李得勤 張述文 段云霞 崔錦

1中國氣象局沈陽大氣環(huán)境研究所，沈陽110016

2蘭州大學(xué)半干旱氣候變化教育部重點(diǎn)實(shí)驗室，蘭州730000

3沈陽市氣象臺，沈陽110168

1 引言

參數(shù)在陸面過程模擬中扮演著重要的角色，陸面模式模擬的效果不僅取決于模型的物理參數(shù)化方案、初始場和驅(qū)動數(shù)據(jù)的質(zhì)量，也取決于參數(shù)的代表性和準(zhǔn)確性（李新等，2007；林朝暉等，2008），對參數(shù)的合理校準(zhǔn)可以改進(jìn)陸面模式的模擬效果（Henderson-Sellers, 1996; Duan et al., 2006；房云龍等，2010；梁曉和戴永久，2010）。陸面模式中的參數(shù)一般分為兩種，一種是具有明確物理涵義的參數(shù)，如地表反照率、植被覆蓋度、地形坡度、葉面指數(shù)、植被高度等，這些參數(shù)一般不取決于模式中所使用的參數(shù)化方案；而另外一種是所謂的“功能性”參數(shù)，如氣孔阻力、空氣動力粗糙度、土壤熱容、土壤水?dāng)U散率和傳導(dǎo)率等，這些參數(shù)很難通過直接測量來獲得（Gupta et al., 1999; Bastidas et al.,2006）。目前，陸面模式對具有實(shí)際物理意義的參數(shù)一般用非均勻區(qū)域（網(wǎng)格）內(nèi)取平均值的辦法來估計，由于不同陸面模式所使用的參數(shù)集（土壤和植被）分辨率不同，加之受到植被和土壤分類方法不同的影響（Gupta et al., 1999），即使是同一種土壤或植被類型，具有相同物理意義的參數(shù)在不同模式中也存在差異，所以，參數(shù)的不確定是造成陸面模式模擬不準(zhǔn)確的主要來源之一。

參數(shù)敏感性分析和優(yōu)化方法，可以先確定關(guān)鍵的敏感參數(shù)，然后通過調(diào)整或訂正這些敏感參數(shù)使模擬結(jié)果更接近于觀測值（Gupta et al., 1999）。早期的參數(shù)優(yōu)化方法通常一次只校正一個參數(shù)，優(yōu)化過程效率較低，并且不考慮參數(shù)間的內(nèi)在聯(lián)系（Bastidas et al., 1999）。多判據(jù)方法的出現(xiàn)部分克服了前期參數(shù)校正中的困難，優(yōu)化方法也可以一次性優(yōu)化多個參數(shù)，大量優(yōu)化算法已陸續(xù)運(yùn)用在陸面模式的參數(shù)優(yōu)化工作中（Chen et al., 2009）；Bastidas et al.（1999）將多判據(jù)廣義敏感性分析方法（Multi-Objective Generalized Sensitivity Analysis, 簡稱MOGSA）應(yīng)用在 BATS（Biosphere-Atmosphere Transfer Scheme）模式的參數(shù)敏感性分析和優(yōu)化中，該方法現(xiàn)為校正陸面模式參數(shù)常用方法之一(Demarty et al., 2004; Bastidas et al., 2006)；Xia et al.（2005）使用貝葉斯隨機(jī)反演（Bayesian Stochastic Inversion, 簡稱BSI）算法探究模式的復(fù)雜度和參數(shù)值域?qū)?shù)估計的影響，發(fā)現(xiàn)模式越復(fù)雜改進(jìn)越明顯，并且模式復(fù)雜度對參數(shù)校正的影響比參數(shù)值域大；Santanello et al.（2007）使用了被動和主動微波遙感土壤濕度的資料，借助于參數(shù)優(yōu)化工具 PEST（Parameter ESTimation model）對Noah模式中土壤成份參數(shù)進(jìn)行優(yōu)化，通過與站點(diǎn)觀測的土壤濕度的對比發(fā)現(xiàn)，參數(shù)優(yōu)化后的Noah模式明顯改進(jìn)了對土壤濕度的模擬。梁曉和戴永久（2010）使用高斯誤差傳播原理分析CoLM（Commom Land Model）模式中參數(shù)的誤差對土壤溫度、土壤濕度、植被蒸散量模擬的誤差貢獻(xiàn)。

隨著土壤濕度觀測手段的提高，通過參數(shù)優(yōu)化方法獲取具有空間一致性和代表性的土壤參數(shù)可能成為現(xiàn)實(shí)（李得勤等, 2012），進(jìn)而校準(zhǔn)陸面模式中的相關(guān)數(shù)據(jù)集和靜態(tài)參數(shù)（Chen et al., 2009；Hogue et al., 2005）。Santanello et al.（2007）認(rèn)為模式中土壤分類大多沒有考慮土壤質(zhì)地的空間異質(zhì)性，經(jīng)過優(yōu)化可以得到具有空間“一致性”和明確“物理意義”的土壤參數(shù)。Scott et al.（2000）使用SCE-UA（Shuffled Complex Evolution Algorithm）方法優(yōu)化了土壤模型參數(shù)后發(fā)現(xiàn)土壤濕度的模擬得到了提高。Mattikalli et al.（1998）使用遙感反演的表層土壤濕度來識別流域尺度土壤質(zhì)地和估計土壤水文參數(shù)，指出微波亮溫和土壤濕度的時間變化與土壤成份和飽和水傳導(dǎo)率之間存在很好的相關(guān)性。除上述單獨(dú)進(jìn)行參數(shù)優(yōu)化外，還有的方法是把參數(shù)優(yōu)化和數(shù)據(jù)同化結(jié)合在一起來提高狀態(tài)量和參數(shù)的估算精度，如Tian et al.（2009, 2010）和Yang et al.（2007）使用雙通道同化系統(tǒng)，第一個通道同化觀測資料的同時進(jìn)行參數(shù)優(yōu)化，優(yōu)化后的參數(shù)在第二個通道中用來再次同化觀測資料；Li and Ren（2011）借助集合卡曼濾波（Ensemble Kalman Filter,簡稱EnKF），使用理想試驗手段驗證了EnKF在同化土壤水勢后，可以有效地改進(jìn)土壤模型的參數(shù)。

