森林立地是指立木生长发育的环境条件。林木生长除其本身的遗传特性和生物学特征影响外，与立地有很大关系。立地质量是影响林地生产能力诸多生境因子（包括生物、气候、土壤等）的综合评价指标。因此，对森林立地的评价是森林经营必不可少的前提条件，是森林经营过程中营林者实施各种经营措施的重要依据^[1]。

在20世纪初，部分美国学者提出以树高生长作为衡量立地质量的指标，这种方法被沿用至今。目前，地位指数法仍然是立地评价的主要方法^[2-3]。Bailey等^[4]于1974年首次提出了通过指定生长方程中一个与形状有关的参数为自由参数，并以此建立差分模型的代数差分法（Algebraic Difference Approach，ADA），并建立了新西兰松木人工林的异形地位指数模型。这种方法不需要事先确定地位指数，基准年龄可变，一些学者通过这种方法建立地位指数模型^[3,5-9]。在此基础上，McDill等^[3]于1992年提出了具有一个速率参数和一个最大树高参数的尺度相容性生长函数，并提出了一种通过指定一部分全局参数（自然条件，如气候等参数，所有林分参数相同）与一部分林分参数（地形地势、土壤条件等参数，林分不同林分参数不同）来拟合模型的方法。随后Amaro等^[6]采用代数差分法（ADA）模拟了桉树（Eucalyptus sp.）优势高生长，通过比较ADA法推导的模型与McDill等^[3]提出的差分方程，为ADA法选择模型形式、全局参数与林分参数的组合提供了依据。

许多研究表明，不同立地条件下的优势木树高生长曲线的形状并非同形，即不同立地条件下优势木的树高生长速度、曲线形状及渐近最大值（潜在树高）均不同^[3-11]。合理的地位指数模型应同时满足多条水平渐近线（树高最大值不同）和多形性（拐点不同）的性质，而使用ADA法所推导的差分地位指数模型，由于只设置一个自由参数，当设置方程中树高最大值参数为自由参数时，曲线满足多条水平渐近线却不能满足多形性；当设置形状参数为自由参数时，曲线满足多形性却不能有多条水平渐近线。也就是说，ADA法只能满足其中的一个条件^[10]。对此，Cieszewski等^[11]于2000年提出了广义代数差分法（GADA）。使用GADA法允许在推导时设置多个自由参数，可以同时设置树高最大值参数与形状参数，因此推导出的模型能同时描述曲线的最大值和形状随立地条件变化的情况，能同时满足多条水平渐近线和多形性两个性质，即异形地位指数曲线^[12]。因此，本研究以长白落叶松（ Larix olgensis ）解析木数据为基础，使用GADA法建立差分地位指数模型，通过拟合与检验筛选出最适用于我国黑龙江省长白落叶松人工林的差分地位指数模型，并将其与通过ADA法建立的模型进行比较分析，进一步验证GADA法所推导异形地位指数模型的优良性能。

长白落叶松又名黄花落叶松，是东北林区的主要造林树种之一，主要分布在完达山、张广才岭和长白山海拔500 ~ 1 800 m湿润山坡及沼泽地区，生长快，材质优良，用途广泛。据第八次全国森林资源清查结果显示：东北地区落叶松人工林面积为205.66万hm²，蓄积为13 364.15万m³，分别占东北三省人工林总面积和总蓄积的35.67%和36.81%^[13]。但是，人工林林分质量差，生产力低，生态功能弱。如何根据长白落叶松人工林树高生长规律，精准评价立地质量，提升人工林质量和木材品质，成为亟待解决的问题。本研究基于优势木的树高生长数据，采用GADA法构建异形地位指数曲线，为长白落叶松人工林立地质量的精准评价和林分质量的精准提升提供依据。

1.   数据与方法

1.1   数据来源

1.1.1   研究区概况

