梁 曉鄭小谷戴永久師春香

1）（國家氣象信息中心，北京100081）2）（北京師范大學(xué)，北京100875）

En KF中誤差協(xié)方差優(yōu)化方法及在資料同化中應(yīng)用

梁 曉1）*鄭小谷2）戴永久2）師春香1）

1）（國家氣象信息中心，北京100081）2）（北京師范大學(xué)，北京100875）

集合卡爾曼濾波（the Ensemble Kalman Filter，簡稱En KF）中將預(yù)報集合的統(tǒng)計協(xié)方差作為預(yù)報誤差協(xié)方差，但該估計可能嚴(yán)重偏離真實的預(yù)報誤差協(xié)方差，影響同化精度。基于極大似然估計理論，發(fā)展了一種優(yōu)化預(yù)報誤差協(xié)方差矩陣的實時膨脹方法，即MLE（the Maximum Likelihood Estimation）方法。利用蒙古國基準(zhǔn)站Delgertsgot（簡稱DGS站）觀測資料，基于En KF方法和 MLE方法，在通用陸面模式（the Common Land Model，簡稱Co LM）中同化了地表溫度和10 cm土壤溫度觀測資料，建立了土壤溫度同化系統(tǒng)。結(jié)果表明：MLE方法對地表溫度和各層土壤溫度（尤其深層土壤溫度）的估計比En KF方法準(zhǔn)確。考慮到淺層和深層土壤溫度的差別，在實施MLE方法時對淺層和深層土壤溫度采用了不同的膨脹因子。對比膨脹因子為單一標(biāo)量時的結(jié)果，多因子膨脹能緩解深層土壤溫度的不合理膨脹，改善同化效果。

數(shù)據(jù)同化；集合卡爾曼濾波；誤差協(xié)方差膨脹

引 言

集合卡爾曼濾波（the Ensemble Kalman Filter，簡稱En KF）是目前常用的一種數(shù)據(jù)同化方法，自Evensen［1］提出和Burgers等［2］發(fā)展后，已應(yīng)用于很多同化研究。理論上，通過將觀測資料融入預(yù)報模型，數(shù)據(jù)同化能夠?qū)崿F(xiàn)模型輸出和觀測資料的最佳融合。但這種“最佳”依賴于預(yù)報（背景）誤差協(xié)方差和觀測誤差協(xié)方差的精度。預(yù)報誤差協(xié)方差和觀測誤差協(xié)方差的估計是否準(zhǔn)確，是決定En KF方法同化精度高低的關(guān)鍵［3］。然而，實際應(yīng)用中很難精確得到這些統(tǒng)計量信息［4］。因此，優(yōu)化誤差協(xié)方差的估計，對En KF方法發(fā)展十分必要［5-6］。

En KF方法中，通常將預(yù)報集合的統(tǒng)計協(xié)方差作為真實的預(yù)報誤差協(xié)方差。但研究表明，由于樣本量有限，該估計中存在的采樣誤差會導(dǎo)致預(yù)報集合的統(tǒng)計協(xié)方差低估真正的預(yù)報誤差協(xié)方差，集合擴散度迅速減小，使得觀測對同化的影響越來越小，甚至完全失效，發(fā)生濾波發(fā)散［7-8］。為此，一種旨在增大集合方差的協(xié)方差膨脹技術(shù)逐漸發(fā)展起來。目前常用的一種方法是給預(yù)報集合成員偏離平均值的偏差乘一個略大于1的常數(shù)［7］，該方法中的膨脹因子通過反復(fù)試驗得到，需多次運行同化系統(tǒng)。對于大型復(fù)雜的實際模型，這種試驗方法的計算代價十分巨大，且膨脹因子在同化時段內(nèi)不變也不合理。因此，有必要發(fā)展無需反復(fù)試驗并能在同化循環(huán)中實時更新膨脹因子的協(xié)方差膨脹技術(shù)。

Dee［9］及其相關(guān)工作［10-11］發(fā)展了誤差協(xié)方差的極大似然估計方法。首先對誤差協(xié)方差矩陣建立參數(shù)化模型，再通過觀測減預(yù)報得到的殘差的極大似然估計來估計模型中的未知參數(shù)。但其工作沒有提出可靠的誤差協(xié)方差矩陣參數(shù)化模型，因此一直未得到普及。近年來，根據(jù)觀測減預(yù)報得到的殘差的統(tǒng)計特征，逐漸提出了一些實時估計膨脹因子的方法［12-14］。然而，這些方法均基于矩估計，而非極大似然估計。