陸面模式中的參數(shù)包括植被和土壤相關(guān)的參數(shù)，植被參數(shù)大多具有季節(jié)性變化特征，而土壤相關(guān)參數(shù)相對比較穩(wěn)定。而在大多數(shù)針對陸面模式參數(shù)的優(yōu)化工作中，常使用地表熱通量（感熱通量、潛熱通量、土壤熱通量）和土壤濕度對模式中所有參數(shù)進(jìn)行優(yōu)化，很少考慮到參數(shù)與觀測變量之間的相關(guān)性。此外，優(yōu)化算法雖然已經(jīng)被用來校正輻射傳輸模式中的部分參數(shù)來提高陸面數(shù)據(jù)同化系統(tǒng)對土壤濕度的模擬（Tian et al., 2009, 2010;Yang et al., 2007），但是仍缺乏真正用來校正和分析陸面模型中土壤濕度相關(guān)參數(shù)的內(nèi)容。加之由于可以同化的物理量相對有限，而陸面模式物理過程的復(fù)雜性和參數(shù)數(shù)量眾多使得陸面模式參數(shù)優(yōu)化工作仍然是一件比較復(fù)雜的工作。本研究將借助SCE-UA方法和觀測系統(tǒng)模擬試驗的手段，探究不同土壤參數(shù)對土壤濕度求解的影響，以及使用土壤濕度來優(yōu)化土壤參數(shù)的可行性和不確定性，進(jìn)而對分析陸面模式誤差來源和將優(yōu)化算法用于改進(jìn)陸面模式中土壤相關(guān)參數(shù)提供一些參考。

2 土壤濕度方程求解和SCE-UA方法

2.1 土壤水分特征曲線和參數(shù)

土壤含水量和水勢之間關(guān)系又叫做土壤水分特征曲線，不同的土壤類型水分傳輸過程不同，土壤水分特征曲線存在明顯的差別（孫菽芬，2005）。現(xiàn)有的描述土壤水分特征曲線的模型有 Brook-Corey模型、Clapp-Hornberger模型、van Genuchten模型、Campbell模型、Mualem模型等。上述模型均為非線性模型，其中van Genuchten模型被認(rèn)為是當(dāng)前描述土壤水分特征曲線最理想的模型，但是由于其高度的非線性，使得求解過程中離散化比較復(fù)雜。相對來說，Clapp-Hornberger模型的應(yīng)用則更廣泛，其具體形式為：

其中，b為與土壤類型有關(guān)的參數(shù)，ψ和θ為土壤水勢和含水量，ψs和θs為飽和水勢和飽和含水量。土壤水導(dǎo)率K為：

其中Ks為飽和土壤水導(dǎo)率。土壤水?dāng)U散率D表示為：

從上面的關(guān)系式可以看出，土壤相關(guān)的參數(shù)主要有：飽和含水量θs、飽和水勢ψs、b和飽和土壤水導(dǎo)率Ks。

目前，陸面模式對上述土壤參數(shù)（b，θs，ψs，Ks）的給出一般采用以下兩種方式：一種是通過土壤成份間接求出，如 CLM（Community Land Model）模式中主要采用這種方法（Oleson et al.,2004, 2010）。另一種方式直接指定每種土壤類型對應(yīng)的土壤性質(zhì)參數(shù)（如b，θs，ψs，Ks），目前大多數(shù)陸面模式中的土壤參數(shù)以這種方式給出。土壤類型一般分表層和深層兩種類型，也有整層土壤采用同種土壤類型。但是受地形數(shù)據(jù)來源和分辨率的不同，即使同一種土壤類型，同一土壤參數(shù)在不同陸面模式中也會出現(xiàn)取值不同現(xiàn)象。

2.2 土壤濕度方程和求解方法

現(xiàn)有陸面模式中對土壤濕度運(yùn)動過程大多采用一維Richards方程來描述，可以寫為土壤含水量、土壤水勢，以及二者混合形式的三種形式。以土壤含水量為例，土壤濕度方程可表示為：

其中，t為時間，z為垂直坐標(biāo)。

土壤濕度方程是一個高度的非線性方程，在一般初始和邊界條件下很難得到解析解，所以常尋求數(shù)值解（Shao and Irannejad, 1999; Xie et al., 1999）。不同陸面模式對其進(jìn)行離散化的方法不盡相同，CLM模式使用增量法（Oleson et al., 2004）求解，澳大利亞聯(lián)邦科學(xué)與工業(yè)研究會（Commonwealth Scientific and Industrial Research Organization，簡稱CSIRO）的CABLE（CSIRO Atmosphere Biosphere Land Exchange）陸面模式用分裂算法（Kowalczyk et al., 2006）求解，中尺度天氣研究與預(yù)報（Weather Research and Forecasting，WRF）模式所采用的Noah陸面模式為全隱式差分格式（Chen and Dudhia,2001）。ECMWF（European Center for Medium range Weather Forecasting）中采用了半隱式的方法求解（Viterbo and Beljaars, 1995），張述文等（2009，2010）對不同的離散化方法求解土壤濕度方程的差別進(jìn)行了對比。

這里采用和Noah陸面模式中相同的全隱式的離散化方法，使用追趕法對一維Richards方程進(jìn)行求解，下邊界條件采用靜力平衡條件，詳細(xì)的離散化過程見 Chen and Dudhia（2001）和張述文等（2010）。

2.3 SCE-UA方法

SCE-UA方法使用種群競爭演變（Competitive Complex Evolution, 簡稱CCE）算法（Duan et al.,1992），將取樣樣本分為若干個種群分區(qū)，每一分區(qū)取樣相互無關(guān)并各自搜索最優(yōu)點(diǎn)，但是種群之間可以相互傳遞搜索得到的信息來更新搜索的分區(qū)樣本，使得所有的參數(shù)同時整體達(dá)到最優(yōu)。大量的研究證實(shí)該方法可以同時優(yōu)化模型中的若干個參數(shù)，成為參數(shù)優(yōu)化的一種有效方法（Duan et al.,1992; Scott et al., 2000; Liu et al., 2005）。