本研究所用数据的调查地点均属黑龙江省长白山系，包括完达山山脉的佳木斯市孟家岗林场、勃利县林业局，以及张广才岭的林口林业局、尚志国有林场管理局、东京城林业局、江山娇林场、宁安市林业局、五常市林业局。各采样点的位置如图1所示。

图  1  长白落叶松人工林数据采样点位置分布

Figure  1.  Location distribution of sampling points of Larix olgensis plantations

下载: 全尺寸图片 幻灯片

完达山位于黑龙江省东部，属长白山山脉北延，是黑龙江省东部主要山地之一。主脉呈西南—东北走向，北抵挠力河，西北与三江平原相接，东南与穆棱—兴凯平原交界，西南接那丹哈达岭，东达乌苏里江，绵延400 km，面积1.2万km²，海拔500 ~ 800 m，主峰神顶山在虎林县与宝清县交界处，海拔831 m。主脊线以南属虎林市，以北在宝清县境内，东北端属饶河县管辖。地形以低山丘陵为主，坡度较为平缓，大多在10° ~ 20°之间，地势东北高，西南低，土壤以暗棕壤为主。气候属中温带大陆性季风气候，年均气温2.5 ℃左右，年均降水量550 mm左右^[14]。

张广才岭是长白山的支脉，位于中国东北地区东部山地北段的中轴部位，大部分在黑龙江省境内，一部分向南伸入到吉林省的敦化市北部，分东西两支，蛟河盆地以西为西老爷岭，蛟河盆地以东为威虎岭。张广才岭是中国东北地区主要山脉之一，主脊南起吉林省敦化市，北接小兴安岭南麓，平均海拔800 m左右，主峰老秃顶子高达1 686.9 m。其中部主要是中山，边缘地形以低山丘陵为主，坡度大多在6° ~ 15°之间，地势由南向北逐渐升高，土壤为暗棕壤。气候属亚寒带大陆性季风气候，年平均气温3.5 ℃左右，年平均降水量550 ~ 750 mm^[15]。

1.1.2   数据获取

本研究所使用的数据为1994—2017年间所调查60株长白落叶松优势木的解析木数据。数据来自1994—1995年黑龙江省孟家岗林场和东京城林业局江山娇林场的9株解析木，2007、2008和2015年孟家岗林场27株解析木，2009年黑龙江重点公益林碳汇项目8株解析木，2016年林口林业局和东京城林业局6株解析木，2017年黑龙江省尚志林场10株解析木。各林场或林业局立地相关因子见表1。

表  1  采样点所在各林场及林业局的立地相关因子

Table  1.  Site quality-related factors of each forest farm and forestry bureau

单位名称 Names of organizations	地形 Terrain	海拔 Altitude/m	坡度 Gradient/(°)	坡向 Slope aspect	坡位 Slope position
孟家岗林场 Mengjiagang Forest Farm	丘陵 Hill	140 ~ 230	0 ~ 15	西南 Southwest	中 Middle
江山娇林场 Jiangshanjiao Forest Farm	山坡 Slope	420 ~ 460	0 ~ 5	北 North	中下 Lower-middle
林口林业局 Linkou Forestry Bureau	山坡 Slope	320 ~ 340	0 ~ 15	西南 Southwest	中下 Lower-middle
东京城林业局 Dongjingcheng Forestry Bureau	山坡 Slope	330 ~ 340	0 ~ 15	西 West	中 Middle
尚志国有林场管理局 Shangzhi State-owned Forest Farm Management Bureau	山坡 Slope	230 ~ 480	0 ~ 15	南 South	中上 Upper-middle
勃利县林业局 Bolixian Forestry Bureau	山坡 Slope	380 ~ 400	5 ~ 10	西北 Northwest	上 Up
五常市林业局 Wuchangshi Forestry Bureau	平地 Flat	270 ~ 280	0 ~ 5		平 Flat
宁安市林业局 Ning’anshi Forestry Bureau	山坡 Slope	350 ~ 360	5 ~ 10	北 North	下 Down

下载: 导出CSV 

| 显示表格

