太陽E10.7 指數(shù)的反演及預報方法

雷蕾 鐘秋珍 王晶晶 師立勤 苗娟

1(中國科學院國家空間科學中心 北京 100190)

2(中國科學院大學 北京 100049)

0 引言

太陽極紫外輻射EUV（Extreme Ultraviolet）是指波長在10～120 nm 范圍的太陽輻射，是中高層大氣基本的物理加熱能量來源，其擾動變化對地球中高層大氣、電離層和空間天氣具有重要影響[1,2]。在極紫外漏洞[3]時期，由于測量技術的限制，地面觀測的太陽黑子數(shù)和10.7 cm 太陽射電流量（F10.7）一直被廣泛用于表征地球大氣接受的太陽極紫外輻射強度。但實際上太陽黑子數(shù)和F10.7對地球大氣沒有任何直接影響[4,5]。Tobiska[6]通過統(tǒng)計1982 年4 月1 日至1983 年8 月9 日的太陽中層大氣探測衛(wèi)星的日海拔衰減數(shù)據(jù)，發(fā)現(xiàn)F10.7的數(shù)值比大氣真實接收的極紫外能量高60%或低50%。為有效提升中高層大氣模型的精度，提出了一種新的太陽輻射指數(shù)E10.7。

E10.7指數(shù)是地球大氣頂部太陽光譜中波長為1～105 nm 的輻射流量E（t）的積分，并以與F10.7指數(shù)一致的單位sfu（1 sfu=10-22W·m-2·Hz-1）給出。表達式為

其單位為 J·cm-2·s-1。利用多項式擬合技術，可以將E(t)轉化為與F10.7同樣的單位，E10.7（單位 sfu）的表達式如下[6]：

對大氣有效加熱效應的太陽輻射主要集中在1～105 nm 的波段，因此E10.7是更有效表征太陽活動對中高層大氣擾動的代理指數(shù)。在針對不同太陽輻射指數(shù)對大氣模型精度的對比分析中，E10.7的應用可有效提升大氣模型的精度[7]。在中國的軌道預報業(yè)務中，后續(xù)也將應用E10.7。

2000 年美國空間環(huán)境技術中心SpaceWX 發(fā)布的Solar2000（S2 K）大氣輻射模型[8,9]中提供了E10.7指數(shù)，可在參考大氣模型、電離層模型以及其他模型中來表征太陽活動的強度。美國空軍高精度大氣拖曳模型（HASDM）的外層溫度的表達式中，就采用E10.7替代F10.7來表征外層溫度的變化[7]。在針對不同太陽輻射指數(shù)對大氣模型精度的對比分析中，E10.7的應用可有效提升大氣模型精度[10]。已有研究利用30.4 nm 和121.6 nm 的萊曼α 射線兩個能道的太陽輻射觀測值嘗試建立了E10.7計算方法[11]。

在中國的軌道預報業(yè)務中，后續(xù)也將應用E10.7指數(shù)。由于目前中國缺少有效的太陽極紫外觀測，無法實時計算得到E10.7，進而提供業(yè)務化應用，因此本文利用TIMED-SEE 提供的0.1～105 nm 的太陽輻射強度數(shù)據(jù)，反演出實時E10.7觀測數(shù)據(jù)，進而開展E10.7的現(xiàn)報和預報研究。

1 E10.7 數(shù)據(jù)與處理

2002 年在軌運行的TIMED 衛(wèi)星搭載的太陽極紫外線實驗儀（SEE）是首次可以直接測量0.1～194 nm 的太陽極紫外輻射強度的星載儀器，時間跨度為2002 年1 月22 日至今[12]。SEE 由分光計和一套光度計組成，用以測量沉積在大氣層的中間層、低熱層以及電離層區(qū)域的太陽軟X 射線、極紫外線和遠紫外線輻射。SEE 上的軟X 光度計系統(tǒng)由9 個涂有薄膜透射濾光片的硅光電二極管構成，能提供分辨率約5 nm 的1～35 nm 內的軟X 光輻照度的測量數(shù)據(jù)。在25～195 nm 內的光輻照度數(shù)據(jù)則是由SEE上的極紫外光柵光譜儀測量，分辨率為0.4 nm[13]。

