基于Richards 方程的冷杉樹高曲線深度神經網絡激活函數*

徐奇剛 雷相東 鄭 宇 胡興國 雷淵才 何 瀟

（1. 中國林業科學研究院資源信息研究所 北京 100091；2. 國家林業和草原局華東調查規劃院 杭州 310000；3. 國家林業和草原局森林經營與生長模擬實驗室 北京 100091；4. 吉林省汪清林業局 汪清 133200）

隨著計算機計算能力的大幅提升以及機器學習算法的突飛猛進，森林生長建模面臨新的機遇和挑戰，既要滿足處理大數據建模的需求，同時對估計精度的要求也越來越高。從“最佳有用”原則出發，機器學習算法，尤其是深度神經網絡算法具有對數據分布不作要求、可處理非線性數據、能處理連續和分類變量、預測精度高、數據適應力強等優勢，在林業建模領域的具有廣泛的應用前景（Humphrieset al.，2018；雷相東，2019），常用機器學習算法，如多元自適應回歸樣條、神經網絡、隨機森林、支持向量機和增強回歸樹等已應用于樹高-胸徑關系（梁瑞婷等，2021；李際平等，1996；陳佳琦等，2020）、單木直徑生長（Ouet al.，2019；Vieiraet al.，2018；馬翔宇等，2009；Reiset al.，2016；浦瑞良等，1999）和枯死模型（Mezeiet al.，2014；Sproullet al.，2015；Oguroet al.，2015；Hasenaueret al.，2001）等森林生長收獲預估領域。

然而，機器學習算法在林業建模應用過程中存在一個潛在缺點，即一個回歸模型的輸出可能不符合生物學邏輯。例如，在一個不被限制輸出的樹高生長神經網絡模型中，若未設置樹高的最大漸進值，隨著年齡增長，模型可能會輸出一個不斷變大偏離常識的樹高值；在林分直徑生長模型中，一個數據清洗不夠徹底、訓練不夠細致的神經網絡模型可能輸出負值。這些問題導致科研工作者在構建類似模型時仍更傾向于采用符合生物學邏輯的傳統回歸方法，如理論生長方程。

神經網絡模型在林業上的應用歷史較長（Hasenaueret al.，2001；Vinícius Oliveira Castroet al.，2013；Tavares Júnioret al.，2019），其處理過程主要通過激活函數實現，激活函數賦予神經網絡非線性的特性，將輸入轉換成輸出（雷相東，2019）。該過程中，最終輸出為最后一層神經元經激活函數轉換得到，故可在輸出層嵌套提前設計好的激活函數，以控制模型輸出使其符合預期，如對多標簽分類問題采用Softmax 函數方式（Liuet al.，2016）。但迄今為止，如何通過激活函數解決模型輸出違背生物學規律的問題仍未得到解決。

鑒于此，本研究提出一個基于理論生長方程（Richards 公式）的激活函數，解決神經網絡算法在森林生長建模時輸出可能不符合生物學規律的問題，并以臭冷杉（Abies nephrolepis）解析木樹高-胸徑數據為例進行驗證。新提出的激活函數可視作傳統方法與機器學習算法的結合，為神經網絡在森林生長建模方面的應用提供一個新的思路和方法。

1 研究區概況

研究區位于吉林省汪清林業局金溝嶺林場，地理位置為129°56′—131°04′E，43°05′—43°40′N。林場森林總面積16 286 hm2，海拔30～1 200 m，全年平均氣溫3.9 ℃左右，屬低山丘陵地貌。森林類型是以紅松（Pinus koraiensis）、云冷杉為優勢樹種的天然針闊葉混交林，主要樹種包括魚鱗云杉（Picea jezoensis）、紅松、楓樺（Betula costata）、臭冷杉、水曲柳（Fraxinus mandschurica）、黃菠蘿（Phellodendron amurense）、白樺（Betula platyphylla）、胡桃楸（Juglans mandshurica）等。

2 數據與方法

2.1 數據

數據源自1988 年8 月吉林省汪清林業局金溝嶺林場臭冷杉解析木數據，樹干解析流程見孟憲宇（2006）。在0 m（0 號盤）、1.3 m（1 號盤）、3.6 m（2 號盤）、5.6 m（3 號盤）等樹高處截取圓盤，查數各圓盤的年輪數。各齡階樹高采用內插法按比例得到，各齡階去皮胸徑通過1.3 m 圓盤（1 號盤）測量得到，應用樹皮系數（本研究臭冷杉樹皮系數為1.072）轉換得到帶皮胸徑。共計96 株伐倒木、458 組樹高-胸徑觀測數據，以8∶2 劃分建模集和測試集，366 組數據用于建模，92 組數據用于檢驗，單木樹高-胸徑統計量見表1。