在SCE-UA方法的具體使用中，根據(jù)待優(yōu)化模型的復(fù)雜程度和參數(shù)的個數(shù)，需要設(shè)置控制概率和一些優(yōu)化問題相關(guān)的參數(shù)。也可根據(jù)待優(yōu)化參數(shù)的個數(shù)，其他參數(shù)采用默認(rèn)設(shè)置：m=(2n+1)；q=(n+1)；y=1；Z=(2n+1)，其中，n為參數(shù)的個數(shù)，m為種群中的樣本數(shù)，q為子種群中的樣本數(shù)，y為每個子種群產(chǎn)生后代的個數(shù)，Z為每個種群和參數(shù)演變的次數(shù)。關(guān)于SCE-UA使用中的具體步驟可以詳細(xì)的參考Duan et al.（1992）文獻(xiàn)。此外，在使用SCEUA算法的過程中，還需要提到以下兩個概念：

（1）目標(biāo)函數(shù)：目標(biāo)函數(shù)用來評價實(shí)際觀測與模型模擬值之間的吻合度，不同的目標(biāo)函數(shù)用來評價模型的不同特征，所以目標(biāo)函數(shù)的選擇也非常重要，Gupta et al.（1999）列出的用于參數(shù)優(yōu)化的目標(biāo)函數(shù)接近 10個，在這里我們采用使用較多的土壤含水量的均方根偏差來定義目標(biāo)函數(shù)F：

其中，I為總的土壤層數(shù)，N為總的積分次數(shù)為第i層、第k個積分時次土壤含水量的模擬值，為在第k個積分時次的土壤含水量的“真值”（在實(shí)際模型參數(shù)優(yōu)化過程中一般取觀測值）。

（2）迭代停止判據(jù)：迭代停止判據(jù)被用來決定迭代過程是否終止，滿足如下三個條件之一：

1）在連續(xù)的 5次迭代中，目標(biāo)函數(shù)值的改進(jìn)程度已經(jīng)不能達(dá)到0.01%，認(rèn)為參數(shù)已經(jīng)收斂到一個很小且平滑區(qū)域上。

2） 搜尋算法已經(jīng)不能通過合理的改變參數(shù)來使得目標(biāo)函數(shù)值進(jìn)一步減小。

3）考慮到計算時間，當(dāng)?shù)螖?shù)達(dá)到超出預(yù)先設(shè)置的次數(shù)以后，迭代停止。在實(shí)際問題中，按照模型的復(fù)雜程度，迭代次數(shù)一般給出一個較大的值，一般情況下在達(dá)到最大值前，優(yōu)化過程已經(jīng)結(jié)束，所以，當(dāng)?shù)螖?shù)越多，也就意味著參數(shù)已接近最優(yōu)了，即基本滿足條件 1），詳細(xì)的介紹參考Duan et al.（1992）。

3 試驗設(shè)計

土壤濕度方程的求解需要給出初始和邊界條件，本研究采用和張述文等（2009）驗證土壤濕度方程求解方法時所使用相同的初始和邊界條件，考慮深度為3 m的土壤層，將其劃分為11層，表1給出了土壤層信息和求解過程中需要的初始和邊界條件信息，這里沒有對所有的土壤類型一一驗證，而只選取了比較典型的壤土（loam）土壤類型作為參考。

表1 土壤層信息、初始條件和邊界條件Table 1 Soil layers, initial and boundary conditions

在2.1節(jié)中給出了土壤濕度方程求解中需要的參數(shù)，以及陸面模式中參數(shù)的確定方法，根據(jù)土壤參數(shù)的確定方法，這里待優(yōu)化模型參數(shù)和對應(yīng)“真值”的給出也分為兩種情況：①以土壤成份的形式：將土壤中粘粒的百分比（Pclay）、粉粒的百分比（Psilt）和沙粒的百分比（Psand）作為待優(yōu)化的參數(shù)，相關(guān)參數(shù)根據(jù)百分比求得（Oleson et al., 2004,2010），由于在土壤質(zhì)地三角圖中每一種土壤類型對應(yīng)的成份存在一定分布區(qū)間，這里隨機(jī)選取壤土土壤區(qū)間中某一固定點(diǎn)對應(yīng)的土壤成份百分比作為 “真值”；②以參數(shù)的形式給出：包括土壤性質(zhì)相關(guān)的參數(shù)b，θs，ψs，Ks，直接根據(jù)Noah模式中壤土土壤對應(yīng)的參數(shù)值作為 “真值”。在具體優(yōu)化過程中，不同參數(shù)的分區(qū)區(qū)間見表 2，這里沒有將土壤層深度和初始條件作為待優(yōu)化的參數(shù)，具體參數(shù)的給出和優(yōu)化過程見第4節(jié)。

表2 土壤參數(shù)分布區(qū)間Table 2 Soil parameter ranges

根據(jù)理想實(shí)驗的設(shè)置方法，根據(jù)上面提到的參數(shù)“真值”的確定方法，將土壤濕度方程在蒸發(fā)的上邊界條件下求解，時間步長為 300 s，總共積分20 d，計算得到“真實(shí)”的土壤濕度廓線。而在優(yōu)化過程中，根據(jù)各個參數(shù)可能的分布區(qū)間，給出初始值，在整個區(qū)間中尋找參數(shù)的最優(yōu)設(shè)置，SCE-UA方法可以同時優(yōu)化若干個參數(shù)。

4 結(jié)果分析和討論

4.1 將土壤成份(Pclay，Psilt，Psand)作為待優(yōu)化參數(shù)

聯(lián)合國糧食和農(nóng)業(yè)組織（Food and AgricultureOrganization, 簡稱 FAO）和美國國家土壤地理（State Soil Geographic, 簡稱STATSGO）根據(jù)土壤中沙粒（sand），粘粒（clay）和粉粒（silt）物質(zhì)含量的不同，將土壤分為 12種類型，每一種土壤的成份均固定。而根據(jù)美國農(nóng)業(yè)部（United States Department of Agriculture, 簡稱USDA）建立的土壤分類標(biāo)準(zhǔn)，土壤類型同樣分為 12種，而每一種土壤成份則規(guī)定在一定的區(qū)間上，12種土壤類型最終可以在土壤質(zhì)地三角圖中表示出來。