在采样点获取数据时，先进行标准地的设置工作。标准地的形状为矩形，大小为20 m × 30 m或20 m × 20 m。对标准地进行设置完成之后，在记录纸上画出标准地的缩略图，同时记录标准地的优势树种、林分起源、坡度、坡向、坡位、地形及海拔等林分因子。标准地内每株胸径大于5 cm树木都要进行实测，进行每木检尺，并通过测得的数据挑选出优势木，将其伐倒进行树干解析。

全部60株长白落叶松人工林优势木树高生长过程如图2所示，建模样本与检验样本的样木调查因子基本统计量见表2。

图  2  长白落叶松人工林优势木树高生长过程曲线

Figure  2.  Growing process curves of tree height of dominant trees for Larix olgensis plantation

下载: 全尺寸图片 幻灯片

表  2  长白落叶松人工林优势木建模数据和检验数据基本统计量

Table  2.  Statistics of fitting data and validation data for dominant trees in Larix olgensis plantations

变量 Variable	建模数据 Fitting data (n = 45)					检验数据 Validation data (n = 15)
	最小值 Min.	最大值 Max.	均值 Mean	标准差 Std.		最小值 Min.	最大值 Max.	均值 Mean	标准差 Std.
年龄/a Age/year	16	53	34.69	10.96		17	47	30.47	10.68
胸径 DBH/cm	11.3	35.7	24.14	4.28		13.8	29.7	21.92	4.96
树高 Tree height/m	11.2	26.8	20.16	3.40		11.0	23.3	18.40	3.53
注：n为样本数。Note: n is the number of samples.

下载: 导出CSV 

| 显示表格

1.2   研究方法

1.2.1   基准模型的选择

在林木和林分生长模型的研究工作中，广泛应用修正Weibull（1978）、Korf（1939）、Richards（1959） 3个理论生长方程^[16-20]。因此，本研究选择这3个方程作为GADA法推导差分地位指数模型的基础方程。

式中：a为树高极限值参数（潜在最大树高），代表一定的立地条件下树高生长的最大值；b为生长率参数（或尺度参数），它影响曲线拐点的位置，在Korf方程中，它是尺度参数，而在修正Weibull 和Richards方程中是与树木增长率有关的参数；c为形状参数，它影响曲线的形状和拐点位置，在Korf方程中，它与植物生长过程中的衰减因子有关，在修正Weibull 和Richards方程中，它与树木生长的同化作用有关^[21]。

1.2.2   差分方程的推导

在Bailey等^[4]于1974提出的代数差分法的基础上，Cieszewski等^[11]于2000年提出了广义代数差分法（GADA）。这种方法允许将生长方程中的一个或多个参数作为SDP（site-dependent parameter，也称“与立地相关参数”或“自由参数”）。这样，推导出来的模型可以同时满足树高生长的多条水平渐近线与多形性的条件，这些条件正是异形地位指数模型应具备的条件^[22]。

通过GADA法推导差分地位指数模型通常包括以下步骤：

（1）选择一个生长方程为基础方程，并根据理论或经验指定其中的两个或多个参数为自由参数。

（2）提出一个与立地质量有关的量X₀，并假设自由参数与它之间有简单的函数关系，在进行这一步时，往往只关注这种函数关系而忽视X₀的具体数值，因为由它所反映出立地条件的不同会通过自由参数体现在方程中。

（3）将上述函数关系代入基础方程，并从中解出X₀的表达式。

（4）将地位指数和基准年龄分别代入表达式中的树高和年龄，并将表达式代入基础方程的X₀，即得到了通过GADA法推导的差分地位指数方程^[23]。

以方程（3）为例，由于参数b代表的树木生长速率是植物固有的属性，与立地关系不紧密，因此将方程（3）中的参数a和c指定为自由参数，并假设参数a与X₀之间有指数关系，参数c与X₀之间有线性关系，即$a = {{\rm{e}}^{{X_0}}}$，c = c₁ + c₂X₀。方程（3）转化为：