表1 單木樹高-胸徑統計量Tab. 1 Statistics of tree height and diameter at breast height

2.2 深度神經網絡模型和新的基于理論生長方程的激活函數

深度神經網絡模型是一個擁有復雜層數和神經元結構的多層感知機模型，由輸入層、多層隱藏層和輸出層組成，層與層的單元全連接，將上一層輸出作為下一層輸入，通過激活函數轉化后繼續作為下一層輸出，逐層向后直至運算到輸出層，用梯度下降來最小化函數近似誤差（Goodfellowet al.，2016）。由于其強大的函數逼近能力，對輸入變量無統計上的分布要求，預測精度高，是現階段應用最廣的人工神經網絡。

對于一個L層的多層神經網絡（Rumelhartet al.，1986；Sibiet al.，2013），令輸入向量為：

輸出向量為：

第l隱藏層神經元的輸出為：

則最后一層神經元的輸出為：

式中：sl為第l層神經元個數；sL-1為第L-1層神經元個數；為第L-1層神經元到第L層神經元的權重為第L層第k個神經元的偏置（bias）；第L層第k個神經元激活函數內部值。

樣本總誤差為：

式中：N為輸入數據集的觀測樣本量；dk(q)為第q個訓練樣本的觀測值；yk(q)為第q個訓練樣本神經網絡輸出值。

反向傳播（BP）算法每次迭代按以下方式對權值和偏置進行更新：

在處理回歸問題時，現階段比較流行的做法是取輸出層（第L層）為線性激活函數，即：

此時，對于單個訓練樣本，BP 算法下輸出層參數的梯度更新公式為：

式中：t為第t輪訓練輪次；E(q)為第q個訓練樣本的損失函數值。

如此處理對輸出層沒有任何限制，但用于森林生長建模時存在缺陷，故本研究對最后一層即輸出層的激活函數進行調整（圖1），引入最常用的理論生長方程Richards 式，作為修正后最后輸出層的激活函數：

圖1 神經網絡模型的激活函數修正Fig. 1 Modified activation function of neural networks

式中：a、b和c為Richards 方程的3 個參數。

在反向傳播算法中，神經網絡模型Richards 激活函數的導數為：

此時，輸出層參數的梯度更新公式變為：

本研究中，Richards 方程的3 個參數a、b和c于R 語言平臺使用nls()函數（Ihakaet al.，1996）進行傳統非線性回歸擬合得到，采用默認的高斯-牛頓迭代算法。最終，a取20.689、b取0.070、c取1.418，即新的激活函數為：

2.3 模型比較與結構設計

除傳統回歸方法外，本研究比較3 種不同激活函數設置的深度神經網絡模型（表2）。

表2 3 個深度神經網絡模型結構與激活函數Tab. 2 Structure and activation functions of three deep neural network models

模型1 為使用Richards 公式作為輸出層（第L層）激活函數的深度神經網絡模型。梯度優化算法選擇“Adam”（Kingmaet al.，2014）。隱藏層激活函數選擇上，由于Sigmoid 函數具有廣泛的飽和性使梯度下降變得困難，整流線性函數ReLU 的行為更接近線性，模型更容易優化，通常情況下表現良好（Goodfellowet al.，2016），故本研究在大多數情況下選用整流線性函數；又考慮到ReLU 函數可能存在Dead ReLU 問題，即對小于0 的輸入的正向輸出和反向梯度更新值均為0，神經元“死亡”，故在第3 層選用Sotfplus 函數。輸出層激活函數使用Richards 公式（式15）。Richards公式本身在x大于0 時才能滿足輸出一定大于0，為使神經網絡模型輸出一定大于1.3，對第L-1 層神經元進行限制，激活函數選擇輸出一定大于0 的Softplus函數，并且通過指定第L層參數use_bias 為False，刪去最后一層的偏置指定kernel_constraint 為NonNeg()，限定非負，使最終模型1 的輸出能一定大于1.3。

對于隱藏層與隱單元數量的最優網絡設計問題，采用反復試錯法‘try and error’（Uzunet al.，2017），觀測模型在測試集上的誤差得到。batch_size 設為64，Epoch 設為1 000 輪。Callbacks 采用ModelCheckpoint策略，save_best_only 設為True，即只保存1 000 輪中損失值最小的模型。

模型2 為使用式（8）作為輸出層（第L層）激活函數的深度神經網絡模型。除第L-1 層使用ReLU 激活函數外，其余結構設計與各項訓練參數與模型1 一致。