為了驗證SCE-UA算法在土壤濕度求解中對參數(shù)優(yōu)化的效果，根據(jù)FAO和STATSGO中土壤分類標(biāo)準(zhǔn)，壤土的成份組成為Pclay=18%，Psilt=39%，Psand=43%，假設(shè)整個土壤柱的成份不隨深度變化，則在優(yōu)化過程中只有兩個參數(shù)（Pclay和Psilt，Psand=1－Pclay－Psilt）。根據(jù)土壤濕度方程的求解，結(jié)合初始條件和邊界條件可以得到土壤濕度廓線，將此作為土壤濕度的“真實(shí)”廓線。在優(yōu)化過程中，Pclay和Psilt的分布區(qū)間均在 0～100%之間（表 1），使用SCE-UA方法在整個區(qū)間上尋找粘粒和粉粒的百分比。

圖1是經(jīng)SCE-UA優(yōu)化100次后Pclay和Psilt分別與目標(biāo)函數(shù)值和迭代次數(shù)的點(diǎn)聚圖，這里只給出了每次迭代中最優(yōu)的種群樣本參數(shù)分布。優(yōu)化中參數(shù)滿足Pclay+Psilt≤100%，所以設(shè)置了當(dāng)不滿足該條件時，目標(biāo)函數(shù)值取 1.0。首先，從參數(shù)的分布和目標(biāo)函數(shù)值間的關(guān)系來看，土壤濕度對Pclay很敏感，隨著優(yōu)化的進(jìn)行，迅速地將參數(shù)Pclay的區(qū)間縮小，達(dá)到 0～40%，而粉粒則沒有粘粒敏感，但是隨著優(yōu)化的進(jìn)行，主要分布范圍為0～90%。其次，圖中看出迭代次數(shù)可以作為參數(shù)優(yōu)化程度的一個判斷標(biāo)準(zhǔn)，當(dāng)?shù)螖?shù)越多，參數(shù)的收斂范圍越小，即更接近于“真值”（黑實(shí)線）。

為了定量地衡量優(yōu)化后參數(shù)相對于“真值”的靠近程度，這里引進(jìn)了參數(shù)偏差絕對百分比（E），即優(yōu)化后參數(shù)與“真值”的誤差相對于“真值”的絕對百分比，具體定義如下：

式中，xe為優(yōu)化后得到的參數(shù)值，xtrue為參數(shù)的“真值”。

Pclay和Psilt作為待優(yōu)化的參數(shù)，當(dāng)?shù)螖?shù)達(dá)到 272次后，得到的參數(shù)估計值分別為Pclay=17.9629%和Psilt=39.1074%，目標(biāo)函數(shù)值為0.0018。根據(jù)式6得到參數(shù)Pclay和Psilt的偏差絕對百分比依次為0.21%和0.28%。所以，將土壤成份作為待優(yōu)化的參數(shù)時，借助于土壤含水量，SCE-UA方法很容易使得參數(shù)達(dá)到最優(yōu)，優(yōu)化效果非常明顯。

圖1 參數(shù)（Pclay和Psilt）與（a, b）目標(biāo)函數(shù)值和（c, d）迭代次數(shù)的點(diǎn)聚圖，黑色直線為“真值”Fig.1 The scatter diagrams of parameters Pclay and Psilt vs.(a, b) the objective function value and (c, d) iteration times, black straight line is “true”

根據(jù) FAO和 STATSGO的全球土壤成份數(shù)據(jù)集，認(rèn)為表層土壤（30 cm以內(nèi)）和深層土壤（30 cm以下）成份不同。在這里分別使用和表示土壤成份在0～31 cm和31～300 cm的粘粒百分比，與之類似，沙粒和粉粒采用同樣的表達(dá)方法。根據(jù)土壤質(zhì)地三角圖，在壤土的區(qū)間內(nèi)隨機(jī)選取兩點(diǎn)，其成份為=18%，=25%，=39%，=45%，=43%和=30%。所以，選取待優(yōu)化的參數(shù)為 4個：，，和。在優(yōu)化過程中，同樣設(shè)置當(dāng)+≥100%或+≥100%時，目標(biāo)函數(shù)值為1.0。圖 2給出了 SCE-UA優(yōu)化以后，總的優(yōu)化次數(shù)為100次，四個參數(shù)的分布與目標(biāo)函數(shù)值和迭代次數(shù)的分布特征?？梢钥吹筋愃朴谥苯訉clay和Psilt作為待優(yōu)化的參數(shù)，使用SCE-UA優(yōu)化后的參數(shù)非常接近于“真值”，當(dāng)?shù)螖?shù)達(dá)到727次后，，，和的優(yōu)化結(jié)果為 17.5476%，24.9736%，40.8750%和 45.0584%，目標(biāo)函數(shù)值為0.001828，根據(jù)式（6）得到優(yōu)化后，，和的偏差絕對百分比依次為 2.51%，4.81%，0.11%和0.13%。可以看出，參數(shù)增多后，保持相同的優(yōu)化次數(shù)的前提下，優(yōu)化后參數(shù)的誤差略有增大，但優(yōu)化效果依然較好。

綜合上面的結(jié)果可以看出，當(dāng)將土壤成份作為待優(yōu)化的參數(shù)時，由于參數(shù)分布區(qū)間比較相近，使用土壤濕度作為判據(jù)，可以得到參數(shù)的最優(yōu)值，且參數(shù)增多的情況下，優(yōu)化效果仍然較好，優(yōu)化后參數(shù)的偏差百分比普遍在 5%以下，結(jié)果說明 SCEUA方法能夠達(dá)到優(yōu)化參數(shù)的目的。同時可以看出SCE-UA也可以作為參數(shù)敏感性判斷的方法，土壤參數(shù)粘粒較粉粒更為敏感，收斂速度更快。