作為極大似然估計方法［9-11］的擴展，Zheng［15］提出了預(yù)報誤差協(xié)方差矩陣的多變量實時膨脹方法，將膨脹因子從常數(shù)標(biāo)量擴展到實時變化的對角矩陣。Liang等［16］基于Zheng［15］的工作，提出了同時膨脹預(yù)報和觀測誤差協(xié)方差矩陣的方法，但上述工作只在簡單的小模型中進行試驗。作為文獻(xiàn)［16］工作的延續(xù)，本文使用陸面過程模型，將Zheng［15］和Liang等［16］的膨脹方法用于陸面同化。本文的膨脹因子已擴展為對角矩陣，即非標(biāo)量。

土壤溫度是陸面過程中的一個重要變量，可以通過模型模擬或儀器觀測得到，但兩種手段各有利弊。站點觀測反映了觀測時間和地點上的真實狀態(tài)，但由于土壤的空間異質(zhì)性，它所能代表的空間范圍有限，無法滿足區(qū)域性研究的要求；遙感觀測可以得到衛(wèi)星過境時刻大范圍的地表溫度信息，但遙感反演過程引入了很大的不確定性，而且只能獲得瞬時信息。陸面過程模型可以給出土壤溫度在時間和空間上的連續(xù)演進，但模型參數(shù)、驅(qū)動場、初始條件等的不完備都會影響模擬精度。如何有效結(jié)合觀測和陸面過程模型模擬得到更真實且時空連續(xù)的信息，尤為重要，因此陸面數(shù)據(jù)同化研究越來越多［17-20］。本文將利用通用陸面模式（the Common Land Model，簡稱Co LM）和兩種同化方法發(fā)展土壤溫度同化方案，改進對土壤溫度廓線的估計。為便于表達(dá)，未使用和使用本文膨脹方法的En KF方法分別簡稱為En KF方法和MLE方法。

1 MLE方法

一個非線性時間離散的預(yù)報系統(tǒng)和線性觀測系統(tǒng)可表示為

其中，i是積分步數(shù)，是真值向量，Mi－1是非線性模型算子，yoi是觀測向量，Hi是將狀態(tài)變量映射到觀測空間的觀測算子，ηi和εi是背景誤差和觀測誤差向量，假定二者相互獨立且時間上不相關(guān)，各自服從平均值為0、協(xié)方差矩陣分別為背景誤差協(xié)方差矩陣Pi和觀測誤差協(xié)方差矩陣Ri的高斯分布。同化目標(biāo)是找到盡可能接近真值的分析狀態(tài)

En KF方法在很多文獻(xiàn)中均有介紹［2，21-22］，不再贅述。En KF方法中通常采用預(yù)報集合的采樣協(xié)方差矩陣Pi（式（3））作為真正的預(yù)報誤差協(xié)方差矩陣

有限的樣本量導(dǎo)致該估計中存在較大的采樣誤差，集合協(xié)方差會低估真正的誤差協(xié)方差，導(dǎo)致濾波發(fā)散［23］。最初，抑制濾波發(fā)散的方法是給預(yù)報集合中每個成員偏離平均值的離差乘略大于1的因子相當(dāng)于將預(yù)報誤差協(xié)方差矩陣Pi膨脹為原來的λ倍。但該膨脹因子λ隨時間不變，且需要多次反復(fù)試驗來確定［7］。

基于觀測減預(yù)報得到的殘差di的極大似然估計，Zheng［15］和 Liang等［16］提出一種隨時間變化的膨脹方法，通過極小化di的似然函數(shù)－2Li（λi）（即對di的概率密度函數(shù)?。?倍的自然對數(shù)）實時估計膨脹因子［λi］：

其中，det表示矩陣的行列式。該方法根據(jù)文獻(xiàn)［9-11］提出的估計誤差協(xié)方差方法發(fā)展而來，膨脹因子［λi］為對角矩陣，不是標(biāo)量。

估計得到［λi］后，預(yù)報誤差協(xié)方差矩陣假定為［λi］Pi［λi］，相當(dāng)于用膨脹后的預(yù)報集合重新計算式（3）中定義的預(yù)報誤差協(xié)方差矩陣：

結(jié)合預(yù)報和觀測信息，可計算出i時刻同化更新后的分析集合

式（7）中，ε′i，j為服從平均值為0、協(xié)方差矩陣為Ri的正態(tài)分布隨機擾動。