从中解出

令方程（5）中h = h₁，t = t₁，其中t₁为指定年龄，h₁为指定年龄下的优势木高。当t₁为基准年龄时，h₁即为地位指数。令方程（4）中h = h₂，t = t₂，其中t₂为预测年龄，h₂为预测年龄下的优势木高。将方程（5）代入方程（4），即得以方程（3）为基础方程，使用GADA法推导的一个差分地位指数模型。它的表达式为：

其中

类似地，在3个生长模型的基础上，使用GADA法推导了6个差分方程（表3）。其中，模型E1与E2由修正Weibull方程导出；模型E3与E4由Korf方程导出；模型E5与E6由Richards方程导出。

表  3  基础方程与差分地位指数方程

Table  3.  Difference site index equations and base equations

基础方程 Basic equation	自由参数 Free parameter	X的解或方程 Solution or equation of X	差分方程 Difference equation	编号 No.
修正Weibull方程 Modified Weibull： $\scriptstyle h = a\left( {1 - {{\rm{e}}^{ - b{t^c}}}} \right)$	$\begin{array}{l}\scriptstyle a = {{\rm{e}}^{{X_0}}}\\\scriptstyle c = {c_1} + {c_2}{X_0}\end{array}$	$\scriptstyle{X_0} = \ln {h_1} - \ln \left( {1 - {{\rm{e}}^{ - b{t_1}^{{c_1} + {c_2}{X_0}}}}} \right)$	$\scriptstyle {h_2} = {{\rm{e}}^{{X_0}}}\left( {1 - {{\rm{e}}^{ - b{t_2}^{{c_1} + {c_2}{X_0}}}}} \right)$	E1*
	$\scriptstyle a = {{\rm{e}}^{{X_0}}}\\c = {c_1} + {c_2}/{X_0}$	$\scriptstyle {X_0} = \ln {h_1} - \ln \left( {1 - {{\rm{e}}^{ - b{t_1}^{{c_1} + {c_2}/{X_0}}}}} \right)$	$\scriptstyle {h_2} = {{\rm{e}}^{{X_0}}}\left( {1 - {{\rm{e}}^{ - b{t_2}^{{c_1} + {c_2}/{X_0}}}}} \right)$	E2*
Korf： $\scriptstyle h = a{{\rm{e}}^{ - b{t^{ - c}}}}$	$\begin{array}{l}\scriptstyle a = {{\rm{e}}^{{X_0}}}\\\scriptstyle c = {c_1} + {c_2}{X_0}\end{array}$	$\scriptstyle {X_0} = \ln {h_1} + b{t^{ - \left( {{c_1} + {c_2}{X_0}} \right)}}$	$\scriptstyle {h_2} = {{\rm{e}}^{{X_0}}}{{\rm{e}}^{ - b{t_2}^{ - \left( {{c_1} + {c_2}{X_0}} \right)}}}$	E3*
	$\begin{array}{l}\scriptstyle a = {{\rm{e}}^{{X_0}}}\\\scriptstyle c = {c_1} + {c_2}/{X_0}\end{array}$	$\scriptstyle {X_0} = \ln {h_1} + b{t^{ - \left( {{c_1} + {c_2}/{X_0}} \right)}}$	$\scriptstyle {h_2} = {{\rm{e}}^{{X_0}}}{{\rm{e}}^{ - b{t_2}^{ - \left( {{c_1} + {c_2}/{X_0}} \right)}}}$	E4*