4.2 將 b，θs，ψs，Ks作為待優(yōu)化參數(shù)

相對于不同于使用土壤成份計算相關(guān)參數(shù)的方法，一些陸面模型中則是借助于靜態(tài)參數(shù)集來給出不同土壤類型的相關(guān)參數(shù)?，F(xiàn)有的大量關(guān)于優(yōu)化陸面模型參數(shù)的工作也大都對土壤相關(guān)的四個變量進(jìn)行了優(yōu)化，但優(yōu)化過程中不考慮觀測物理量（感熱通量，潛熱通量或土壤濕度）與這些參數(shù)之間的關(guān)系，而使用所有的觀測物理量對模型中幾乎所有的參數(shù)同時進(jìn)行優(yōu)化。由于土壤相關(guān)參數(shù)對土壤濕度廓線的模擬影響較大，與已往工作不同，本文就使用土壤濕度廓線觀測對描述土壤性質(zhì)的四個參數(shù)（b，θs，ψs和Ks）進(jìn)行優(yōu)化，考察SCE-UA對參數(shù)訂正的效果和準(zhǔn)確度。

土壤含水量“真值”同樣用 loam對應(yīng)的參數(shù)設(shè)置模擬得到（b=5.25,θs=0.439 cm3/cm3,ψs=－355.0 cm，Ks=0.000338 cm/s），初始和邊界條件見表1，時間步長同樣選取300 s，共積分20 d。

表2給出參數(shù)b，θs，ψs和Ks的范圍，可以看出四個參數(shù)值的分布區(qū)間差異較大，尤其是土壤水傳導(dǎo)率Ks的分布區(qū)間跨度較大，最小值和最大值相差3個量級，且其值本身較小，由于其直接控制著土壤中水分的傳導(dǎo)和擴(kuò)散過程，不同土壤類型對應(yīng)的土壤水傳導(dǎo)率Ks的差別較大（根據(jù) USGS的參數(shù)集，如沙壤土（loamy sand）和粘土（clay）的水傳導(dǎo)率相差兩個量級，所以對該參數(shù)的優(yōu)化難度較大，當(dāng)同時對所有的參數(shù)進(jìn)行優(yōu)化時，對其他參數(shù)能否達(dá)到最優(yōu)的影響究竟如何也需要進(jìn)一步驗證。

在 4.1節(jié)中土壤成份作為參數(shù)時，經(jīng)過優(yōu)化100次以后參數(shù)基本接近“真值”，而在這里優(yōu)化100次后，優(yōu)化效果并不明顯，因此將優(yōu)化次數(shù)增加到300次。同時，從圖1和圖2可以發(fā)現(xiàn)，參數(shù)的最優(yōu)值更容易從迭代次數(shù)上明顯地看出，因此在圖3中給出了參數(shù)b，θs，ψs和Ks與迭代次數(shù)的分布圖?？梢钥闯觯c對土壤成份優(yōu)化結(jié)果不同，只有參數(shù)b能快速并最終與“真值”吻合，其他個參數(shù)則出現(xiàn)多個次優(yōu)點(diǎn)，最終結(jié)果也并不能和“真值”完全一致。

當(dāng)總的優(yōu)化次數(shù)為 300次時，迭代次數(shù)超出500次的次優(yōu)點(diǎn)出現(xiàn)了9個，只選擇迭代次數(shù)多且保證目標(biāo)函數(shù)值最小的次優(yōu)點(diǎn)作為最終優(yōu)化的結(jié)果。當(dāng)?shù)螖?shù)為 830次以后，目標(biāo)函數(shù)值達(dá)到0.000003，而當(dāng)?shù)螖?shù)為749次以后，目標(biāo)函數(shù)值達(dá)到0.000002，所以，將后者的參數(shù)作為最終結(jié)果。最終，優(yōu)化得到b，θs，ψs和Ks的值分別為5.2499，0.4166，－467.2969和 0.000167，對應(yīng)的偏差絕對百分比為0，5.1%，31.63%和50.59%?？梢钥闯?，當(dāng)參數(shù)的不確定分布區(qū)間較大時，優(yōu)化以后的結(jié)果略差，并且不同的參數(shù)優(yōu)化效果不同。

SCE-UA也可以對參數(shù)進(jìn)行敏感性分析：當(dāng)參數(shù)敏感時容易得到該參數(shù)的最優(yōu)值，而參數(shù)不敏感時，則易陷入次優(yōu)值，難找到最優(yōu)值。Bastidas et al.（2006）和Roundy（2009）分別使用不同站點(diǎn)觀測資料和不同參數(shù)優(yōu)化算法對Noah模式中不同參數(shù)的敏感性進(jìn)行了檢驗，檢驗結(jié)果發(fā)現(xiàn)四個參數(shù)敏感性的次序依次為Ks＞b＞θs＞ψs。而在本文試驗中，由于Ks的取值范圍較大，取值比較小，達(dá)到“真值”比較困難，而b易得到最優(yōu)值，從四個參數(shù)偏差絕對百分比上來看，參數(shù)敏感性的順序依次為b＞θs＞ψs＞Ks，這和 Roundy（2009）的結(jié)果不完全相同。θs，ψs，Ks三個參數(shù)具有一個明顯的“敏感區(qū)間”，如果將迭代次數(shù)和目標(biāo)函數(shù)值依然作為敏感區(qū)間的判據(jù)依據(jù)，則loam土壤θs的敏感區(qū)間為[0.4166, 0.4569]，ψs的敏感區(qū)間為[－467.2969，－287.7977]，Ks的敏感區(qū)間為[0.000167，0.000580]。值得注意的是，參數(shù)的優(yōu)化不僅得到了參數(shù)的敏感性區(qū)間，同時得到了參數(shù)的組合信息，只有合適的組合才能得到最小的目標(biāo)函數(shù)值（詳見4.3節(jié)）。

圖2 參數(shù)（，，和）與目標(biāo)函數(shù)值（上）和迭代次數(shù)（下）的點(diǎn)聚圖，黑色直線為“真值”Fig.2 The scatter diagrams of parameters, , and vs.the objective function value (top) and iteration times (bottom), black straight line is “true”