本文將使用同化結(jié)果的均方根誤差來評價同化效果。通常來說，均方根誤差越小，同化結(jié)果越好。

2 試驗設(shè)計

本文選取國際協(xié)同強化觀測期計劃（Coordinated Energy and Water Cycle Observations Project，簡稱CEOP）蒙古國基準(zhǔn)站Delgertsgot（簡稱DGS站），使用 En KF和 MLE兩種同化方法，在Co LM［24-25］中以12 h 1次的頻率同化DGS站地表溫度和10 cm土壤溫度觀測資料，比較兩種方法對土壤溫度廓線的估計效果。同化系統(tǒng)狀態(tài)向量包括地表溫度、10層土壤溫度和10層土壤濕度。本文還分別使用經(jīng)驗算法得到的葉面積指數(shù)和MODIS LAI產(chǎn)品，用來分析葉面積指數(shù)精度對土壤溫度的影響。

2.1 研究區(qū)概況及觀測數(shù)據(jù)

蒙古國基準(zhǔn)站位于蒙古高原，地表覆蓋為均質(zhì)矮草，觀測數(shù)據(jù)較好［18，26］。其中 DGS站（46.1°N，106.4°E）于2002年10月1日—2003年9月30日進行了地表溫度和土壤溫度及氣象要素的連續(xù)觀測，本文選取該站開展同化試驗。土壤溫度的觀測深度為3，10，40 cm和100 cm，相應(yīng)還有土壤濕度觀測，30 min 1次。在該站點已開展了一些陸面同化研究［17-18］，本文在選取站點觀測資料和研究時段等方面均予以考慮。本研究中，在Co LM中同化了DGS站每日08：00（世界時，下同）和20：00的地表溫度和10 cm土壤溫度的觀測數(shù)據(jù)，并使用逐時連續(xù)的各層土壤溫度和土壤濕度觀測數(shù)據(jù)來驗證同化結(jié)果。

Co LM運行所需的大氣驅(qū)動數(shù)據(jù)中，除入射太陽短波輻射和大氣長波輻射沒有觀測資料外，其他變量均來自自動氣象站資料。日本氣象廳利用其全球數(shù)據(jù)同化系統(tǒng)為CEOP的基準(zhǔn)站提供了始于2002年10月1日的時間分辨率為1 h的大氣數(shù)據(jù)集，Huang等［17］從該數(shù)據(jù)集獲取了蒙古國基準(zhǔn)站的太陽短波輻射和大氣長波輻射，本研究所需的輻射信息來源于此。

為分析葉面積指數(shù)（leaf area index，簡稱LAI）對土壤溫度的影響，本研究還使用MODIS LAI產(chǎn)品代替原模式經(jīng)驗算法計算的葉面積指數(shù)，比較使用經(jīng)驗葉面積指數(shù)和MODIS LAI產(chǎn)品得到的結(jié)果。MODIS LAI產(chǎn)品提供了全球的1 km分辨率的葉面積指數(shù)（每8 d 1次），本文從中提取DGS站的2002年10月1日—2003年9月30日的葉面積指數(shù)和質(zhì)量控制說明數(shù)據(jù)。

2.2 MLE方法的兩種情況

本文使用的膨脹因子［λi］是對角矩陣，而不再局限為標(biāo)量。當(dāng)［λi］的主對角線上含有不同數(shù)值時，即可對不同的狀態(tài)變量使用不同系數(shù)進行不同程度的膨脹；而當(dāng)［λi］的主對角線上均為相同數(shù)值時，膨脹因子［λi］便相當(dāng)于標(biāo)量，即對所有的狀態(tài)變量都使用同一系數(shù)進行相同程度的膨脹。

