4個竹種高生長模型的比較

凡美玲，方水元，馮俊嬌，張春霞，梁大洪，時培建

(1.南京林業(yè)大學南方現(xiàn)代林業(yè)協(xié)同創(chuàng)新中心，竹類研究所，江蘇 南京 210037；2.貴州省桐梓縣林業(yè)局，貴州 桐梓 563200)

竹子屬于禾本科竹亞科植物，竹林為我國南方重要的森林資源之一，其生長不同于木本植物，生長速度快、種群密度高使得竹子得以在我國南方廣泛分布[1]。用合理的模型對植物生長數(shù)據(jù)進行擬合，能夠模擬植物生長的動態(tài)變化，這將有益于探索和比較生物個體的生長模式[2]。上世紀以來，相繼許多用于描述溫度變化對外溫生物發(fā)育速率的影響非線性模型被構建，許多模型被運用于擬合動植物的動態(tài)生長。澳大利亞學者Ratkowsky等[3]提出了平方根模型，此模型最初被用于描述細菌的種群增長速率，后被擴展到對節(jié)肢動物溫發(fā)育速率的研究中。Shi等發(fā)現(xiàn)去平方根模型和除平方根的Ratkowsky模型對昆蟲的溫發(fā)育數(shù)據(jù)的擬合效果較好[4]，該模型是否能夠用于植物生長發(fā)育有待研究。Yin等[5]首次提出beta模型用于描述溫度對農(nóng)作物的影響，Yin等[6]通過對beta模型的積分提出了一個新模型(BSE)，該模型能預測植物達到最大生物量的時間，為農(nóng)作物的生長培育提供理論依據(jù)。Shi等[7]用BSE模型與其它4個模型對6種農(nóng)作物的干重進行擬合對比，發(fā)現(xiàn)BSE模型擬合效果優(yōu)于其它4個模型，該模型能適應于描述多種農(nóng)作物的生長。Shi等用BSE模型[2]和LRF模型[8]對6種作物的干重、4個竹種高生長、2種樹木原始年輪數(shù)據(jù)和16個動物的鮮重數(shù)據(jù)進行擬合比較，擬合隨著時間推移各個動植物生物量(或者體型大小)指標變化的數(shù)據(jù)，發(fā)現(xiàn)2種增長方程與實際上動植物的觀測數(shù)據(jù)吻合，且可以通過這2個模型預測動植物達到最大生物量和最高高度的時間，揭示了這2種模型功效的普適性[2]。

唐紅等[9]2014年對湖南省植物園的圣音竹(Phyllostachysedulisf.tubaeformis)高生長研究發(fā)現(xiàn)：用Logistic方程能很好的擬合幼竹平均高生長數(shù)據(jù)，晝夜生長量可達16.48 cm。周芳純[1]對毛竹高生長與出土時間關系研究發(fā)現(xiàn)：兩者關系呈“S”形曲線，可用Logistic方程擬合。曾慶南等[10]對江西林業(yè)科學院竹種園花毛竹(Ph.edulis‘Tao Kiang’)的高生長節(jié)律研究表明Logistic方程同樣能夠較好地擬合竹高生長過程。竹類植物在幼竹高生長時期的生長速度迅速，不同竹種每天的高生長量不同，通過對4個幼竹高生長速度的比較，以期為4個竹種的栽培管理提供理論依據(jù)。以往用Logistic方程對竹子的高生長擬合較多，本研究旨在比較不同增長方程對竹子高生長數(shù)據(jù)的擬合效果的優(yōu)劣，用模型擬合生物個體間不同的生長模式，并對不同模型的適用性和擬合效果進行比較，探討適用于擬合竹類植物高生長的最優(yōu)方程。

1 材料與方法

1.1 研究區(qū)概況

研究地點位于南京林業(yè)大學南京本部校園內(nèi)(32°04′N，118°48′E)，海拔約25 m，屬于亞熱帶季風氣候，研究地年平均氣溫為(15.6±0.9)℃，年平均降雨量為(1 045±242.7)mm左右，該地季節(jié)分明，每年6-7月的梅雨季節(jié)約22 d(1986-2016年氣候數(shù)據(jù)來源于中國氣象數(shù)據(jù)網(wǎng))。

1.2 研究材料

研究試驗材料為4個竹種的竹筍-幼竹，從竹筍露出地面1～2 cm為基準，開始對竹筍高度進行跟蹤調(diào)查，于2016年春季選擇4個竹種(每個竹種4株)進行跟蹤調(diào)查16株幼竹高生長。竹種分別為紅殼竹(Ph.iridescens)、美竹(Ph.mannii)、唐竹(Sinobambusatootsik)、斑苦竹(Pleioblastusmaculatus)，發(fā)筍期均為春季。