圖3 參數(shù)b，θs，ψs和Ks與迭代次數(shù)的點(diǎn)聚圖，黑色直線為“真值”Fig.3 The scatter diagrams of parameters b, θs, ψs, and Ks vs.iteration times, black straight line is “true”

為了驗證參數(shù)的優(yōu)化區(qū)間對優(yōu)化效果的影響，縮小參數(shù)優(yōu)化區(qū)間，將四個參數(shù)b，θs，ψs和Ks優(yōu)化的區(qū)間依次定義為[4.0，7.0]，[0.4，0.5]，[－500，－200]和[0.0001，0.001]，優(yōu)化次數(shù)同樣選取了300次。迭代次數(shù)超過500次的次優(yōu)點(diǎn)共出現(xiàn)29個，且目標(biāo)函數(shù)值均為0.000002。圖4給出四個參數(shù)的最優(yōu)值對應(yīng)的迭代次數(shù)，類似于圖 3，只有參數(shù)b能達(dá)到最優(yōu)并接近于真值，而其它三個參數(shù)則出現(xiàn)更多的次優(yōu)點(diǎn)。和上面的試驗一樣，選取迭代次數(shù)最多的點(diǎn)作為最終優(yōu)化結(jié)果，迭代次數(shù)最多達(dá)到735次，優(yōu)化后的參數(shù)b，θs，ψs和Ks值分別為5.2500，0.4324，－384.3251和 0.000276，偏差絕對百分比依次為0，1.5%，8.26%和18.34%。所以，通過縮小參數(shù)的搜索區(qū)域，可以在一定程度上提高參數(shù)優(yōu)化的程度。

為了看出次優(yōu)點(diǎn)與真實(shí)土壤含水量廓線的差別，圖5分別給出用次優(yōu)點(diǎn)參數(shù)和“真值”進(jìn)行模擬得出的第一天（黑色）、第10天（紅色）和第20天（藍(lán)色）土壤含水量廓線，作為對比也給出由未訂正參數(shù)模擬得出的廓線?？梢钥闯?，雖然次優(yōu)點(diǎn)的參數(shù)值和“真值”存在一定的差別，但是計算得到的土壤含水量廓線和使用“真值”的土壤含水量廓線基本吻合。這說明基于土壤含水量的均方根誤差來定義目標(biāo)泛函（5）式很難在眾多的次優(yōu)點(diǎn)中準(zhǔn)確找到參數(shù)的“真值”，隱含著問題的超定性，需要添加其他附加約束條件。

圖4 同圖3，但縮小了參數(shù)范圍Fig.4 Same as Fig.3, but with reduced parameter ranges

圖5 分別使用參數(shù) “真值”、優(yōu)化參數(shù)和初始參數(shù)設(shè)置模擬的第1天（黑色）、第10天（紅色）和第20天（藍(lán)色）土壤含水量廓線Fig.5 Soil moisture profiles simulated respectively with the inputs of true parameters, optimal parameters, and priori parameters with no calibration on the first day (black), tenth day (red), twentieth day (blue)

從上面結(jié)果看出，SCE-UA方法可以用于對土壤濕度模型中的參數(shù)優(yōu)化和參數(shù)敏感性檢驗，優(yōu)化程度依賴于模型、待優(yōu)化參數(shù)的性質(zhì)和優(yōu)化次數(shù)。當(dāng)待優(yōu)化參數(shù)分布區(qū)間較大且不敏感時，很難搜索到參數(shù)的“真值”。通過增加優(yōu)化次數(shù)在一定程度上改善了優(yōu)化的效果?？傮w而言，參數(shù)優(yōu)化算法對土壤參數(shù)的優(yōu)化效果較明顯。

4.3 參數(shù)個數(shù)和參數(shù)組合對參數(shù)優(yōu)化的影響

從4.2節(jié)中的優(yōu)化結(jié)果來看，同時對b，θs，ψs和Ks四個參數(shù)優(yōu)化時，部分參數(shù)無法完全與“真值”吻合，出現(xiàn)超定情況。考慮到優(yōu)化參數(shù)的個數(shù)和參數(shù)之間的相關(guān)性對優(yōu)化效果的影響，這里通過減少優(yōu)化參數(shù)的個數(shù)和使用不同的參數(shù)組合來對模型參數(shù)進(jìn)行優(yōu)化。試驗發(fā)現(xiàn)：在b，θs，ψs和Ks四個參數(shù)中，選取其中三個組合進(jìn)行優(yōu)化，只要每個組合包含b均能收斂到“真值”（圖略）；當(dāng)把參數(shù)個數(shù)進(jìn)一步減少到2個后，不論何種參數(shù)組合，優(yōu)化效果均較好（圖略）。

從4.2節(jié)結(jié)果中可知b是一個較敏感的參數(shù)，在不同的試驗中優(yōu)化效果均較好。在下面優(yōu)化試驗不包含b，只選θs，ψs和Ks作為待優(yōu)化參數(shù)，試驗結(jié)果見圖 6。迭代超過 300次的次優(yōu)點(diǎn)竟出現(xiàn) 19個，其中有9個次優(yōu)點(diǎn)目標(biāo)函數(shù)值達(dá)到最小。按照上面確定最優(yōu)點(diǎn)的判斷方法，最大迭代次數(shù)529且目標(biāo)函數(shù)值最小時，參數(shù)θs，ψs和Ks的值分別為0.4552，－293.3358和0.000552，對應(yīng)的偏差絕對百分比依次為3.69%，17.37%和63.31%。對比發(fā)現(xiàn)上述結(jié)果類似于4.2節(jié)中將b，θs，ψs和Ks同時作為待優(yōu)化參數(shù)的試驗，并且優(yōu)化參數(shù)后對應(yīng)的土壤含水量廓線和“真實(shí)”廓線吻合得很好（圖略），這也與 4.2節(jié)中結(jié)果類似。同時發(fā)現(xiàn)迭代次數(shù)為431次的次優(yōu)點(diǎn)更接近于真值，θs，ψs和Ks取值分別為0.4374，－362.0123和0.000321，對應(yīng)的偏差絕對百分比依次為0.36%，1.98%和5.03%。以上結(jié)果說明，隨意減少參數(shù)個數(shù)不一定提高參數(shù)的優(yōu)化結(jié)果，同時說明最優(yōu)點(diǎn)的最終選擇還與參數(shù)敏感性和最優(yōu)點(diǎn)判據(jù)有關(guān)。