Co LM中地表溫度和深層土壤溫度的性質(zhì)差異大，如地表溫度隨時間變化快且變化幅度大，而深層土壤溫度隨時間變化平緩。如果對整個土壤溫度廓線均采用同一膨脹因子，可能無法兼顧淺層和深層土壤溫度膨脹的合理性，最終影響土壤溫度廓線的同化效果。因此，本文對淺層和深層土壤溫度采用不同的膨脹因子，即膨脹因子矩陣［λi］的主對角線上含有兩個待估計的標(biāo)量膨脹因子λi，1和λi，2。其中［λi］主對角線上的前3個元素值均為λi，1，用于膨脹地表溫度和上兩層土壤溫度（2.8 cm以上）；后8個元素值均為λi，2，用于膨脹下8層土壤溫度（2.8 cm以下）；與土壤濕度對應(yīng)的10個元素值則為1.0，即 不 膨 脹 土 壤 濕 度。由 式 （5）可 知，－2Li（λi）實際上是關(guān)于Hi［λi］的函數(shù)，即無論膨脹因子矩陣［λi］中包含多少未知參數(shù)，通過極小化－2Li（λi）只能估計出Hi［λi］能包含的那些參數(shù)。如果被同化的觀測信息只有地表溫度，則Hi［λi］只包含用于膨脹地表溫度和上兩層土壤溫度的參數(shù)λi，1，即 MLE方法通過極小化－2Li（λi）只能估計出λi，1，無法估計出λi，2。為此，本文在同化系統(tǒng)中又引入10 cm土壤溫度觀測信息，這樣Hi［λi］中就同時包含了λi，1和λi，2，從而保證 MLE方法可以同時估計得到這兩個參數(shù)。

為探討對淺層和深層土壤溫度進行不同程度膨脹的必要性，本文將比較淺層和深層土壤溫度使用同一系數(shù)和不同系數(shù)進行膨脹時的結(jié)果。為便于表達(dá)，淺層和深層土壤溫度使用同一膨脹系數(shù)和不同膨脹系數(shù)的MLE方法分別簡稱為MLE1方法和MLE2方法。

3 同化結(jié)果

3.1 確定觀測算子和觀測誤差

本研究同化的觀測變量為地表溫度和10 cm土壤溫度，狀態(tài)變量為地表溫度、10層土壤溫度和10層土壤濕度。

本研究同化的地表溫度觀測是站點觀測數(shù)據(jù)，觀測變量與模型中的地表溫度變量一致，故地表的觀測算子為1。本研究還同化了站點觀測的10 cm土壤溫度，而這一深度無法與Co LM中定義的任一層土壤節(jié)點的深度直接對應(yīng)。因此，利用相鄰兩層模擬深度上的土壤溫度經(jīng)線性內(nèi)插得到10 cm土壤溫度，該插值算法即為10 cm深度觀測算子：

其中，y1和y2分別為地表和10 cm土壤溫度的模擬值，d3和d4分別為Co LM中第3層和第4層土壤節(jié)點的深度（6.2 cm和11.9 cm），x0為Co LM 的地表溫度模擬值，x3和x4分別為Co LM中第3層和第4層的土壤溫度模擬值。該觀測算子將Co LM輸出的模型狀態(tài)變量與觀測變量聯(lián)系起來。

在DGS站，地表溫度使用4000.4G紅外溫度傳感器觀測，土壤溫度使用溫度計觀測。查閱相關(guān)儀器手冊，最終確定地表溫度和10 cm土壤溫度的觀測誤差分別為1.0 K和0.5 K。

3.2 同化結(jié)果

使用DGS站2002年10月1日—2003年9月30日的觀測數(shù)據(jù)進行試驗。Co LM從2002年10月1日積分到2003年5月29日，作為spin-up使模型達(dá)到平衡，2003年5月30日同化開始，積分步長為1 h。基于En KF方法和MLE方法，Co LM每12 h同化1次地表溫度和10 cm土壤溫度觀測數(shù)據(jù)（08：00和20：00），集合樣本量為100。選用研究時段為2003年9月1—30日不同深度的土壤溫度來驗證同化效果。

淺層土壤溫度的時間變化較快且變幅較大，很難區(qū)分其在30 d內(nèi)的觀測場、模擬場和同化場。因此，對地表溫度和3 cm土壤溫度在30 d內(nèi)的結(jié)果進行平均，得到其日變化（圖1）。En KF方法和MLE方法對地表溫度（圖1a）和3 cm土壤溫度（圖1b）的估計都比模擬值準(zhǔn)確，尤其在09：00—23：00時段。與模擬相比，En KF方法在01：00—07：00有所低估，MLE方法在此時段內(nèi)的結(jié)果則可以與模擬場相當(dāng)甚至更優(yōu)。

圖1 2003年9月1—30日地表溫度（a）和3 cm土壤溫度（b）觀測場、模擬場和同化場平均日變化Fig.1 The diurnal variation of observed，simulated and assimilated soil temperature at 0 cm（a）and 3 cm（b）averaged from 1 Sep to 30 Sep in 2003