1.3 研究方法

研究調(diào)查日期從2016年4月4日開始，直至6月25日結束調(diào)查。調(diào)查中，選用卷尺、可伸縮魚竿測量幼竹高度(每次測定幼竹全高，即為每個幼竹從地面到幼竹的最高處)，調(diào)查初期為每天調(diào)查1次，后期直至幼竹基本完成高生長后，調(diào)查周期為2-3 d，調(diào)查時間均為中午12點，直至幼竹高度值穩(wěn)定，結束調(diào)查。

研究共選擇3個非線性模型對4個竹種高生長數(shù)據(jù)進行擬合，每個竹種均選擇4個個體，每個模型都對16個個體的高生長數(shù)據(jù)進行擬合，3個模型的參數(shù)均為3個，比較3個模型對16個個體高生長數(shù)據(jù)的擬合優(yōu)度及模型參數(shù)的線性近似表現(xiàn)，3個模型如下所示。

1.3.1von Bertalanffy模型 von Bertalanffy模型[11]首次由von Bertalanffy提出，2001年West等[12]提出了一個新的個體增長模型(OGM)。2014年Shi等將OGM模型用于擬合3種淡水魚身體長度的變化，指出不同的物種應該有不同的新陳代謝率指數(shù)取值[13]。研究使用的von Bertalanffy模型，其具體形式如下：

h=H·(1-exp(-3/4·K(t-t0)))4/3

此處，h表示在時間t處生物體一維幾何參數(shù)的大小(長度或者高度)。k為常數(shù)；t0亦是一個常量，表示是當竹子高度為0時的起始時間；參數(shù)H表示漸近線趨于最大值的一維幾何參數(shù)(即最大長度或者最大高度)。

1.3.2BSM模型 Yin等[5]提出一個beta模型，最初是為了描述溫度對玉米發(fā)育速率的影響。Yin等[6]在beta模型的基礎上建立了1個3個參數(shù)的個體生長模型BSE，該模型生物量為0時對應的時間被定義為0。該模型具體形式如下：

h=copt·t·(2tmax-topt-t)/(2tmax-topt)(t/topt)topt/(tmax-topt)

此處，h表示在時間t處生物體一維幾何參數(shù)的大小(長度或者高度)。模型共有3個待擬合參數(shù)，tmax表示幼竹生長到最高高度的時間，topt表示最適幼竹生長的時間(生長速率最大時)，copt表示幼竹最大生長速率的值。

1.3.3Logistic模型 Logistic模型描述動物和植物的生物量和身體長度數(shù)據(jù)運用最廣泛的模型，曲線的對稱性為該模型最典型的特點。該模型具體形式如下：

h=H/(1+exp(b-rt))

此處，共有3個參數(shù)，h表示隨時間t變化的幼竹高度，H代表待擬合的幼竹最高高度，r代表幼竹的瞬時增長率，b代表一個常數(shù)。

1.4 參數(shù)擬合方法及模型評估的指標

研究使用R語言自編程序?qū)τ字窀呱L進行擬合，使用進化遺傳算法[14]來執(zhí)行對模型參數(shù)的估計，該方法是針對于解決求一個函數(shù)最大值的優(yōu)化問題。研究使用決定系數(shù)(R2)和曲率評判模型的優(yōu)劣，R2越大擬合優(yōu)度越好[15]。在非線性回歸中曲率分為均方根內(nèi)稟曲率和均方根參數(shù)效果曲率(下文將省略“均方根”)，其中內(nèi)稟曲率是受到模型本身制約的無法通過對模型參數(shù)的的再參數(shù)化來減小，而參數(shù)效果曲率則可以通過對模型參數(shù)的再參數(shù)化來實現(xiàn)變化，參數(shù)效果曲率相對于臨界曲率越小，表明模型參數(shù)近似線性的表現(xiàn)越好，模型越合理[16-17]。

2 結果與分析

2.1 3個模型的擬合參數(shù)比較