圖6 同圖4，但是不包括參數(shù)bFig.6 Same as Fig.4, but not including b

5 結(jié)論和討論

陸面模式中包含大量土壤和植被性質(zhì)的參數(shù)，參數(shù)的不確定性是陸面模式誤差的主要來源之一。地表通量（感熱通量和潛熱通量）和土壤濕度被廣泛用來對模式參數(shù)進(jìn)行敏感性分析和校準(zhǔn)。本文使用參數(shù)優(yōu)化算法SCE-UA結(jié)合一維土壤濕度方程的求解來檢驗參數(shù)優(yōu)化算法對土壤濕度方程中相關(guān)參數(shù)的優(yōu)化效果，初步得到以下結(jié)論：

（1）類似于CLM模式中，土壤濕度方程相關(guān)的參數(shù)可以通過土壤成份計算得到，當(dāng)將土壤成份作為模式待優(yōu)化的參數(shù)時，SCE-UA算法的優(yōu)化效果很好，并且土壤分層增多時參數(shù)個數(shù)增多對優(yōu)化效果的影響并不大。

（2）當(dāng)把土壤濕度方程中參數(shù)直接作為待優(yōu)化參數(shù)時，由于不同參數(shù)的分布區(qū)間差別較大，特別是水導(dǎo)率Ks的值比較小，且分布區(qū)間跨度較大，使用SCE-UA方法優(yōu)化以后，除b外其余參數(shù)均存在次優(yōu)點(diǎn)，優(yōu)化過程無法與“真值”吻合，但可通過縮小參數(shù)分布范圍提高優(yōu)化效果。

（3）SCE-UA算法也可以作為參數(shù)敏感性分析的工具，對于較敏感的參數(shù)，較容易找到最接近“真值”的最優(yōu)點(diǎn)。而對于不敏感的參數(shù)和分布區(qū)間跨度較大的參數(shù)容易陷入局部最優(yōu)。但是從土壤含水量廓線上難于分辨，說明問題可能存在超定性。

（4）最后，考察減少參數(shù)個數(shù)和改變參數(shù)組合對參數(shù)優(yōu)化的影響，發(fā)現(xiàn)隨意減少參數(shù)個數(shù)不一定提高優(yōu)化效果，參數(shù)組合反應(yīng)了不同參數(shù)之間某種內(nèi)在聯(lián)系，通過減少不敏感參數(shù)個數(shù)以及合理的參數(shù)設(shè)置可以提高參數(shù)優(yōu)化的效果。

本研究用觀測系統(tǒng)模擬試驗手段，目的是消除未知的模式誤差來源，達(dá)到對SCE-UA算法公正合理評價的目的。實(shí)際陸面模型參數(shù)個數(shù)可能更多，具體的判據(jù)也可能不止一個，也會涉及到不同判據(jù)的歸一化問題，不確定來源較多，難于分辨和定量化，更需要更多高質(zhì)量的觀測數(shù)據(jù)，所有這些都值得進(jìn)一步研究。

致謝感謝北京師范大學(xué)全球變化與地球系統(tǒng)科學(xué)研究院段青云教授提供SCE-UA方法代碼。感謝兩位審稿專家提出的寶貴意見。

（References）

Bastidas L A, Gupta H V, Sorooshian S, et al.1999.Sensitivity analysis of a land surface scheme using multicriteria methods [J].J.Geophys.Res.,104 (D16): 19481–19490.

Bastidas L A, Hogue T S, Sorooshian S, et al.2006.Parameter sensitivity analysis for different complexity land surface models using multicriteria methods [J].J.Geophys.Res., 111: D20101.

Chen F, Dudhia J.2001.Coupling an advanced land surface-hydrology model with the Penn State-NCAR MM5 Modeling System [J].Mon.Wea.Rev., 129 (4): 569–585.

Chen Wen, Zhu Deqin, Liu Huizhi, et al.2009.Land-air interaction over arid/semi-arid areas in China and its impact on the East Asian summer monsoon.Part I: Calibration of the land surface model (BATS) using multicriteria methods [J].Advances in Atmospheric Sciences, 26 (6):1088–1098.

Demarty J, Ottlé C, Braud I, et al.2004.Using a multiobjective approach to retrieve information on surface properties used in a SVAT model [J].J.Hydrol., 287 (1–4): 214–236.

Duan Q Y, Sorooshian S, Gupta V K.1992.Effective and efficient global optimization for conceptual rainfall-runoff models [J].Water Resour.Res.,28 (4): 1015–1031.

Duan Q, Schaake J, Andréassian V, et al.2006.Model parameter estimation experiment (MOPEX): An overview of science strategy and major results from the second and third workshops [J].J.Hydrol., 320 (1–2):3–17.

房云龍, 孫菽芬, 李倩, 等.2010.干旱區(qū)陸面過程模型參數(shù)優(yōu)化和地氣相互作用特征的模擬研究 [J].大氣科學(xué), 34 (2): 290–306.Fang Yunlong, Sun Shufen, Li Qian, et al.2010.The optimization of parameters of land surface model in arid region and the simulation of land–atmosphere interation [J].Chinese Journal of Atmospheric Sciences(in Chinese), 34 (2): 290–306.

Gupta V K, Bastidas L A, Sorooshian S, et al.1999.Parameter estimation of a land surface scheme using multicriteria methods [J].J.Geophys.Res.,104 (D16): 19491–19503.

Henderson-Sellers A.1996.Soil moisture: A critical focus for global change studies [J].Global and Planetary Change, 13 (1-4): 3–9.

Hogue T S, Bastidas L, Gupta H.et al.2005.Evaluation and transferability of the Noah land surface model in semiarid environments [J].Journal of Hydrometeorology, 6 (1): 68–84.

