基于立地类型哑变量的马尾松优势木枝下高模型研究

童建明

(国家林业和草原局中南调查规划设计院，长沙 410014)

枝下高(Heighttocrownbase)一般指的是直立树干上第一根活树枝到地面的高度，它是单木树冠特征的一个重要指标[1]，不仅能反映树木的生长活力和立地生产力，还是林分内竞争水平等的外在直观体现[2]。大量研究表明，枝下高常被作为输入变量应用于冠形模型、冠幅模型以及生物量模型等，此外也可用于计算单株或林分生长收获模型的关键预测变量，如冠长、冠长率和其它等。但是，在实际生产经营活动中对立木枝下高的获取较为困难，尤其是在林况通透性差、郁闭度高的林分中，其测量误差较大，造成了大量人力和物质资源的浪费。同时，因为工作人员的自我意识，在测量时对第一活枝的判断会存在偏差，从而影响数据的准确性；另外，同一树种由于立地类型的不同，其枝下高也会存在较大的差异。因此在枝下高的广义模型中，仅仅以树高和胸径作为自变量建立枝下高的单一模型，不足以反映所有可能的枝下高与相关测树因子的关系，而基于立地类型哑变量为枝下高模型的研究，尝试提供了一种有效的思路与途径。

哑变量(dummyvariable)，也叫做指示变量、假变量。在哑变量基础上的回归分析方法，可参照李希非等[3]的探索。目前大部分的探索成果表明，哑变量是处理分类变量和定性因子经常使用的一种方法，统计学中的每一种量化方法全部都有涉及哑变量处理的问题[4-6]。同时哑变量模型广泛运用于不同的建模实践与回归分析中[7]，该模型能有效表示调查中的定性因子，模型内既能体现复杂的林分结构，也能反映海拔、坡度等立地因子的不同水平。

1 研究区概况

湖南省位于中国中南部，24°38′—30°08′N，108°47′—114°15′E，在长江中游，省会长沙。陆地总面积大概是3 174.35万km2，其中：山地占51%、盆地占7%、平原占13%、丘陵占29%。东部、西部和南部的地形很高，形成了一个东北有开口的马蹄形形状。属亚热带季风潮湿的大陆气候，年均气温16～18 ℃度至摄氏度，年日照时间1 300～1 800 h；年均降水量为1 200～1 700 mm，降雨量丰富，是我们国家降雨量较多的一个省份[8]；土壤重点是黄色土壤与红色土壤，冲积土壤与紫色土壤其次。

马尾松(Pinusmassoniana)，是松科树种，广泛分布于我国中南部，北自河南、鲁南，南自广东、广西，东自湖南、台湾，西自四川中部、贵州[9]，是中国南部主要的木材品种，具有很高的经济价值。

2 研究方法

2.1 数据来源

数据采集于2020年9—11月，在湖南省怀化、邵阳和益阳三地共调查了42块马尾松人工林临时样地，测量样地内胸径大于5 cm的每株立木。样地调查内容主要涉及树高、枝下高、胸径和其它测树因子，以及坡面、海拔、坡度、土壤类型、土壤厚度和其它立地因素。优势树种(组)由实测数据计算，选取优势木高(HD)等变量指标。

2.2 优势木选取

以样地内所有马尾松胸高断面积(ΣBAi)与所有树种胸高断面积(BA)之和为商，计算出各样地马尾松的组成系数(XSi)，也就是XSi=BAi/BA。为划分优势树种(组)，若ΣXSi≥0.65，则认为该临时样地优势树种(组)为马尾松，样地保留，否则舍弃。

在优势树种(组)为马尾松的样地中，在优势树种(组)中选取3棵优势树(包括次优势树)，分别统计树高和胸径，取其平均值，得到样地平均优势木[10]树高和胸径，根据这两个数值在样地内选取与之最接近的一株立木，该立木数据即为构建模型的原始数据。所有样地构建模型数据，详见表1。

表1 建模数据基本统计量项目枝下高(HCB)/m树高(HD)/m胸径(D)/cm高径比(RHD)平均值 6.915.617.10.92最小值 18.5110.69最大值 13.52224.41.29标准差 3.23.23.50.14变异系数0.50.20.20.15

2.3 模型选择

本研究中，选择5种常见的枝下高模型[11-14](表2)，作为研究枝下高的基础模型，讨论立地类型哑变量马尾松优势木枝下高模型的最佳参数形式。

表2 枝下高基础模型模型 表达式M1HCB= HD(1+exp(X))M2HCB= HD(1+exp(X))M3HCB= HD6(1+exp(X))M4HCB= HD×(1-exp(X))M5HCB= HD×(a+exp(X)) 注: 式中,X=b0+b1∗D;其中 a ,b0,b1为模型参数。

2.4 模型及参数检验

模型用确定系数(R2)、均方根误差(RMSE)和平均系统误差(MSE)3个标准进行评估和比较。其中R2反映模型的适用性，R2越接近1，模型就越准确；RMSE值越小，模型就越准确；MSE是反映一定范围(3%或5%)内拟合效果的关键指标，接近0时，可以得到最佳效果。具体的公式为：