模型3 為使用式（8）作為輸出層（第L層）激活函數的深度神經網絡模型。由于模型2 的激活函數行為接近線性，故加入模型3 作為對比，其結構設計與各項訓練參數與模型2 一致。與模型2 不同的是，模型3 最后一層隱藏層（第5 層）的激活函數選擇Tanh函數，第4 層選擇Sigmoid 函數，以此觀察隱藏層中帶有常見非線性激活函數的深度神經網絡模型表現。

深度神經網絡模型的建立和實現采用Python 語言的keras（Gulliet al.，2017）和tansorflow（Abadiet al.，2016）模塊完成。

2.4 模型評價

選取決定系數（R2）、均方根誤差（root mean squared error，RMSE）和平均絕對誤差（mean absolute error，MAE）作為模型評價指標：

式中：n為樣本數；Hi為第i株單木的樹高觀測值；為第i株單木的樹高預測值；為樣本內單木樹高觀測值的平均值。

3 結果與分析

3.1 基于Richards 激活函數的深度神經網絡模型

4 種模型訓練集和測試集精度檢驗結果如表3 所示。傳統模型、使用Richards 激活函數的模型1、與隱藏層使用Tanh 和Sigmoid 函數的模型3 精度表現非常接近。模型1 表現最好，從訓練集來看，R2比傳統模型提升0.051%，RMSE 和MAE 比傳統模型分別下降0.322%和0.117%；從測試集來看，R2比傳統模型提升0.175%，RMSE 和MAE 比傳統模型分別下降2.282%和4.011%。模型2 訓練集R2比傳統模型降低0.140%，RMSE 和MAE 比傳統模型分別提高0.876%和0.036%；測試集R2比傳統模型降低0.598%，RMSE和MAE 比傳統模型分別提高4.591%和2.378%。模型3 訓練集R2比傳統模型訓練集提升0.422%，RMSE和MAE 比傳統模型分別降低2.681%和2.242%；測試集R2比傳統模型降低0.476%，RMSE 和MAE 比傳統模型分別提高3.529%和0.185%。

表3 4 種模型精度檢驗結果Tab. 3 Accuracy test statistics based on 4 models

3.2 模擬環境下4 種模型的輸出表現

對于已建立好的臭冷杉樹高-胸徑傳統模型和3種深度神經網絡模型，本研究模擬0～100 cm 胸徑區間觀察4 種模型的輸出表現。輸入數據為服從[0,100)均勻分布的1 000 組胸徑數據，輸出樹高如圖2 所示。使用Richards 激活函數的深度神經網絡模擬圖與使用傳統方法表現類似，均有明顯的樹高最大漸進值20.69 m。模型2 輸出樹高隨胸徑增長呈線性遞增關系，但沒有樹高最大漸進值，不符合生物學邏輯。模型3 雖然表現出樹高在20.0 m 達到最大值，但曲線不光滑，模擬的胸徑超出訓練數據胸徑最大值（29.38 cm）后模型預測樹高達到最大值后保持恒定。

圖2 4 種模型在胸徑0～100 cm 情況下的輸出樹高Fig. 2 Predicted tree heights from 4 models with simulated diameter from 0 to 100 cm

4 討論

本研究提出一個基于Richards 方程的樹高-胸徑深度神經網絡激活函數，從模型檢驗結果來看，使用Richards 激活函數的深度神經網絡模型相比傳統非線性回歸模型（直接使用Richards 公式擬合）精度略有提高，測試集R2提升0.175%，RMSE 和MAE 分別降低2.282%和4.011%，總體上看神經網絡模型與傳統非線性回歸模型結果接近，與以往研究結果類似（徐奇剛等，2019；?z?el?ket al.，2010；2017；Vahedi，2016）。Shen 等（2020）采用人工神經網絡方法建立廣東省楊樹（Populus）人工林樹高-胸徑模型，與傳統非線性回歸方法進行比較，結果顯示神經網絡樹高-胸徑模型的R2增加10.3%，RMSE 和MAE 分別減少12%和13.5%，神經網絡方法擁有更高的泛化潛力。本研究中模型1 與傳統方法構建的Richards 模型最關鍵的區別在于其在胸徑與樹高的映射中間增加了一層映射關系，即最后一層神經元值與胸徑之間的映射，通過基于Richards 公式的激活函數估計樹高，這樣處理可使模型表達出更精細的信息。使用Richards 激活函數的模型1 相比使用傳統激活函數的模型3 精度也有所提升，測試集R2提升0.653%，RMSE 和MAE 分別降低5.613%與4.189%，這可能是因為Richards 的a、b和c參數值與Richards 本身的模型形式作為先驗信息加入神經網絡訓練過程，提高了模型精度。