圖2為2003年9月1—30日10 cm和40 cm土壤溫度的觀測場、模擬場和同化場。在10 cm深度（圖2a），En KF方法和MLE方法均比模擬更準(zhǔn)確地估計出了土壤溫度日變化的趨勢和幅度，但MLE方法相比En KF方法的優(yōu)勢不大。在40 cm深度（圖2b），與觀測場相比，模擬場明顯低估，而兩種同化方法均有效改善了這種低估，且MLE方法同化得到的土壤溫度比En KF方法更接近觀測。100 cm的結(jié)果與40 cm相似（圖略）。

圖2 2003年9月1—30日10 cm（a）和40 cm（b）土壤溫度觀測場、模擬場和同化場變化Fig.2 The observed，simulated and assimilated soil temperature at 10 cm（a）and 40 cm（b）from 1 Sep to 30 Sep in 2003

表1為各層土壤溫度模擬場和同化場在2003年9月1—30日平均的均方根誤差。En KF方法和MLE方法同化得到的均方根誤差均低于模擬結(jié)果的均方根誤差，而且膨脹因子為對角矩陣的MLE2方法比En KF方法更好。表1還表明，同化對土壤溫度估計效果的改善在10 cm，40 cm和100 cm深度比0 cm和3 cm深度明顯。淺層的土壤溫度變化較快且變幅較大，因此，雖然同化時刻的同化結(jié)果會接近觀測值，但同化時刻后，在持續(xù)較長的預(yù)報時間內(nèi)，變化很快的淺層土壤溫度會以較快的速度向模型模擬結(jié)果靠攏，導(dǎo)致上一步同化對淺層土壤溫度的改善在下一步同化到來時影響不大。這樣，同化時刻對淺層土壤溫度造成的影響就很難在整個預(yù)報過程中得到有效保持和傳遞。為進一步分析該問題，本文單獨針對2003年9月1—30日每日08：00和20：00的同化時刻，統(tǒng)計了膨脹因子為對角矩陣的MLE2方法同化得到的地表溫度、3 cm和10 cm土壤溫度在同化更新前的背景場均方根誤差和同化更新后的分析場均方根誤差?？偟膩碚f，時間平均的地表溫度、3 cm和10 cm土壤溫度同化場均方根誤差（分別為0.402，1.920 K和0.157 K）均顯著低于背景場均方根誤差（分別為3.028，3.509 K和1.110 K），這表明同化更新后的同化場比背景場更加接近觀測。這與前面的分析一致：在每日兩次的同化時刻，同化的確對淺層土壤溫度產(chǎn)生了較大影響，但這種影響很難在整個預(yù)報過程中得到有效保持。

考慮到淺層和深層土壤溫度的性質(zhì)不同，本文在實施MLE方法時，對淺層和深層土壤溫度采用了不同的膨脹因子。從2003年9月1—30日平均結(jié)果看，MLE1方法（淺層和深層土壤溫度使用同一系數(shù)進行膨脹）中地表溫度和10層土壤溫度的膨脹因子值均為2.2，而MLE2方法（淺層和深層土壤溫度使用不同系數(shù)進行膨脹）中地表溫度和上兩層（0.7～2.8 cm）土壤溫度的膨脹因子值為3.2，下8層（6.2～286.5 cm）土壤溫度的膨脹因子值為1.3?？梢?，MLE1方法和MLE2方法都通過膨脹因子而對預(yù)報誤差協(xié)方差矩陣的離散度發(fā)揮了作用。但MLE1方法沒有對淺層和深層土壤溫度區(qū)別膨脹，導(dǎo)致深層土壤溫度膨脹過度；而MLE2方法對淺層和深層土壤溫度采用了不同的膨脹因子，使得深層土壤溫度的膨脹相對合理。表1中分別列出了MLE1和MLE2方法得到的土壤溫度的均方根誤差。可以看到，MLE1方法對各層土壤溫度的估計均沒有MLE2方法準(zhǔn)確，尤其對于深層土壤溫度（40 cm和100 cm），MLE1方法同化得到的均方根誤差比MLE2方法高50%以上，甚至比En KF方法同化得到的均方根誤差還高。圖3為En KF方法、MLE2方法和MLE1方法同化得到的40 cm和100 cm土壤溫度。由圖3可見，由于MLE1方法對整個土壤溫度廓線都采用同一因子進行膨脹，將適用于變化很快的淺層土壤溫度的膨脹因子，也用于膨脹變化較慢的深層土壤溫度，造成深層土壤溫度被不合理膨脹，進而影響了深層的同化效果。

表1 2003年9月1—30日土壤溫度模擬場和同化場平均的均方根誤差（單位：K）Table 1 Root mean square error of the simulated and assimilated soil temperature from 1 Sep to 30 Sep in 2003（unit：K）