Kowalczyk E A, Shao Y P, Law R M, et al.2006.The CSIRO atmosphere biosphere land exchange (CABLE) model for use in climate models and as an offline model [R].CRIRO marine and atmospheric research paper 013, 20–24.

Li Chao, Ren Li.2011.Estimation of unsaturated soil hydraulic parameters using the ensemble Kalman filter [J].Vadose Zone Journal, 10 (4):1205–1227.

李得勤, 段云霞, 張述文.2012.土壤濕度觀測、模擬和估算研究 [J].地球科學(xué)進(jìn)展, 27 (4): 424–434.Li Deqin, Duan Yunxia, Zhang Shuwen.2012.Soil moisture measurement and simulation: A review [J].Advances in Earth Science (in Chinese), 27 (4): 424–434.

李新, 黃春林, 車濤, 等.2007.中國陸面數(shù)據(jù)同化系統(tǒng)研究的進(jìn)展與前瞻 [J].自然科學(xué)進(jìn)展, 17 (2): 163–173.Li Xin, Huang Chunlin, Che Tao, et al.2007.Development of a Chinese land data assimilation system:Its process and prospects [J].Process in Natural Science (in Chinese), 17(2): 881–892.

林朝暉, 劉輝志, 謝正輝, 等.2008.陸面水文過程研究進(jìn)展 [J].大氣科學(xué), 32 (4): 935–949.Lin Zhaohui, Liu Huizhi, Xie Zhenghui, et al.2008.Recent progress in the land-surface and hydrological process studies [J].Chinese Journal of Atmospheric Sciences (in Chinese), 32(4):935–949.

梁曉, 戴永久.2010.陸面模式中土壤和植被經(jīng)驗參數(shù)隨機(jī)誤差的傳播研究 [J].大氣科學(xué), 34 (2): 457–470.Liang Xiao, Dai Yongjiu.2010.Soil and plant parameters-related stochastic uncertainty propagation in the common land model [J].Chinese Journal of Atmospheric Sciences (in Chinese), 34 (2): 457–470.

Liu Y Q, Gupta H V, Sorooshian S, et al.2005.Constraining land surface and atmospheric parameters of a locally coupled model using observational data [J].Journal of Hydrometeorology, 6 (2): 156–172.

Mattikalli N M, Engman E T, Jackson T J, et al.1998.Microwave remote sensing of temporal variations of brightness temperature and near-surface soil water content during a watershed-scale field experiment, and its application to the estimation of soil physical properties [J].Water Resour.Res., 34 (9): 2289–2299.

Oleson K W, Dai Y J, Bonan G, et al.2004.Technical description of the Community Land Model (CLM) [R].NCAR Tech.Note NCAR/TN-461_STR, 113–123.

Oleson K W, David M L, Bonan G B, et al.2010.Technical description of Version 4.0 of the Community Land Model (CLM) [R].NCAR Tech.Note NCAR/TN-461_STR, 136–160.

Roundy J K.2009.Uncertainty Analysis for Land Surface Model Predictions: Application to the Simple Biosphere 3 and Noah Models at Tropical and Semiarid Locations [M].USA: Utah State University, 28–54.

Santanello J A Jr, Peters-Lidard C D, Garcia M E, et al.2007.Using remotely-sensed estimates of soil moisture to infer soil texture and hydraulic properties across a semi-arid watershed [J].Remote Sens.Environ., 110 (1): 79–97.

Scott R L, Shuttleworth W J, Keefer T O, et al.2000.Modeling multiyear observations of soil moisture recharge in the semiarid American Southwest [J].Water Resour.Res., 36 (8): 2233–2247.

Shao Y P, Irannejad P.1999.On the choice of soil hydraulic models in land-surface schemes [J].Boundary-Layer Meteor., 90 (1): 83–115.

孫菽芬.2005.陸面過程的物理、生化機(jī)理和參數(shù)化模型 [M].北京: 氣象出版社, 29–46.Sun Shufen.2005.Physics, Biochemistry and Parameterization of Land Surface Model [M].Beijing: China Meteorological Press (in Chinese), 29–46.

Tian X D, Xie Z H, Dai A G, et al.2009.A dual-pass variational data assimilation framework for estimating soil moisture profiles from AMSR-E microwave brightness temperature [J].J.Geophys.Res., 114:D16102.

Tian X D, Xie Z H, Dai A G, et al.2010.A microwave land data assimilation system: Scheme and preliminary evaluation over China [J].J.Geophys.Res., 115: D21113.

Viterbo P, Beljaars A C M.1995.An improved land surface parameterization scheme in the ECMWF model and its validation [J].J.Climate, 8 (11): 2716–2748.

Xia Y L, Yang Z L, Stoffa P L, et al.2005.Optimal parameter anduncertainty estimation of a land surface model: Sensitivity to parameter ranges and model complexities [J].Advances in Atmospheric Sciences,22 (1): 142–157.

Xie Z H, Luo Z D, Zeng Q C, et al.1999.Numerical simulation of moisture content and flux for the unsaturated soil water flow problem [J].Progress in Natural Science, 9 (9): 679–686.

Yang K, Watanabe T, Koike T, et al.2007.Auto-calibration system developed to assimilate AMSR-E data into a land surface model for estimating soil moisture and the surface energy budget [J].J.Meteor.Soc.Japan, 85A: 229–242.

張述文, 李得勤, 邱崇踐.2009.三類陸面模式模擬土壤濕度廓線的對比研究 [J].高原氣象, 28 (5): 988–996.Zhang Shuwen, Li Deqin, Qiu Chongjian.2009.A comparative study of the three land surface models in simulating the soil moisture profile [J].Plateau Meteorology (in Chinese),28 (5): 988–996.

張述文, 李得勤, 劉彥華, 等.2010.土壤濕度方程求解方法的比較 [J].蘭州大學(xué)學(xué)報 (自然科學(xué)版), 46(4): 46–53.Zhang Shuwen, Li Deqin,Liu Yanhua, et al.2010.A comparative study of the difference schemes for solving soil moisture equation [J].Journal of Lanzhou University(Natural Sciences) (in Chinese), 46 (4): 46–53.