基于PROSAIL 模型的玉米冠層參數敏感性分析

李 青

（天津市測繪院有限公司，天津 300381）

植被葉面積指數等特征參數是表征植被生長發育狀態的重要量化指標，被廣泛應用于農業管理、環境檢測、生態保護和氣候變化等領域。植被冠層輻射傳輸物理模型建立了植被特征參數與植被光譜之間的數學關系，能夠通過遙感手段觀測光譜，進而反演植被特征參數。PROSAIL 模型是一種基于輻射傳輸理論建立的植被冠層反射率模擬模型，其輸入參數包括葉片結構和生理特征、冠層結構、土壤反射率、天空光散射比和觀測幾何等。為分析不同參數對應的光譜特征，并探究不同輸入參數在不同波段對模型輸出結果的影響，本文以玉米冠層為研究對象，利用OAT 和EFAST 敏感性分析方法，分析各輸入參數對冠層反射率變化的貢獻。

1 模型與方法

1.1 PROSAIL 模型

PROSAIL 模型是建立在PROSPECT 葉片模型和SAILH 冠層模型之上的。在葉片尺度上，PROSPECT提供了400～2 500 nm 光譜上葉片的方向半球反射率和透過率，并將其作為SAILH 的輸入參數。PROSPECT 模型以Allen 的平板模型（plate model）為基礎，將葉片描述為一個或多個疊加起來的表面粗糙的各向同性平板[1]。這個模型假設葉片是由N 層同性層疊加起來的，每兩層之間被空氣間隔，葉片結構參數N 就是用來描述葉片中的平均空氣/細胞間隔的。在葉子內部，光線被認為是各向同性的[2]。PROSPECT 模型的輸入參數包括：葉片結構參數，葉綠素含量、類胡蘿卜素含量、褐色素含量、葉片含水量和葉片干物質含量等。

SAILH 模型是Kuusk 在SAIL（Scattering by Arbitrarily Inclined Leaves）模型的基礎上引入熱效應得到的改進模型。SAIL 模型的數學形式是一個基于輻射傳輸方程的4 方程9 參數的線性微分方程組，是Verhoef[3]和Bunnik 根據Suits 模型擴充得到的。由于輻射傳輸方程微分-積分的形式無法求得其嚴格的解析解，應用比較廣泛的解法是KM 理論，即用一組微分方程近似輻射傳輸方程，只考慮下行漫散射和上行漫散射。Suits模型在此基礎上加入了觀測方向的輻亮度[4]。SAIL 模型則以Suits 為基礎，加入葉片傾角和分布函數，可以計算任意天空散射光比例、太陽高度、觀測方向下的冠層反射率[3]。由于SAIL 模型不能很好地對熱點效應進行模擬，因此Kuusk[5]在SAIL 模型的基礎上加入熱點效應從而發展出了SAILH 模型，引入熱點效應因子，熱點效應因子為葉片長度和冠層高度的函數。SAILH模型的輸入參數包括：葉面積指數、平均葉傾角（葉傾角分布函數）、熱點參數、土壤反射率、天空光散射比、太陽天頂角、觀測天頂角和相對方位角等。

1.2 模型參數設置

本文以玉米冠層為研究對象，設定PROSAIL 模型輸入參數見表1。根據玉米的特征，玉米葉片大多平展且狹長，在生長的過程中葉片長度與植株高度的比值從接近1 逐漸減小到0.25 左右。因此將平均葉傾角（ALA）的默認值設為50°，取值范圍為30～60°；熱點參數（hspot）默認值設為0.25，取值范圍為0.25～1。除此之外，由于褐色素含量與葉綠素含量相沖突，且其含量不易測量，因此將褐色素含量設為0，即衰老的植株不予考慮。

表1 PROSAIL 模型輸入參數及取值范圍

1.3 敏感性分析方法

敏感性分析可以確定輸入參數對模型模擬結果影響的大小，分析方法主要分為局部敏感性分析和全局敏感性分析。

局部敏感性分析是研究單個輸入參數在局部范圍內變化時對模型輸出響應的影響，局部敏感性分析方法簡單、計算量小，但其缺點在于無法充分描述模型參數的空間分布形態，并且忽略了參數之間的相互作用[6]。本文采用一次一個變量法（OAT）計算每個參數在其默認值附近進行變化導致的模型輸出的變化率，其變化率的絕對值代表了該參數的敏感性。由于之后的全局敏感性分析會得出個輸入參數在各波段的一階敏感性，因此這一部分不做定量分析，僅定性展示各輸入參數的局部敏感性。