2.5 哑变量模型构建

立地类型对林木枝下高影响较显著，为探索这种显著性，首先将初始立地类型作为哑变量应用在基础模型中拟合，再采用R语言k均值算法(k-means)对初始立地类型分级处理，其分类数标准为聚类精度≥0.99[15]，基于聚类结果最终构建含聚类后的立地类型哑变量模型。

本文在构建马尾松优势木枝下高模型时，在确保预测准确性的基础上，尽可能简化了模型，建立了具有立地类型哑变量的优势木枝下高模型，以确定模型的最佳模式。

那么基于模型M1的哑变量模型的形式可表示为：

式中:Zi为哑变量，bi为对应的具体参数或局部参数。

按照立地类型的分类、聚类，立地类型涉及6个哑变量；即Ⅰ，Ⅱ，Ⅲ，Ⅳ，Ⅴ和Ⅵ。如果立地类型是Z1的时候，取Z1=1；那么Z2，Z3，…，Z6均为0，以此类推。

2.5.1 初始立地类型哑变量

对于初始立地类型的划分，通过数量化方法I得出6个立地因子(海拔、坡度、坡向、坡位、土层厚度、土壤类型)对枝下高的显著影响，并根据《湖南省森林资源规划设计调查技术规程(2013年版)》对所有立地因子进行组合，每个立地因子水平组合为一个立地类型[16-17]，42块样地共划分为21个初始立地类型。

2.5.2 聚类立地类型哑变量

基于最优基础模型考虑21个初始立地类型哑变量拟合的得分值，为了研究方便以及考虑模型精度，采用R语言k均值算法(k-means)对各因子分级处理，允许合并后因子水平信息损失≤1%。

3 结果与分析

3.1 基础模型拟合与评价

基于样地的平均优势木数据，对基本模型的枝下高模型进行拟合，结果详见表3。

表3 基础模型拟合结果模型R2RMSEMSE/%M10.320 42.738 3-2.594 0M20.320 22.738 8-2.632 5M30.319 82.739 6-2.678 2M40.321 52.736 2-2.430 5M5不收敛——

结果显示，模型M1—M4的确定系数都较低。其中，模型M4的确定系数最大，为0.321 5；模型M3的确定系数最小，为0.319 8；模型M5拟合结果不收敛。考虑到4种基础模型(M1—M4)拟合的确定系数差值均小于0.001 7，并且该模型具有固有的可变性。所以，四个基本模型都作为构建初始立地类型哑变量枝下高的模型，其后再根据评价指标选出最优模型。

3.2 哑变量模型拟合与评价

利用Forstat软件中非线性混合效应模块，以初始立地类型为哑变量加在4种候选枝下高模型不同参数上拟合；依据R2，AIC和BIC进行评价。再分析哑变量在不同模型上的拟合效果，选出最优模型来构建含聚类后立地类型哑变量的候选模型。

3.2.1 初始立地类型哑变量拟合结果

根据数量化方法I结果得出6个立地因子对枝下高均显著，不同因子水平组合为21个初始立地类型。利用初始立地类型作为哑变量加在4种基础模型不同参数上拟合，所有模型拟合指标结果详见表4。

表4 初始立地类型哑变量拟合结果模型R2AICBICM1-b00.773 8198.135 4221.114 9M1-b10.772 4313.831 5 336.733 4 M2-b00.773 8183.764 2210.505 6M2-b10.770 2301.255 2324.234 7M3-b00.773 8147.386 9170.366 4M3-b10.767 7265.059 7288.039 2M4-b00.773 8237.653 1260.632 6M4-b10.771 8354.972 9377.952 4

由表4可知，加入初始立地类型哑变量后，模型确定系数从0.319 8～0.321 5提高到0.767 7～0.773 8，且收敛度均较好，4种模型的确定系数差值均小于0.006 1，都是哑变量加在模型参数b0上最大，均为0.773 8。另外哑变量加在同一模型的不同参数上时，其确定系数差值很小。

基于4种模型的拟合结果，综合考虑评价指标，最终选择模型M1作为后期研究的候选模型。同时，由于初始立地类型哑变量加在模型M1上的两个不同参数(b0和b1)得到拟合的确定系数很接近，分别为0.773 8和0.772 4，考虑到模型内在的相容性，将分别以模型M1两个不同参数的拟合结果来构建含聚类后立地类型哑变量模型，最终选择最优模型形式。

3.2.2 聚类立地类型哑变量拟合结果

根据42块样地初始划分的21个立地类型，分别以模型M1两个参数形式拟合的初始立地类型得分值聚类，聚类后分类数对应的立地编号及数量如表5。

表5 立地类型聚类模型M1-b0模型M1-b1立地类型编号数量立地类型编号数量Ⅰ16Ⅰ17Ⅱ4Ⅱ6Ⅲ12Ⅲ10Ⅳ5Ⅳ4Ⅴ1Ⅴ1Ⅵ4Ⅵ4

