本发明涉及一种林分胸径分布的估算方法，特别是一种毛竹林分胸径分布的估算方法。
背景技术：
胸径结构是林分结构最重要的指标，毛竹林分具有独特的生物学特征，为研究林分胸径结构分布提供了一个很好的试验样本，毛竹林分胸径分布对毛竹林生物量估算及其经营活动起着非常重要的作用，不仅可预估毛竹林各径阶株数、生长量测定、产量预估、枯损量计算等，还能对毛竹林抚育间伐设计及其碳库估算提供科学依据。
目前对毛竹林分胸径结构研究最多的是林分胸径分布函数的构建、比较与选择，而通过构建林分胸径分布通用模型来估算毛竹林分胸径分布还尚待深入研究，毛竹林分胸径分布通用模型的构建可直接估算任意毛竹林分的胸径分布，不仅对生产经营有重要指导意义，而且能减少林分调查所消耗的人力、物力与财力。因此，亟需建立一种毛竹林分胸径分布估算的通用模型。
技术实现要素：
本发明的目的在于，提供一种毛竹林分胸径分布的估算方法。它具有可以估算任意毛竹林分胸径分布的优点。
本发明的技术方案：毛竹林分胸径分布的估算方法，其特征在于，包括以下步骤：
步骤一，收集待估年份的n个毛竹固定样地连续清查数据；
步骤二，用3参数weibull概率密度函数分别拟合步骤一中每个毛竹样地连续清查数据并用柯尔莫哥洛夫检验法对每个样地的weibull概率密度函数进行拟合优度检验；
步骤三，根据柯尔莫哥洛夫拟合优度检验结果，计算能用3参数weibull概率密度函数描述毛竹样地胸径分布的样地数占总毛竹样地数的比例，得毛竹林分服从weibull分布；
步骤四，从n个毛竹样地中取m个样地为建模样本，其余的n-m个样地为检验样本；
步骤五，设样地平均胸径dg与样地平均胸径的平方dg2都服从gamma分布，联接函数均为对数函数，以待估年份毛竹样地weibull概率密度函数的位置参数a、尺度参数b与形状参数c为解释变量，毛竹样地的平均胸径dg与dg2为因变量，用广义线性模型建立样地平均胸径与位置参数a、尺度参数b、形状参数c的函数关系式dg＝exp(ω+βa+γb+ηc)，即为g函数，其中ω，β，γ，η为g函数参数，用广义线性模型建立样地平均胸径的平方dg2与位置参数a、尺度参数b的关系式dg2＝exp(ε+θa+μb)，即为h函数，其中ε，θ，μ为h函数参数；
步骤六，用matlab软件的glmfit函数分别拟合g函数的参数ω，β，γ，η与h函数的参数ε，θ，μ，拟合函数的数据采用建模样地连续清查数据；
步骤七，用matlab软件的glmfit函数分别计算g函数与h函数的拟合误差w1与w2，计算拟合误差的数据采用建模样地连续清查数据，再查χ2分布表得
步骤八，对m个建模样本与建模样本的weibull概率密度函数参数a，b，c，用非线性最小二乘拟合g函数，得建模精度，用
步骤九，对m个建模样本与建模样本的weibull概率密度函数参数a，b，用非线性最小二乘拟合h函数，得建模精度，用
步骤十，建立n个毛竹样地的weibull概率密度函数尺度参数b与形状参数c的关系式c＝φbυ，即f函数，其中φ，υ为f函数参数，经非线性最小二乘法拟合得参数φ，υ与r2的值；
步骤十一，把步骤十中得到的参数φ，υ带入f函数得p函数；
步骤十二，建立待估年份毛竹样地胸径分布的通用模型，即由l函数，k函数，p函数与3参数weibull概率密度函数组成的方程组。
前述的毛竹林分胸径分布的估算方法中，所述步骤四中毛竹样地数n为≥170个，建模样本数m为≥120个，检验样本数n-m为≥50个。
前述的毛竹林分胸径分布的估算方法中，对于给定任意一个毛竹样地或毛竹林分，根据毛竹林分胸径分布的通用模型，由该毛竹林分的平均胸径可以得出该毛竹林分的weibull概率密度函数3参数，进而获得该毛竹林分某胸径的概率值，从而估算出该毛竹林分的胸径分布。
与现有技术相比，本发明通过收集待估年份的n个毛竹固定样地连续清查数据，用3参数weibull概率密度函数分别拟合每个毛竹样地连续清查数据并用柯尔莫哥洛夫检验法对每个样地的weibull概率密度函数进行拟合优度检验；根据柯尔莫哥洛夫拟合优度检验结果，计算能用3参数weibull概率密度函数描述毛竹样地胸径分布的样地数占总毛竹样地数的比例；设样地平均胸径dg与样地平均胸径的平方dg2都服从gamma分布，联接函数均为对数函数，以待估年份建模样地weibull概率密度函数的位置参数a、尺度参数b与形状参数c为解释变量，建模样地的平均胸径dg与dg2为因变量，用广义线性模型建立样地平均胸径与位置参数a、尺度参数b、形状参数c的函数关系式dg＝exp(ω+βa+γb+ηc)，即为g函数，用广义线性模型建立样地平均胸径的平方dg2与位置参数a、尺度参数b的关系式dg2＝exp(ε+θa+μb)，即为h函数；用matlab软件的glmfit函数分别拟合g函数的参数ω，β，γ，η与h函数的参数ε，θ，μ，拟合函数的数据采用建模样地连续清查数据；将参数ω，β，γ，η带入g函数得l函数；将参数ε，θ，μ带入h函数得k函数；建立n个毛竹样地的weibull概率密度函数尺度参数b与形状参数c的关系式c＝φbυ，即f函数，经非线性最小二乘法拟合得参数φ，υ与r2的值；将参数φ，υ带入f函数得p函数；建立待估年份毛竹样地胸径分布的通用模型，即由l函数，k函数，p函数与3参数weibull概率密度函数组成的方程组；l函数、k函数与p函数联合在一起，在只知道林分平均胸径的条件下可以解出任何一个林分的weibull概率密度函数3参数，从而可以建立待估年份林分胸径分布的通用模型。综上，本发明具有可以估算任意毛竹林分胸径分布的优点。
具体实施方式
下面结合实施例对本发明作进一步的说明，但并不作为对本发明限制的依据。
实施例：毛竹林分胸径分布的估算方法，其特征在于，按以下步骤进行：
步骤一，收集待估年份的n个毛竹固定样地连续清查数据；
步骤二，用3参数weibull概率密度函数分别拟合步骤一中每个毛竹样地连续清查数据并用柯尔莫哥洛夫检验法对每个样地的weibull概率密度函数进行拟合优度检验；
步骤三，根据柯尔莫哥洛夫拟合优度检验结果，计算能用3参数weibull概率密度函数描述毛竹样地胸径分布的样地数占总毛竹样地数的比例，以95％为最低值，超过95％说明毛竹林分服从weibull分布；
步骤四，从n个毛竹样地中取m个样地为建模样本，其余的n-m个样地为检验样本；
步骤五，设样地平均胸径dg与样地平均胸径的平方dg2都服从gamma分布，联接函数均为对数函数，以待估年份毛竹样地weibull概率密度函数的位置参数a、尺度参数b与形状参数c为解释变量，毛竹样地的平均胸径dg与dg2为因变量，用广义线性模型建立样地平均胸径与位置参数a、尺度参数b、形状参数c的函数关系式dg＝exp(ω+βa+γb+ηc)，即为g函数，其中ω，β，γ，η为g函数参数，用广义线性模型建立样地平均胸径的平方dg2与位置参数a、尺度参数b的关系式dg2＝exp(ε+θa+μb)，即为h函数，其中ε，θ，μ为h函数参数；
步骤六，用matlab软件的glmfit函数分别拟合g函数的参数ω，β，γ，η与h函数的参数ε，θ，μ，拟合函数的数据采用建模样地连续清查数据；
步骤七，用matlab软件的glmfit函数分别计算g函数与h函数的拟合误差w1与w2，计算拟合误差的数据采用建模样地连续清查数据，再查χ2分布表得
步骤八，对m个建模样本与建模样本的weibull概率密度函数参数a，b，c，用非线性最小二乘拟合g函数，得建模精度，用
步骤九，对m个建模样本与建模样本的weibull概率密度函数参数a，b，用非线性最小二乘拟合h函数，得建模精度，用
步骤十，建立n个毛竹样地的weibull概率密度函数尺度参数b与形状参数c的关系式c＝φbυ，即f函数，其中φ，υ为f函数参数，经非线性最小二乘法拟合得参数φ，υ与r2的值，如果r2的值大于0.9，则说明f函数精度很高；
步骤十一，把步骤十中得到的参数φ，υ带入f函数得p函数；
步骤十二，建立待估年份毛竹样地胸径分布的通用模型，即由l函数，k函数，p函数与3参数weibull概率密度函数组成的方程组。
所述步骤四中毛竹样地数n为≥170个，建模样本数m为≥120个，检验样本数n-m为≥50个，建模样本与检验样本均是大样本数据(样本数大于50即为大样本)，符合统计要求。
对于给定任意一个毛竹样地或毛竹林分，根据毛竹林分胸径分布的通用模型，由该毛竹林分的平均胸径可以得出该毛竹林分的weibull概率密度函数3参数，进而获得该毛竹林分某胸径的概率值，从而估算出该毛竹林分的胸径分布。
如建立2014年毛竹样地胸径分布通用模型：
步骤一，收集待估年份的n个毛竹固定样地连续清查数据；浙江省于1979年建立了森林资源连续清查体系，以5年为一个复查周期。共设置固定样地4250个，样点格网为4km×6km，样地形状为正方形，边长28.28m，面积800m2。本研究利用2009年与2014年两期森林连续清查数据，选择两期毛竹重复样地172个。2009年毛竹样地的基本情况为：每个样地毛竹株数22～897不等，样地内毛竹胸径为2～15cm，年龄1～4度以上(当年生竹记为1度竹；2～3年生竹记为2度竹，依此类推)；2014年毛竹样地的基本情况为：每个样地毛竹株数37～940不等，样地内毛竹胸径为2～15cm，年龄1～4度以上(当年生竹记为1度竹；2～3年生竹记为2度竹，依此类推)。样地的调查因子包括土层厚度、坡向(1：北，2：东北，3：东，4：东南，5：南，6：西南，7：西，8：西北，9：无坡向即山顶)、坡位(1：脊部，2：上部，3：中部，4：下部，5：谷部，6：平地)、坡度、海拔(10m～1200m)、样地平均胸径、样地株数及样地内每株毛竹的胸径与年龄等；本发明利用2014年172个毛竹样地的连续清查数据构建毛竹林分胸径分布估算的通用模型；
步骤二，用3参数weibull概率密度函数分别拟合步骤一中每个毛竹样地连续清查数据并用柯尔莫哥洛夫检验法对每个样地的weibull概率密度函数进行拟合优度检验，结果见表1；
步骤三，根据柯尔莫哥洛夫拟合优度检验结果，计算能用3参数weibull概率密度函数描述毛竹样地胸径分布的样地数占总毛竹样地数的比例，以95％为最低值，超过95％说明毛竹林分服从weibull分布；2009年、2014年两期连续清查数据的毛竹重复样地172个，其中2014年服从weibull分布的毛竹样地数占总样地数的比例为0.9826，而2009年的比例为0.9709，因毛竹独特的生物学特性与采伐措施，使毛竹林分胸径结构的周期为5年，故对毛竹林分而言，其林分胸径分布是稳定的，服从weibull分布的毛竹样地数占总样地数的比例低于0.95是不会出现的，因此，验证得2014年毛竹林分服从weibull分布；
步骤四，从n个毛竹样地中取m个样地为建模样本，其余的n-m个样地为检验样本；这里以122个样地为建模样本，其余的50个样地为检验样本；
步骤五，设样地平均胸径dg与样地平均胸径的平方dg2都服从gamma分布，联接函数均为对数函数，以待估年份毛竹样地weibull概率密度函数的位置参数a、尺度参数b与形状参数c为解释变量，毛竹样地的平均胸径dg与dg2为因变量，用广义线性模型建立样地平均胸径与位置参数a、尺度参数b、形状参数c的函数关系式dg＝exp(ω+βa+γb+ηc)，即为g函数，其中ω，β，γ，η为g函数参数，用广义线性模型建立样地平均胸径的平方dg2与位置参数a、尺度参数b的关系式dg2＝exp(ε+θa+μb)，即为h函数，其中ε，θ，μ为h函数参数；
步骤六，用matlab软件的glmfit函数分别拟合g函数的参数ω，β，γ，η与h函数的参数ε，θ，μ，拟合函数的数据采用建模样地连续清查数据；
步骤七，用matlab软件的glmfit函数分别计算g函数与h函数的拟合误差w1与w2，计算拟合误差的数据采用建模样地连续清查数据，再查χ2分布表得
步骤八，对m个建模样本与建模样本的weibull概率密度函数参数a，b，c，用非线性最小二乘拟合g函数，得建模精度，用
步骤九，对m个建模样本与建模样本的weibull概率密度函数参数a，b，用非线性最小二乘拟合h函数，得建模精度，用
步骤十，建立n个毛竹样地的weibull概率密度函数尺度参数b与形状参数c的关系式c＝φbυ，即f函数，其中φ，υ为f函数参数，经非线性最小二乘法拟合得参数φ，υ与r2的值，如果r2的值大于0.9，则说明f函数精度很高；经多次尝试得172个样地形状参数b与尺度参数c有如下关系：c＝φbυ，φ，υ为模型参数，经非线性最小二乘拟合得(φ，υ)＝(0.989575,0.933805)，r2＝0.9830，可见精度很高；
步骤十一，把步骤十中得到的参数φ，υ带入f函数得p函数；
步骤十二，建立待估年份毛竹样地胸径分布的通用模型，即由l函数，k函数，p函数与3参数weibull概率密度函数组成的方程组。
这里的l函数为dg＝exp(1.0904+0.1131a+0.1130b+0.0001c)，k函数为dg2＝exp(2.1829+0.2261a+0.2261b)，p函数为c＝0.989575b0.933805；l函数、k函数与p函数联合在一起，在只知道林分平均胸径的条件下可以解出任何一个林分的weibull概率密度函数3参数，从而可以建立2014年林分胸径分布的通用模型。表1weibull概率密度函数参数与检验
表1中，a，b，c为weibull概率密度函数
的位置参数，尺度参数，形状参数，x为自变量，在本发明中为毛竹林分中毛竹的某个胸径；dn为经验分布函数与理论分布函数偏差的最大值，dna为柯尔莫哥洛夫检验临界值。当dn<dna时，说明该毛竹样地服从weibull分布，此时h＝0，当dn>dna时，说明该毛竹样地不服从weibull分布，此时h＝1。从表1可知，服从weibull分布的毛竹样地数占总样地数的比例为0.9826。