向南海

（旌德縣蔡家橋林場(chǎng)，安徽 宣城 242600）

淺析理查德方程擬合樹木高生長(zhǎng)模型及參數(shù)求定

向南海

（旌德縣蔡家橋林場(chǎng)，安徽 宣城 242600）

在林學(xué)方面，描述樹木及林分生長(zhǎng)過程時(shí)，理查德方程是近代應(yīng)用最為廣泛、適應(yīng)性較強(qiáng)的一類生長(zhǎng)曲線方程。筆者就速生樹種杉木進(jìn)入中齡林后用理查德方程擬合樹木高生長(zhǎng)模型并求定參數(shù)，在此與大家做一交流。

理查德方程；參數(shù)

1 樹木生產(chǎn)量測(cè)定意義和樹木生長(zhǎng)方程概念

1.1 樹木生長(zhǎng)量測(cè)定意義

找出樹木生長(zhǎng)規(guī)律測(cè)定樹木生長(zhǎng)量，在森林經(jīng)營管理上有著很重要的意義，它既反映立地條件的好壞，又可以作為判斷營林效果和森林生產(chǎn)能力以及確定年伐量和主伐年齡的重要依據(jù)。

在現(xiàn)代林業(yè)特別是集體林權(quán)制度深化改革中，森林資源作為資產(chǎn)納入資本運(yùn)作，給林業(yè)發(fā)展增添了新的活力，且森林資源資產(chǎn)評(píng)估工作已經(jīng)作為必不可少的工作環(huán)節(jié)納入多個(gè)經(jīng)濟(jì)領(lǐng)域。在一個(gè)森林經(jīng)理期內(nèi)，測(cè)定林木生長(zhǎng)量的預(yù)期收獲值直接影響森林資產(chǎn)價(jià)值評(píng)估結(jié)果，給林業(yè)資本運(yùn)作也帶來一定的牽引作用。

1.2 樹木生長(zhǎng)方程的基本概念

樹木生長(zhǎng)方程是描述某樹種（組）各調(diào)查因子總生長(zhǎng)量y（t）隨年齡（t）生長(zhǎng)變化規(guī)律的數(shù)學(xué)模型。它是該樹種某調(diào)查因子的平均生長(zhǎng)過程，也就是在均值意義上的生長(zhǎng)方程。主要特點(diǎn)如下：

當(dāng)t＝0時(shí)y（t）＝0。此條件稱之為樹木生長(zhǎng)方程應(yīng)滿足的初始條件。

y（t）存在一條漸進(jìn)線y（t）＝A，A是該樹木生長(zhǎng)極大值。

樹木的生長(zhǎng)是不可逆的，使得y（t）是關(guān)于年齡（t）的單調(diào)非減函數(shù)。

y（t）是關(guān)于t的連續(xù)且光滑的函數(shù)曲線。

2 理查德方程產(chǎn)生背景與應(yīng)用范圍

2.1 理查德方程產(chǎn)生背景

理查德（Richards）方程是基于著名的Bertalanffy生長(zhǎng)理論發(fā)展而來。方程的基本形式為：y=A（1±Bexp-kt）1/（1-m）（當(dāng)m＞1時(shí)取“+”號(hào)，0＜m＜1時(shí)取“-”號(hào)）。該方程的參數(shù)或參數(shù)的組合有明顯的生物學(xué)意義，而且隨著方程中參數(shù)m的取值不同，它構(gòu)成了不同的生物生長(zhǎng)或種群動(dòng)態(tài)方程。因此，理查德方程有著廣泛的適用性，邏輯性強(qiáng)，參數(shù)可由獨(dú)立的試驗(yàn)加以驗(yàn)證，即參數(shù)可作出生物學(xué)解釋，從理論上對(duì)未來的趨勢(shì)進(jìn)行預(yù)測(cè)。其產(chǎn)生背景簡(jiǎn)述如下：

Bertalanffy通過分析動(dòng)物的生長(zhǎng)，發(fā)現(xiàn)在動(dòng)物生長(zhǎng)期間，動(dòng)物的體重增長(zhǎng)速率為同化速率與消耗速率之差，而后兩者分別和同化器官的大小以及動(dòng)物體重成比例，即：

式中：F—同化器官重，W—體重，Ra—同化速率，Rt—消耗速率。

由相對(duì)生長(zhǎng)關(guān)系，有F=rWm，因此：

式中：η=ar