Richards： $\scriptstyle h = a{\left( {1 - {{\rm{e}}^{ - bt}}} \right)^c}$	$\begin{array}{l}\scriptstyle a = {{\rm{e}}^{{X_0}}}\\\scriptstyle c = {c_1} + {c_2}{X_0}\end{array}$	$\scriptstyle {X_0} = \frac{ \scriptstyle{\ln {h_1} - {c_1}\ln \left( {1 - { {\rm{e} }^{ - b{t_1} } } } \right)} }{ \scriptstyle{1 + {c_2}\ln \left( {1 - { {\rm{e} }^{ - b{t_1} } } } \right)} }$	$\scriptstyle {h_2} = {{\rm{e}}^{{X_0}}}{\left( {1 - {{\rm{e}}^{ - b{t_2}}}} \right)^{\left( {{c_1} + {c_2}{X_0}} \right)}}$	E5
	$\begin{array}{l}\scriptstyle a = {{\rm{e}}^{{X_0}}}\\\scriptstyle c = {c_1} + {c_2}/{X_0}\end{array}$	$\begin{array}{l}\scriptstyle {X_0} = \frac{\scriptstyle{1}}{\scriptstyle{2}}\left[ {F + \sqrt { {F^2} - 4{c_2}\ln \left( {1 - { {\rm{e} }^{ - b{t_1} } } } \right)} } \right]\\\scriptstyle F = \ln {h_1} - {c_1}\ln \left( {1 - { {\rm{e} }^{ - b{t_1} } } } \right)\end{array}$	$\scriptstyle {h_2} = {{\rm{e}}^{{X_0}}}{\left( {1 - {{\rm{e}}^{ - b{t_2}}}} \right)^{\left( {{c_1} + {c_2}/{X_0}} \right)}}$	E6
	$\scriptstyle a = {X_0}$	$\scriptstyle {X_0} = {h_1}{\left( {1 - {{\rm{e}}^{ - b{t_1}}}} \right)^{ - c}}$	$\scriptstyle {h_2} = {h_1}\frac{ { \scriptstyle{ {\left( {1 - { {\rm{e} }^{ - b{t_2} } } } \right)}^c} } }{ { { \scriptstyle{\left( {1 - { {\rm{e} }^{ - b{t_1} } } } \right)}^c} } }$	E0
注：*代表模型E1 ~ E4在推导过程中无法解出X0的显式解，在拟合过程中使用迭代法求出其数值解来代替。Notes: * represents the numerical solutions were used instead, since the explicit solutions of X0 in model E1−E4 cannot be deduced mathematically.

下载: 导出CSV 

| 显示表格

同时，为了比较GADA法与ADA法的差异，以Richards方程为基础，使用ADA法推导出了一个差分模型E0（表3）。

1.2.3   模型的拟合与检验

本研究使用R3.5.1软件来完成模型参数的估计和拟合统计量的计算。

为了全面验证模型的拟合优度和预测能力，本研究采用拟合样本计算决定系数（R²），用均方根误差（RMSE）来评价模型的拟合优度；采用独立样本计算模拟效率（EF），用平均绝对误差（MAE）来评价模型的预测效果。决定系数（R²）和均方根误差（RMSE）用于反映模型对于建模样本的拟合优度，决定系数越接近于1、均方根误差越小，说明模型越好；模拟效率（EF）和平均绝对误差（MAE）用于反映模型对独立样本的预测能力，其中模拟效率越接近于1、平均绝对误差越小，说明模型预测效果越好。各指标计算公式如下：

式中：${h_i}$为优势高实测值，${{\hat h}_i}$为优势高估计值，$\bar h$为优势高实测值的平均值，n为样本个数。

1.2.4   GADA法优势的验证

相比于ADA法，使用GADA法推导的方程能同时具有多条水平渐近线和多形性两项性质。为了进一步验证这一点，本研究在选出最优模型后，通过绘制地位指数曲线簇、比较两组曲线的最大值参数a与形状参数c，求出并比较两组曲线的拐点位置、拐点处斜率（连年生长量最大值）、拐点处函数值（拐点处优势木高），来进一步说明两种方法的差异。若一组曲线满足多条水平渐近线，在不同的地位指数下的曲线会有不同的最大值，即不同的参数a；若一组曲线满足多形性，则不同地位指数下的曲线应具有不同的拐点位置，同时，曲线之间也不应有简单的比例关系。

2.   结果与分析

2.1   模型拟合结果