TIMED 衛(wèi)星SEE 儀器提供0～105 nm 的太陽輻射積分值，為了區(qū)分Tobiska 定義的E10.7，以下稱為ETS。ETS由戈達德飛行中心的協(xié)調數(shù)據(jù)分析網站（CDA Web）提供*https:--cdaweb.sci.gsfc.nasa.gov-index.html-，時間分辨率約為97 min，單位W·m-2，為準實時數(shù)據(jù)，數(shù)據(jù)延遲時間為4 天。為了獲取實時E10.7數(shù)據(jù)，本文利用歷史數(shù)據(jù)建立ETS和E10.7之間的經驗關系式，并據(jù)此反演得到近實時的E10.7，以通過預報得到未來E10.7預報值，滿足大氣密度模型和軌道預報的需求。

使用的E10.7數(shù)據(jù)是美國海洋氣象局（NOAA）發(fā)布的觀測值*https:--cdaweb.sci.gsfc.nasa.gov，時間分辨率為1天，單位為sfu。歷史數(shù)據(jù)的時間跨度為2002 年2 月7日至2006 年6 月30 日。圖1 給出了2002 年2 月7 日至2006 年6 月30 日的ETS和E10.7的對比，可以看出二者的波動趨勢相似，并且ETS數(shù)據(jù)在時間分辨率上占絕對優(yōu)勢。

從圖1 還可以看出，2002 年中段的ETS部分數(shù)據(jù)明顯偏低，且在ETS數(shù)據(jù)的某些波峰上方有少許偏高的點（紅色點）。前者是由SEE 的某些濾光器異常等導致，后者是由觀測儀器在X 級太陽耀斑發(fā)生期間響應延遲產生[13]。表1 列出了2002 年2 月7 日至2006年6 月30 日發(fā)生的X 級耀斑。本文在后續(xù)計算中根據(jù)表1 中X 級太陽耀斑發(fā)生的峰值時刻剔除這些無效數(shù)據(jù)。

表1 2002 年2 月7 日至2006 年6 月30 日X 級耀斑Table 1 X class solar flares from 7 February 2002 to 30 June 2006

圖1 2002 年2 月7 日至2006 年6 月30 日的ETS和E10.7（紅色點為異常數(shù)據(jù)）Fig.1 ETS and E10.7 from 7 February 2002 to 30 June 2006 (The red points are abnormal)

為了分析E10.7和ETS之間的關系，將分辨率97 min的ETS數(shù)據(jù)進行線性平均得到日均值，缺失的ETS和E10.7數(shù)據(jù)用相鄰兩個數(shù)據(jù)的均值插入補足。通過統(tǒng)計計算可知，E10.7和ETS日均值的線性相關系數(shù)為0.986。為了確定E10.7和ETS的關系，采用最小二乘法對E10.7和ETS的日均值進行線性擬合，可以得到E10.7與ETS的對應關系為

擬合方程的均方根誤差為5.724，結果如圖2 所示。

圖2 2002 年2 月7 日至2006 年6 月30 日NOAA發(fā)布的E10.7與ETS 日均值關系。R 為二者的線性相關系數(shù)，RMSE 為線性擬合均方根誤差，紅色實線為線性擬合結果Fig.2 Relationship between E10.7 and ETS daily average values from 7 February 2002 to 30 June 2006.R is the linear correlation coefficient,RMSE is the root mean square error of linear fitting,and the red solid line is the result of linear fitting

由于TIMED 任務期間存在觀測儀器退化等原因，長時間連續(xù)觀測的ETS數(shù)據(jù)會受到影響。為了減小這些干擾，采用滑動擬合構建E10.7和ETS日均值的線性關系，即

其中，a0和a1分別為線性擬合關系式的常數(shù)項和一次項系數(shù)。滑動窗口寬度為2 個太陽活動自轉周，即54 天。