理查德（Richards）認(rèn)為，m的取值為2/3，從而得到著名的理查德生長(zhǎng)方程：

2.2 理查德方程在林業(yè)中的應(yīng)用

（1）在胸徑和斷面積生長(zhǎng)中的應(yīng)用：林木胸徑和斷面積生長(zhǎng)曲線滿足t=t0（生長(zhǎng)至1.3 m 所需的年齡），y＝0的初始條件。

（2）在樹木高生長(zhǎng)中的應(yīng)用滿足當(dāng)t=0時(shí)、y=0的初始條件。

因此，通過引入?yún)?shù)使對(duì)樹木生長(zhǎng)具有廣泛的適應(yīng)能力，主要應(yīng)用在以下幾個(gè)方面：1）與林分密度有關(guān)的林分?jǐn)嗝娣e生長(zhǎng)模型；2）樹高生長(zhǎng)及地位指數(shù)方面研究；3）林分直徑結(jié)構(gòu)研究方面。

3 理查德方程應(yīng)用實(shí)例

3.1 擬合對(duì)象

擬合對(duì)象為杉木中齡林高生長(zhǎng)過程。林業(yè)調(diào)查中，按齡組、齡級(jí)劃分條件定義杉木中齡林生長(zhǎng)階段為林齡10至20年，5年1個(gè)齡級(jí)跨Ⅲ和Ⅳ二個(gè)齡級(jí)。

3.2 擬合模型

將方程兩邊各取對(duì)數(shù)得：

式中：t為樹木年齡觀測(cè)值，H為樹高，a、k、c為參數(shù)，其中k為模擬生長(zhǎng)規(guī)律時(shí)制定的修正值，一般在1.0以內(nèi)。

3.3 參數(shù)求定

3.3.1 參數(shù)k值求定

利用計(jì)算機(jī)應(yīng)用程序求定參數(shù)k值，應(yīng)用程序環(huán)境為：Microsoft Excel 2010。

第一步，列出杉木中齡林年齡與樹高觀測(cè)值于Excel表中。在實(shí)際調(diào)查中，可根據(jù)不同坡向坡位等不同立地條件調(diào)查多組數(shù)據(jù)，這里僅列出三組調(diào)查數(shù)據(jù)如下：

年齡( t ) 1 0 1 1 1 2 1 3 1 4 1 5 1 6 1 7 1 8 1 9 2 0樹高( h ) 6 . 5 7 . 6 8 . 2 9 . 0 9 . 9 1 1 . 0 1 1 . 7 1 2 . 3 1 3 . 2 1 4 . 0 1 4 . 6年齡( t ) 1 0 1 1 1 2 1 3 1 4 1 5 1 6 1 7 1 8 1 9 2 0樹高( h ) 6 . 4 7 . 5 8 . 0 9 . 1 1 0 . 0 1 0 . 6 1 1 . 3 1 2 . 4 1 3 . 0 1 3 . 8 1 4 . 7年齡( t ) 1 0 1 1 1 2 1 3 1 4 1 5 1 6 1 7 1 8 1 9 2 0樹高( h ) 6 . 8 7 . 3 8 . 0 9 . 1 9 . 5 9 . 9 1 0 . 6 1 1 . 4 1 1 . 5 1 1 . 9 1 3 . 0

第二步，繪制樹高曲線圖。選擇年齡與樹高觀測(cè)值數(shù)據(jù)區(qū)域，點(diǎn)擊插入“散點(diǎn)圖”，然后添加趨勢(shì)線。趨勢(shì)線格式選擇“多項(xiàng)式”，并勾選“顯示公式”和“顯示R平方值”，見下圖1。圖2為樹高曲線圖繪制結(jié)果。

圖2 中，公式y(tǒng)=-0.008 2x2+0.987 4x-2.467為Excel自動(dòng)擬合的樹高與年齡的回歸方程，R2值是回歸方程擬合優(yōu)度，取值范圍是[0，1]，越接近1擬合程度越好，反之越差。

第三步，分別取三個(gè)年齡觀測(cè)值t1、t2、t3，將它們代入回歸方程中，求出三個(gè)相對(duì)應(yīng)的理論高H1、H2、H3。在這里，t1、t2、t3分別取值為10、15、20，即杉木中齡林的下限、中限、上限三個(gè)年齡值。理論高求算結(jié)果為6.6、10.5、14。

第四步，根據(jù)算出的H1、H2、H3，求觀測(cè)比。

第五步，由（5）式可得：