研究共用3種模型分別對4個竹種不同個體的幼竹高生長數(shù)據(jù)進行擬合，計算結果(R2)如表1所示，從表1中的擬合結果可以看出：3個模型分別對16株幼竹個體高生長數(shù)據(jù)擬合效果都很好，R2最小也可以達到0.981 66，為von Bertalanffy模型對紅殼竹4號個體的擬合結果，而擬合效果最好的出現(xiàn)在BSM模型對斑苦竹1號個體的擬合效果，R2最大值可達0.999 28。對擬合結果R2分析發(fā)現(xiàn)，Logstic模型擬合后的R2的計算值有6個數(shù)值達0.998 00，BSM模型有3個數(shù)值達0.998 00，但von Bertalanffy模型擬合后R2最大為0.996 52，沒有數(shù)值達0.998 00。不同模型對不同竹種的擬合結果也不同，Logstic模型和BSM模型對斑苦竹擬合效果最佳，但von Bertalanffy模型對美竹(Ph.mannii)擬合效果較差，該模型2014年被應用于描述淡魚體長變化，研究將其應用于植物，效果欠缺，有待改善模型結構及參數(shù)使其適應性更廣泛。綜上所述，Logistic模型對16株幼竹個體高生長數(shù)據(jù)擬合效果最優(yōu)，BSM模型擬合效果較好，von Bertalanffy模型的擬合結果較差， Logistic模型對于諸多竹種幼竹高生長的擬合效果最好，可以被廣泛運用于動植物個體增長研究。

表1 3個模型對4個竹種16株個體高生長數(shù)據(jù)集的R2擬合結果Tab.1 Coefficient of determination (R2) obtained from using the three models to fit the height growth data sets of sixteen individual shoots of four bamboo species

說明：編號Ⅰ、Ⅱ、Ⅲ、Ⅳ分別代表紅殼竹、美竹、唐竹和斑苦竹，編號1、2、3、4代表每個竹種的不同個體。
The number of Ⅰ、Ⅱ、Ⅲ、Ⅳ representPh.iridescens、Ph.Mannii、Sinobambusatootsik、Pleioblastusmaculatus，and the number of 1、2、3、4 represent different individuals of each bamboo species.

2.2 3個模型參數(shù)曲率的比較

表2為3個模型對16個幼竹高生長數(shù)據(jù)擬合后的曲率對比，內(nèi)稟曲率和參數(shù)效果曲率用于評估模型的非線性表現(xiàn)和參數(shù)結構的合理性。此處，曲率指標共有3個：內(nèi)稟曲率(RMS IC)、參數(shù)效果曲率(RMS PEC)和臨界曲率(CC)，當模型擬合后的內(nèi)稟曲率和參數(shù)效果曲率都小于臨界曲率，表明該模型近似線性的表現(xiàn)越強，參數(shù)結構更為合理。由表2的擬合結果可以看出：3個模型Logstic模型擬合結果中，內(nèi)稟曲率和參數(shù)效果曲率均小于臨界曲率，表明該模型參數(shù)的線性近似表現(xiàn)較好，模型參數(shù)設置合理；對于BSM模型，內(nèi)稟曲率和參數(shù)效果曲率均小于臨界曲率，與Logstic模型擬合結果相比，BSE模型擬合得到的內(nèi)稟曲率和參數(shù)效果曲率均大于Logstic模型，表明Logstic模型的擬合結果優(yōu)于BSM模型；von Bertalanffy模型擬合結果中，內(nèi)稟曲率均小于臨界曲率，但有4個參數(shù)效果曲率大于臨界曲率，參數(shù)效果曲率均>BSM模型>Logstic模型。綜上所述，Logistic模型對16株幼竹個體高生長數(shù)據(jù)擬合優(yōu)度最高，參數(shù)設置更為合理，其次為BSM模型，von Bertalanffy模型的擬合優(yōu)度最低[17]。

2.3 3個模型的擬合曲線比較

圖1為BSM模型、von Bertalanffy模型和Logistic模型3種不同模型對斑苦竹中1號幼竹高生長數(shù)據(jù)的擬合結果(3個模型對該個體擬合效果最好)。研究中的斑苦竹種植于南京林業(yè)大學竹種園內(nèi)，1號斑苦竹于2016年4月29日破土而出，直至2016年6月25日完成其高生長，最高達12 m。開始生長速度較慢，最快速度出現(xiàn)于6月6日，此時增長速率達到72 cm·d-1，隨后生長速度減慢，其生長曲線呈典型的“慢-快-慢”趨勢。

表2 3個模型對16個個體高生長數(shù)據(jù)集擬合后得到的曲率比較Tab.2 Curvatures obtained from using three models to fit the height growth data sets of sixteen bamboo shoots