圖3 給出了利用式（3）計算所得的E10.7數(shù)據(jù)（紅色實線）與滑動擬合得到的E10.7數(shù)據(jù)（藍色實線，以下稱為T10.7）的對比。后者的均方根誤差為5.445，低于前者。滑動擬合方法可以排除長時間觀測儀器退化對數(shù)據(jù)的影響，且計算所需數(shù)據(jù)明顯更少，擬合的誤差也更低，因此在實際應用中具有明顯的優(yōu)勢。這里建模所用數(shù)據(jù)就是由滑動擬合得到的T10.7。

圖3 2002 年2 月7 日至2006 年6 月30 日NOAA發(fā)布的E10.7 日均值、利用式（3）計算的E10.7 數(shù)據(jù)及滑動擬合得到的E10.7 數(shù)據(jù)Fig.3 E10.7 of NOAA calculaed by Eq.(3),and obtained by sliding window linear fitness from 7 February 2002 to 30 June 2006

通過以上方法，得到了與E10.7單位同為sfu的T10.7歷史值，數(shù)值保留一位小數(shù)。這使得后續(xù)在各模型中能夠更方便地直接應用T10.7數(shù)據(jù)[6,9]。

2 E10.7 中期預報模型構建

對于E10.7數(shù)據(jù)這類時間序列，通常能夠通過統(tǒng)計方法計算出其相關關系，并利用這個關系來預測未來的數(shù)據(jù)。常用的預報模型主要包括自回歸（Autoregressive，AR）模型、滑動平均（Moving Average，MA）模型以及自回歸-滑動平均（Autoregressive Moving Average，ARMA）模型。ARMA 模型的適用性更廣，但是其建模計算非常復雜；AR 模型建模相對簡單，并且一個高階的AR 模型可以近似于一個ARMA 模型。因此，本文選擇利用高階自回歸模型對E10.7做中期預報。

AR 模型是利用同一變數(shù)，例如X的往期即Xt-p至Xt-1來預測本期Xt的表現(xiàn)，并假設其為一線性關系，表達式為

其中，p為AR 模型的階數(shù)；?t為白噪聲序列，視為預報誤差。模型的參數(shù)估計采用伯格（Burg）算法[14,15]。

為了確定自回歸預報模型的參數(shù)和階數(shù)，利用2002 年2 月7 日至2006 年6 月30 日的T10.7數(shù)據(jù)進行自相關和偏自相關分析。圖4 給出了分析結果，可以看出T10.7的自相關系數(shù)隨天數(shù)的增加逐漸減小，并在天數(shù)為600 的時候最趨于0；而T10.7的偏自相關系數(shù)在52 天時最趨于0，且其后天數(shù)的偏自相關系數(shù)均在正負0.05 以內。根據(jù)自相關系數(shù)的分析，T10.7數(shù)據(jù)的自相關性隨著相隔天數(shù)的增加下降，600 天時相關系數(shù)接近0，即沒有相關性。因此，采用連續(xù)600 天的觀測數(shù)據(jù)建立自回歸模型的樣本數(shù)據(jù)，即用連續(xù)600 天的觀測數(shù)據(jù)擬合式（5）中的模型參數(shù)。根據(jù)偏自相關系數(shù)的分析，與當前數(shù)據(jù)存在獨立相關性的數(shù)據(jù)主要在相隔52 天以內的數(shù)據(jù)，并且52 天接近兩個太陽自轉周的時間。因此，預報E10.7自回歸模型的階數(shù)選用52 階。

圖4 2002 年2 月7 日至2006 年6 月30 日的T10.7 的自相關（a）和偏自相關（b）結果（紅色點為相關系數(shù)為0 的點）Fig.4 Autocorrelation function (a) and the partial autocorrelation function (b) of T10.7 from 7 February 2002 to 30 June 2006 (Red dot is the point with the correlation coefficient of 0)

實際應用中，通常用E10.7的81 天中心平均值來表示近3 個太陽自轉周極紫外輻射的平均狀態(tài)。本文將過去40 天、未來40 天和當日E10.7的日均值定義為E10.7的81 天中心平均值（E81C），表達式為