因此，觀測(cè)比λ可改寫為：

已知t1＝10、t2＝15、t3＝20，代入（7）式如下：

已知λ＝0.617，令e-10k=x，這里x＞0，且x＜1，則：

第六步，在Excel中規(guī)劃求解（15）式中x值。方法如下：

首先，裝載規(guī)劃求解加載宏。依次點(diǎn)擊“文件－選項(xiàng)－加載宏－管理Excel加載項(xiàng)轉(zhuǎn)到……”，將分析工具庫和規(guī)劃求解加載項(xiàng)一并勾選，見圖3，點(diǎn)擊“確定”；然后，點(diǎn)擊“數(shù)據(jù)”菜單“規(guī)劃求解”按鈕，彈出對(duì)話框進(jìn)行設(shè)置，見圖4。

設(shè)置目標(biāo)單元格是要顯示（15）式中左邊公式內(nèi)容所在的單元格，目標(biāo)值選定（15）式中等號(hào)右邊值，即1。在這里我們已經(jīng)于規(guī)劃求解前在$T$8單元格中預(yù)先用函數(shù)寫下了公式“=POWER（1+$T $9，0.617）*（1-$T$9）+POWER（$T$9，3/2）”，式中$T $9是（12）式中變量x引用的單元格絕對(duì)位置，也是要規(guī)劃求解的值，并添加約束規(guī)則2項(xiàng)：$T$9＜=1，$T$9＞=0。勾選“使無約束變量為非負(fù)數(shù)”，選擇求解方法為“非線性GRG”。接著點(diǎn)擊“選項(xiàng)”按鈕，見圖5、圖6。在所有方法選項(xiàng)卡中，如果想要看到迭代計(jì)算結(jié)果就勾選“顯示迭代結(jié)果”，要直接看到計(jì)算結(jié)果就取消勾選，其他可不用填寫；在非線性GRG選項(xiàng)卡中派生選擇“中心”，并勾選“使用多初始點(diǎn)”，確定退出。

回到圖4界面，點(diǎn)擊“求解”按鈕，得到全局最優(yōu)解，見圖7對(duì)話框。

關(guān)閉對(duì)話框，回到Excel工作表中，見圖8。在$T $9單元格中已經(jīng)看到我們要的x變量結(jié)果：x= 0.62224955。如果最優(yōu)解得到的是1或0，重復(fù)求解過程便可。

第七步，將x值代入定義的e-10k=x式中，求得：k=0.047 441 406。

3.3.2 參數(shù)a、c值求定

將年齡觀測(cè)值t1、t2、t3和對(duì)應(yīng)的用回歸方程求解的樹高理論值H1、H2、H3以及參數(shù)k值一并代入式（7）、（8）、（9）中，在Excel中用函數(shù)輕松求解，得到：a=29.918 477 22，c=1.552 514 895。

3.4 杉木中齡林高生長(zhǎng)方程模型擬合結(jié)果

將上述參數(shù)求定值代入理查德方程模型，擬合結(jié)果為：

將年齡觀測(cè)值代入公式驗(yàn)證：當(dāng)t1=10時(shí)，解得H1＝6.6。

當(dāng)t2=15時(shí)，解得H2＝10.5。

當(dāng)t3=20時(shí)，解得H3＝14.0。驗(yàn)證通過。

4 理查德方程應(yīng)用注意事項(xiàng)

一般的直線方程或曲線方程有a、b兩個(gè)待定參數(shù)，而理查德方程 H=a(1-e-kt)c中有三個(gè)待定參數(shù)a、k、c，問題顯然復(fù)雜。我們?cè)谏鲜鰧?shí)例中，設(shè)定k值為模擬生長(zhǎng)規(guī)律時(shí)制定的修正值，單獨(dú)求算出，然后再代入方程求算出其他兩個(gè)參數(shù)，這樣便于求解。在這里，需要特別注意的是：k值模擬樹木生長(zhǎng)規(guī)律求算時(shí)，要盡量滿足樹木生長(zhǎng)函數(shù)曲線要求，縮短擬合區(qū)間范圍，以各齡組分段擬合為最優(yōu)。

S711

A

2095-0152（2015）03-0071-04

2015-03-08

2015-05-07

向南海(1967- )，男，林業(yè)工程師，主要從事林業(yè)調(diào)查規(guī)劃設(shè)計(jì)和森林資源資產(chǎn)評(píng)估咨詢工作。E-mail：376141298@qq.com

楊婷婷）