研究共選擇了3個非線性模型對同一竹種不同個體高生長數(shù)據(jù)進行擬合，3個模型都設置了3個待擬合參數(shù)(BSM模型中tmin設定為0)。von Bertalanffy模型是一個通用型模型，可以用于諸多生物的生長，研究該模型新陳代謝率指數(shù)選用的為3/4，但von Bertalanffy模型對幼竹高度的擬合曲線與實際數(shù)據(jù)有很大出入(圖1c和圖1d)，該擬合曲線在前10 d左右零生長，隨后生長速度驟增，但在曲線末端未見曲線平緩，有繼續(xù)增加的趨勢，這與實際觀測數(shù)據(jù)不符，幼竹在完成高生長后，高度會趨于一個穩(wěn)定值，其下方對幼竹生長速率的擬合曲線同樣可以看出，該幼竹前期生長速率為0，隨后在10 d左右達到最大生長速率，達到最大生長速率共用24 d，其后幼竹繼續(xù)生長，生長速率線始終沒有趨于0，且整條生長速率曲線不對稱，擬合后決定系數(shù)計算結果為0.988 45，該模型對其幼竹高生長擬合效果差。BSM模型最大的特點是能夠包容對稱和不對稱曲線，能準確估計植物到達最大生物量的時間。圖1a和圖1b展現(xiàn)了BSM模型對斑苦竹中1號幼竹高生長數(shù)據(jù)及生長速率的擬合結果，該擬合結果表現(xiàn)為前期生長緩慢，隨后進入高生長盛期，曲線最后趨于平緩，幼竹高度達到最大高度后稍有下降(幼竹開始伸枝展葉)，最后穩(wěn)定于一個值，擬合曲線與幼竹實際高度重合率高。從該模型對其生長速率擬合曲線可以發(fā)現(xiàn)：高生長階段達到最大生長速率的時間為距離出筍35 d，52 d后結束高生長，該生長速率曲線左側寬于右側表明，竹子前期到達最大生長速率的時間略多于后期生長的時間，擬合后決定系數(shù)計算結果可達0.999 28，該模型對其擬合效果較優(yōu)。Logistic模型被廣泛應用于描述動植物的生物量，從圖1e和圖1f可以看出，幼竹高生長呈“慢-快-慢”趨勢，擬合曲線與幼竹高生長數(shù)據(jù)高度重合，與BSM模型不同的是，該模型對生長速率的擬合曲線更平緩，且曲線兩側對稱，幼竹最大生長速率時間為距出筍33 d，與實際數(shù)據(jù)相符，擬合后決定系數(shù)計算結果能夠達到0.999 01，綜上所述，Logistic模型擬合結果優(yōu)于BSM模型優(yōu)于von Bertalanffy模型，在本研究中，3個模型對不同幼竹個體高生長數(shù)據(jù)的擬合結果表明，Logistic模型為最優(yōu)模型。

圖1 用3種模型對斑苦竹(Pleioblastus maculatus)1號個體高生長數(shù)據(jù)的擬合結果Fig.1 Fitted results of the height growth for No.1 bamboo of Pleioblastus maculatus by using three models說明：a、b為BSM模型；c、d為von Bertalanffy模型；e、f為Logistic模型。a、b are the BSM model; c、d are the von Bertalanffy model; and e、f are Logistic model

3 討論

竹子每年的發(fā)筍時間主要受當年的溫度和降雨量等因素影響，按發(fā)筍期發(fā)筍數(shù)量可分為初期、中期和末期，末期筍常因營養(yǎng)不足或蟲害退筍[18-20]。實驗選擇是4個竹種發(fā)筍期差異主要因為竹種差異造成，選擇發(fā)筍前期到中期的筍為研究對象，調(diào)查發(fā)現(xiàn)紅殼竹與美竹發(fā)筍日期與高生長期相近，于4月初開始發(fā)筍，6月伊始結束高生長，唐竹和斑苦竹發(fā)筍期略晚，分別為4月中期和末期，高生長期約為40-55 d。