全局敏感性分析則能夠分析整個參數空間對模型結果的影響，考慮參數之間的相互作用。常用的全局敏感性分析方法包括回歸分析法、Morris 篩選法、Sobol’方法和EFAST 方法等[7]。PROSAIL 模型數據間呈現非線性關系且各輸入參數間相互作用，不宜采用回歸分析法。Morris 篩選法可以簡單快速地篩選出敏感性高的參數，但無法量化影響，Sobol 方法可以定量分析各參數敏感性大小，但對于復雜模型計算成本高。因此本文選用能定量分析且計算量較小的EFAST 方法進行敏感性分析。

EFAST 方法是由Saltelli 等結合Sobol’法和傅里葉振幅靈敏度檢驗法（Fourier amplitude sensitivity test，FAST）的優點提出的基于方差分解的全局敏感性分析方法[7]，可以計算得到每個參數及參數相互作用對模型結果的方差的貢獻率。該方法的核心思想是通過搜索函數對模型的每個參數引入一個共同獨立參數的函數，并給每個輸入參數定義一個整數頻率，最后使模型成為獨立參數的周期函數，對模型的輸出結果進行傅里葉變換，獲得傅里葉級數的頻譜曲線，通過該頻譜曲線分別得到參數的敏感性。

式中：Vi是參數xi引起的方差，Vij、Vijm、…、V1，2，…，k是各參數相互作用貢獻的方差。那么，通過歸一化處理，參數xi的一階敏感性指數可定義如下

一階敏感性指數反應的是參數對模型輸出總方差的直接貢獻率。同理，參數xi的二階及三階敏感性指數可定義為

參數xi對模型輸出總方差的直接貢獻和通過參數間相互作用對模型輸出方差的貢獻可以用總敏感性指數STi來表示，即各階敏感性指數之和

該方法采樣的總個數為

式中：Nr是搜索曲線的重復次數，可取2；M 是干涉因子，一般取4 或更高；ωmax為最大采樣頻率；n 為參數個數[7]。

EFAST 敏感性分析借助專業敏感性分析軟件Simlab（Version 2.2.1）[8]完成。方法實現過程包括樣本生成、模型選擇、Monete Carlo 模擬及模擬結果的不確定性和敏感性分析。具體步驟如下所示。

1）對PROSAIL 模型的所有輸入變量在取值范圍內按均勻分布采樣；

為確保在行洪安全的前提下，充分合理地開發利用沂沭河汛期洪水資源，實現洪水資源化，臨沂市通過新建和改建河道蓄水工程，對河道來水進行節節攔蓄、梯級開發，汛期調蓄洪水，非汛期攔蓄來水，從而有效緩解流域內工農業用水矛盾，改善區域水環境，實現區域水循環。目前，已在祊河、沂河沿城區段建成7座橡膠壩和8座攔河閘，在臨沂城形成了88.5 km長、面積達48.5 km2的水面，一次性蓄水可達1.62億m3。

2）使用MATLAB 語言編寫的程序讀入樣本，運行PROSAIL 模型進行模擬；

3）將PROSAIL 模型輸出的反射率按Simlab 軟件樣本輸出格式將寫入文件；

4）在Simlab 中選擇結果文件，執行Monete Carlo模擬和EFAST 敏感性分析。

2 結果與分析

2.1 局部敏感性分許

局部敏感性分析采用的方法是OTA，對PROSAIL模型的葉綠素含量（Cab）、類胡蘿卜素含量（Car）、干物質含量（Cw）、葉片水分含量（Cm）、葉片結構參數（N）、葉面積指數（LAI）、熱點參數（hspot）和平均葉傾角（ALA）8個冠層參數進行敏感性分析，每次只改變一個輸入參數的值，其他參數取表1 中的默認值，觀察輸入參數的改變對冠層反射光譜在400～2 500 nm 波段的影響，結果如圖1 所示。

圖1 冠層參數局部敏感性OTA 分析結果