每个样本的初始立地类型被转换成相应的类，并作为哑变量添加到模型M1的相应参数中进行拟合，也就是说模型M1的参数b0上的聚类立地类型哑变量仅被添加到参数b0中进行拟合，而参数b1上的聚类立地类型哑变量仅被添加到参数b1中进行拟合，分析并比较两种拟合结果。

基于聚类后立地类型哑变量，获得在模型M1的两个不同参数下模拟的RMSE，R2，MSE，AIC和BIC，详见表6。

表6 模型不同参数拟合结果模型R2RMSEMSE/ %AICBICM1-b00.769 21.595 8-1.031 3184.413 1196.855 8M1-b10.767 01.603 6-1.140 5211.969 1224.411 9

结果表明，在聚类后添加立地类型哑变量后，模型的拟合结果要比基本模型好得多，但相比初始立地类型哑变量拟合的确定系数有所降低。其中，哑变量加在模型M1参数b0上的R2最大，为0.769 2；RMSE和MSE最小，分别为1.595 8和-1.0313；AIC和BIC值也优于参数b1上的拟合结果。因此，聚类立地类型哑变量中，模型M1参数b0上的形式拟合最优。

3.2.3 模型参数估计

基础模型参数值和立地类型哑变量模型参数估计，详见表7和表8。

表7 基础模型参数值参数 估计值渐近标准差参数下限参数上限b00.564 30.596 4-0.641 11.769 6b1-0.019 80.032 4-0.085 20.045 6

在模型M1中，基础模型和哑变量模型的全部参数均不错。立地类型哑变量模型中哑变量的估计值有很大不同，即此项探索中，聚类后的立地类型哑变量模型较好地体现了不同立地类型之间的不同，将模型形式进行了简化。另外，哑变量模型中的形式参数都具有较好的稳定性。

表8 立地类型哑变量模型参数估计参数估计值SDtPⅠb00.491 40.539 10.911 50.368 2Ⅱb0-0.852 70.519 1-1.642 80.109 4Ⅲb0-0.395 00.570 0-0.692 90.492 9Ⅳb0-0.164 90.457 8-0.360 30.720 8Ⅴb01.748 30.979 41.785 00.082 9Ⅵb01.439 60.659 22.183 90.035 8b10.007 60.028 20.269 70.789 0

3.2.4 最优模型残差图

按照拟合的成果，将选择一种最理想的模型M1最好的哑变量参数形式。用枝下高实测值减去枝下高预测值得到残差，残差图以枝下高预测值为横轴，残差为纵轴进行描绘，具体结果见模型M1基础的残差图(图1)和聚类后立地类型哑变量最优参数形式的残差图(图2)。

图1 基础模型残差图

图2 聚类立地类型哑变量残差图

图1和图2直观地反映了哑变量模型方法的优势，以聚类立地类型哑变量进行估计后的枝下高残差图要优于基础模型的拟合；另外，哑变量模型的枝下高估计值与实测值相差程度都较小，也就是说，各对应的立地类型哑变量模型得到的估计值与基础模型估计值差异较小，这验证了模型之间的相容性。

3.3 模型检验

利用建模样本对构建的枝下高模型进行检验。由于调查的样地数量有限，在全部样本数用于建模的情况下，本文模型检验采取分径阶检验的方式，即将样本数以胸径分组，最终分别以径阶12、14、16、18、20和22～24六组数据检验模型，检验指标选择确定系数(R2)、均方根误差(RMSE)和平均系统误差(MSE)3个评价指标进行评价和比较。结果见表9。

表9 模型分径阶拟合结果径阶/cm样本数R2RMSEMSE/%1250.999 00.060 2-0.221 5 1480.694 71.886 50.134 7 1690.960 20.509 1-0.124 9 1870.855 51.301 9-1.521 4 2080.751 71.046 7-0.318 0 22～2450.669 92.491 2-3.045 9

从表9中可知，建模样本检验中，各径阶拟合结果均较好，但差异较大。在各径阶上，聚类后立地类型哑变量模型确定系数均优于基础模型，均方根误差(RMSE)和平均系统误差(MSE)也均表明良好。检验结果较好地显示了模型M1所建的马尾松枝下高模型具有良好的全面切合性能。

4 结论与讨论

因为此项探索的数据收集仅限于湖南三个地方，所以模型的应用范围具有局限性。在优势树的选择中，选择优势树种(组)中的三棵优势木，没有与其他选取方法数据作对比，在以后的研究中可尝试以两种以上不同优势木来分析，如最高优势木(优势树种组中最高的一株立木)，平均优势木等。除此之外，关于立地类型哑变量，立地类型哑变量模型法的应用受到限制，因为在此项探索中全部的样地初始的立地类型并不包括湖南地区所有的立地类型。此篇论文运用R语言的k-means算法将初始立地类型聚类成若干类作为哑变量，简化哑变量类型是为了说明此方法的应用，因为不包括这一区域的全部立地类型，其实用性需进一步探讨。研究考虑了立地类型哑变量，为区域性马尾松人工林枝下高模型的研究提供了一种有效的思路与途径。