研究使用的von Bertalanffy模型，參數(shù)=3/4，但總體擬合效果差，同時也嘗試參數(shù)=2/3，擬合優(yōu)度低于前者，2014年Shi等[2]用該模型和OGM模型同時對3種淡魚體長數(shù)據(jù)進行擬合，他們發(fā)現(xiàn)當OGM模型中的參數(shù)α越接近1，模型的擬合效果越好，與本研究的結果一致。本研究中，使用的BSM模型中限制δ=1、tmin為0，但并不會影響該模型描述數(shù)據(jù)的有效性，結果表明BSM模型對4種幼竹高生長數(shù)據(jù)擬合效果較好，對幼竹生長到最大生物量的時間估計與實際相符，Shi等[2]用2種生長模型分別對動植物生長進行擬合，作者指出：如果可以預先確定動植物增長的起始時間，那對于BSM模型的參數(shù)δ可以不限制為1，可以放寬δ作為一個參數(shù)，這樣將幫助對比動植物種間增長模式的差異。對于常用的生長模型Logistic模型，擬合結果同前人研究一致，對本研究中的4種竹種個體高生長數(shù)據(jù)擬合效果好，Shi等13用4種模型同時對相同數(shù)據(jù)擬合比較發(fā)現(xiàn)，Logistic模型對于描述動植物生物量和身體長度都是通用的，且擬合優(yōu)度都優(yōu)于其它模型，是最適合描述動植物的生長模型。Richards[21]提出一個廣義的Logstic模型，即Richards模型，該模型更加靈活，不再限制生長速率曲線的對稱性，且該模型對于在其最小值到最大值之間的任何值都能有其最大斜率，但該模型也存在一些爭議，Ratkowsky指出該模型對系數(shù)估計值不準確[21]，Birch[23]也指出該模型不能反映一些蔬菜和牧草的生長，如黑麥草(LoliumperenneL.)在低密度時就能達到其最大凈產(chǎn)量，與實際不符。研究表明：許多非線性模型不僅適用于溫度對動物生長速度和內(nèi)稟增長率的影響，對于植物的生長也是適用的。

最初一些非線性模型被提出用于描述溫度對昆蟲發(fā)育速率的影響，如表現(xiàn)型模型[24]、Taylor模型[25]和Ratkowsky模型[4]等等，這些模型可以幫助預測一些重要參數(shù)，對農(nóng)業(yè)生產(chǎn)及科學研究有著重要意義，但因為模型結構設置及參數(shù)選擇不同，對不同數(shù)據(jù)的擬合優(yōu)度不同，選擇一個通用性的最優(yōu)模型具有重要意義。評價模型擬合優(yōu)度的好壞可以從幾個方面評估，模型結構的復雜程度、模型適用性以及模型預測曲線與實際是否相符等[2]。隨后，王如松等將數(shù)學模型用于描述溫度變化對昆蟲內(nèi)稟增長率的影響，擬合結果顯示相關系數(shù)高，擬合結果較優(yōu)，指出可以利用數(shù)學模型研究一些與溫度相關的現(xiàn)象[26]。2016年Shi等提出BSM的通用表達式(即生長發(fā)育初始時間不為0時的BSM模型)能夠擬合植物的動態(tài)生長，并能幫助預測植物生長到最大生長量的時間，這為森林資源的合理利用提供了參考價值[7]。Shi等[13]將1種個體增長模型運用于作物的生長，結果表明Logistic模型能夠更好地擬合動植物生長的數(shù)據(jù)，這與本研究結果一致。本研究中，BSM模型與Logistic模型對幼竹生長到最大生物量時間的預測相近，與實際數(shù)據(jù)吻合，而von Bertalanffy模型預測幼竹最大生物量的時間提前10 d，正如圖1d所展示的生長速率曲線左偏。同樣，Logistic模型擬合的生長曲線同前人研究一樣，曲線兩側關于對稱軸對稱，對數(shù)據(jù)缺乏靈活性[27]，但總體Logistic模型適合描述植物生長過程中生物量的變化趨勢。

研究使用3個非線性模型對不同竹種幼竹高生長數(shù)據(jù)進行擬合，將以往描述溫度對昆蟲發(fā)育速率影響的模型用于擬合植物的生長，將自變量溫度調(diào)換為時間，結果表明，生長模型也能很好地描述植物的生長[2]。調(diào)查發(fā)現(xiàn)，不同竹種前期與盛期筍生長為成竹的概率較大，末期筍易發(fā)生退筍現(xiàn)象，所以在選擇留筍養(yǎng)竹時，應選擇初期和盛期萌發(fā)的筍為宜。利用3個模型對不同幼竹高生長數(shù)據(jù)進行擬合，4個竹種高生長均符合“S”型生長，符合“慢-快-慢”是生長規(guī)律，這為4個竹種的撫育管理提供了重要的理論依據(jù)。本研究發(fā)現(xiàn)，Logistic模型對動植物生長的擬合效果最優(yōu)，BSM模型與von Bertalanffy模型有待改進，這些生長模型的適用性有待進一步探索，生物有機體生長模式的比較也將成為一個研究熱點。