使用长白落叶松人工林优势木的建模样本（表2）对各个差分地位指数方程（表3）进行拟合，结果见表4。注意到模型E6的参数c₁在5%水平下差异不显著，将参数c₁固定为0重新拟合模型E6，得到模型E7，拟合与检验结果见表4。将模型E6与E7的残差进行方差分析的结果显示两模型差异不显著（P = 0.871），计算两模型的赤池系数（AIC），结果为AIC_E6 = 3 385.388，AIC_E7 = 3 383.414，差异不显著，且E7要略优于E6。可见E7不仅模型形式比E6简单，拟合效果也比E6要好。因此以E7代替E6参与比较。

表  4  长白落叶松人工林差分地位指数模型的参数估计值及拟合优度

Table  4.  Parameter estimates and goodness of fit of the difference site index equations for Larix olgensis plantations

编号 No.	参数 Parameter	估计值 Estimating value	标准误 Std. err.	P值 P value	R2	RMSE	EF/%	MAE
E1	b	0.017 0	0.000 4	< 0.001	0.987 27	0.758 7	97.22	0.888 5
	c1	0.544 6	0.236 1	0.021
	c2	0.215 5	0.069 8	0.002
E2	b	0.017 0	0.000 3	< 0.001	0.987 27	0.759 0	97.22	0.888 7
	c1	1.915 6	0.271 4	< 0.001
	c2	− 2.175 8	0.909 1	0.017
E3	b	6.965 7	0.085 0	< 0.001	0.986 81	0.772 5	97.04	0.908 4
	c1	− 2.100 8	0.018 5	< 0.001
	c2	0.606 0	0.005 8	< 0.001
E4	b	6.957 8	0.084 4	< 0.001	0.986 77	0.775 4	97.04	0.909 9
	c1	1.167 5	0.021 6	< 0.001
	c2	− 2.920 8	0.015 6	< 0.001
E5	b	0.046 8	0.000 8	< 0.001	0.987 43	0.748 2	97.04	0.904 9
	c1	2.663 7	0.304 9	< 0.001
	c2	− 0.378 6	0.089 8	< 0.001
E6	b	0.047 0	0.001 1	< 0.001	0.987 43	0.749 2	97.04	0.904 0
	c1	− 0.000 9	0.005 5	0.871
	c2	4.683 6	0.066 8	< 0.001
E7	b	0.046 8	0.000 8	< 0.001	0.987 41	0.749 3	97.04	0.904 0
	c2	4.675 4	0.044 1	< 0.001
E0	b	0.045 7	0.000 8	< 0.001	0.987 20	0.763 3	97.17	0.891 7
	c	1.377 6	0.015 3	< 0.001

下载: 导出CSV 

| 显示表格

由表4可知：所有模型参数估计值都很小，说明所有的模型参数都具有良好的稳定性；决定系数（R²）均达到0.986以上，均方根误差（RMSE）均在0.78以下，说明所有模型对于建模数据的拟合效果较好；模拟效率（EF）均达到97%以上，平均绝对误差（MAE）在0.88 ~ 0.91之间，说明所有模型都具有良好的预测能力。其中E5和E7决定系数相对较高，分别为0.987 43和0.987 41，说明基于Richards方程推导出的差分地位指数模型能较好地解释不同立地条件下的黑龙江地区长白落叶松人工林优势木的高生长过程；E5、E7均方根误差相对较小，说明基于Richards方程推导的差分地位指数模型能更为精确地拟合优势木的高生长过程；同时，E1、E2的模拟效率和平均绝对误差均相对较小，说明基于修正Weibull方程推导的差分地位指数模型具有较好的预测能力。

综合考虑各项指标，挑选模型E1、E2和E7为较优模型，进行进一步的检验。

2.2   模型检验