在傳統非線性回歸建模過程中，Richards 等理論生長方程的各參數均具有生物學意義，且在樹高-胸徑模型中表現優異（Xuet al.，2022）。Lei 等（1997）論證樹高-胸徑關系時指出，在樹高-胸徑基礎模型選擇期間，應同時考慮數據相關及合理的生物學標準。一個合理的模型形式應該具有單調增量（dH/dD＞0）、S拐點（存在二階導為0）和最大漸進值。對于樹高-胸徑模型，Richards 公式中a代表樹高生長的最大漸進值，b代表相對生長速率參數。各樹種的參數值相對接近，林業工作者能夠依據工作建模經驗選擇合適的參數初值，使模型更加容易收斂，同時也能使模型輸出保證在符合生物學邏輯的合理范圍內。而基于機器學習算法的模型在建模數據處理不夠細致、超參數調整粗糙的情況下，可能會出現樹高無限增大甚至出現負值的現象。本研究模擬結果顯示，使用傳統神經網絡激活函數的模型2 和模型3 的輸出存在樹高沒有最大值以及模型在訓練集數據區間外預測值失真的情況，由于樹高生長存在極限，且樹高與胸徑并不是線性關系。而使用Richards 激活函數的模擬結果與傳統方法建立模型的模擬結果類似，避免樹高隨胸徑的線性增長，更符合生物學規律。

需要指出的是，在傳統回歸中，模型漸進值的確定取決于數據和基礎模型本身，如對于同一組觀測，舒馬赫模型預測的漸進值往往會比其他模型更大（Lee，2000），又如若用于擬合的數據本身樹高偏低，得到的樹高最大漸進值較低會導致出現與實際樹高最大值不相符合的情況（區域內的樹高最大值出現在樣地之外）。本研究中，Richards 公式作為激活函數加入神經網絡輸出層，公式中代表樹高最大漸進值的參數a和代表相對生長速率的參數b是來自預先擬合好的非線性回歸模型的參數結果，若數據本身覆蓋范圍不夠大，模型的泛化性能可能受到影響。預先設置好的參數值在模型訓練后不會被改變，神經網絡通過梯度更新，將輸入層變量與輸出層變量的函數關系均映射在網絡本身的權值中，并不會影響輸出層激活函數。與此同時，a、b和c參數也因為反向傳播算法影響每次更新的梯度變化值，對模型訓練產生直接影響；如果a、b和c參數設置不合理，會對模型訓練產生負面效果。

神經網絡模型是數據驅動模型，在輸出層加入Richards 激活函數并不能保證模型輸出一定存在漸近線，因為當訓練數據存在“倒長”，即數據中胸徑較大時對應的樹高小于胸徑較小時對應的樹高，那么神經網絡模型會被影響，以致最后一層的神經元值與胸徑在該區間內呈負相關，模型則會失去漸近線，同時也會造成模型在該區間預測時曲線并非單調遞增，但Richards 函數本身形式保證模型輸出一定存在一個合理的最大值。同理，由于神經網絡的權值更新由建模數據驅動，雖然基于Richards 公式激活函數的神經網絡模型不能保證在胸徑為0 時樹高一定等于1.3，但對第L-1 層的神經元進行合理限制后可保證模型輸出一定大于1.3。在神經網絡模型中，若試圖使dH/dD一定大于0，則需要施加限制使網絡中的權值W、偏置b以及組成網絡各激活函數的導數一定大于0，這種情況下，BP 算法擬合模型極其容易梯度消失，難以收斂，因此本研究并未加入這方面的限制。

本研究也比較了林業中其他常用的理論生長方程，如邏輯斯蒂和考爾夫公式，擬合效果均不如Richards 公式，故選擇Richards 公式作為輸出層的激活函數參考形式。雖然本研究最終確定的激活函數為Richards 方程，但從原理上講，基于理論生長方程修正神經網絡激活函數從而保證模型輸出更具有生物合理性這一方法具有普適性。在未來森林生長建模工作中，可對深度神經網絡模型的激活函數進行建模者本身的經驗選擇，以得到最佳可用的模型。

5 結論

與傳統回歸相比，機器學習算法具有對數據分布不作要求、可處理非線性數據、能處理連續和分類變量、預測精度高、數據適應力強等優勢，但存在輸出可能會偏離生物學規律的問題。本研究提出一個基于Richards 公式的激活函數控制神經網絡模型輸出，并得到一個輸出符合生物學邏輯的樹高-胸徑神經網絡模型。使用Richards 公式的激活函數具有如下優點：1） 輸出一定存在一個合理的最大值；2） 配合合理的神經網絡結構可使輸出一定大于1.3；3） 將傳統回歸方法擬合得到的參數作為神經網絡模型輸入，能使神經網絡的訓練得到先驗知識。