其中，t為距離預報當天的天數(shù)。

由于T10.7觀測數(shù)據(jù)有4 天的延遲時間，若采取此數(shù)據(jù)做預報，在t=-3時，就要開始對T10.7進行預報。對于計算未來第27 天的E81C，需要至未來第67 天的E10.7。因此，需要預報共71 天（t∈[-3,67]）的T10.7。其中前4 天的值為延遲時間內的預報值，5 至31 天的值為未來27 天的預測值，71 天的預報值和實測值用于滑動計算未來27 天的E81C預報值。的預報模型構建流程如下。

第1步利用最小二乘法將過去600 天的E10.7和ETS歷史觀測數(shù)據(jù)進行滑動線性擬合，滑動窗口寬度為54，得到共600 天的T10.7數(shù)據(jù)。

第2步將第1 步所得共600 天的T10.7數(shù)據(jù)作為訓練數(shù)據(jù)，帶入式（5）的52 階自回歸模型，擬合得到自回歸模型的參數(shù)。

第3步把前52 天的T10.7數(shù)據(jù)輸入自回歸模型，計算得到當天的E10.7預報值；以前51天T10.7數(shù)據(jù)和當天E10.7預報值作為輸入，計算得到未來1天E10.7的預報值。以此類推，計算71 天的E10.7預報值。由于T10.7觀測數(shù)據(jù)有4 天的延遲時間，實際應用中，第5 至31 天的值為未來27 天的E10.7預報值。

第4步利用第3 步所得的71天E10.7預報值和E10.7歷史實測值，通過式（6）計算未來27 天的E81C預報值。

3 E10.7 預報模型評估

為了評估預報模型的效果，這里以2005 年1 月1 日至2006 年6 月30 日期間的數(shù)據(jù)為測試集，對E10.7指數(shù)進行預報測試。

評估預報結果主要采用逐日相對誤差、相對誤差、平均逐日相對誤差和平均相對誤差。其中，逐日相對誤差（）為單次預報的未來第t天的相對誤差，相對誤差（Emr）為單次預報未來27 天的相對誤差的平均值，平均逐日相對誤差（）為測試集的的平均值，平均相對誤差（Emr,av）為全部測試集的Emr的平均值。這4 個量的表達式為

其中，N為測試集樣本數(shù)，t為預測未來的天數(shù)，i為測試集中的單個樣本，ot,i為E10.7的實測日值，ft,i為E10.7的預報日值。

利用E10.7自回歸預報模型預報2005 年5 月20 日未來67 天的E10.7和未來27 天的E81C，如圖5所示。T10.7觀測數(shù)據(jù)有4 天的延遲時間，因此利用式（5）的52 階自回歸模型，共預報71 天（橫坐標-3 至67）的T10.7，其中前4 天（橫坐標-3 至0）的值為延遲時間內的預報值，5 至31 天（橫坐標1 至27）的值為未來27 天的預測值。從圖5 可以看出，模型對5 月20 日之后兩個太陽自轉周內的E10.7的周期變化有較好的預報效果，且未來27 天的E81C預報值與實測值也有較好吻合。模型27天E10.7預報值的平均逐日相對誤差為3.05%，模型27 天的E81C預報值平均逐日相對誤差為0.62%。

利用E10.7自回歸預報模型預報2005 年5 月1 日未來67 天的E10.7和未來27 天的E81C，如圖6所示。與圖5 相比，本次預報E10.7和E81C的平均逐日相對誤差分別增大了10.4% 和7.39%。從圖6 可以看出，52 階自回歸預報模型可以很好地確定太陽活動指數(shù)的27 天周期，但是預報值與實測值存在一定差距，預報值的周期波動幅度隨預報時間增長而減小。因此，當未來27 天的E10.7比前一個周期有增幅時，52 階自回歸預報模型的預報誤差較大。