模型E1、E2和E7的残差图如图3所示。从3个模型的残差图来看，残差点较为均匀地分布在x轴两侧，大部分的点集中在x轴附近，在预估值小于15时有少量点分布较为分散。这就说明，3个模型对于长白落叶松较低的优势木高拟合残差绝对值均比较大，而随着优势木高预估值的增大，残差值的分布范围逐渐减小，在预估值较大时3个模型预估精度差别很小；但模型E7在预估值较低时的残差点分布相对集中，说明模型E7对于优势木幼龄时期高生长的预估能力更强。

图  3  差分地位指数模型E1、E2及E7的残差图

Figure  3.  Residuals against predicted dominant tree height for the difference site index models E1, E2 and E7

下载: 全尺寸图片 幻灯片

结合样本点，绘制3个模型基准年龄为25年（t_b）、地位指数级距为2 m、地位指数范围从12 m至22 m的地位指数曲线簇（图4）。从图4可以看出，3个模型都满足了多形性以及多条水平渐近线的条件，均能较好地预测优势木15年以上的树高生长情况；但对于15年以下的预测，模型E7可以更明显地区分不同立地条件下优势高的生长过程，这与残差图比较的结果一致。因此，选择E7为本研究所建立的最优模型，它的表达式为：

图  4  差分地位指数模型E1、E2及E7所绘制的地位指数曲线簇

Figure  4.  Site index curves generated with the difference site index models E1, E2 and E7

下载: 全尺寸图片 幻灯片

其中

式中：h₂是待预测年龄t₂时的树高，h₁是已知年龄t₁时的树高。当t₁为基准年龄t_b时（如25年），相应的h₁即为地位指数SI，此时该模型可以估计一定地位指数下长白落叶松人工林优势木的树高生长过程；当t₂为基准年龄t_b时，相应的h₁即为地位指数SI，此时该模型可在已知长白落叶松人工林在某一年龄下的优势木高时估计该林分的地位指数。

2.3   基于GADA法推导的模型与ADA法模型的比较

为检验GADA法与ADA法的差异，将上述最优模型E7与通过ADA法推导的E0进行比较分析。首先根据E7和E0绘制地位指数为12、17、22 m时两模型的曲线簇和拐点连线（图5）。从图5可以看出：两个模型在平均立地条件下优势木的生长过程基本一致；当年龄小于25年时，模型E7在地位指数为最小值处的曲线在模型E0的下方，而地位指数为最大值处的曲线在模型E0的上方；当年龄大于25年时情况相反。从散点图中观察到，在年龄较小时，不同立地条件下长白落叶松人工林优势木高生长之间差距较大，这说明模型E7能更好地区分年龄较小时不同立地条件下优势木的生长过程。而当年龄较大时，由于解析木数据量较少，不能完全从图中得出不同立地条件下优势木的生长情况。

图  5  差分地位指数模型E0和E7的对比

Figure  5.  Comparison of difference site index models E0 and E7

下载: 全尺寸图片 幻灯片

为了进一步分析模型E0与E7的差异，计算了两个模型再参数化的参数在不同地位指数下的值、拐点位置和拐点处的优势木高（表5）。由表5的结果可以看出：通过ADA法建立的差分地位指数模型E0，在不同地位指数下曲线的极值之间只有简单的比例关系，它们的拐点位置完全相同，而拐点处的斜率和预估值同样成比例，满足了多条水平渐近线的性质，然而没有满足多形性，这是因为再参数化的过程中，仅再参数化了代表树木生长极大值的a，而决定曲线形状的参数b和c保持固定。通过GADA法建立的差分地位指数模型E7，随着地位指数的提高，相邻曲线之间极值的差距逐渐减小，拐点位置明显左移，相邻曲线拐点处斜率之间的差距逐渐增大，预估值之间的差距逐渐减小，同时满足了多条水平渐近线和多形性的性质，可以更好地预测不同立地条件下长白落叶松人工林优势木的生长过程。

表  5  差分地位指数模型E0与E7在不同地位指数下再参数化参数、拐点、连年生长量最大值和拐点处预估值对比

Table  5.  Comparison of reparameterized parameters, inflections, maximum height growth rates and fitted values at inflection of difference site index models E0 and E7

