用SASmixed过程拟合林分线性差分生长模型.docVIP

用SASmixed过程拟合林分线性差分生长模型 　　【摘 要】本研究的目的在于研究如何用SAS的proc mixed过程拟合线性代数差分模型。所用数据来源于148个集约经营火炬松人工林。直接拟合了一个胸高断面积的收获模型，而非代数差分生长模型。模型拟合过程如下：i）.同时确定随林分变化的参数和最优拟合的方差结构模型；ii）.依据AIC、BIC和极大似然比检验化简期望模型；iii）.用代数差分法将拟合的收获模型转化为代数差分生长模型。 　　【关键词】线性代数差分模型；mixed过程；模型筛选；林分生长与收获预估 　　0.前言 　　在林分生长与收获预估的模型中，差分生长模型得到了广泛的应用。线性差分模型基本上为Schumacher模型的变型，广泛应用于林分蓄积、胸高断面积的建模，以及单位面积株数和优势木树高生长模型。 　　差分生长模型拟合方法有“直接最小二乘估计法”和“分类变量回归法”[1]。一般认为后???可以获得近似无偏的估计，而前者则导致检验统计量如RMSE的失真[2]。传统上差分生长模型的拟合主要是直接拟合差分生长模型，然后根据拟合统计量如RMSE、R2等确定最优拟合模型。与传统方法不同，本文直接拟合生长模型，在获得参数估计值后，再用代数差分法导出相应的代数差分生长模型。这样做的优越之处在于非常便于对期望模型和方差结构模型进行筛选。更为重要的是，可以通过模型拟合识别最适合的随林分变化参数。 　　本文详细讨论了如何用“分类变量回归法”和SAS的mixed过程拟合代数差分生长模型，可简述如下：i）.直接以生长收获模型为对象，同时确定一个随林分变化的参数和最优拟合方差结构模型；ii）. 保持方差结构模型不变，根据拟合统计量逐步化简期望模型；iii）.在确定最优拟合的期望模型后，运用代数差分法导出相对应的代数差分生长模型。所有拟合与筛选均用SAS的mixed过程完成，并给出了详细的SAS代码和代码解释。 　　1.方法与材料 　　1.1数据 　　数据来源于148个集约经营的火炬松实验人工林逐年观测的固定样地数据（样地约0.152公顷）。SAS的数据集basal的内容如表（1）。 　　表1 模型拟合的基本数据结构 　　Table 1 data structure for mode fitting 　　age=林分年龄；fert=经营措施，分别取值为H=施加除草剂以控制竞争植物、F=施肥以增加土壤肥力、HF=除草剂和施肥并用、C=对照；code=样地代码，每个样地有一个唯一代码；logba=样地胸高断面积的自然对数值；logdh=样地优势木树高的自然对数值；iage=林分年龄的倒数；logtpa=样地株数的自然对数值。共有148个样地数据，1491个记录。 　　1.2数学模型 　　所考虑的胸高断面积收获的数学模型为 　　E（lny）=α+αlnN+αlnH+α+α （1） 　　这里lny、lnH和lnN分别对应logba、logdh和logtpa；1/t对应iage，其余为模型参数。通过不同假设，以上模型可以导出很多广泛应用的差分生长模型。例如假设α0为随林分变化的参数，运用代数差分法可以导出Pienaar和Shiver （1986）的胸高断面积模型[3]；假设α4=0且α0为随林分变化的参数则可导出Forss等人（1996）提出的差分生长模型[4]；假设α4=0和α2=0且α3为随林分变化的参数，则可以导出Souter （1986）的模型[5]。 　　1.3参数估计方法 　　从以上代数差分法的讨论中，总是假设模型（1）中有一个参数随林分变化而变化，因而使用“分类变量回归法”以便考虑这一特点。如果用SAS的reg过程，必须构造特殊的数据结构。用SAS的reg过程拟合线性差分生长模型有两点不足之处：1）.需要构造复杂的数据文件。以本研究的为例，数据文件中需要增加额外的148个变量；2）.reg过程无法考虑重复观测数据的自相关性和异质方差结构。由于所考虑的胸高断面积模型均为线性模型，而且数据为典型的重复观测数据，因而模型的拟合与筛选均用SAS的proc mixed过程。 　　1.4模型筛选和筛选指标选用 　　模型筛选包括期望模型和最优拟合的自相关与异质方差结构模型的筛选。Ngo和Brand（1997）推荐了两种模型选择方法[6]。一种方法就是首先列出所有可用的数学期望模型和所有可能的方差结构模型，拟合二者间的所有组合，然后根据模型拟合指标选择最优拟合模型。另一种方法为Wolfinger和Diggle提出的方法，即首先考虑最复杂的期望模型，并保持期望模型不变，然后选择最优拟合的方差结构模型。一旦选出最优拟合的方差模型，保持方差模型不变，再逐步简化期望模型。模型筛选指标较多，但常用的指标为极大似然比检验