圖5 2005 年5 月20 日預報與觀測的E10.7和E81C。t 為距2005 年5 月20 日的天數(shù)，負值為距2005 年5 月20 日的歷史天數(shù)Fig.5 E10.7 and E81C that are predicted and observed on 20 May 2005.Positive values in abscissa are the numbers of days after 20 May 2005,the negative values in abscissa are the numbers of days before 20 May 2005

圖6 2005 年5 月1 日預報與觀測的E10.7和E81C。t 為距離2005 年5 月1 日的天數(shù)，負值為距2005 年5 月1 日的歷史天數(shù)Fig.6 E10.7 and E81C that are predicted and observed on 1 May 2005.Positive values in abscissa are the numbers of days after 1 May 2005,the negative values in abscissa are the numbers of days before 20 May 2005

為衡量自回歸模型的整體預報效果，對2005 年1 月1 日至2006 年6 月30 日期間的預報結果進行整體評估。圖7 為此期間未來27天E10.7（黑色菱形）和E81C（紅色三角）預報結果的平均逐日相對誤差。由于T10.7的實測數(shù)據(jù)存在4 天延時，圖中0 天是提前4 天的預報結果，-3 天是提前1 天的預報結果。從圖7 可看出，模型E10.7預報值的誤差在前7 天內（即提前10 天的預報）較快增長，平均相對誤差由3.5%（-3 天）增長到8%（第7 天），7 天至24 天的平均相對誤差基本維持在8%附近，24 天后又出現(xiàn)增長。E81C預報值的平均相對誤差隨預報天數(shù)的增加而緩慢增大，從2.5%（-3 天）增加到4.5%（第27 天）。

圖7 模型預測未來27天E10.7和E81C 的平均逐日相對誤差Fig.7 Average daily relative error of E10.7 and E81C that are predicted in the next 27 days

表2 給出了預報結果的誤差統(tǒng)計。在整個測試時間段，E10.7預報值27 天平均逐日相對誤差最小值為2.4%，最大值為15.08%，平均值為7.83%；E81C預報值27 天平均逐日相對誤差最小值為0.29%，最大值為10.10%，平均值為3.63%。對于未來大氣模型的應用來說，E81C預報值平均3.63%的誤差可以較好地滿足應用需求。

表2 測試集的E10.7和E81C 預報誤差Table 2 Predicted error of E10.7 and E81C of all test sets

4 討論與結論

針對高層大氣密度預報和軌道預報業(yè)務中對新型太陽紫外輻射指數(shù)E10.7的需求，利用TIMED/SEE觀測到的高時間分辨率的太陽輻射強度數(shù)據(jù)，基于最小二乘法建立了反演E10.7指數(shù)的方法。在分析反演數(shù)據(jù)自相關和偏自相關的基礎上，構建了一個52 階的自回歸預報模型，可以實現(xiàn)E10.7的中期預報。初步結論如下。

TIMED/SEE 觀測儀提供的ETS實測值具有高時間分辨率、延遲時間短和易獲得的優(yōu)勢，利用最小二乘法擬合可反演出準實時的E10.7，均方根誤差為5.445 sfu，能夠基本滿足實際應用。利用該方法可以為中國相關用戶提供準實時的E10.7指數(shù)。

利用52 階自回歸模型對E10.7的中期預報效果尚好，預報誤差隨預報時間的增加而增大，未來27 天的預報值平均相對誤差為7.83%。利用同樣方法也開展了對E10.7的81 天滑動平均值E81C的未來27 天預報試驗。試驗結果表明，與E10.7類似，預報平均相對誤差隨預報天數(shù)的增加而增大，但誤差較低，未來27 天的預報值平均相對誤差僅為3.63%。

隨著航天事業(yè)發(fā)展，中國衛(wèi)星軌道預報業(yè)務對軌道大氣密度的精確度要求也越來越高，將需要逐步在大氣模型中采用E10.7指數(shù)替代F10.7指數(shù)，提高大氣密度計算的準確性。本文建立的E10.7反演方法和自回歸預報模型可用于未來E10.7的現(xiàn)報和預報。

致謝ETS實測值由TIMED/SEE 觀測儀提供。