地位指数 Site index/m	E0					E7
	a	拐点/a Inflection/ year	连年高生长量最大值 Annual maximum height growth rate/m	拐点处优势木高 Dominant tree height at inflection/m		a	c	拐点/a Inflection/ year	连年高生长量最大值 Annual maximum height growth rate/m	拐点处优势木高 Dominant tree height at inflection/m
12	20.37	7.01	0.571	3.344		21.20	1.53	9.10	0.565	4.190
14	23.77	7.01	0.666	4.075		24.16	1.47	8.20	0.662	4.512
16	27.16	7.01	0.761	4.580		27.09	1.42	7.45	0.761	4.789
18	30.56	7.01	0.857	5.153		30.00	1.37	6.80	0.863	5.025
20	33.96	7.01	0.952	5.725		32.89	1.34	6.23	0.966	5.223
22	37.35	7.01	1.047	6.298		35.76	1.31	5.72	1.073	5.386

下载: 导出CSV 

| 显示表格

3.   结论与讨论

本研究以修正Weibull、Korf和Richards方程为基础方程，使用GADA法推导了6个差分地位指数模型，基于黑龙江省长白落叶松人工林解析木调查数据，拟合出模型中的参数，选出了拟合效果较好的3个模型E1、E2和E7。绘制了3个较优模型的残差图和地位指数曲线簇进行比较，挑选出以Richards生长方程为基础，自由参数为$a = {{\rm{e}}^{{X_0}}}$、c = c₂/X₀所推导的E7为本研究所建立的最优模型。它的决定系数为0.987 41，均方根误差为0.749 3，模拟效率为97.04%，平均绝对误差为0.904 0。此外，本研究通过对比模型E7与E0的图形与曲线性质，证明了GADA法相比于ADA法的优势。结果表明：

（1）在推导地位指数模型时，根据GADA法，指定多个参数为自由参数所推导出的差分模型不仅具有良好的拟合效果，也能同时符合多条水平渐近线与多形性的性质，而ADA法只能满足其中的一个条件。

（2）最优模型为${h_2} =$$\exp ( {X_0} \; - \; 0.046 \; 8{t_2}^{ - 4.675 \; 4/{X_0}} )$，${X_0} = \dfrac{1}{2} \left[ {\ln {h_1} + }\right.$$\left. { \sqrt {{{\left( {\ln {h_1}} \right)}^2} - 18.701 \; 6\ln \left( {1 - {{\rm{e}}^{ - 0.046 \; 8{t_1}}}} \right)} } \right]$，具有良好的拟合效果与预测能力，可以作为黑龙江省长白落叶松人工林地位指数预估的依据。

（3）根据模型E7，随着地位指数的提升，优势木高生长曲线的极值逐渐增大，拐点位置逐渐左移，拐点处斜率增大的速率逐渐增加，而模型极值和拐点处的预估值增大的速率逐渐降低。这说明较好的立地条件下长白落叶松人工林优势木树高连年生长量会更快地增长，并更早达到最大值，但达到最大值之后，其增长的速率会更快地下降。也就是说，较好的立地条件对长白落叶松人工林年龄较小时优势木高生长的促进作用逐渐增强，但对于年龄较大时优势木高生长的促进作用会逐渐减弱。这可能由于立地条件较差时，限制优势木高生长的因素主要来自于其所处的环境条件；而立地条件较好时，限制优势木高生长的因素主要由其自身的性质产生。

本研究所选取的采样点散布于黑龙江境内长白落叶松的分布区，包含了不同的立地类型，在空间分布上具有较强的代表性；而所选取的解析木年龄包含了人工林经营收获的不同阶段，在时间上同样具有较强的代表性。总的来说，使用GADA法推导，利用黑龙江省长白落叶松人工林优势木解析木数据拟合的差分地位指数模型具有良好的拟合效果与预测能力，为黑龙江省长白落叶松人工林立地质量精准评价提供了基础。