由圖1（a）（b）可知，對于葉片中的色素來說，Cab和Car對冠層反射光譜的影響主要集中在可見光波段，其中對Cab最敏感的波段為550～670 nm，對Car最敏感的波長為525 nm。Cab對反射率的影響程度遠遠大于Car。隨著色素含量的增加，冠層的反射率減小，且減小的幅度越來越小。這可能是由于色素對可見光特定波段的強吸收造成的。由圖1（c）可知，水分子在1 200、1 400、1 950 nm 處存在吸收谷，因此Cw在1 200～2 500 nm 的波段范圍內對植被反射率均有較大影響，反射率隨著葉片水分的增加而減小。對于冠層葉面積指數來說，當LAI＞2 時，冠層反射率隨著LAI 的增大而增大；當LAI≤2 時，冠層反射率在400～700 nm 和1 250～2 500 nm 處有明顯的抬高，這是因為地表植被稀疏造成土壤背景裸露，模擬得到的反射光譜中包含了土壤的光譜特性。ALA 對反射率的影響主要集中在700～1 400 nm 波段，葉傾角越大（葉片越平）反射率越小。除此之外，干物質含量（Cm）、葉片結構參數（N）和熱點因子（hspot）沒有表現出在特定波段的敏感性，對反射率的影像不大。葉片結構參數和熱點因子與冠層反射率成正相關，干物質含量則相反。

針對觀測幾何參數，進行如下分析：①保持觀測天頂角和相對方位角均為0°，改變太陽高度角，發現隨著太陽高度角的減小，植被冠層的反射率增加。分析原因主要是由于觀測角度和太陽角度逐漸靠近，并且植被陰影逐漸減少。②保持太陽高度角為30°，改變觀測天頂角和相對方位角，計算觀測半球上的冠層反射率，得到結果如圖2 所示。從圖中可以發現，整體上隨著觀測天頂角的減小和相對方位角的增大，冠層反射率減小，但相對方位角為0°時除外，在這個角度可以觀察到明顯的熱點效應。

圖2 觀測半球的反射率（黑圈為太陽入射位置）

2.2 全部敏感性分析

對于表1 中的所有冠層生理結構參數進行EFAST全局敏感性分析，對400～2 000 nm 的波段范圍結果進行輸入輸出均勻采樣，得33 組的敏感性分析結果，同時計算得到了各輸入參數的一階敏感性（局部敏感性）和全局敏感性。實驗結果表明，各輸入參數在各波段的一階敏感性和全局敏感性趨勢基本一致。全局敏感性較強的輸入參數包括：LAI、Cab、Cw和ALA。在400～500 nm、750～950 nm 和1 900～2 000 nm 波段，LAI 對反射率方差的貢獻最大，全局敏感性指數均超過0.5。Cab在500～750 nm 波段占據絕對優勢，對反射率方差貢獻高達0.7。而Cw則是在1 150～1 850 nm 成為敏感性最高的輸入參數，全局敏感性大多超過0.6。

針對能夠反應植被光譜特征的藍光、綠光、紅光、近紅外波段，對輸入參數的敏感性進行了排序和具體分析。如圖3 所示，在藍光波段和近紅外波段LAI 敏感性最高，分別超過0.8 和0.4；在綠光和紅光波段，Cab的敏感性最高，能夠解釋超過一半的反射率方差變化，這是由于葉綠素吸收紅光反射綠光的特性造成的。

圖3 植被光譜特征波段EFAST 敏感性分析結果

3 結論

本文對PROSPECT 模型的輸入參數進行了定性的局部敏感性分析和定量的EFAST 全局敏感性分析，得到以下結論：

全局敏感性較強的輸入參數包括：冠層葉面積指數LAI、葉片葉綠素含量Cab、葉片水分含量Cw和平均葉傾角ALA。

1）對于葉片葉綠素含量Cab，其對冠層反射光譜的影響主要集中在綠光和紅光波段，能夠解釋超過一半的反射率方差變化。隨著葉綠素含量的增加，冠層的反射率減小。

2）水分子在1 200、1 400、1 950 nm 處存在吸收谷，葉片水分含量Cw在1 150～1 850 nm 波段為敏感性最高的輸入參數，反射率隨著葉片水分的增加而減小。

3）對于冠層葉面積指數來說，在400～500 nm、750～950 nm 和1 900～2 000 nm 波段，LAI 對反射率方差的貢獻最大，全局敏感性指數均超過0.5。當LAI＞2時，冠層反射率隨著LAI 的增大而增大；當LAI≤2時，地表植被稀疏造成土壤背景裸露，模擬得到的反射光譜中包含了土壤的光譜特性，導致冠層反射率在400～700 nm 和1 250～2 500 nm 處有明顯的抬高。

4）隨著太陽高度角的減小，植被冠層的反射率增加；保持太陽高度角為30°，觀測半球上的冠層反射率整體上隨著觀測天頂角的減小和相對方位角的增大，冠層反射率減小。相對方位角為0°時可以觀察到明顯的熱點效應。