圖3 MLE1方法和MLE2方法分別同化得到2003年9月1—30日40 cm（a）和100 cm（b）土壤溫度Fig.3 Soil temperature at 40 cm（a）and 100 cm（b）assimilated by MLE1 and MLE2 from 1 Sep to 30 Sep in 2003

圖4為2003年9月1—30日3 cm近地表土壤濕度的觀測場、模擬場和同化場。與觀測相比，Co LM的模擬普遍高估了土壤濕度，而En KF方法和MLE方法均對此有所緩解，且MLE方法的改善更明顯。這表明，雖然同化過程中沒有引入土壤濕度的觀測信息，同化過程仍可以通過誤差協(xié)方差矩陣和模式傳輸而影響土壤濕度。

圖4 2003年9月1—30日3 cm土壤濕度觀測場、模擬場和同化場Fig.4 The observed，simulated and assimilated soil moisture at 3 cm from 1 Sep to 30 Sep in 2003

4 結(jié)論和討論

本文基于極大似然估計理論，發(fā)展了一種優(yōu)化預(yù)報誤差協(xié)方差的實時膨脹方法——MLE方法。當(dāng)然，除了協(xié)方差膨脹方法外，國內(nèi)外在預(yù)報誤差協(xié)方差的分析和優(yōu)化方法上還有很多其他研究［27-30］，這些超出了本文目前的研究范圍，有待于下一步的跟蹤和關(guān)注。本研究表明：

1）總的來說，MLE方法對地表溫度和各層土壤溫度的估計都比En KF方法準(zhǔn)確，對深層土壤溫度（40 cm和100 cm）尤其如此。

2）考慮到淺層和深層土壤溫度的差別，本文對淺層和深層土壤溫度采用不同的膨脹因子，即MLE2方法。與膨脹因子為單一標(biāo)量值時的MLE1方法比較表明，由于對整個土壤溫度廓線采用相同的膨脹因子，MLE1方法對各層（尤其深層）土壤溫度的估計沒有MLE2方法準(zhǔn)確。這表明，多因子膨脹可以緩解深層土壤溫度的不合理膨脹問題。

3）為分析葉面積指數(shù)對土壤溫度的影響，本文在Co LM中引入MODIS LAI產(chǎn)品代替經(jīng)驗算法計算的葉面積指數(shù)。結(jié)果表明，與使用經(jīng)驗葉面積指數(shù)相比，使用MODIS LAI可以改善Co LM對地表溫度和3 cm土壤溫度的模擬，對地表溫度和各層土壤溫度的同化也更加準(zhǔn)確。

本文在估算模型誤差協(xié)方差矩陣時使用對角矩陣（協(xié)方差項為0），沒有考慮不同土壤層變量之間的相關(guān)。Jin等［31］研究表明，當(dāng)合理給定模型誤差協(xié)方差矩陣中的協(xié)方差項時，同化系統(tǒng)能將表層優(yōu)化后的信息快速傳遞給深層，從而快速顯著改善整個土壤廓線狀態(tài)變量的估計。但其工作是單獨優(yōu)化土壤溫度或濕度，未涉及土壤溫度與濕度之間的誤差協(xié)相關(guān)。在同步優(yōu)化土壤溫度和土壤濕度時，如何確定土壤溫度和濕度之間的協(xié)相關(guān)，值得研究。

本文膨脹因子隨時間變化，但在空間上不變，這實際上不合理，尤其是觀測資料在空間上分布不均勻時。適用于觀測資料濃密區(qū)的膨脹因子，若持續(xù)用于所有地區(qū)，可能導(dǎo)致觀測資料稀疏區(qū)過度膨脹。為對淺層和深層土壤溫度進行不同程度的膨脹，本文引入了10 cm土壤溫度的觀測信息。如果隨時空變化的膨脹方法得以實現(xiàn)，則僅同化地表溫度就能針對不同深度的土壤溫度估計出不同的膨脹因子。Anderson［32］提出的方法中膨脹因子時空上均可變化，但其工作只考慮了觀測誤差空間獨立的情況，不利于遙感資料在同化中的應(yīng)用。因此，有必要發(fā)展適用于空間相關(guān)觀測系統(tǒng)的、時空變化的膨脹方法。

［1］ Evensen G.Sequential data assimilation with a nonlinear quasi-geostrophic model using Monte-Carlo methods to forecast error statistics.J Geophys Res，1994，99（C5）：10143-10162.

