非线性混合效应模型及其在林业上应用
|
摘要
非线性混合效应模型（Nonlinear Mixed Effects Models，简称NLMEMs）是依据回归函数依赖于固定效应参数和随机效应参数的非线性关系而建立的。它是分析纵向数据、多水平数据及重复调查数据的近代统计学方法，既能反映总体的平均变化趋势，又能刻画个体差异的统计分析工具。近几年，NLMEMs受到越来越多的学者关注，并且被引入到多门学科，如医学、工学、农学及林学等。经过30多年的发展历程，学者已提出了单水平和嵌套多水平NLMEMs的多种参数估计方法，典型的计算软件有SAS和S-Plus。但在实际应用中发现2个问题，一是主流软件SAS和S-Plus计算经常不收敛，尤其是模型中待估参数较多时最为明显，二是现有的混合模型没有包括所有类型的随机效应的组合（例如交互作用，这在林业中是常用的），因而限制了模型的应用。本文的目的在于（1）提出一种收敛性好的计算方法来分析单水平和嵌套多水平NLMEMs；（2）提出一个包括所有随机效应类型的NLMEMs的统一表达式，并且给出一种参数计算方法；（3）完成上述2个算法的程序代码并在ForStat上实现；（4）应用我们的程序解决一个林学上的实际问题，该问题用已有的程序是无法解决的。
本研究实现了上述4个目的。具体内容如下：1）根据一阶条件期望线性化—数学期望极大方法（First-order conditional expectation linearization–expectation maximation，简称FOCE-EM）的理论推导出计算单水平和嵌套多水平NLMEMs的计算公式、设计了计算流程。2）本研究提出了一种正态NLMEMs的标准表达式，它包括了正态NLMEMs的所有随机效应类型，给出该模型相应的一种参数估计方法，即线性逼近—逐步2次规划算法。3）发现并通过实例说明SAS中mixed模块不能保证随机效应参数方差为非负定矩阵，所以在本文提出采用线性逼近—逐步2次规划算法，同时给出十种方差类型满足正定或半正定的条件，因而算法可以保证不出现类似SAS中的错误。4）模型的标准表达式可以处理固定效应和随机效应参数分级（即数量化问题），并且指出了S-Plus中nlme函数分级差法计算有缺陷，而线性逼近—逐步2次规划算法克服了此问题。5）首次利用带有交互作用的两因素（林分密度和地位级指数）NLMEMs分析了落叶松树高—直径模型。在此基础上进一步分析了随机效应与海拔的关系。
通过本项研究，可以得出以下主要结论：1）本研究提出了一种正态NLMEMs的标准表达式，它包含多种类型的非线性混合效应模型（随机效应参数服从正态分布），例如单水平NLMEMs、逐级嵌套多水平NLMEMs、只含主效应的多因素NLMEMs、包括主效应和交互效应的NLMEMs以及某几种类型组合的一般性NLMEMs等。模型中固定效应参数和随机效应参数可以考虑分级（即数量化）。同时还把正态NLMEMs的标准表达式推广到参数方差与某些因素（称为组变量）有关的NLMEM（s考虑组变量的NLMEMs）。因此模型比传统的NLMEMs表达式更为一般化，具有更广的用途。2）在计算单水平或逐级嵌套的多水平NLMEMs时，FOCE-EM算法与SAS及S-Plus提供的Lindstrom andBates（LB）算法计算精度非常接近，数值实例表明，至少四位小数相同。但FOCE-EM算法从理论上保证了线性步计算收敛，从而使得该算法计算收敛性明显好于LB算法。3）本研究给出线性逼近—逐步2次规划算法计算一般类型的NLMEMs。为保证方差非负定采用逐步2次规划是必要的。该算法同样也能计算考虑组变量的NLMEMs。4）从计算速度上讲，FOCE-EM算法比线性逼近—逐步2次规划算法要快，因此建议在计算单水平和逐级嵌套多水平NLMEMs时，使用FOCE-EM算法，而计算其它类型的NLMEMs时，使用线性逼近—逐步2次规划算法。5)通过两因素NLMEMs对落叶松树高—直径模型研究得出，考虑林分密度与地位级的交互作用能明显提高模型的预测精度，而且把海拔高度作为组变量时，还可以进一步提高模型预测精度。
Nonlinear mixed effects models (NLMEMs) are built on the regression function, whichrelying on the nonlinear relationship of fixed effects parameters and random effects parameters.They are modern statistical methods to analyze longitudinal data, multilevel data and repeatedsurvey data. They also can be considered as powerful statistical tools, which not only reflectthe overall variation trend of population, but also describe individual differences. In recentlyyears, NLMEMs are getting much attention from the scholars in various disciplines, such asmedicine, engineering, agriculture and forestry and so on. After30years’ development,several scholars have proposed the single level NLMEMs and the multi-level NLMEMs, andtypical statistical software such as SAS and S-Plus were developed to estimate them. Inpractice, we found two main problems in existing NLMEMs and estimation methods. One isthat the main software such as SAS and S-Plus may not always reach convergence, especiallywhen the model contains many estimable parameters. The other is that the current NLMEMsmay not cover all types of random effects combinations (such as interaction that is commonlyused in forestry), which limits application of NLMEMs. The objectives of this paper aretherefore to (1) propose an algorithm with good convergence to estimate single level andmulti-level NLMEMs,(2) put forward a unified expression that includes all types of randomeffects of NLMEMs and corresponding parameter estimation method,(3) complete and realizethe coding for the above two algorithm in ForStat software,(4) and use our programs toresolve a practical problem in forestry, which is not resolved yet.
The study has successfully achieved the above four objectives. Specific contents of thepaper are as follows:1) Based on the first-order conditional expectation linearization–expectation maximation (FOCE-EM), it draws calculation formulas and processes of singlelevel and multi-level NLMEMs.2) The study proposed a normal standard expression forNLNEMs, which contains all types of random effects, and put forward the relevant parameterestimation method, naming linearization approximation-sequential quadratic programming algorithm.3) We found the mixed procedure in SAS may not guarantee the variance ofrandom effects parameters as non-negative definite, so we used one specific case to verify thisdrawback. Therefore, we proposed the linearization approximation-sequential quadraticprogramming algorithm to estimate the parameters in the NLMEMs standard expression, andwe present ten common variance types of constraint conditions to guarantee each variancepositive definite or half positive definite.4) The standard expression of this model can dealwith classification (namely quantification problem) of both fixed effects and random effectsparameters. We also pointed out some defects of classification differential method in nlmefunction of S-Plus, and put forward the solutions such as linearization approximation-sequential quadratic programming algorithm to overcome this drawback5) It is the first timethat we used two factors (stand density and site index class) NLMEMs within interactions toanalyze the height–diameter model for Larix olgensis. We also further analyzed therelationship of random effects and altitude based on the two factors NLMEMs.
Based on the analysis and results from the study, the following conclusions can be drawn:1) the proposed normal NLMEMs standard expression contain various types of NLMEMs(random effects is normal distribution), such as single level NLMEMs, grading multi-levelNLMEMs, multi-factor NLMEMs with main effects, multi-factor NLMEMs with main effectsand interaction effects, and general NLMEMs that combined several types of random effects..We also spread the normal NLMEMs standard expression to those NLMEMs that thevariances of random effects are related to some factors (also called as group variable). Ascompared with the traditional models, the proposed models are more general and useful.2)When calculating single level or grading multi-level NLMEMs, the FOCE-EM algorithm hasa similar computational accuracy with the LB algorithm of SAS and S-Plus. However,FOCE-EM algorithm in theory ensures the linear step convergence, therefore, the convergenceof this algorithm is significantly better than LB algorithm.3) The study proposed linearizationapproximation-sequential quadratic programming algorithm to estimate the general types ofNLMEMs. It is necessary to choose this algorithm to guarantee the random effects variance asnon-negative definite.4) The computational speed of FOCE-EM algorithm is more quickly than linearization approximation-sequential quadratic programming algorithm. Therefore, it isbetter to use the FOCE-EM algorithm to analyze single level and multi-level NLMEMs. Forother types of NLMEMs, linearization approximation-sequential quadratic programmingalgorithm would be better.5) Through using the two factor NLMEMs to analyze the height-diameter model for Larix olgensis, it can be drawn that the prediction precision can beobviously increased by considering the interactions of stand density and site index class, and ifthe altitude was considered as group variable in the mixed model, the prediction precisionwould be further improved.
本研究实现了上述4个目的。具体内容如下：1）根据一阶条件期望线性化—数学期望极大方法（First-order conditional expectation linearization–expectation maximation，简称FOCE-EM）的理论推导出计算单水平和嵌套多水平NLMEMs的计算公式、设计了计算流程。2）本研究提出了一种正态NLMEMs的标准表达式，它包括了正态NLMEMs的所有随机效应类型，给出该模型相应的一种参数估计方法，即线性逼近—逐步2次规划算法。3）发现并通过实例说明SAS中mixed模块不能保证随机效应参数方差为非负定矩阵，所以在本文提出采用线性逼近—逐步2次规划算法，同时给出十种方差类型满足正定或半正定的条件，因而算法可以保证不出现类似SAS中的错误。4）模型的标准表达式可以处理固定效应和随机效应参数分级（即数量化问题），并且指出了S-Plus中nlme函数分级差法计算有缺陷，而线性逼近—逐步2次规划算法克服了此问题。5）首次利用带有交互作用的两因素（林分密度和地位级指数）NLMEMs分析了落叶松树高—直径模型。在此基础上进一步分析了随机效应与海拔的关系。
通过本项研究，可以得出以下主要结论：1）本研究提出了一种正态NLMEMs的标准表达式，它包含多种类型的非线性混合效应模型（随机效应参数服从正态分布），例如单水平NLMEMs、逐级嵌套多水平NLMEMs、只含主效应的多因素NLMEMs、包括主效应和交互效应的NLMEMs以及某几种类型组合的一般性NLMEMs等。模型中固定效应参数和随机效应参数可以考虑分级（即数量化）。同时还把正态NLMEMs的标准表达式推广到参数方差与某些因素（称为组变量）有关的NLMEM（s考虑组变量的NLMEMs）。因此模型比传统的NLMEMs表达式更为一般化，具有更广的用途。2）在计算单水平或逐级嵌套的多水平NLMEMs时，FOCE-EM算法与SAS及S-Plus提供的Lindstrom andBates（LB）算法计算精度非常接近，数值实例表明，至少四位小数相同。但FOCE-EM算法从理论上保证了线性步计算收敛，从而使得该算法计算收敛性明显好于LB算法。3）本研究给出线性逼近—逐步2次规划算法计算一般类型的NLMEMs。为保证方差非负定采用逐步2次规划是必要的。该算法同样也能计算考虑组变量的NLMEMs。4）从计算速度上讲，FOCE-EM算法比线性逼近—逐步2次规划算法要快，因此建议在计算单水平和逐级嵌套多水平NLMEMs时，使用FOCE-EM算法，而计算其它类型的NLMEMs时，使用线性逼近—逐步2次规划算法。5)通过两因素NLMEMs对落叶松树高—直径模型研究得出，考虑林分密度与地位级的交互作用能明显提高模型的预测精度，而且把海拔高度作为组变量时，还可以进一步提高模型预测精度。
Nonlinear mixed effects models (NLMEMs) are built on the regression function, whichrelying on the nonlinear relationship of fixed effects parameters and random effects parameters.They are modern statistical methods to analyze longitudinal data, multilevel data and repeatedsurvey data. They also can be considered as powerful statistical tools, which not only reflectthe overall variation trend of population, but also describe individual differences. In recentlyyears, NLMEMs are getting much attention from the scholars in various disciplines, such asmedicine, engineering, agriculture and forestry and so on. After30years’ development,several scholars have proposed the single level NLMEMs and the multi-level NLMEMs, andtypical statistical software such as SAS and S-Plus were developed to estimate them. Inpractice, we found two main problems in existing NLMEMs and estimation methods. One isthat the main software such as SAS and S-Plus may not always reach convergence, especiallywhen the model contains many estimable parameters. The other is that the current NLMEMsmay not cover all types of random effects combinations (such as interaction that is commonlyused in forestry), which limits application of NLMEMs. The objectives of this paper aretherefore to (1) propose an algorithm with good convergence to estimate single level andmulti-level NLMEMs,(2) put forward a unified expression that includes all types of randomeffects of NLMEMs and corresponding parameter estimation method,(3) complete and realizethe coding for the above two algorithm in ForStat software,(4) and use our programs toresolve a practical problem in forestry, which is not resolved yet.
The study has successfully achieved the above four objectives. Specific contents of thepaper are as follows:1) Based on the first-order conditional expectation linearization–expectation maximation (FOCE-EM), it draws calculation formulas and processes of singlelevel and multi-level NLMEMs.2) The study proposed a normal standard expression forNLNEMs, which contains all types of random effects, and put forward the relevant parameterestimation method, naming linearization approximation-sequential quadratic programming algorithm.3) We found the mixed procedure in SAS may not guarantee the variance ofrandom effects parameters as non-negative definite, so we used one specific case to verify thisdrawback. Therefore, we proposed the linearization approximation-sequential quadraticprogramming algorithm to estimate the parameters in the NLMEMs standard expression, andwe present ten common variance types of constraint conditions to guarantee each variancepositive definite or half positive definite.4) The standard expression of this model can dealwith classification (namely quantification problem) of both fixed effects and random effectsparameters. We also pointed out some defects of classification differential method in nlmefunction of S-Plus, and put forward the solutions such as linearization approximation-sequential quadratic programming algorithm to overcome this drawback5) It is the first timethat we used two factors (stand density and site index class) NLMEMs within interactions toanalyze the height–diameter model for Larix olgensis. We also further analyzed therelationship of random effects and altitude based on the two factors NLMEMs.
Based on the analysis and results from the study, the following conclusions can be drawn:1) the proposed normal NLMEMs standard expression contain various types of NLMEMs(random effects is normal distribution), such as single level NLMEMs, grading multi-levelNLMEMs, multi-factor NLMEMs with main effects, multi-factor NLMEMs with main effectsand interaction effects, and general NLMEMs that combined several types of random effects..We also spread the normal NLMEMs standard expression to those NLMEMs that thevariances of random effects are related to some factors (also called as group variable). Ascompared with the traditional models, the proposed models are more general and useful.2)When calculating single level or grading multi-level NLMEMs, the FOCE-EM algorithm hasa similar computational accuracy with the LB algorithm of SAS and S-Plus. However,FOCE-EM algorithm in theory ensures the linear step convergence, therefore, the convergenceof this algorithm is significantly better than LB algorithm.3) The study proposed linearizationapproximation-sequential quadratic programming algorithm to estimate the general types ofNLMEMs. It is necessary to choose this algorithm to guarantee the random effects variance asnon-negative definite.4) The computational speed of FOCE-EM algorithm is more quickly than linearization approximation-sequential quadratic programming algorithm. Therefore, it isbetter to use the FOCE-EM algorithm to analyze single level and multi-level NLMEMs. Forother types of NLMEMs, linearization approximation-sequential quadratic programmingalgorithm would be better.5) Through using the two factor NLMEMs to analyze the height-diameter model for Larix olgensis, it can be drawn that the prediction precision can beobviously increased by considering the interactions of stand density and site index class, and ifthe altitude was considered as group variable in the mixed model, the prediction precisionwould be further improved.
引文
Adame P, Rio M D, Canellas I. A Mixed Nonlinear Height-Diameter Model for Pyrenean oak(Quercuspyrenaica Willd.). Forest Ecology and Management,2008a,256:88~98.
Bates D M, Watts D G. Nonlinear Regression Analysis and Its Applications. New York: Wiley,1988.
Bates D, Maechler M, Dai B. Lme4: Linear Mixed-effects Models Using S4Classes. R package version0.999375-22. http:lme4.r-forge.r-project.org,2008.
Baum L E, Petrie T, Soules G, Weiss N. A Maximization Technique Occurring in the Statistical Analysis ofProbabilistic Functions of Markov Chains. The Annals of Mathematical Statistics,1970,41,164~171.
Beal S L, Sheiner L B. NONMEM User's Guide. University of California, San Francisco, NONMEMProject Group,1992.
Beal S, Sheiner L B. Estimating Population Kinetics. Critical Reviews in Biomedical Engineering.1982,8:94~128.
Beal S, Sheiner L. The NONMEM System. American Statistician,1980,34:118~119.
Becker M P, Yang I, Lange K. EM Algorithms without Missing Data. Statistical Methods in MedicalResearch,1997,6:38~54.
Berkey C S, Laird N M. Nonlinear Growth Curve Analyses: Estimating Population Parameters. Annals ofHuman Biology,1986,13:111~128.
Booth J G, Hobert J P. Maximizing Generalized Linear Mixed Model Likelihoods with an Automated MonteCarlo EM algorithm. Journal of the Royal Statistical Society, Series B (Statistical Methodology),1999,61:265~285.
Breslow N E, Clayton D G. Approximate Inference in Generalized Linear Mixed Models. Journal of theAmerican Statistical Associaation,1993,88:9~25.
Breslow N E, Lin X. Bias Correction in Generalized Linear Mixed Models with a Single Components ofDispersion. Biometrika,1995,82:81~91.
Browne W J, and Draper D. A comparison of Bayesian and likelihood-based methods for fitting multilevelmodels. Bayesian Analysis,2006,1:473~550.
Budhathoki C B, Lynch T B, Guldin J M. February28-March4. Individual Tree Growth Models for NaturalEven-Aged Shortleaf (Pinus echinata Mill.). In:Proceedings of the13th Biennial SouthernSilvicultural Conference,Memphis,TN,2005.
Budhathoki C B, Lynch T B, Guldin J M. Nonlinear mixed modeling of basal area growth for shortleaf pine.Forest Ecology and Management,2008,255,3440~3446.
Budhathoki, Chakra B, Lynch T B, Guldin J M. Nonlinear Mixed Modeling of Basal Area Growth forShortleaf Pine. Forest Science,2008,255:3440~3446.
Calama R, Montero G. Interregional Nonlinear Height-Diameter Model with Random Coefficients for StonePine in Spain.Canadian Journal of Forest Research,2004,34:150~163.
Calegario N, Daniels R F, Maestri R, Neiva R. Modeling Dominant Height Growth Based on NonlinearMixed-effects Model:A Clonal Eucalyptus Plantation Case Study. Forest Ecology and Management,2005,204:11~21.
Clutter J L. The development of compatible analytic models for growth and yield of Loblolly pine. Ph. D.dissertation, Duke University,1961.
Cudeck R, Du Toit S H C. Nonlinear Multilevel Models for Repeated Measures Data. In S. P.Resie andN.Dunn(Eds.),Multilevel Modeling:Methodogical Advances, Issues and Applications (pp.1-24).Mahwah, NJ:Erlbaum,2003.
Davidian M, Gallant A R. The Nonlinear Mixed Effects Model with a Smooth Random EffectsDensity.Biometrika,1993,80:475~488.
Davidian M, Gallant, A R. Smooth Nonparametric Maximum Likelihood Estimation for PopulationPharmacokinetics, with Application to Quinidine. Journal of Pharmacokinetics and Biopharmaceutics,1992,20:529~556.
Davidian M, Giltiman D M. Nonlinear models for repeated measurement data: an overview and update.Journal of Agricultural, Biological, and Environmental Statistics,2003,8:387~419.
Davidian M, Giltinan D M. Nonlinear Models for Repeated Measurement Data. Chapman&Hall, New York,1995.
Davis P J, Rabinowitz P. Methods of Numerical Integration,2nd ed., Academic press, New York.1984.
Dempster A P, Laird N, Rubin D B. Maximum Likelihood from Incomplete Data Via the EM Algorithm.Journal of the Royal Statistical Society, Series B (Statistical Methodology),1977,39:1~38.
Dorado F C, Ulises D A, Marcos B A, Marina S R, Klaus V G. A Generalized Height-Diameter ModelIncluding Random Components for Radiata Pine Plantations in Northwestern Spain. Forest Ecologyand Management,2006,229:202~213.
Draper N R, Smith H. Applied Regression Analysis,2nd edition.1981, New York: Wiley.
Eerik inen K. Stem volume models with random coefficients for Pinus kesiya in Tanzania, Zambia, andZimbabwe. Canadian journal of forest research,2001,31:879~888.
Fang Z, Bailey R L. Nonlinear Mixed Effects Modeling for Slash Pine Dominant Height Growth FollowingIntensive Silvicultural Treatments. Forest Science,2001a,47:287~300.
Fehrmann L, Lehtonen A, Kleinn C, Tomppo R. Comparison of Linear and Mixed-Effect RegressionModels and a K-nearest Neighbor Approach for Estimation of Single-tree Biomass. Canadian journalof forest research,2008,38:1~9.
Garber S M, Maguire D A. Modeling Stem Taper of Three Central Oregon Species Using Nonlinear MixedEffects Models and Autoregressive Error Structures. Forest Ecology and Management,2003,179:507~522.
Geyer C J, Thompson E A. Constrained Monte Carlo Maximum Likelihood for Dependent Data. Journal ofthe Royal Statistical Society, Series B (Statistical Methodology),1992,54:657~699.
Gregoire T G, Schabenberger O, Barrett J P. Linear modeling of irregularly spaced, unbalanced, longitudinaldata from permanent plot measurements. Canadian journal of forest research,1995,25:137~156.
Gregorie T G, Schabenberger O. A Nonlinear Mixed-effects Model to Predict Cumulative Bole Volume ofStanding Trees.Journal of Apply Statistics,1996,23:257~271.
Hall D B, Bailey R L. Modeling and Prediction of Forest Growth Variables Based on Multilevel NonlinearMixed Models. Forest Science,2001,47(3):311~321.
Hall D B, Clutter M. Multivariate multilevel nonlinear mixed effects models for timber yield predictions.Biometrics,2004,60:16~24.
Harring J R. Nonlinear mixed effects mixture model:A model for clustering nonlinear longitudinal profiles.The university of Minnesota PHD thesis,2005,6~26.
Hartley H O. Maximum likelihood estimation from incomplete data. Biometrics,1958,14,174~194.
Harville D A. Maximum likelihood approaches to variance component estimation and to related problems.Journal of the American Statistical Association,1997,72:320~340.
Henderson, C R. Applications of Linear Models in Animal Breeding. University of Guelph,1984.
Huang S, Morgan D J, Klappstein G, Heidt J, Yang Y, Greidanus G. A Growth and Yield ProjectionSystem(GYPSY)for Natural and Regenerated Lodgepole Pine Stands in Alberta. Alberta SustainableResource Development Tech. Rep.Publ. T/485.193P,2001.
Jin Z. Some statistical methods for analysis of Nonlinear mixed effects models. The Columbia universityPHD thesis,1998,12~32.
Jordan L, Daniels R F, Clark A, He R. Multilevel Nonlinear Mixed-Effects Models for the Modeling ofEarlywood and Latewood Microfibril Angle. Forest Science,2005,51(4):357~371.
Ke C, Wang Y. Semiparametric Nonlinear Mixed-Effects Models and Their Applications. Journal of theAmerican Statistical Association,2001,96,1272~1281.
Kuhn E, Lavielle M. Maximum Likelihood Estimation in Nonlinear Mixed Effects Models. ComputationalStatistric and Data Analysis,2005,49:1020~1038.
Laird N M, Ware J H. Random effects models for longitudinal data. Biometrics,1982,38:963~974.
Laird N, Lange, N. Stram, D. Maximum Likelihood Computations with Repeated Measures: Application ofthe EM algorithm. Journal of the American Statistical Association,1987,82:97~105.
LANDAU, S.and EVERITT, B.S. A Handbook of Statistical Analyses using SPSS. Boca Raton: Chapman&Hall/Crc,2004.
Lappi J, Bailey R L. A Height Prediction Model with Random Stand and Tree Parameters:An Alternative toTraditional Site Methods. Forest Science,1988,34:907~927.
Lappi J, Bailey R L. A Height Prediction Model with Random Stand and Tree Parameters: An Alternative toTraditional Site Methods. Forest Science,1991,34:907~927.
Lappi J. A Longitudinal Analysis of Height/Diameter Curves. Forest Science,1997,43:555~570.
Leites L P, Robinson A P. Improving Taper Equations of Loblolly Pine with Crown Dimensions in aMixed-Effects Modeling Framework. Forest Science,2004,50:204~212.
Lewis J, Richard F D, Clark A, He R. Multievel Nonlinear Mixed-effects Models for the Modeling ofEarlywood and Latewood Microfibril Angle. Forest Science,2005,51:357~371.
Li R, Weiskittel A, Kershaw Jr J. Modeling annualized occurrence, frequency, andcomposition of ingrowthusing mixed-effects zeroiflated models and permanent plots in the acadian Forest region of NorthAmerica. Canadian journal of forest research,2011,41:2077~2089.
Lin X, Breslow N E. Bias Correction in Generalized Linear Mixed Models with Multiple Components ofDispersion. Journal of the American Statistical Association,1996,91,1007~1016.
Lindstrom M J, Bates D M. Nonlinear Mixed Effects Models for Repeated Measures Data. Biometrics,1990,46:673~687.
Lindstrom M J, Bates, D M. Newton-Raphson and EM Algorithms for Linear Mixed-Effects Models forRepeated-Measures Data (corr:94v89p1572). Journal of the American Statistical Association,1988,83:1014~1022.
Littell R C, Milliken G A, Stroup W W, Wolfinger R D, Schabenberber, O. SAS for Mixed Models.2nd ed.,SAS Institute Inc, Cary, NC,2006.
Littell R C, Milliken G A, Stroup W W, Wolfinger R D. SAS System for Mixed Models. SAS Institute Inc.,Cary, NC,1996.
Little R J A, Rubin D B. Statistical Analysis with Missing Data (2nd Ed). John Wiley and Sons, Inc.,Hoboken, New Jersey,2002.
Lynch T B, Holley A G, Stevenson D J. A random-parameter height-diameter model for Cherry-bark Oak.Southern journal of applied forest,2005,29:22~26.
Mahadev S, John P. Height-Diameter Equations for Boreal Tree Species in Ontario Using a Mixed-EffectsModeling Approach. Forest Ecology and Management,2007,249:187~198.
Mallet A, Mentre F, Steimer J L, Lokiec F. Nonparameteric Maximum Likelihood Estimation for PopulationPharmacokinetics, with Application to Cyclosporine. Journal of Pharmacokinetics andBiopharmaceutics,1988,16:311~327.
McCulloch C E. Maximum Likelihood Algorithms for Generalized Linear Mixed Models. Journal of theAmerican Statistical Association,1997,92:162~170.
Meht talo L. A Longitudinal Height-Diameter Model for Norway Spruce in Finland. Canadian Journal ofForest Research,2004,34:131~140.
Meng S X, Huang S. Improved Calibration of Nonlinear Mixed-Effects Models Demonstrated on a HeightGrowth Function. Forest Science,2009,55(3):239~248.
Meng X L, Pedlow S. EM: a bibliographic review with missing articles, proc. Statist. Computing Sec.Journal of the American Statistical Association,1992,24~27.
Nanos N, Rafael C, Gregorio M, Luis G. Geostatistical Prediction of Height-Diameter Models. ForestEcology and Management,2004,195:221~235.
Nerman O, Olsson M. Fitting Phase-type distributions via the EM algorithm. Scandinavian Journal ofStatistics.1996.23(4):419~441
Orchard T, Woodbury M A. A missing information principle:theory and applications. Proceedings of theSixth Berkeley Symposium on Mathematical Statistics and Probability,1972,1:697~715.
Pandhard X, Samson A. Extension of the SAEM algorithm for nonlinear mixed models with2levels ofrandom effects. Biostatistics,2009,10(1):121~135.
Pinheiro J C, Bates D M. Approximations to the Loglikelihood Function in the Nonlinear Mixed EffectsModel. Journal of Computational and Graphical Statistics,1995,4:12~35.
Pinherio J C, Bates, D M. Mixed-Effects Models in S and S-PLUS. Spring-Verlag, New York, NY,2000.
Rao C R. Linear Statistical Inference and its Application, John Wiley and Sons,1973.
Reader A J. EM algorithm system modeling by image-space techniques for PET reconstruction. NuclearScience.2003.50(5):1392-1397.
Sakamoto Y, Ishiguro M, Kitagawa G. Akaike Information Criterion Statistics, Reidel, Dordrecht,Holland,1986.
Schwarz G. Estimating the dimension of a model. Annals of Statistics,1978.6:461–464.
Sharma M J. Parton.Height-Diameter Equations for Boreal Tree Species in Ontario Using a mixed-effectsModeling Approach. Forest Ecology and Management,2007,249:187~198.
Sheiner L B, Beal S L. Evaluation of Methods for Estimating Population Pharmacokinetic Parameters. I.Michaelis-Menton Model Routine Clinical pharmacokinetic Data. Journal of Pharmacokinetics andBiopharmaceutics,1980,8:553~571.
Shih M C. Estimation in Nonlinear mixed effects models: parametric and nonparametric approaches.Stanford University PHD thesis.1999,20~21.
Skaug H, Fournier D. Random effects in AD Model Builder(ADMB-RE) User Guide.Otter Resrarch Ltd,Canada,2006.
Strathe A B, Danfaer A, S rensen, and Kebreab E. A Multilevel Nonlinear Mixed-Effects Approach toModel Growth in Pigs. Journal of Animal Science,2010,88:638~649.
Sundberg R. Maximum likelihood theory for incomplete data from an exponential family. ScandinavianJournal of Statistics,1974,1:49~58.
Swallow W H, Monahan J F. Monte Carlo Comparison of AVOVA, MIVQUE, REML, and ML estimatorsof variance components. Technometrics,1984,28:47~57.
Tassisa G, Burkhart H E. An Application of Mixed Effects Analysis to Modeling Thinning Effects on StermProfile of Loblolly Pine. Forest Ecology and Management,1998,103:87~101.
The MIXED/NLMIXED Procedure. SAS/STAT User’s Guide,SAS On-Line Document. SAS Institute Inc.1999.
Train K E. EM algorithms for nonparameteric estimation of mixing distributions. Journal of ChoiceModelling.2008.1(1):40~69.
Wang Y, Valerie M L, Thomas G B. Modelling and Prediction of Dominant Height and Site Index ofEucalyputs globules Plantations Using a Nonlinear Mixed-Effects Model Approach. Canadian Journalof Forest Research,2007,37:1390~1403
Wlofinger R. Laplace’s Approximation for Nonlinear Mixed Models. Biometrika,1993,80:791~795.
Wolfinger R D, Lin X. Two Taylor-Series Approximation Methods for Nonlinear Mixed Models.Computational Statistric and Data Analysis,1997,25:465~490.
Wolfinger R, Tobias R, Sall J. Computing Gaussian likelihoods and their derivatives for general linearmixed models. Journal on Scientific Computing,1994,15:1294~1310.
Wu C F J. On the convergence properties of the EM algorithm. Annals of Statistics,1983,11:95~103.
Yang Y, Huang S. Comparison of different methods for fitting nonlinear mixed forest models and formaking predictions. Canadian journal of forest research,2011,41(8):1671~1686
Zhang P. Nonparametric Importance Sampling. Journal of the American Statistical Association,1996,91:1245~1253.
Zhang Y J, Borders B E. Using a System Mixed-Effects Modeling Method to Estimate Tree CompartmentBiomass for Intensively Managed Loblolly Pines—an Allometric Approach. Forest Ecology andManagement,2004,194:145~157.
Zhao D, Wilson M, Borders, B E. Modeling response curves and testing treatment effects in repeatedmeasures experiments: A multilevel nonlinear mixed-effects model approach. Canadian journal offorest research,2005,35:122~132.
陈宝林.最优化理论与算法.北京:清华大学出版社,2002.
杜纪山,唐守正,王洪良.天然林区小班森林资源数据的更新模型.林业科学,2000,36(2):26~32
符利勇,李永慈,李春明,唐守正.基于两水平非线性混合模型对杉木林优势高生长量研究.林业科学研究,2011,24(6):720~726.
符利勇，李永慈，李春明，唐守正.利用2种非线性混合效应模型（2水平）对杉木林胸径生长量的分析.林业科学,2012,48(5),36~43.
李春明，唐守正.基于非线性混合模型的落叶松云纱林分断面积模型.林业科学，2010,46（7）：106~113.
李春明，张会儒．利用非线性混合模型模拟杉木林优势木平均高．林业科学,2010,46(3):89~95
李永慈，唐守正．用Mixed和Nlmixed过程建立混合生长模型．林业科学研究,2004,17(3):297~283
唐守正,李勇.生物数学模型的统计学基础.北京:科学出版社,2002.
唐守正，郎奎建，李海奎.统计和生物数学模型计算.北京:科学出版社，2008
Bates D M, Watts D G. Nonlinear Regression Analysis and Its Applications. New York: Wiley,1988.
Bates D, Maechler M, Dai B. Lme4: Linear Mixed-effects Models Using S4Classes. R package version0.999375-22. http:lme4.r-forge.r-project.org,2008.
Baum L E, Petrie T, Soules G, Weiss N. A Maximization Technique Occurring in the Statistical Analysis ofProbabilistic Functions of Markov Chains. The Annals of Mathematical Statistics,1970,41,164~171.
Beal S L, Sheiner L B. NONMEM User's Guide. University of California, San Francisco, NONMEMProject Group,1992.
Beal S, Sheiner L B. Estimating Population Kinetics. Critical Reviews in Biomedical Engineering.1982,8:94~128.
Beal S, Sheiner L. The NONMEM System. American Statistician,1980,34:118~119.
Becker M P, Yang I, Lange K. EM Algorithms without Missing Data. Statistical Methods in MedicalResearch,1997,6:38~54.
Berkey C S, Laird N M. Nonlinear Growth Curve Analyses: Estimating Population Parameters. Annals ofHuman Biology,1986,13:111~128.
Booth J G, Hobert J P. Maximizing Generalized Linear Mixed Model Likelihoods with an Automated MonteCarlo EM algorithm. Journal of the Royal Statistical Society, Series B (Statistical Methodology),1999,61:265~285.
Breslow N E, Clayton D G. Approximate Inference in Generalized Linear Mixed Models. Journal of theAmerican Statistical Associaation,1993,88:9~25.
Breslow N E, Lin X. Bias Correction in Generalized Linear Mixed Models with a Single Components ofDispersion. Biometrika,1995,82:81~91.
Browne W J, and Draper D. A comparison of Bayesian and likelihood-based methods for fitting multilevelmodels. Bayesian Analysis,2006,1:473~550.
Budhathoki C B, Lynch T B, Guldin J M. February28-March4. Individual Tree Growth Models for NaturalEven-Aged Shortleaf (Pinus echinata Mill.). In:Proceedings of the13th Biennial SouthernSilvicultural Conference,Memphis,TN,2005.
Budhathoki C B, Lynch T B, Guldin J M. Nonlinear mixed modeling of basal area growth for shortleaf pine.Forest Ecology and Management,2008,255,3440~3446.
Budhathoki, Chakra B, Lynch T B, Guldin J M. Nonlinear Mixed Modeling of Basal Area Growth forShortleaf Pine. Forest Science,2008,255:3440~3446.
Calama R, Montero G. Interregional Nonlinear Height-Diameter Model with Random Coefficients for StonePine in Spain.Canadian Journal of Forest Research,2004,34:150~163.
Calegario N, Daniels R F, Maestri R, Neiva R. Modeling Dominant Height Growth Based on NonlinearMixed-effects Model:A Clonal Eucalyptus Plantation Case Study. Forest Ecology and Management,2005,204:11~21.
Clutter J L. The development of compatible analytic models for growth and yield of Loblolly pine. Ph. D.dissertation, Duke University,1961.
Cudeck R, Du Toit S H C. Nonlinear Multilevel Models for Repeated Measures Data. In S. P.Resie andN.Dunn(Eds.),Multilevel Modeling:Methodogical Advances, Issues and Applications (pp.1-24).Mahwah, NJ:Erlbaum,2003.
Davidian M, Gallant A R. The Nonlinear Mixed Effects Model with a Smooth Random EffectsDensity.Biometrika,1993,80:475~488.
Davidian M, Gallant, A R. Smooth Nonparametric Maximum Likelihood Estimation for PopulationPharmacokinetics, with Application to Quinidine. Journal of Pharmacokinetics and Biopharmaceutics,1992,20:529~556.
Davidian M, Giltiman D M. Nonlinear models for repeated measurement data: an overview and update.Journal of Agricultural, Biological, and Environmental Statistics,2003,8:387~419.
Davidian M, Giltinan D M. Nonlinear Models for Repeated Measurement Data. Chapman&Hall, New York,1995.
Davis P J, Rabinowitz P. Methods of Numerical Integration,2nd ed., Academic press, New York.1984.
Dempster A P, Laird N, Rubin D B. Maximum Likelihood from Incomplete Data Via the EM Algorithm.Journal of the Royal Statistical Society, Series B (Statistical Methodology),1977,39:1~38.
Dorado F C, Ulises D A, Marcos B A, Marina S R, Klaus V G. A Generalized Height-Diameter ModelIncluding Random Components for Radiata Pine Plantations in Northwestern Spain. Forest Ecologyand Management,2006,229:202~213.
Draper N R, Smith H. Applied Regression Analysis,2nd edition.1981, New York: Wiley.
Eerik inen K. Stem volume models with random coefficients for Pinus kesiya in Tanzania, Zambia, andZimbabwe. Canadian journal of forest research,2001,31:879~888.
Fang Z, Bailey R L. Nonlinear Mixed Effects Modeling for Slash Pine Dominant Height Growth FollowingIntensive Silvicultural Treatments. Forest Science,2001a,47:287~300.
Fehrmann L, Lehtonen A, Kleinn C, Tomppo R. Comparison of Linear and Mixed-Effect RegressionModels and a K-nearest Neighbor Approach for Estimation of Single-tree Biomass. Canadian journalof forest research,2008,38:1~9.
Garber S M, Maguire D A. Modeling Stem Taper of Three Central Oregon Species Using Nonlinear MixedEffects Models and Autoregressive Error Structures. Forest Ecology and Management,2003,179:507~522.
Geyer C J, Thompson E A. Constrained Monte Carlo Maximum Likelihood for Dependent Data. Journal ofthe Royal Statistical Society, Series B (Statistical Methodology),1992,54:657~699.
Gregoire T G, Schabenberger O, Barrett J P. Linear modeling of irregularly spaced, unbalanced, longitudinaldata from permanent plot measurements. Canadian journal of forest research,1995,25:137~156.
Gregorie T G, Schabenberger O. A Nonlinear Mixed-effects Model to Predict Cumulative Bole Volume ofStanding Trees.Journal of Apply Statistics,1996,23:257~271.
Hall D B, Bailey R L. Modeling and Prediction of Forest Growth Variables Based on Multilevel NonlinearMixed Models. Forest Science,2001,47(3):311~321.
Hall D B, Clutter M. Multivariate multilevel nonlinear mixed effects models for timber yield predictions.Biometrics,2004,60:16~24.
Harring J R. Nonlinear mixed effects mixture model:A model for clustering nonlinear longitudinal profiles.The university of Minnesota PHD thesis,2005,6~26.
Hartley H O. Maximum likelihood estimation from incomplete data. Biometrics,1958,14,174~194.
Harville D A. Maximum likelihood approaches to variance component estimation and to related problems.Journal of the American Statistical Association,1997,72:320~340.
Henderson, C R. Applications of Linear Models in Animal Breeding. University of Guelph,1984.
Huang S, Morgan D J, Klappstein G, Heidt J, Yang Y, Greidanus G. A Growth and Yield ProjectionSystem(GYPSY)for Natural and Regenerated Lodgepole Pine Stands in Alberta. Alberta SustainableResource Development Tech. Rep.Publ. T/485.193P,2001.
Jin Z. Some statistical methods for analysis of Nonlinear mixed effects models. The Columbia universityPHD thesis,1998,12~32.
Jordan L, Daniels R F, Clark A, He R. Multilevel Nonlinear Mixed-Effects Models for the Modeling ofEarlywood and Latewood Microfibril Angle. Forest Science,2005,51(4):357~371.
Ke C, Wang Y. Semiparametric Nonlinear Mixed-Effects Models and Their Applications. Journal of theAmerican Statistical Association,2001,96,1272~1281.
Kuhn E, Lavielle M. Maximum Likelihood Estimation in Nonlinear Mixed Effects Models. ComputationalStatistric and Data Analysis,2005,49:1020~1038.
Laird N M, Ware J H. Random effects models for longitudinal data. Biometrics,1982,38:963~974.
Laird N, Lange, N. Stram, D. Maximum Likelihood Computations with Repeated Measures: Application ofthe EM algorithm. Journal of the American Statistical Association,1987,82:97~105.
LANDAU, S.and EVERITT, B.S. A Handbook of Statistical Analyses using SPSS. Boca Raton: Chapman&Hall/Crc,2004.
Lappi J, Bailey R L. A Height Prediction Model with Random Stand and Tree Parameters:An Alternative toTraditional Site Methods. Forest Science,1988,34:907~927.
Lappi J, Bailey R L. A Height Prediction Model with Random Stand and Tree Parameters: An Alternative toTraditional Site Methods. Forest Science,1991,34:907~927.
Lappi J. A Longitudinal Analysis of Height/Diameter Curves. Forest Science,1997,43:555~570.
Leites L P, Robinson A P. Improving Taper Equations of Loblolly Pine with Crown Dimensions in aMixed-Effects Modeling Framework. Forest Science,2004,50:204~212.
Lewis J, Richard F D, Clark A, He R. Multievel Nonlinear Mixed-effects Models for the Modeling ofEarlywood and Latewood Microfibril Angle. Forest Science,2005,51:357~371.
Li R, Weiskittel A, Kershaw Jr J. Modeling annualized occurrence, frequency, andcomposition of ingrowthusing mixed-effects zeroiflated models and permanent plots in the acadian Forest region of NorthAmerica. Canadian journal of forest research,2011,41:2077~2089.
Lin X, Breslow N E. Bias Correction in Generalized Linear Mixed Models with Multiple Components ofDispersion. Journal of the American Statistical Association,1996,91,1007~1016.
Lindstrom M J, Bates D M. Nonlinear Mixed Effects Models for Repeated Measures Data. Biometrics,1990,46:673~687.
Lindstrom M J, Bates, D M. Newton-Raphson and EM Algorithms for Linear Mixed-Effects Models forRepeated-Measures Data (corr:94v89p1572). Journal of the American Statistical Association,1988,83:1014~1022.
Littell R C, Milliken G A, Stroup W W, Wolfinger R D, Schabenberber, O. SAS for Mixed Models.2nd ed.,SAS Institute Inc, Cary, NC,2006.
Littell R C, Milliken G A, Stroup W W, Wolfinger R D. SAS System for Mixed Models. SAS Institute Inc.,Cary, NC,1996.
Little R J A, Rubin D B. Statistical Analysis with Missing Data (2nd Ed). John Wiley and Sons, Inc.,Hoboken, New Jersey,2002.
Lynch T B, Holley A G, Stevenson D J. A random-parameter height-diameter model for Cherry-bark Oak.Southern journal of applied forest,2005,29:22~26.
Mahadev S, John P. Height-Diameter Equations for Boreal Tree Species in Ontario Using a Mixed-EffectsModeling Approach. Forest Ecology and Management,2007,249:187~198.
Mallet A, Mentre F, Steimer J L, Lokiec F. Nonparameteric Maximum Likelihood Estimation for PopulationPharmacokinetics, with Application to Cyclosporine. Journal of Pharmacokinetics andBiopharmaceutics,1988,16:311~327.
McCulloch C E. Maximum Likelihood Algorithms for Generalized Linear Mixed Models. Journal of theAmerican Statistical Association,1997,92:162~170.
Meht talo L. A Longitudinal Height-Diameter Model for Norway Spruce in Finland. Canadian Journal ofForest Research,2004,34:131~140.
Meng S X, Huang S. Improved Calibration of Nonlinear Mixed-Effects Models Demonstrated on a HeightGrowth Function. Forest Science,2009,55(3):239~248.
Meng X L, Pedlow S. EM: a bibliographic review with missing articles, proc. Statist. Computing Sec.Journal of the American Statistical Association,1992,24~27.
Nanos N, Rafael C, Gregorio M, Luis G. Geostatistical Prediction of Height-Diameter Models. ForestEcology and Management,2004,195:221~235.
Nerman O, Olsson M. Fitting Phase-type distributions via the EM algorithm. Scandinavian Journal ofStatistics.1996.23(4):419~441
Orchard T, Woodbury M A. A missing information principle:theory and applications. Proceedings of theSixth Berkeley Symposium on Mathematical Statistics and Probability,1972,1:697~715.
Pandhard X, Samson A. Extension of the SAEM algorithm for nonlinear mixed models with2levels ofrandom effects. Biostatistics,2009,10(1):121~135.
Pinheiro J C, Bates D M. Approximations to the Loglikelihood Function in the Nonlinear Mixed EffectsModel. Journal of Computational and Graphical Statistics,1995,4:12~35.
Pinherio J C, Bates, D M. Mixed-Effects Models in S and S-PLUS. Spring-Verlag, New York, NY,2000.
Rao C R. Linear Statistical Inference and its Application, John Wiley and Sons,1973.
Reader A J. EM algorithm system modeling by image-space techniques for PET reconstruction. NuclearScience.2003.50(5):1392-1397.
Sakamoto Y, Ishiguro M, Kitagawa G. Akaike Information Criterion Statistics, Reidel, Dordrecht,Holland,1986.
Schwarz G. Estimating the dimension of a model. Annals of Statistics,1978.6:461–464.
Sharma M J. Parton.Height-Diameter Equations for Boreal Tree Species in Ontario Using a mixed-effectsModeling Approach. Forest Ecology and Management,2007,249:187~198.
Sheiner L B, Beal S L. Evaluation of Methods for Estimating Population Pharmacokinetic Parameters. I.Michaelis-Menton Model Routine Clinical pharmacokinetic Data. Journal of Pharmacokinetics andBiopharmaceutics,1980,8:553~571.
Shih M C. Estimation in Nonlinear mixed effects models: parametric and nonparametric approaches.Stanford University PHD thesis.1999,20~21.
Skaug H, Fournier D. Random effects in AD Model Builder(ADMB-RE) User Guide.Otter Resrarch Ltd,Canada,2006.
Strathe A B, Danfaer A, S rensen, and Kebreab E. A Multilevel Nonlinear Mixed-Effects Approach toModel Growth in Pigs. Journal of Animal Science,2010,88:638~649.
Sundberg R. Maximum likelihood theory for incomplete data from an exponential family. ScandinavianJournal of Statistics,1974,1:49~58.
Swallow W H, Monahan J F. Monte Carlo Comparison of AVOVA, MIVQUE, REML, and ML estimatorsof variance components. Technometrics,1984,28:47~57.
Tassisa G, Burkhart H E. An Application of Mixed Effects Analysis to Modeling Thinning Effects on StermProfile of Loblolly Pine. Forest Ecology and Management,1998,103:87~101.
The MIXED/NLMIXED Procedure. SAS/STAT User’s Guide,SAS On-Line Document. SAS Institute Inc.1999.
Train K E. EM algorithms for nonparameteric estimation of mixing distributions. Journal of ChoiceModelling.2008.1(1):40~69.
Wang Y, Valerie M L, Thomas G B. Modelling and Prediction of Dominant Height and Site Index ofEucalyputs globules Plantations Using a Nonlinear Mixed-Effects Model Approach. Canadian Journalof Forest Research,2007,37:1390~1403
Wlofinger R. Laplace’s Approximation for Nonlinear Mixed Models. Biometrika,1993,80:791~795.
Wolfinger R D, Lin X. Two Taylor-Series Approximation Methods for Nonlinear Mixed Models.Computational Statistric and Data Analysis,1997,25:465~490.
Wolfinger R, Tobias R, Sall J. Computing Gaussian likelihoods and their derivatives for general linearmixed models. Journal on Scientific Computing,1994,15:1294~1310.
Wu C F J. On the convergence properties of the EM algorithm. Annals of Statistics,1983,11:95~103.
Yang Y, Huang S. Comparison of different methods for fitting nonlinear mixed forest models and formaking predictions. Canadian journal of forest research,2011,41(8):1671~1686
Zhang P. Nonparametric Importance Sampling. Journal of the American Statistical Association,1996,91:1245~1253.
Zhang Y J, Borders B E. Using a System Mixed-Effects Modeling Method to Estimate Tree CompartmentBiomass for Intensively Managed Loblolly Pines—an Allometric Approach. Forest Ecology andManagement,2004,194:145~157.
Zhao D, Wilson M, Borders, B E. Modeling response curves and testing treatment effects in repeatedmeasures experiments: A multilevel nonlinear mixed-effects model approach. Canadian journal offorest research,2005,35:122~132.
陈宝林.最优化理论与算法.北京:清华大学出版社,2002.
杜纪山,唐守正,王洪良.天然林区小班森林资源数据的更新模型.林业科学,2000,36(2):26~32
符利勇,李永慈,李春明,唐守正.基于两水平非线性混合模型对杉木林优势高生长量研究.林业科学研究,2011,24(6):720~726.
符利勇，李永慈，李春明，唐守正.利用2种非线性混合效应模型（2水平）对杉木林胸径生长量的分析.林业科学,2012,48(5),36~43.
李春明，唐守正.基于非线性混合模型的落叶松云纱林分断面积模型.林业科学，2010,46（7）：106~113.
李春明，张会儒．利用非线性混合模型模拟杉木林优势木平均高．林业科学,2010,46(3):89~95
李永慈，唐守正．用Mixed和Nlmixed过程建立混合生长模型．林业科学研究,2004,17(3):297~283
唐守正,李勇.生物数学模型的统计学基础.北京:科学出版社,2002.
唐守正，郎奎建，李海奎.统计和生物数学模型计算.北京:科学出版社，2008