［2］ Burgers G，Leeuwen P J V，Evensen G.Analysis scheme in the Ensemble Kalman Filter.Mon Wea Rev，1998，126：1719-1724.

［3］ Evensen G.The ensemble Kalman filter：Theoretical formulation and practical implementation.Ocean Dynam，2003，53（4）：343-367.

［4］ Senegas J，Wackernagel H H，Rosenthal W，et al.Error covariance modeling in sequential data assimilation.Stoch Env Res Risk A，2001，15：65-86.

［5］ Daley R.Atmospheric Data Analysis.Cambridge：Cambridge University Press，1991.

［6］ Kalnay E.Atmospheric Modeling，Data Assimilation，and Predictability.Cambridge：Cambridge University Press，2002.

［7］ Anderson J L，Anderson S L.A Monte Carlo implementation of the non-linear filtering problem to produce ensemble assimilations and forecasts.Mon Wea Rev，1999，127：2741-2758.

［8］ Constantinescu E M，Sandu A，Chai T，et al.Ensemble-based chemical data assimilation.I：General approach.Q J R Meteorol Soc，2007，133：1229-1243.

［9］ Dee D P.On-line estimation of error covariance parameters for atmospheric data assimilation.Mon Wea Rev，1995，123（4）：1128-1145.

［10］ Dee D P，da Silva A M.Maximum-likelihood estimation of forecast and observation error covariance parameters.PartⅠ：Methodology.Mon Wea Rev，1999，127：1822-1834.

［11］ Dee D P，Gaspari G，Redder C，et al.Maximum-likelihood estimation of forecast and observation error covariance parameters.PartⅡ：Applications.Mon Wea Rev，1999，127：1835-1849.

［12］ Wang X，Bishop C H.A comparison of breeding and ensemble transform Kalman filter ensemble forecast schemes.J Atmos Sci，2003，60：1140-1158.

［13］ Li H，Kalnay E，Miyoshi T.Simultaneous estimation of covariance inflation and observation errors within an ensemble Kalman filter.Q J R Meteorol Soc，2009，135（639）：523-533.

［14］ Miyoshi T.The gaussian approach to adaptive covariance inflation and its implementation with the Local Ensemble Transform Kalman Filter.Mon Wea Rev，2011，139：1519-1535.

［15］ Zheng X G.An adaptive estimation of forecast error covariance parameters for Kalman filtering data assimilation.Adv Atmos Sci，2009，26：154-160.

［16］ Liang X，Zheng X G，Zhang S P，et al.Maximum likelihood es-timation of inflation factors on error covariance matrices for ensemble Kalman filter assimilation.Q J R Meteorol Soc，2011，138：263-273，doi：10.1002／qj.912.

［17］ Huang C L，Li X，Lu L.Retrieving soil temperature profile by assimilating MODIS LST products with ensemble Kalman filter.Rem Sens Environ，2008，112：1320-1336.

［18］ Yang K，Koike T，Kaihotsu I，et al.Validation of a dual-pass microwave land data assimilation system for estimating surface soil moisture in semiarid regions.J Hydrometeorology，2009，10（3）：780-793.

［19］ 楊曉峰，陸其峰，楊忠東.基于AMSR-E土壤濕度產(chǎn)品的LIS同化試驗.應(yīng)用氣象學(xué)報，2013，24（4）：435-445.

［20］ 吳統(tǒng)文，宋連春，劉向文，等.國家氣候中心短期氣候預(yù)測模式系統(tǒng)業(yè)務(wù)化進展.應(yīng)用氣象學(xué)報，2013，24（5）：533-543.

［21］ 王莉，黃嘉佑.Kalman濾波的試驗應(yīng)用研究.應(yīng)用氣象學(xué)報，1999，10（3）：276-282.

［22］ 趙曉琳，朱國富，李澤椿.基于TIGGE資料識別適應(yīng)性觀測敏感區(qū)的應(yīng)用研究.應(yīng)用氣象學(xué)報，2010，21（4）：405-415.

［23］ Whitaker J S，Hamill T H.Ensemble data assimilation without perturbed observations.Mon Wea Rev，2002，130：1913-1924.

［24］ Dai Y J，Zeng X B，Dickinson R E，et al.The common land model.Bull Amer Meteor Soc，2003，84：1013-1023.

［25］ 孟春雷，張朝林.路面氣象數(shù)值預(yù)報模型及性能檢驗.應(yīng)用氣象學(xué)報，2012，23（4）：451-458.

［26］ Koike T.Coordinated Enhanced Observing Period （CEOP）-An initial step for integrated global water cycle observation.World Meteorological Organization Bulletin，2004，53（2）：115-121.

［27］ 龔建東，趙剛.全球資料同化中誤差協(xié)方差三維結(jié)構(gòu)的準(zhǔn)確估計與應(yīng)用：背景誤差協(xié)方差調(diào)整與數(shù)值試驗分析.氣象學(xué)報，2006，64（6）：669-682.

［28］ 曹小群，黃思訓(xùn)，張衛(wèi)民，等.區(qū)域三維變分同化中背景誤差協(xié)方差的模擬.氣象科學(xué)，2008，28（1）：8-14.

［29］ 馬旭林，莊照榮，薛紀(jì)善，等.GRAPES非靜力數(shù)值預(yù)報模式的三維變分資料同化系統(tǒng)的發(fā)展.氣象學(xué)報，2009，67（1）：50-60.

［30］ 王曼，李華宏，段旭，等.WRF模式三維變分中背景誤差協(xié)方差估計.應(yīng)用氣象學(xué)報，2011，22（4）：482-492.

［31］ Jin R，Li X.Improving the estimation of hydrothermal state variables in the active layer of frozen ground by assimilating in situ observations and SSM／I data.Sci China Ser D：Earth Sci，2009，39（9）：1220-1231.

［32］ Anderson J L.Spatially and temporally varying adaptive covariance inflation for ensemble filters.Tellus，2009，61：72-83.

A Method of Improving Error Covariances in En KF and Its Application to Data Assimilation

Liang Xiao1）Zheng Xiaogu2）Dai Yongjiu2）Shi Chunxiang1）

1）（National Meteorological Information Center，Beijing100081）
2）（Beijing Normal University，Beijing100875）

In the ensemble Kalman filter（En KF），the forecast error covariance matrix is estimated as the sampling covariance matrix of the forecast ensemble.However，previous studies suggest that the sampling error resulting from finite-size ensembles may make such estimations far from the true forecast error covariance，and finally degrade the performance of En KF.A common way to address this problem is covariance inflation with a time-constant inflation factor.A time-dependent infiation approach on forecast error covariance matrix（i.e.，MLE method）is developed based on the maximum likelihood estimation theory，so as to improve estimates of forecast error covariances.At Delgertsgot（DGS）Station in the Mongolian Plateau reference site，point observations of ground temperature and soil temperature at the depth of 10 cm are assimilated into the Common Land Model（Co LM）with a frequency of every 12 hours，using two assimilation algorithms（En KF method and MLE method），in order to test the effectivity of MLE in practical assimilation.In this way，a soil temperature assimilation system is constructed on the point scale.

Results indicate that MLE method performs better than En KF method for ground temperature and soil temperatures at most depths（especially for soil temperatures at deeper depths）.Moreover，considering differences between soil temperatures at shallower depths and those at deeper depths，different inflation factors are adopted to them when implementing MLE method.Compared to results of MLE using a single scalar inflation factor，it shows that multiple-factor inflation is able to alleviate the unreasonable inflation of soil temperatures at deeper depths and therefore get better assimilation results.In addition，the leaf area index（LAI）in the Co LM is updated dynamically by MODIS LAI products，and results derived using MODIS LAI are compared to those derived using LAI computed by experiential formula，so as to study the effect of the LAI accuracy on simulated and assimilated soil temperatures.It shows that using MODIS LAI can get better simulation of soil temperature at depths of 0 cm and 3 cm，as well as more accurate assimilation of soil temperature at most depths.

The inflation factor is set to be variable in time，but constant in space.However，variables such as soil temperature and soil moisture behave quite differently at shallow surfaces and deep depths，and observations may be unevenly distributed in space in regional assimilation researches.Therefore，it is necessary to adopt different inflation factors to different variables or to the same variable at different locations.In the future，it is necessary to develop a time-and-space dependent inflation method and test its capability in real applications.

data assimilation；En KF；error covariance inflation

梁曉，鄭小谷，戴永久，等.En KF中誤差協(xié)方差優(yōu)化方法及在資料同化中應(yīng)用.應(yīng)用氣象學(xué)報，2014，25（4）：397-405.

2014-01-06收到，2014-05-05收到再改稿。

公益性行業(yè)（氣象）科研專項（GYHY201206013，GYHY201306045），國家國際科技合作專項（2011DFG23150）

*email：liangx＠cma.gov.cn