"Forestry, Faculty of"@en . "DSpace"@en . "UBCV"@en . "The\u00CC\u0081rien, Guillaume"@en . "2011-02-10T20:20:05Z"@en . "1990"@en . "Doctor of Philosophy - PhD"@en . "University of British Columbia"@en . "Permanent sample plots have become the main source of information for estimating models which quantify the dynamic processes of a forest. Fitted models allow for projecting inventories, used to determine timber production and many forest management decisions. The quality of these models is largely dependent on the quality of the information provided by the permanent sample plots. However, the pool of information contained in recent permanent sample plots is limited. Efficient estimation techniques must use all the information available from such plots.\r\nCurrent estimation techniques can be improved. Existing techniques employed in forestry have failed to recognize the random nature of the individual model characterizing each plot. On the other hand, techniques designed for remeasured entities in other scientific fields do not address particular forestry situations such as the small number of remeasurements or the irregularity of remeasurements. A framework for estimating forestry growth models which recognizes the individuality of each plot and special forestry situations is presented in this dissertation.\r\nThe proposed framework is a two-stage estimation technique, in which the growth rate of a permanent sample plot is considered analogous to the interest rate on a bank account. The first stage estimates the growth rate after removing the time effect. The second stage, based on Von Bertalanffy's growth curve, relates growth rate to site index and volume at the beginning of the growing season. The proposed predictor of future growth rates, the \"weighted predictor,\" is a weighted average between the growth rate observed on a plot and the growth rate predicted from the second-stage model. The weighted predictor is then used to compound the current volume of a plot. An estimate of the variance of the prediction can also be computed."@en . "https://circle.library.ubc.ca/rest/handle/2429/31159?expand=metadata"@en . "G R O W T H P R E D I C T I O N OF R E C E N T P E R M A N E N T S A M P L E P L O T S F O R F O R E S T I N V E N T O R Y P R O J E C T I O N By Guillaume Therien B. Sc. App. (Foresterie) Universite Laval A THESIS SUBMITTED IN PARTIAL F U L F I L L M E N T OF T H E REQUIREMENTS FOR T H E DEGREE OF D O C T O R OF PHILOSOPHY in T H E FACULTY OF GRADUATE STUDIES FORESTRY We accept this thesis as conforming to the required standard T H E UNIVERSITY OF BRITISH COLUMBIA April 1990 \u00C2\u00A9 Guillaume Therien, 1990 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of F O R b S T ^ V The University of British Columbia Vancouver, Canada Date \\fie. 21s, mo DE-6 (2/88) A b s t r a c t Permanent sample plots have become the main source of information for estimating mod-els which quantify the dynamic processes of a forest. Fitted models allow for projecting inventories, used to determine timber production and many forest management decisions. The quality of these models is largely dependent on the quality of the information pro-vided by the permanent sample plots. However, the pool of information contained in recent permanent sample plots is limited. Efficient estimation techniques must use all the information available from such plots. Current estimation techniques can be improved. Existing techniques employed in forestry have failed to recognize the random nature of the individual model characterizing each plot. On the other hand, techniques designed for remeasured entities in other scientific fields do not address particular forestry situations such as the small number of remeasurements or the irregularity of remeasurements. A framework for estimating forestry growth models which recognizes the individuality of each plot and special forestry situations is presented in this dissertation. The proposed framework is a two-stage estimation technique, in which the growth rate of a permanent sample plot is considered analogous to the interest rate on a bank account. The first stage estimates the growth rate after removing the time effect. The second stage, based on Von Bertalanffy's growth curve, relates growth rate to site index and volume at the beginning of the growing season. The proposed predictor of future growth rates, the \"weighted predictor,\" is a weighted average between the growth rate observed on a plot and the growth rate predicted from the second-stage model. The ii weighted predictor is then used to compound the current volume of a plot. An estimate of the variance of the prediction can also be computed. m Table of Contents Abstract ii List of Tables vii List of Figures viii Acknowledgement ix 1 Introduction 1 2 Literature Review 7 2.1 Forestry Growth and Yield Literature 7 2.2 The Econometrics Literature 14 2.2.1 Estimated Generalized Least Squares 14 2.2.2 A l l Parameters Are Constant 16 2.2.3 Only the Intercepts Vary 17 2.2.3.1 The General Model 17 2.2.3.2 The Intercepts Are Fixed 18 2.2.3.3 The Intercepts Are Random 19 2.2.3.4 Recapitulation 19 2.2.4 A l l Parameters Vary 20 2.2.4.1 The General Model 20 2.2.4.2 The Parameters Are Fixed 21 2.2.4.3 The Parameters Are Random 22 iv 2.2.4.4 Recapitulation 25 2.2.5 Summary of the Econometrics Literature 26 2.3 The Growth Curve Literature 26 2.4 Conditional Densities and Expectations 28 2.5 Overview 31 3 Statistical Methods 33 3.1 The Timber Capital 33 3.1.1 The Interest Rate Analogy 34 3.1.2 Uncertainty and Growth Rate 37 3.2 The Statistical Model 38 3.2.1 Notational Conventions 39 3.2.2 First-Stage Model 40 3.2.2.1 Model Definition 40 3.2.2.2 Error Components 42 3.2.2.3 Assumptions 43 3.2.2.4 Matrix Notation ' . . 47 3.2.3 Second-Stage Model 49 3.2.3.1 Background 50 3.2.3.2 Model Definition 52 3.2.3.3 Assumptions 53 3.2.3.4 Matrix Notation 54 3.2.4 Combined Model 55 3.2.4.1 Assumption . . . 55 3.2.4.2 Matrix Notation 55 3.2.5 Summary 56 v 3.3 Parameter Estimation 56 3.3.1 Estimation When Covariance Matrices Are Known 57 3.3.1.1 Parameter Estimators 57 3.3.1.2 Estimation of the Variance of the Predictor 61 3.3.2 Covariance Matrix Estimators 64 3.3.2.1 First-Stage Estimators 64 3.3.2.2 Second-Stage Estimators 67 3.3.3 The Best Linear Unbiased Predictor and Its Variance 69 3.4 Model Predictions 69 3.4.1 Volume Predictor 70 3.4.2 Confidence Intervals Estimator 72 4 Numerical Example 75 4.1 Description of the Data Base 75 4.2 Goodness of Fit Criteria 79 4.3 Results 81 5 Discussion 89 5.1 The Basic Concept 89 5.2 The Statistical Model 91 5.3 Results 96 6 Conclusions 99 References Cited 102 v i Lis t of Tables 4.1 Time Intervals 77 4.2 Frequency of Volumes by Initial Age and Projection Time Class 78 4.3 Frequency of Volumes by Initial Age and Total Lag Class 79 4.4 Estimation of Observed Growth Rates 82 4.5 Overall Results from the First- and the Second-Stage Estimators and from the Weighted Predictor 84 4.6 Bias and (Relative Bias) by Initial Age and Projection Time Class . . . . 86 4.7 Bias and (Relative Bias) by Initial Age and Total Lag Class 86 4.8 Root Mean Squared Error and (Relative Root Mean Squared Error) by Initial Age and Projection Time Class 87 4.9 Root Mean Squared Error and (Relative Root Mean Squared Error) by Initial Age and Total Lag Class 87 4.10 Mean Absolute Deviation and (Relative Mean Absolute Deviation) by Ini-tial Age and Projection Time Class 88 4.11 Mean Absolute Deviation and (Relative Mean Absolute Deviation) by Ini-tial Age and Total Lag Class 88 vii Lis t of Figures 3.1 Growth Rate over Time for a Sigmoidal Yield Curve 41 4.2 Distribution of Site Index 76 4.3 Individual Plot Volume versus Age 76 viii Acknowledgement I would like to thank my supervisor Dr. A. Kozak as well as the rest of my supervisory committee, Dr. P. Marshall, Dr. K. lies, Dr. P. de Jong, and Dr. D. Tait for their direction and guidance. Their contributions to my thesis and forest biometrics education were enormous, as was their patience. I extend my appreciation to Dr. N. E. Heckman and Dr. J. G. Cragg for their short but invaluable contributions to my work. Funding for this research was supplied by the University of British Columbia through numerous scholarships, fellowships, teaching assistantships, research assistantships, etc. Special thanks to MacMillan Bloedel Ltd. for allowing me access to the best data base in the country. Finally, I would like to thank my peers of the first hour: Mrs Karen Jahraus, Dr. Margaret Penner, Dr. Valerie LeMay, and Dr. Jim Thrower, for their constant encour-agement and support. They initiated an instructive forum where ideas and opinions were debated freely. Many thanks to Valerie who also commented on zillions of earlier drafts, and to Margaret whose own acknowledgement \"inspired\" this one. ix Chapter 1 Introduction The main purpose of forest management is to maintain or increase the yield of the different forest resources. This goal requires a knowledge of current inventory and growth potential. Growth potential is the change in the current inventory that would occur without harvesting. However, estimating growth potential is a difficult task because the forest is a dynamic system. Quantifying the dynamic processes in a forest is essential for forecasting the future inventory. One way of quantifying the timber aspect of the dynamic processes of a forest is to use permanent sample plots. Monitoring the same piece of land over time represents a logical way to assess timber growth. In British Columbia, permanent sample plots were first established for growth and yield purposes in the late twenties and the early thirties. They became more popular in the sixties and seventies because of an increased interest in predicting yield from managed stands (Marshall and Jahraus 1987). Today, many agencies have a network of recent permanent sample plots, which have had only a few measurements. Timber dynamics can be represented by a mathematical model that ideally has the following properties: representativity, tractability and simplicity. This model is called a growth model. The model should be representative of the biological behaviour and of the stochastic nature of this behaviour; it should be tractable so that estimation is possible, and simple in its parameterization. 1 Chapter 1. Introduction 2 Modeling growth of permanent sample plots allows for the projection of sample plot inventories. Projections of temporary sample plot inventories can also be obtained. This information will be needed to estimate the timber production of a management unit and to determine an appropriate level of harvesting for this management unit. Most techniques used today to quantify inventory changes do not take advantage of all the information available from permanent sample plots. Techniques especially designed for remeasured entities are required. As most forest agencies have a large proportion of recent permanent sample plots, these techniques must be compatible with sample plots that have had only a few measurements. The objective of this dissertation is to define a modeling framework that will make use of all the information available from recent permanent sample plots in order to project an inventory with precision and accuracy. Estimates of the variability of the projections is also of special interest. For simplicity, only linear models are considered. The general linear statistical model of permanent sample plot collections can be ex-pressed by: Yij \u00E2\u0080\u0094 floij + PlijXuj + . . . + 0Kij^Kij + e;j where Yij is the dependent variable measured on individual i [i = 1, . . . ,N) at time j (j = 1,... , T); Xpij is the variable p (p = 1,..., K) measured on individual i at time j; 0pij is the parameter associated with X^j, fioij is the intercept, and e^ - is the error term on individual i at time j. Judge et al. (1985, p. 515-516) classified the parameter estimation techniques of such models in 5 broad categories: 1. all parameters are constant over individuals and time; 2. the intercept varies over individuals only, the other parameters are constant over individuals and time; Chapter 1. Introduction 3 3. the intercept varies over i n d i v i d u a l s a n d t ime, the other parameters are constant over i n d i v i d u a l s and t ime; 4. a l l parameters vary over i n d i v i d u a l s on ly ; 5. a l l parameters vary over i n d i v i d u a l s a n d t ime . In forestry, most techniques developed for f i t t ing models o n permanent sample plots have assumed constant parameters over t i m e a n d i n d i v i d u a l s . T h i s fails to recognize the i n d i v i d u a l i t y of each sample plot . O t h e r scientific fields dea l ing w i t h collections of remea-sured i nd iv idua l s have used models w i t h at least one parameter va ry ing over i nd iv idua l s . These models can \" i n d i v i d u a l i z e \" the general m o d e l . However , recent permanent sample plots have some pa r t i cu la r problems that models used i n other fields do not address. These inc lude : 1. the plots are not the same age at the first measurement; 2. the t ime in terva l between measurements is not a lways the same w i t h i n a plot; 3. plots are not remeasured at the same t ime; 4. the number of remeasurements are s m a l l . A s a general rule, specia l a t ten t ion must be g iven to the error s t ructure of the m o d e l and the assumpt ions about the covar iance m a t r i x of the r a n d o m terms i f efficient es t ima-tors are required. T h e u s u a l assumpt ions required for s imple or m u l t i p l e l inear regression are frequently not v a l i d w i t h remeasured i n d i v i d u a l s . Because each plot is remeasured at different t imes, the plo ts are often referred to as a time-series. Observat ions f rom a time-series are often correlated; th is type of corre la t ion is ca l led serial correla t ion. O b -servations measured at the same t ime can also be correlated; this type of corre la t ion is Chapter 1. Introduction 4 called contemporaneous correlation. Both types of correlation provide information about the remeasured entities needed to estimate model parameters efficiently. In this dissertation, a new framework for projecting recent permanent sample plot inventories is proposed. This new framework addresses the complexities of remeasured entities particular to forestry. The variable of interest is growth rate, that is, the change in volume per unit of initial volume.1 For short projections, an intuitive technique to predict future volumes on a given plot is to compound the current volume with the growth rate recently observed on this plot. Another technique, requiring a collection of permanent sample plots, is to regress the observed growth rates on stand attributes. The volume is then updated with the growth rate predicted from the regression. Another option, the one favoured in this thesis, is to use an average of both techniques, with each technique's growth rate estimates weighted proportionally to its reliability. The framework employed a two-stage model. In the first stage, average recent growth rates are obtained for individual plots. In the second stage, model parameters are es-timated assuming growth rates to be a function of site index and initial volume. The predictor for a sample plot is a weighted average between the growth rate estimated in the first stage and the growth rate predicted from the second-stage model. Statistically, the framework can be explained in the following way. A general model exists, explaining the growth rate from two stand attributes: site index and initial volume. However, each individual permanent sample plot follows its own model; these individual models vary randomly around the general model. The individual model can be expressed as the general model plus an error term specific to a particular plot. The actual growth 1This variable is also sometimes referred to as relative growth rate. In this thesis the term growth rate will be used. Chapter 1. Introduction 5 rate observations also vary, but around the individual model. The framework can there-fore be considered a hierarchical model: the observations are random variations of the general model's random variation. Repeated observations on the same individual and an estimate of the general model can yield a prediction of the random error term specific to a particular plot. Adding this term to the estimated general model gives a prediction of the individual model. The individual model is \"predicted\" and not \"estimated\" because the individual model is a random variation of the general model. The predicted individual model is then used to provide a prediction of the growth rate used to update the current volume. There are many reasons to prefer an estimate of the individual growth rate rather than a global estimate. The individual estimate can be used as a decision aid for plot remeasurement scheduling. It can also be employed in sampling with partial replacement to update the plots that have not been remeasured. A third possible use of individual growth rates is the ordination of the sample plots in groups to study the causes of the differences among groups. Most techniques used in forestry have failed to recognize the stochastic nature of this individual relationship. They also have failed to take into account that many permanent sample plots in North America are less than 25 years of age and have been measured less than 5 times. The proposed framework will provide a better understanding of the information given by permanent sample plots, especially recent ones, and will result in more efficient fitting techniques. The dissertation has been organized according to the following format: a brief presen-tation of the past work in growth and yield and in other fields using data sets containing remeasured individuals is discussed in Chapter 2; the statistical model, the various as-sumptions and the estimators are discussed in Chapter 3; a numerical example and results Chapter 1. Introduction 6 are presented in Chapter 4; these results, along with the statistical model, are discussed in Chapter 5, and, finally, conclusions are drawn in the last chapter. Chapter 2 Literature Review Permanent sample p lo t da ta bases are a co l lec t ion of plots remeasured a number of t imes. Remeasured da t a bases are not un ique to forestry. Research scientists i n business, economics, b io logy and medic ine have deve loped techniques cus tomized to their specific problems i n v o l v i n g entities w i t h repeated measurements over t ime . A n overview of these techniques is presented i n the fol lowing pages. A br ief h i s tory of forestry g rowth a n d y i e l d can be found i n the first sect ion. T h e second section groups the techniques developed i n econometr ics . M u c h a t ten t ion has been given i n business and economics research to da ta bases s imi la r to pe rmanent sample p lo t da ta bases. Research on b io log ica l g rowth curves is reviewed i n the t h i r d section. Because trees are b io logica l organisms, general knowledge of b io log ica l g r o w t h can help unders tand tree or s tand g rowth . F i n a l l y , the last section is a br ief presenta t ion of cond i t i ona l densities and expectat ions. Useful m a t r i x invers ion rules are also given. These basic rules w i l l help i n unders tand ing the s ta t i s t ica l m o d e l suggested i n C h a p t e r 3 to project inventory. 2.1 Forestry Growth and Yield Literature P r e d i c t i n g y i e ld , a long w i t h de te rmin ing the current inventory, is a c o m m o n forest m a n -agement ac t iv i ty . A v a i l a b i l i t y of re l iable i n fo rma t ion , c o m p u t i n g technology a n d infor-m a t i o n requirements have inf luenced the q u a l i t y of the in format ion desired f rom forest inventory project ions. Pe rmanen t sample p lo t d a t a bases have become the m a i n source 7 Chapter 2. Literature Review 8 of information for these projections. A historical perspective of the uses of permanent sample plot data as a source of information is presented in the following pages. Yield is highly correlated with a few well-known variables: species, age, ecological association, climate, edaphic conditions, among others. Nineteenth century German foresters used the variables age, species, and site class to create their normal yield curves. \"Normal\" in this case merely meant that the variability in stocking was omitted by using only fully-stocked, even-aged, monospecific temporary sample plots. At any point in time, these curves represented the upper bound or potential yield of a stand. The interpretation of predictions from these tables was left to the manager's judgment and experience. To make these tables applicable to understocked forests, Gerhardt (1930, in Davis and Johnson 1987; p. 106) developed a simple correction factor. Growth was estimated from the percentage of stocking compared to that of a fully-stocked forest, the predicted growth for a fully-stocked forest and the species. With the correction factor, it was assumed that all stands were moving toward full stocking at a rate depending on the species and the stocking. Schumacher (1939) published a simple empirical yield function where the logarithm of volume was a function of site index (S) and the reciprocal of age (1/A): InV = 30+ Pi{l/A) + 82{S) + 33(S/A) No statistical assumptions were made to develop this model; no statistical inference could be drawn from this model. The parameters were fitted using least squared error as the criterion of optimality. Schumacher's model was among the first of the mathematical expressions for yield. Buckman (1962) recognized that serially correlated data from permanent plots were providing misleading significance levels and confidence intervals when analyzed as inde-pendent observations. Chapter 2. Literature Review 9 Clutter (1963) recognized the need for compatible growth and yield models, where predicted growth is the derivative of predicted yield. He justified the use of volume on a logarithmic scale to fit his growth and yield equations by assuming a better compatibil-ity in general with \"the statistical assumptions customarily made in regression analysis (linearity, normality, additivity and homogeneity of variance).\" Leak (1966) developed a technique for providing a large-sample estimator of the over-all slope of a random coefficient model. He first assumed a linear relationship between the dependent and the independent variable for each sample plot, with the error terms inde-pendently distributed as N(0,crf), where af is the variance of an error term observed on sample plot i. These individual regressions were solved using the least squares technique. On each plot, he assumed the same number of observations and a similar distribution of the independent variable, X, around the average, X (i.e., the corrected sums of squares of the variable X were approximately equal). The expected value of any individual slope was the overall slope. This was estimated as the simple average of the individual slopes. Assuming that a? was equal among plots, the standard error of the estimated overall slope was the standard error of the mean slope. Curtis (1967) developed a method to estimate the gross yield of Douglas-fir. He recognized that the usual assumption of independence was probably violated when fit-ting equations to permanent sample plots. To handle the problem, he first fitted linear growth equations on each plot without regard to autocorrelation using the least squares technique. Then, the coefficient of determination was computed for u^_i and ui, the con-secutive estimated residuals on a plot. If the coefficient of determination r2 was small, autocorrelation was considered negligible; otherwise, the analysis was repeated using a single randomly selected point per plot. Swindel (1968) showed the theoretical bounds on the bias of the variance of parameter estimates when serially correlated data were assumed independent for multiple linear Chapter 2. Literature Review 10 models estimated with least squares. Sullivan and Reynolds (1976), using Swindel's work, provided an example for data with two measurements. Sullivan and Clutter (1972) developed growth and yield models using permanent sample plots with two measurements per plot. They assumed common parameters over plots and time, homogeneous variance within a measurement period, serial correlation p between the first and the second measurements, and no contemporaneous correlation, in that the volumes measured in a given year were uncorrelated. If = (Vn V{2) is the vector containing the first and second measurements of volume on plot i, then, for any i 7^ j , V; and Vj are uncorrelated and where p is the correlation between the first and the second measurement. They derived the maximum likelihood estimates of multiple linear models by an iterative procedure. Seegrist and Arner (1978) and Amer and Seegrist (1979, 1980) extended this method to Ferguson and Leech (1978), fitted a two-stage model drawing on the theory relating to random coefficient modeling. In the first stage, they fitted a quadratic model to each plot with the reciprocal of age as the independent variable using ordinary least squares. They assumed that the coefficients were different among plots: where V{j is the volume measured on plot i at time j, /3,-j. is the coefficient k on plot i, Aij is the age of plot i at time j, and e^ - is an error term. Using matrix notation, V is the column-vector of the logarithmic volumes; 6 is the vector of parameters to be estimated; X is the matrix of independent variables, and e is the vector of error terms. The first 1In econometrics nomenclature, this type of estimation is called classical pooling. allow multiple and varying number of remeasurements.1 ln V{j = Ao + ftiUMy) + ft2(l/Ay )2 + Chapter 2. Literature Review 11 stage can be written V = X8 + e. With the assumptions that e is distributed with mean vector zero and variance ft and is uncorrelated with /? or X , the least squares estimate of 8, 8, is unbiased with variance They then considered 8 as a random variable that could be predicted from certain independent variables (second-stage models). The set of independent variables to predict 8 do not need to be related to the set of independent variables used to predict V . They fitted three different models: to predict the intercept, the coefficient of the linear term and the coefficient of the quadratic term. The second-stage model structure was where d^i is the parameter associated with the variable I for the coefficient k\ Ziki is the corresponding explanatory variable on plot i, and u;*. is the error term. Using matrix notation, the second stage can be written Assuming that u is distributed with mean vector zero and variance S and is uncorrelated with 8 and Z, the random variable 8 is distributed with mean vector Z9 and variance S . Because 8 was considered a random variable, the distribution of 8 in the first stage was actually the conditional distribution of 8 given 8 with [x'rr 'x]- 1. a = ze + u. E\j3\8] = 8, and vai0\8) = [X'f2 _ 1 X]-- l Combining the conditional distribution of 8 given 8 with the distribution of 8, the unconditional distribution of 8 will have ze, and Chapter 2. Literature Review 12 var(/3) = [ X ' f r ' X ] \" 1 + E. Because Vt and \u00C2\u00A3 are unknown, they must be replaced by estimates. In the first-stage model, Ferguson and Leech assumed homogeneous variance within plots, no serial correlation, and heterogeneous variance among plots. In the second-stage, they assumed a homogeneous covariance matrix for the coefficients and independence between sets of coefficients measured on different plots. Wi th these assumptions they developed a feasible generalized least squares (FGLS) algorithm. They also investigated more restrictive assumptions about the error structure. Commenting on this study, Davis and West (1981) claimed that no proper confidence intervals nor test of significance can be computed with this approach. Ferguson and Leech (1981) replied that asymptotically normal statistics have proved adequate for most practical purposes using large samples. Garcia (1983) fitted a height growth model with a stochastic differential equation. Height growth was assumed to be a function of height (if): dH/dt = f(H) The function of height was assumed to be (^if* 2 \u00E2\u0080\u0094 d>3H, where i, d>2, and 3 are the function parameters. This function integrates to the well-known Chapman-Richards function (Chapman 1961, Richards 1959). The model included an error term to account for the random environmental variation in growth and the measurement error. The environmental error was assumed to be a Wiener process. This means that the variation in H accumulated over a short time interval is normally distributed with zero mean and a variance increasing with the interval length, and that errors for nonoverlapping time intervals are independent. The measurement error was assumed to be a function of observed height. The plots were assumed statistically independent. The parameters were simultaneously estimated with an iterative maximum likelihood procedure. Some parameters were global, that is identical across plots; some were local, that is site specific. Chapter 2. Literature Review 13 Gregoire (1987) applied error component modeling to permanent sample plot data. The basic model was Yij \u00E2\u0080\u0094 ^ ^dpXpij -f- Ujj p = e; + ij + Vij with XUJ defined as 1. The error was divided into three components: one for plot-dependent variation (e,), one for time-dependent variation (rj) and one, (V;J), for the variation unexplained by plot-dependent or time-dependent effects; was the dependent variable measured on plot i at time j, and X^j was the variable p measured on plot i at time j. The expected value of each component of the error term was assumed to be zero. The error term components were also assumed to be uncorrelated with the independent variables and uncorrelated with each other. Alternative error covariance structures were analyzed. When compared to the least squares technique, with the undivided error terms assumed independently and identically distributed, the results were not conclusive. In general, the error component models performed better than the least squares for the likelihood criterion but worse for the minimum prediction error criterion. Gregoire also pointed out the need to find good statistics to compare alternative models. Biging (1985) derived a site index model using random coefficient modeling theory in order to account for between-tree differences in individual tree height growth. The error structure was identical to Swamy's (1970) model, and the parameters were estimated using Maddala's (1977) estimators. Finally, Tait et al. (1988) developed a simple, two-equation growth model where the growth of the average tree was a function of the average tree volume, site quality and the number of trees per hectare, and the mortality was a function of an index of the rate of area expansion of a tree. They simultaneously estimated the best model parameters (i.e., the global parameters) with the best initial values for each plot (the local parameters). Chapter 2. Literature Review 14 This type of analysis is close to an analysis of covariance where each plot is considered a different treatment with two covariates, the initial stand density and volume. This brief historical report illustrates the change of attitude toward the information pool contained in permanent sample plot data bases. It is now recognized that an ob-servation at a point in time on a given individual provides partial information, on one hand, on the same individual at other points in time and, on the other hand, on other individuals at the same point in time. Efficient parameter estimation techniques must exploit this information. 2.2 The Econometrics Literature In econometrics literature, a data base that provides repeated measurement for each of a number of entities is called a pooled cross-sectional and time-series data base (Diel-man 1983). The term panel data is also often used. This type of data base is very common in econometrics. The abundance of information has led econometricians to de-velop techniques that take advantage of this information to explain and predict economic phenomena. Because of the similarities with forestry permanent sample plots, there has been a growing interest in forestry to adapt techniques developed in econometrics (e.g., Furnival and Wilson 1971, Ferguson and Leech 1978, Gregoire 1987, LeMay 1988). The literature about pooled cross-sectional and time-series data is extensive and only a brief overview of the major developments will be presented here. 2.2.1 Estimated Generalized Least Squares When a simple or a multiple linear model is fitted to predict or to explain the outcome of a dependent variable, the most common assumptions about the error term are zero expectation and spherical disturbance, (i.e., the covariance matrix of the error term is Chapter 2. Literature Review 15 0, Hm P[\a - a | > (] = 0. The EGLS technique is the general approach for obtaining the BLUE of the param-eters of any linear model. The OLS technique is a special case, particularly restrictive, of the EGLS technique. With pooled cross-sectional and time-series data, the assumed structure of 0 will have a supreme importance on the structure of the covariance matrix of the error terms. This is why the EGLS technique for pooled cross-sectional and time-series data has been organized in three different sections corresponding to the general assumptions about the structure of the model parameters: all parameters are constant, only the intercepts vary, and all parameters vary. 2.2.2 All Parameters Are Constant Each individual measured over time can be represented by its own model Yi = Xifr + U i Y{1 Xi = 1 \u00E2\u0080\u00A2 XKH 0iO u; = U i l YIT 1 X L I T . \u00E2\u0080\u00A2 X K i T 0iK MiT Chapter 2. Literature Review 17 When it is assumed that all individuals follow a similar behaviour, the model pa-rameters are assumed to be equal over individuals. The variation between individuals is then explained by the stochastic error term only. This is called classical pooling (CP). If the error terms on the observations from a single individual were uncorrelated, the observations from this individual would oscillate around their expected values. However, because of serial correlation, a given individual can perform constantly above or below the expected behaviour. Classical pooling is used when there is no reason to believe that the different behaviour between individuals can be explained other than by serial correlation and chance. 2.2.3 O n l y the Intercepts V a r y Sometimes, it is possible to explain partially the difference between individuals by their initial states. In these cases, it is assumed that the partial derivatives, that is the slopes of each independent variable, are constant over individuals and time, and that the intercepts will vary over individual and time. In this section, the general model describing such a situation is presented. Alternative approaches to parameter estimation depend on the assumptions about the nature of the intercepts. The intercepts can be assumed to be fixed or random. Both alternatives lead to different estimation techniques that are discussed separately. Finally, the model and the estimation techniques are summarized in a short recapitulation. 2.2.3.1 T h e G e n e r a l M o d e l In the most general case, the error term of the general linear model can be divided into three components: one to account for the difference between individuals, one to account for the difference between initial time periods, and a third for the variables unaccounted Chapter 2. Literature Review 18 for in the model that are not individual-dependent or time-dependent. This gives K Yij = B0 + ]T BpXpij + e,- + rj + Vjj P = i where Yij is measurement j on plot i , /30 is the average intercept, Xpij is the variable p measured on plot i at time j, ej is the individual error term, ij is the temporal error term, and Vij is the unexplained variation. The individual intercept is defined by combining ei and Tj with the average intercept. These different intercepts will likely arise from omission of important unobservable or unmeasurable variables in the model (Dielman 1989, p. 49). 2.2.3.2 The Intercepts Are Fixed In this case, the individual error component and the temporal error component are as-sumed to be fixed with mean zero (i.e., ^ e i = 0 a r m 2~Z r j = 0) a i m the third unexplained component is random with mean zero and variance a2. The estimation technique is called analysis of covariance ( A N C O V A ) or least squares with dummy variables. Mund-lak (1961) was among the first to use the A N C O V A for analyzing pooled cross-sectional and time-series data. There are a few disadvantages of A N C O V A . Since the different intercepts are consid-ered fixed, they are using up as many degrees of freedom. Also, the different components of the error term represent some ignorance; the model does not explain the differences be-tween initial states. Maddala (1971) suggested that this specific ignorance was not really different from the general ignorance un and should be considered in a similar way. This would also decrease the number of fixed parameters to be estimated, increasing therefore the number of degrees of freedom. This led to the assumption of random intercepts. Chapter 2. Literature Review 19 2.2.3.3 The Intercepts Are Random In this case, it is assumed that the various components of the error term are random effects with mean zero and specific variance: E[ei] = E[xA = E[v{j} = 0 m = *j E[vl] = 4 . It is also usually assumed that e;, ij and V J J are mutually independent. The estimation technique based on these assumptions is called error component modeling ( E C M ) . Wallace and Hussain (1969) suggested estimating the variance components using the residuals obtained from the OLS estimate of the individual models. These estimators are consistent when the independent variables are assumed nonstochastic. Other estimators of the error components have been developed (e.g., Amemiya 1971, Swamy and Arora 1972, Rao 1972, Fuller and Battese 1974). A l l these estimators have similar asymptotic properties. Consequently, the choice of an estimator cannot be based on its asymptotic properties. One must consider the simplicity of the estimator, finite sample results and Monte Carlo studies of small sample performance (Baltagi 1981). The finite sample properties of some estimators have been investigated analytically by Swamy and Arora (1972) , Swamy and Mehta (1979) and Taylor (1980). Monte Carlo studies were published by Arora (1973), Maddala and Mount (1973) and Baltagi (1981). Swamy and Mehta (1973) analyzed the error component model in a Bayesian framework. 2.2.3.4 Recapitulation The difference between individuals can sometimes be explained by different initial states (or model intercepts). These different intercepts come from partitioning the error term Chapter 2. Literature Review 20 into generally three components: an individual error term, a temporal error term, and an unexplained error term, and from the combination of the common intercept with the individual and the temporal error terms. If the individual and the temporal error terms are assumed to be fixed, then the ANCOVA is the appropriate EGLS technique. If the individual and temporal error terms are assumed to be random, then ECM must be used. 2.2.4 All Parameters Vary-Some time-series observations of different individuals are highly correlated with the same set of independent variables but with slopes varying over individuals. Not only are their initial states different but also the increment in the dependent variable for one unit of change in one of the independent variables. This will occur if an important variable for all the individuals is omitted from the model or if a variable important for some individuals is omitted. This section follows the same format as the previous section. The general model describing the situation is presented first. The parameter estimation technique when the parameters are assumed to be fixed comes next, followed by the technique for random parameters. Last, the model and the techniques are briefly recapitulated. 2.2.4.1 The General Model The individual model is Yi = Xidi + Ui where Yn 1 0i = Ui = U i l YiT 1 XliT \u00E2\u0080\u00A2 \u00E2\u0080\u00A2 xKiT BiK Chapter 2. Literature Review 21 Then, the full model can be written: Y = XB + u where ' Pi Y = X = P = u = XJV _0N UJV 2.2.4.2 The Parameters Are Fixed With fixed parameters, efficient parameter estimates can be obtained using OLS if there is no contemporaneous correlation. When contemporaneous correlation is present, Zellner (1962) developed a technique using this supplementary information. He called the tech-nique seemingly unrelated regressions (SUR). Estimators obtained with this technique are more efficient than simply applying least squares to each individual independently. The main assumptions about the error structure are: \u00C2\u00A3[U;J] = 0, -Efujuf] = O^IT, and X,- is a fixed matrix in repeated samples. TV is the number of individuals, K is the number of independent variables in the model, and T is the number of observations on each individual. Based on the assumptions about the error terms, the covariance matrix of the error terms can be written il = S I T where S = O i l &12 \u00C2\u00B0~21 \u00C2\u00B0~22 CT2N 0~N1 0~N2 \u00E2\u0080\u00A2\u00E2\u0080\u00A2\u00E2\u0080\u00A2 \u00C2\u00B0~NN When the contemporaneous correlation between individuals i and I is equal to zero, an = 0 for i / /, the estimates are easily shown to be identical to the OLS estimates Chapter 2. Literature Review 22 computed independently on each individual. The SUR estimator will also be the same as the OLS estimator if the X; are identical for all individuals. If u is normally distributed, and hmT_ > 0 0 [XTi _ 1 X] _ 1 is a finite nonsingular matrix, the EGLS estimate is a minimum variance unbiased estimator, and is asymptotically efficient. To estimate ft, Zellner (1962) suggested using the residuals from the regression of Y{ on either the independent variables in the ith equation alone or all the independent variables in the system. Writing the residuals from either method as hij, the elements of S (or ft) are estimated by o-n = uijUij/{T - K). Zellner (1962) also suggested iterating the estimator by recomputing the residuals using 8 to get a new estimate of ft, and so on until convergence occurred. Dhrymes (1971) showed that the iterative estimator is equivalent to the maximum likelihood estimator. Parks (1967) and Kmenta and Gilbert (1970) extended the theory to include the case where the disturbances are both serially and contemporaneously correlated. Schmidt (1977) examined the SUR model with missing observations. Swamy and Mehta (1975) also studied SUR models with missing observations but in a Bayesian framework. Phillips (1985) derived the exact finite sample properties of the SUR estimator. Zellner (1987, pp. 240-256), Srivastava (1973), and Zellner and Vandaele (1975) analyzed the SUR model in a Bayesian framework. 2.2.4.3 The Parameters Are Random The sets of different parameters can also be assumed to be random. The EGLS tech-nique for such models has been called random coefficient modelling (RCM; Swamy 1970, 1971, 1974). Because the model parameters link the same independent variables to the dependent variable, they could very well be related. If it can be assumed that many Chapter 2. Literature Review 23 independent unobservable or unmeasurable factors can influence a model parameter, this model parameter can be seen as a random variable. These excluded variables can be assumed to counteract each other, so that their overall effect is negligible. Consequently, the expected value of the model parameter is the mean parameter, and its variance is the variance of the aggregate effect of these excluded variables. The parameter 8{ can be rewritten Pi = P + e; where 0 = Po PK and e i 0 Replacing this new expression for Pi in the original individual model yields Y i = Xi{p + e;) + Ui = Xij3 + X ; e ; + U i = Xi/3 + Wi The error term w; is then the sum of two components, X i e i and Ui. The full model can be written Y = X ^ + w where (Xjei + (Xjv-Civ Swamy (1970, 1971) made the following assumptions: the number of individuals and the number of observations on each individual are larger than the number of indepen-dent variables; the U, are independently and identically distributed with -E[uij] = 0 and .E[uiU^] = oo. Carter and Yang (1986) showed that this result was also valid when N \u00E2\u0080\u0094> oo or T \u00E2\u0080\u0094> oo if an = a2 (i.e., if all the variances of the disturbance terms are equal). Amemiya (1978) investigated the case where the individual parameters are a linear function of some explanatory variables Pi = Z;0 + e i . Swamy (1974) extended the theory to include the case with contemporaneous correlation and serial correlation. Finally, there are also techniques to estimate models with different sets of parameters varying across individuals and time. These models were not considered relevant to this dissertation and are not discussed (for information concerning these models, see Judge et. al. 1985; p. 545-550). 2.2.4.4 Recapitulation The difference between individuals can sometimes be explained by different sets of model parameters. These parameters can be fixed in which case the SUR technique is the ap-propriate EGLS technique if contemporaneous correlation is present. If the parameters Chapter 2. Literature Review 2 6 are f ixed, a n d there is no contemporaneous corre la t ion, then the O L S technique is ap-propr ia te . T h e parameters can also be assumed to vary r a n d o m l y a r o u n d a c o m m o n average. In th is case, R C M is the appropr ia te E G L S technique. 2.2.5 S u m m a r y of the Econometr i c s L i t era ture T h e best, l inear , unb iased est imate ( B L U E ) of the parameters of a l inear m o d e l is given i n its general f o r m by the es t imated generalized least squares ( E G L S ) est imate. W h e n the error terms are assumed to be independent ly and iden t ica l ly d i s t r ibu ted , the E G L S est imate is equivalent to the o rd ina ry least squares ( O L S ) est imate. However , this rarely is the case w i t h models f i t ted on poo led cross-sectional and time-series da ta . T h e struc-ture of the covar iance m a t r i x of the error t e rm is h igh ly inf luenced by the a s sumpt ion about the i n d i v i d u a l m o d e l parameters. There are three general cases: the i n d i v i d u a l models share the same parameters , the i n d i v i d u a l models have different intercepts bu t c o m m o n slopes, or the i n d i v i d u a l models have different sets of parameters altogether. These different parameters can be assumed to be f ixed or r a n d o m . C l a s s i c a l p o o l i n g ( C P ) is the appropr i a t e technique i f a l l sets of parameters are assumed to be iden t i ca l across i n d i v i d u a l s a n d t ime. A n a l y s i s of covariance ( A N C O V A ) is used w h e n the inter-cepts are a s sumed to be different but f ixed, whi le error component m o d e l l i n g ( E C M ) is app l i cab le for different r a n d o m intercepts. T h e seemingly unre la ted regression ( S U R ) technique is used to es t imate different sets of f ixed parameters, whi le e s t ima t ing different sets of r a n d o m parameters requires us ing r a n d o m coefficient mode l l i ng ( R C M ) . 2.3 T h e G r o w t h C u r v e L i t era ture I n d i v i d u a l m o n i t o r i n g is a c o m m o n p r o b l e m i n b io log ica l research. Therefore, various models a n d techniques have been invest igated to analyze g rowth curves. In the b io log ica l Chapter 2. Literature Review 27 l i terature , time-series da t a are usua l ly cal led l o n g i t u d i n a l da ta or repeated measurement data . Po thof f and R o y (1964) suggested the mu l t i va r i a t e analysis of variance ( M A N O V A ) m o d e l to est imate g r o w t h curves. T h e M A N O V A m o d e l requires that a l l i n d i v i d u a l s be measured at the same age. T h i s m o d e l a l lowed confidence in te rva l es t imat ion and hypothesis tes t ing. R a o (1965) extended Po thof f a n d R o y ' s mode l , t ransforming i t to a r a n d o m coefficient m o d e l where the m a t r i x of independent variables were iden t i ca l across i nd iv idua l s (i.e., X; = X for any i). G r i z z l e a n d A l l e n (1969) used an analysis of covariance technique where the weight ing was ob ta ined from a subset of covariates. Geisser (1970) p rov ided a Bayes i an analysis of the M A N O V A mode l . Fearn (1975) app l i ed the general Bayes i an l inear m o d e l proposed by L i n d l e y and S m i t h (1972) to g rowth curves. H e also s tud ied predic t ions given a sample f rom this m o d e l . R a o (1975) gave an e m p i r i c a l Bayes so lu t ion to a r a n d o m coefficient m o d e l where the parameters of the i n d i v i d u a l models were assumed r a n d o m l y d r a w n f rom a c o m m o n d i s t r i bu t i on whose expected value was the average of the i n d i v i d u a l parameters . L a i r d and W a r e (1982), cons ider ing the i n d i v i d u a l parameters as miss ing data , used the E M a l g o r i t h m (Demps te r et al. 1977) to est imate t h e m . W i t h the E M a lgo r i t hm, slow convergence or convergence to a l oca l rather t h a n g loba l m a x i m u m are major concerns. Berkey (1982) extended R a o ' s (1975) technique to nonl inear models . H u i and Berger (1983) also used an e m p i r i c a l Bayes i an approach to g rowth curves. T h e i r technique cou ld be cal led an e m p i r i c a l Bayes es t imat ion of a r a n d o m coefficient mode l . R a o (1987) discussed the p r o b l e m of p red ic t ing future observat ions w i t h g rowth models . H e ana lyzed paramet r ic and nonparamet r i c models . Fo r the pa rame t r i c models , he considered b o t h Bayes i an and e m p i r i c a l Bayes i an techniques to deal w i t h u n k n o w n parameters . Bayes ian and empi r i ca l Bayes i an techniques have p redomina ted b io logica l g rowth curve l i terature . W i t h o u t s ta r t ing a ph i l o soph i ca l debate on how they differ and w h i c h Chapter 2. Literature Review 28 one, if any, is more appropriate, one can notice that both are partially based on some auxiliary information. This auxiliary information, particular to an individual, individu-alizes a general model and its prediction of future observations. The prediction of future observations is dependent on the outcome of this individual-specific auxiliary information and on an estimate of the general model. 2.4 Conditional Densities and Expectations Many techniques for analyzing remeasured individuals and for predicting future observa-tions described so far assume that future observations are conditional on past, random, observations. The modeling process, therefore, involves conditional densities and expec-tations. A conditional density is a probability density function (pdf) where the outcome of a random variable X depends on the outcome of another random variable Y. The variable of interest is thus X given Y, or X | Y, and its pdf is noted p(X \Y). Most of the theory concerning conditional pdf's was developed by Bayes (1763), and a brief summary is presented here. Some basic matrix inversion rules are useful for deriving conditional densities and expectations. These rules can be found in certain texts (see Smith 1973, for instance). Rule (1) (D + E F E ' ) - 1 = D 1 - D 1 E ( E ' D ~ 1 E + F ^ ^ E ' D \" 1 Rule (2) ( D + B ) - 1 = D 1 - D ^ D 1 + B 1 ) - ^ 1 Rule (3) (D + B ) - ^ = I - ( D + B) - 1 D Conditional densities can be defined for univariate or multivariate variables. The multivariate case is more general. Let X and Y be two random vectors of n elements with pdf p(X) and p(Y), respectively. Their joint pdf is p(X,Y). Using Bayes' theorem on conditional probabilities, the joint pdf can be rewritten: p(X,Y) = p(X|Y)p(Y) Chapter 2. Literature Review 29 p ( Y | X ) p ( X ) From this theorem, the conditional pdf of X given Y can be defined as: p ( Y | X ) p ( X ) P ( X | Y ) = P(Y) When Y is given, the denominator p (Y) contains no unknown. It can be seen as a proportionality constant, so that the multiple integral of p ( X | Y ) over Xi, ... ,Xn yields one. This constant will be equal to the inverse of the integral of the numerator over the domain of X : p(Y) = l / / / . . . / P(Y\X)p(X)dX1dX2...dXn For instance, if Y given X is distributed as a i V ( X , f i ) , p ( Y | X ) = 1 (2TT)\" / 2 \ f l \ ^ V - - ( Y - x y n - ^ Y - x ) and X is distributed as a N(p, A ) , 1 P(X) = -exp -(X-p)'A-\X-p) (2TT)\"/ 2 I A I1/2 the unconditional variable Y will be distributed as a N(p,A + ft) P ( Y ) = (2TT)\"/ 2 I A + fl\V: ;exp hY-p)'(A + ftr^Y-p) To show this, rather than using the integral of the product of p ( Y | X ) and p(X), an easy alternative can be used. The variable Y given X can be seen as the sum of the variable X , normally distributed, plus an error term, e, uncorrelated with X , normally distributed too, with expected value zero and variance ft. Obviously, Y will be normally distributed since the sum of normal variables yields a normal variable. The expected value and the variance of this normal variable can be found easily: Chapter 2. Literature Review 30 Y = X + e, E[Y] = E[X] + E[e] = fi, var(Y) = var(X) + var(e) = A + ft. Using the definition of a conditional density, it can be shown that the variable X given Y will also be normally distributed: P[X|Y] = where (27r)\u00E2\u0084\u00A2|A | 1 /2|f2 | l /2 e X P -H(Y - xyn-^Y - x) + (x - ^A-^ X - M)> (i^ /'lW'\"* [-H(Y - ^ )'(A + \")\"1(Y - /*)}] A + n i 1 / 2 X exp \u00E2\u0080\u00A2-{(Y - xyn-^Y - x) + (x - U)'A-\X - p) -(Y-py(A + ny\Y-p)}} 1 (2TT)\"/2 I (A-1 + 0 - 1 ) - 1 ! 1 / 2 exp [-1 [(X - E[X\Y})'(A-1 + f l ^ f X - \u00C2\u00A3[X| Y]) \u00C2\u00A3[x|Y] = (A^ + n-^^n^Y + ^ ^ + n-1)-^-1/*-The expected value of the conditional pdf X given Y is a weighted average between the expected values of the variables Y given X and X, where the weights are inversely proportional to their respective variance. The variance of X given Y is the inverse of the sum of the inverses of the respective variances. The expected value i\u00C2\u00A3[X|Y] can be seen to be a random variable since it depends on the random outcome of Y. It is the sum of a constant, (A-1 + fi-1)-1 A_1/i, and a random vector normally distributed, Y, multiplied by a vector of constants. Therefore, it will be normally distributed, with the Chapter 2. Literature Review 31 following expected value and variance: E[E[X\Y}} E[X]=p va r (\u00C2\u00A3 [X |Y] ) v a r ( X ) - \u00C2\u00A3 [ v a r ( X | Y ) ] = A - ( A - 1 + ft-1)- 1 If it is assumed that A and fi are two positive definite matrices, then it can be shown that var(X) is always larger than va,i(E[X | Y]) . A matrix A is larger than another matrix B if the matrix A \u00E2\u0080\u0094 B is positive definite. Using matrix inversion rules (2) and (3) defined above, the expected value of the conditional density is equivalent to As demonstrated, statistical computations with conditional variables follow the same general rules as with unconditional variables, despite the more elaborate algebra. It is important to be famihar with the algebra of conditional variables to understand models making use of past observations to predict future observations. 2.5 Overview Permanent sample plot attributes can be assumed to be explained by some general bio-logical phenomenon. Observations of these attributes on a given plot at a point in time will also be influenced by the individual character of the plot and time-related variables, such as climate. Repeated measurements on a group of plots provide information on both the influence of the individual character and the time period effect. This suggests that permanent sample plot observations can be explained by a general model plus random terms specific to the individual and the time period. Since future observations must be predicted on a given plot, a prediction of the individual character of the plot as well as an estimate of the general model are needed E[X\Y] = A ( A + + ft(A + ft)~V-Chapter 2. Literature Review 32 to individualize the volume prediction. This approach recognizes the random nature of the relationship between observations coming from a sample plot. Future observations depend on the general model, and the outcome of the individual plot character. The individual plot character can be predicted from past observations. Therefore, predictions of future observations are dependent on past observations. Many approaches can be used to predict future observations based on past observa-tions. The approach suggested in this dissertation was inspired by the empirical Bayes models described by Rao (1975, 1987) and Hui and Berger (1983), the random coefficient model (Swamy 1970, 1971, 1974; Amemyia 1978; Ferguson and Leech 1978; Biging 1985), as well as the theory related to error component modeling. Chapter 3 Statistical Methods Future volumes can be expressed as the current volume plus a volume increment. The analogy with the interest on the money in a bank account is often used to explain forestry ideas such as volume growth and allowable annual cut to non-foresters. Going from interest to interest rate to biological growth rate are steps easy to understand, and the mathematics relevant to interest rates can simply be adapted to a forestry context. Consequently, the parameter of interest in the proposed model is the growth rate. A biologically-reasonable model for predicting volume growth rate is defined and estimated. This predicted volume growth rate is used to compound the current volume of a plot, the same way money can be compounded with a certain interest rate. This chapter is divided into four parts. In the first section, the similarities between biological growth rate and interest rate are analyzed and explained. In the second sec-tion, a two-stage statistical model based on interest rate mathematics and fundamental principles of biological growth is suggested to predict volume growth rate. Parameter es-timation is discussed and estimators are proposed in the third section. Finally, algorithms for compounding the current volume of a sample plot and for computing a confidence interval around the predicted volumes are presented in the last section. 3.1 The Timber Capital The trees on a permanent sample plot can be compared to a timber capital destined for yielding interest. The mathematics relevant to interest rates are simple and well-known; 33 Chapter 3. Statistical Methods 34 converting them to growth rates is straightforward as it is shown in section 3.1.1 The main difference between money and trees is the uncertainty coming with observations on biological organisms such as trees. This uncertainty must be recognized and incorpo-rated into the mathematical expression of the compounded capital. In section 3.1.2, the problem of uncertainty attached to volume observations is addressed. 3.1.1 The Interest Rate Analogy The growth in volume for a tree in a particular year is the amount of volume the tree gained between the beginning and the end of the growing season. When the growth is divided by the initial volume, it can be seen as a growth rate (i.e., the amount of volume per unit of initial volume gained in a year). The growth rate is the annual \"interest rate\" on the initial volume of a particular tree in a given year. This image can be extended to a plot or a stand without problem. In this section, the mathematics relevant to interest rates and their similarities with growth rates are explained. The interest rate at a particular bank will likely fluctuate over time. It will also likely differ between banks; two different banks can have different interest rates. In the same way, the growth rate of a tree (or a plot or stand) will be specific to a given individual and time. Mathematically, this gives V - y + 1 ) = ^ ( 1 + ^ . ) (3.1) where Vij = volume measured on individual i at time j, Pij = annual growth rate on individual i at time j, and (1 + pij) = annual multiplicative rate. In forestry, permanent sample plots are commonly used to monitor growth. However, annual measurements on permanent sample plots are rarely available. The sample plots Chapter 3. Statistical Methods 35 are remeasured at a regular frequency, usually between three and ten years. When the time interval between measurements j and j + 1 is more than one growing season, periodic growth is then observed. In this case, the ratio between the volume at the end of the period and the initial volume is the periodic multiplicative rate. The periodic multiplicative rate is the product of the annual multiplicative rates. Since the annual multiplicative rates are not observed, the periodic multiplicative rate can be considered as the geometric mean of the annual multiplicative rates raised to a power, Lj, equal to the time lag between measurements j and j + 1: Vi{j+1) = V k ( l ( 3 . 2 ) where Pij = [ H \u00C2\u00A7 = 1 ) ( 1 + P i j)}^ - 1. The periodic multiplicative rate is the percentage of initial volume obtained after a certain period of time. Therefore, the notation used for annual remeasurements is somewhat inadequate for periodic remeasurements since the subscript j in refers to the period between measurements j and j + 1 while in Vij it refers to the instantaneous moment at measurement j. Model (3.2) can be rewritten using the subscript p to refer to the time interval between measurements j and j + 1 with the initial volume on plot i and the volume at the end of the period: = V\u00C2\u00B0(l + pip)L\u00C2\u00BB. (3-3) With this notation, the growth rate is located in space and time. The subscript i indicates a geographic space, (delimited by some longitude and latitude) and the subscript p is temporal space (delimited by an initial year and a final year). This approach to growth modeling could be taken for other stand attributes, such as diameter, height, basal area or number of stems per hectare. For illustrative purposes, the approach is shown for plot volume. For the remainder of this thesis, the term volume will refer to plot volume. Chapter 3. Statistical Methods 36 The periodic multiplicative rate in model (3.3) can also be expressed using the Napierian base (e = 2.7182818 ...) K = V\u00C2\u00B0 x exp[8tpLf (3.4) where 9ip = ln( l + pip). The parameter &ip is easily isolated 1 0* = \u00E2\u0080\u0094 I n V-1 ' P vip (3.5) Using the logarithm rules and the definition of a derivative, it can be shown that, as the time period goes to zero, 6{p becomes the instantaneous growth rate lim 9j\u00E2\u0080\u009E = lim *0 Lr In V 1 v\u00C2\u00B0 v ip \u00E2\u0080\u0094 lim hm At\u00E2\u0080\u0094o d\n[V] dt 1 dV V~dt Aln[V] At (3.6) where t is time. Similarly, 0ipLp can be shown to be the integral of with respect to time between the initial and the final measurements. Time can be expressed on a relative scale. The initial measurement is taken at time 0, and the final one is taken at time Lp. This yields dipdt (3.7) which is the periodic multiplicative rate expressed in a logarithmic scale. The concept of growth as the derivative of volume with respect to time (or volume as the integral of growth) was recognized in the early sixties (Buckman 1962, Clutter Chapter 3. Statistical Methods 37 1963). This is a simplification of reality; it assumes that growth is a constant process, which is not exactly true since growth stops in winter. A discrete analog would be to take the annual interest rate analogy. This approach has some shortcomings too; there is no growth assumed before the end of the growing season, similar to the annual interest on a bank account that is paid once a year. Another discrete approach would be to follow the modern bank system and use daily interest rate, but this would demand a tremendous data set. Assuming that growth is a continuous process leads to relatively easy mathematical computations while being close to reality. 3.1.2 Uncertainty and Growth Rate The deterministic model (3.4) could be used if the true initial and final volumes were observed. However, these quantities will likely be measured with some error. The ob-served volumes are random variables. The variability in observed volumes will introduce some variation into the observed growth rate. Because the growth rate is essentially the ratio of two variables, it becomes a random variable itself. A stochastic term must be incorporated in model (3.4) to account for this variability. A priori assumptions are then made about this error term in order to draw inferences using statistical theory after mathematically fitting the model. This stochastic term can be included in different ways. Three common ones are iVi/VP) = exp[0ipJLp] + u t p (3.8) iVi/VP) = exp[(0ipLp) + uip] (3.9) {V\u00C2\u00A3/V?) = exp[(0ip + uip)Lp}. (3.10) The choice of a model depends on its mathematical tractability and its biological interpretation. Model (3.8) is a nonlinear model and is particularly difficult to fit to recent permanent sample plots. Model (3.9) is linear after a logarithmic transformation. Chapter 3. Statistical Methods 38 It adds an error term to the periodic growth rate. The error term will undoubtedly depend on the length of the time period even if not explicitly expressed in the model. Model (3.10) is also log-linear. The error term is included with the instantaneous growth rate. Its dependency on the time lag is therefore implicitly expressed in the model. This model has the advantages of being simple to fit and of corresponding to the interpretation of a growth rate whose variation is dependent on the period length. For these reasons, model (3.10) was considered the best approach. Bank account and trees have this point in common: they become bigger as time passes. However, while knowing the exact amount of money in a bank account is no problem, direct observations of plot volume is impossible. An error term can be included in the mathematical expression of periodic growth rate to account for this uncertainty. Model (3.10) is the basic expression of growth rate with uncertainty. Within this basic expression, it is possible to assume that the growth rate, &{v, can be explained by some biological phenomena. A relationship explaining #;p from some variables can be substi-tuted in the basic expression. This framework leads to the statistical model to predict future growth rates, which is discussed in the next section. 3.2 The Statistical Model The framework proposed for estimating the statistical model is a two-stage procedure. The first stage is the basic expression describing the mathematics of growth rate, and the second stage is a biological explanation of the true growth rate. The framework was designed specifically for the irregular remeasurements of perma-nent sample plots. These irregular remeasurements, along with overlapping time periods will create notation difficulty. Conventions must be defined to avoid confused notation. Once the notation is defined, each stage of the procedure can be discussed. The first Chapter 3. Statistical Methods 39 stage deals with the problem of a small number of observations on each plot and the error components. The second stage makes use of Von Bertalanffy's growth curve to estimate the growth rate from site index and volume. Finally, both stages are joined in the combined model. 3.2.1 Notational Conventions It will be assumed that there are N plots covering M intervals. Some time intervals may be overlapping. For instance, if some plots were measured in 1970, 1975 and 1980, and the other plots were measured in 1973, 1978, and 1983, there would be four time intervals: (1970,1975], (1973,1978], (1975,1980], and (1978,1983]1. Not every plot will have an observation in each time interval. Furthermore, the number of plots measured each year will vary with the money budgeted annually for the remeasurement of permanent sample plots. Consequently, the number of remeasurements will differ among plots and among time periods. It will be assumed that there are T; growth observations on plot i [i = 1, . . . , N), Sp observations in time interval p ( p = l , . . . , M ) , and W observations in total where N M W = Y,Ti = Y,SP; with W > N and Sp > 1. i = i P = i Because time intervals may be overlapping, consecutive measurements on a plot do not lead to observations in consecutive time intervals. For the example mentioned above, the first set of plots would have observations in time intervals (1970,1975] and (1975,1980], while the second set of plots would have observations in time intervals (1973, 1978] and (1978, 1983]. The time interval numbering and the plot numbering are arbitrary. They are only coordinates to link an observation with geographic and time locations. Also, it is assumed that all plots used in the analysis come from the same geographic region, and 1This notation supposes that the plots were measured at the end of the growing season. Chapter 3. Statistical Methods 40 are dominated by the same species. No notation is used to identify these two qualitative variables. 3.2.2 F ir s t -S tage M o d e l The first-stage model is the basic expression for periodic growth rate, with uncertainty included with the instantaneous growth rate. Recent permanent sample plots only have a few growth observations. Using an approximation, the proposed first-stage model can be estimated from plots that have had only two growth observations. This approximation is explained in the model definition section. Next, the error term is decomposed to take advantage of the measurements of many individuals in a given time interval. Assumptions about the error term and their justifications follow. Finally, the matrix notation of the model is given to simplify the notation and derivation of the parameter estimators. 3.2.2.1 M o d e l Def in i t ion The observed periodic multiplicative rate is the ratio Vip/Vip. The average annual mul-tiplicative rate observed during this time period is the observed periodic multiplicative rate raised to the power l/Lp. Because the model (3.10) is in exponential form, it will be more convenient to work with the natural logarithm of this ratio. So let It will be assumed that the variable Y{p represents an estimate of the instantaneous growth rate, 8ip, since Lp is usually small. It is the sum of the true instantaneous growth rate plus an error term: (3.11) (3.12) Chapter 3. Statistical Methods 41 p l o t oge Figure 3.1: Growth Rate over Time for a Sigmoidal Yield Curve which corresponds to the logarithmic transformation of model (3.10). Recent permanent sample plots cover only a few time intervals, and if these time intervals are short, the whole period covered by the sample plot can still be considered close to zero, especially when, it is compared to the lifespan of the stand. Therefore^ the instantaneous growth rate could also be estimated from the ratio between the last and the first measurements on a sample plot. This corresponds to assuming a constant growth rate on a plot. Hence, the growth rates estimated in the intervals can be considered as repeated observations of the same parameter, the constant growth rate. If the yield is assumed to follow a sigmoidal shape over time, the growth rate in the early ages declines rapidly to level off as the stand gets older (Figure 3.1). Therefore, the assumption of constant growth rate seems reasonable for a large part of the plot's life. Obviously, this assumption is a compromise to approximate the trend of growth rate over time for a sigmoidal yield curve and to allow for parameter estimation. Chapter 3. Statistical Methods 42 This approximation will simplify the model, but also bring some limitations. It likely could cause some bias. Assuming a constant growth rate on a plot will tend to un-derestimate the growth rate at the first measurement and overestimate it at the last measurement. This bias should be negligible, or at least reasonably small, if the whole period covered by the plot is short enough. The period of time where the approximation is reasonable depends on the species and the age of the plot. The bias should be more pronounced on young plots. Also, the rate of change in the instantaneous growth rate is species-dependent; fast growing species will level off faster than slow growing species. Considering these limitations, the approximation can still prove useful in many cases such as older plots and fast growing species. Model (3.12) can be approximated by This model can be estimated with as few as two growth observations on a perma-nent sample plot. That is, the model requires a minimum of three measurements on a permanent sample plot since the number of measurements is equal to the number of growth observations plus one. Also, measurements on different individuals in the same time period will allow decomposing u; p , the error term. 3.2.2.2 Error Components Following the principles of error component modeling, when measured on a given individ-ual (i.e., the individual is fixed), the error term u; p can be divided into two components: the error due to the time-dependent variables unaccounted for in the model and a term for the unaccounted variables that are not time-dependent. This gives Yip \u00E2\u0080\u0094 Oi + U;. (3.13) u i p = r p + Vj. ip (3.14) where Chapter 3. Statistical Methods 43 rp = error term explained by the time-dependent variables unaccounted for in the model. It remains constant for all individuals in a time period but varies among time periods. V i p = unexplained residual noise. It accounts for the observed difference between indi-viduals with the same true growth rate in a given time period. Some unobservable or unmeasurable variables will have a similar effect on geograph-ically close plots in a time period (Gregoire 1987). For example, climatic variables, and major insect and disease epidemics will uniformly affect a small region in a given time period. Total number of days in the growing season, rainfall, quantity and quality of sunlight, etc. should be roughly similar for all plots in a time period. However, it would be impossible to include all these variables in the model. Their cumulative effect will vary from time period to time period. With a two-way lay-out such as a panel data set, it becomes possible to estimate the variation of the random effect caused by these variables. It is therefore possible to distinguish the time-dependent ignorance from the unexplained ignorance. The error term in model (3.13) was decomposed to represent this knowledge about the stochastic variation in the model. 3.2.2.3 Assumptions Combining equations (3.13) and (3.14), the model for estimating the growth rate on plot i becomes: Yip = 6i + Tp + vip. (3.15) Assumption 1. E[TP] \u00E2\u0080\u0094 E\v{p\ = 0 Assumption 2. E [ T p i q ] ' = 62 iip = q = Spq iipf^q Chapter 3. Statistical Methods 44 Assumption 3. E[vipVjq] = cr2 if i = j and p = q = 0 otherwise Assumption 4. E[TpViq] \u00E2\u0080\u0094 0 Assumption 5. rp and v;p are normally distributed. These assumptions imply that u;p is normally distributed and E[uip] =0 E[uipujq] = o2 + 52 if i = j and p = q = S2 ii i ^ j and p = q = Spg Hpf^q Assumption 1 means that the expected value of the observed growth rate on a plot is assumed to be the true growth rate for this plot. Since Y\u00C2\u00B1v is the logarithm of a ratio of random variables raised to a certain power, this assumption will make sense if and only if the expected value of the logarithm of the plot volume measured on plot i at time j is the logarithm of the true volume E[Yip] = 6i t=^> E[ln Vn] = l n where Vij is the true volume. The error of the logarithm of the observed volume has often been assumed to have zero mean in growth and yield literature (Schumacher 1939, Clutter 1963, Curtis 1967, Ferguson and Leech 1978, Gregoire 1987), and the same assumption is made here. The bias created by the approximation of a constant growth rate is assumed to be negligible. Assumption 2 concerns the time-dependent error terms affecting the growth rate. The effects of these error terms should be correlated with each other since the time Chapter 3. Statistical Methods 45 periods are overlapping. Following Garcia (1983), it can be assumed that the correlation is proportional to the overlap between two time periods. Nonoverlapping time periods are considered independent; the weather varies randomly from one period to another. Overlapping time periods share the same climatic conditions for a certain number of years which should indicate the presence of correlation. This assumption also implies that consecutive growth observations on a given individual are considered independent because they are always taken in nonoverlapping time periods. Let Cpq be the number of overlapping years between periods p and q, then the covari-ance between rp and rq, 8m, can be assumed to be 2 C p q ( L p + Lq) *\u00C2\u00BB ~ 5 2LpLq ( 3 - 1 6 ) where Lp and L q are the respective lengths of periods p and q. To explain this estimator, let r p be divided into two parts: r p l and r p 2 where r p i is the temporal error term for the first Lp \u00E2\u0080\u0094 Cm years, and r p 2 is the temporal error term for the remaining Cm years. Similarly, let rg be divided into r g l and r g 2 where r g l is the temporal error term for the first Cpq years, and r g 2 is the temporal error term for the remaining Lq \u00E2\u0080\u0094 Cpq years. Both r p 2 and r g i are error terms for the same time period. Within a time period, the variance is assumed constant. This means that the variance of r p 2 is Cpq82/Lp when estimated from time interval p. Similarly, the variance of the error explained by the same time interval but estimated from time interval q, iqi, would be Cpq82/Lq. Also, the correlation between r p 2 and r g l is really the variance of the time-dependent error term for the Cpq overlapping years. Therefore, Cpq82/Lp and Cpq82/Lq are two different estimates of the correlation between r p 2 and r g i . Averaging out both estimates will yield a single estimate of the correlation between the two error terms, providing in the same time an estimate of the variance of the time-dependent error term for these Cm overlapping years. Chapter 3. Statistical Methods 46 Mathematically, this yields E[(iPi + rp2 ) ( r qi + r g 2)] E[ipliql] + E[Tpliq2\ + E[ip2iql} + \u00C2\u00A3 [ r p 2 r 9 2 ] \u00C2\u00A3 [ r p 2 r g i ] \{E[rl2] + E[?ql]) -(var(r p 2) + var(r g l)) 2 Lp Lq = 82 (3.17) 2LpLq Two nonoverlapping time periods have Cpq \u00E2\u0080\u0094 0 satisfying the assumption of indepen-dence; two completely overlapping time periods have Lp = Lq \u00E2\u0080\u0094 C^ which would yield a covariance of 82. Assuming a known form for 8pq is convenient since it decreases the number of unknown parameters to be estimated in the model. Assumption 3 implies that the unexplained error terms are independently distributed with common variance within a time period. This error term represents the stochastic variation between individuals observed in given time period after correction for the true growth rate. It accounts for the general ignorance in the model. For various reasons, some plots can be more sensitive than others in unusual but not extreme climatic conditions. For example, rainfalls slightly below the average can have serious effects on some plots and leave other plots almost unaffected. Other special conditions can prevail on a particular stand in a time period. A local insect infestation, for example, could affect only a small number of plots. For extreme conditions, it can effect all plots similarly. For instance, a major insect infestation would probably have the same kind of impact on all plots. The variance of this error term is assumed to differ Chapter 3. Statistical Methods 47 among time periods. The variation among plots within a time period will differ with the magnitude of the effect of climatic conditions on each plot. Some climatic conditions will cause a wide range of impact on individual plot growth rates while other conditions will affect all plots in a similar fashion. This error term can be seen as the interaction between time periods and plots, and it is assumed that the variance of this interaction varies with time periods. It is also assumed that the interaction between a given plot and a time period is uncorrelated with the interaction between the same plot and another time period, or with the interaction between another plot and the same time period. Assumpt ion 4 implies that the time-dependent error term is uncorrelated with the unexplained error term. The individual plot behaviour in certain climatic conditions is unrelated with the average behaviour of all plots in the same conditions. Assumpt ion 5 was justified by Aitchison and Brown (1957; p. 1-2) who mentioned that a variable which is the product of numerous independent random variables tends to be lognormally distributed around its true value. If the growth rate at any time is considered a random proportion of the initial volume, and if the growth rates are considered independent for a given individual, then the volume tends to be lognormally distributed. Therefore, the observed growth rate tends to be normally distributed. 3.2.2.4 M a t r i x Nota t ion Model (3.15) and its assumptions can be expressed with matrices. To define the model using matrix notation, the observations are sorted initially by time periods and secondly by plots. The first observation is then the observation measured in period 1 on the plot with the smallest plot number from all plots measured in period 1, and is denoted Y\u00E2\u0080\u009E\. Chapter 3. Statistical Methods 48 The last observation is the observation measured in period M on the plot with the largest plot number from all plots measured in period M ( Y i M ) . There are W observations in the matrix containing the growth observations. The first-stage model is conditional on fl. The second-stage model will express fl; as a function of site index and initial plot volume. This approach will make the estimation of the covariance matrix of the random terms easier. The first-stage model is: Y = F9 + Gr + v (3.18) with F = G = Gi V . l Y = 9 = r = v = YiM 8N TM Fa . . . F\u00E2\u0080\u009E Gp = F, Fiw G p i G pW where Fik \u00E2\u0080\u0094 1 if the kth element of Y was measured on plot i = 0 otherwise and Gpk = 1 if the kth element of Y was measured in time period p = 0 otherwise. Chapter 3. Statistical Methods 49 T h e m a t r i x of the s u m of b o t h error components w o u l d be u = Gr + v . T h e various assumpt ions about the error terms lead to the fo l lowing matrices E[w'] = S = d iag where E[rr'] A = S2A A Ml AM2 1 A12 . . . A i M A 2 i 1 . . . A 2 M (3.19) (3.20) a n d _ Cvq(Lv + Lq) m ~ 2LpLq \u00C2\u00A3 [ u u ' ] = t o G ' + S = n . (3.21) T h i s concludes the first-stage mode l . In this m o d e l , a constant g rowth rate is assumed for each plot . Since the i n d i v i d u a l is f ixed, the first-stage m o d e l yields an est imate of the app rox ima ted constant g rowth rate of the p lo t . T h i s est imate comes f rom the i n d i v i d u a l g rowth observations and can be cal led the observed ins tantaneous g rowth rate. U s i n g a group of i n d i v i d u a l s , the g rowth rate can also be pred ic ted f rom independent variables. T h i s is the second-stage m o d e l . 3.2.3 Second-Stage Model T h e second-stage m o d e l provides a b io log ica l re la t ionsh ip between the true instantaneous g rowth rate and the site i n d e x and the vo lume. T h i s re la t ionship is der ived f rom V o n Bertalanffy 's g rowth curve. T h e b a c k g r o u n d exp l a in ing the re la t ionship is presented first. T h e s ta t i s t ica l m o d e l and its assumpt ions fol low. F i n a l l y , the m a t r i x n o t a t i o n of the m o d e l is g iven. Chapter 3. Statistical Methods 50 3.2.3.1 Background It is reasonable to assume that the individual growth rate can be explained by some variables plus an individual-specific random term. For example, it can be assumed that the individual growth rate is a function of site index (SI) and initial volume (V). Any other factor that could have an influence on the individual growth rate can be included in a stochastic term. This function can be assumed to be identical for all individual growth rates, with fixed parameters and an additive error term Fitting the relationship expressed in (3.22) would provide an estimate of the growth rate of a plot as explained by its site index and its initial volume. Also, the stochastic nature of the individual growth rate is recognized by the additive error term. The structure of the deterministic part of equation (3.22) can be derived from the biological interpretation leading to Von Bertalanffy's growth curve. The fundamental assumption underlying this curve is as follows. The potential growth is the difference between anabolic gains and catabolic losses. Anabolic gains are assumed to be propor-tional to the absorptive surface. The absorptive surface is assumed to be proportional to the volume raised to the power 2/3. Catabolic losses are assumed to be proportional to the volume. This gives The potential growth is the production of surplus metabolic products (Pienaar and Turn-bull 1973), which is the gain in biomass. The gain in biomass is distributed among various components of the trees: fine roots, large roots, stem, branches, leaves and bark. There-fore, the actual volume growth depends on the percentage of the total gain in biomass going to stem biomass. This percentage, \u00C2\u00AB 3 , is an efficiency coefficient relating volume 6i=f(SIi,Vi) + ei. (3.22) Potential growth = K^V2^3 \u00E2\u0080\u0094 K2V. (3.23) Chapter 3. Statistical Methods 51 production to production of surplus metabolic products: Actual growth = K 3 ( K I V 2 ^ 3 \u00E2\u0080\u0094 K2V) \u00E2\u0080\u0094 = K 3 K I V 2 / 3 - K3K2V. (3.24) Keyes and Grier (1981) showed that the efficiency coefficient was correlated with site quality. Richer sites have a better availability of nutrients to put on stem growth. Poorer sites are poorer in nutrients, and more energy must be dedicated to fine roots production in order to get the nutrients. Kurz and Kimmins (1987) similarly found that the percent-age of total biomass going into below-ground production is inversely proportional to site quality. As site quality decreases, a larger proportion of total biomass goes below ground, and consequently, a smaller proportion goes above ground. If the percentage of above-ground biomass is proportional to site quality, it can be expected that the percentage of biomass going into stem biomass will be correlated to site quality as well. Axelsson and Axelsson (1986) showed that the efficiency coefficient increases with fertilization. Richer sites do not need as developed of fine root system as the poorer sites, and are there-fore more efficient at transforming their production of surplus metabolic products into growth. Site index can be assumed to be highly correlated with site quality; therefore, it is reasonable to assume that the efficiency coefficient is proportional to site index. Let K 3 be K4SI. Equation (3.24) can be rewritten ~ = K A S I K I V 2 / 3 - K A S I K 2 V ot = fa SI V2/Z - 02 SI V. (3.25) The growth rate is the growth divided by the volume. Dividing both sides of equation (3.25) by V yields Chapter 3. Statistical Methods 52 3.2.3.2 Model Definition The left-hand side of model (3.26) is the instantaneous growth rate, which is approx-imated by fl;. Consequently, model (3.26) represents the deterministic part of model (3.22). The complete model (3.22) can be rewritten: SI-9 i = h y ^ i - h M i + * - (3-27) The parameter fl; is assumed to be constant during the time periods covered by plot i. However, this does not correspond to the concept of instantaneous growth rate. This concept can be retrieved in the second-stage model if fl; is assumed to be measured at a single point, within the interval between the first measurement and the last measurement on plot i. The volume at this point can be used as the initial volume. This definition of initial volume will allow for simulating an instantaneous observation of the growth rate. Hui and Berger (1983) adopted a comparable strategy in a similar situation, mentioning that it should bring some robustness to the model. The constant growth rate is an average of the instantaneous (but unknown) growth rates. Therefore, its value is the value of an instantaneous growth rate located between the first and the last measurements. This suggests that an \"average\" initial volume would lead to a more appropriate definition of initial volume for an instantaneous growth rate. There are many possible definitions of this \"average\" initial volume, and since the real growth rate is unknown, there seems to be no definite rule to define this initial volume. One possible estimator is to take the mean volume, that is the arithmetic average of the observed volumes on a plot. V\u00C2\u00B0 = (3.28) where Vij is the volume observed on plot i at measurement j. When the growth rate is as-sumed to be estimated at this point, it better corresponds to the concept of instantaneous growth rate represented by fl;. Chapter 3. Statistical Methods 53 3.2.3.3 Assumptions Once an appropriate definition of the initial volume is accepted, stating the assumptions about the error terms is the next step. These assumptions are: Assumption 6. Sli and ( V ^ 0 ) 1 / 3 are fixed and measured without error Assumption 7. Efe] = 0 Assumption8. E[eiej] = 2/(V?)1/3 if i = j = 0 otherwise Assumption 9. e; is normally distributed. Assumption 6 is a common assumption in regression. As a consequence, the error term e; is assumed to be uncorrelated with both independent variables, Sli and V\u00C2\u00AE. Assumption 7 implies that the expected growth rate on a plot is assumed to be a function of site index and initial volume. Assumption 8 has two implications. First, it implies that the growth rates are inde-pendent between plots. Knowing the true growth rate on a plot gives no idea about the true growth rate on another plot. The eight assumption also states that the variance of the individual growth rates is inversely proportional to initial volume raised to a power of 1/3. The variance of the individual growth rate should decrease as the plot gets older. Variation in growth rate is more important when the growth rate is large (i.e., when the plots are young). If it is assumed that the variance is proportional to the expected individual growth rate var(0{) oc &i Chapter 3. Statistical Methods 54 (i.e., i f the variance expressed as a p r o p o r t i o n of the expected g rowth rate is assumed to be constant) , i t w i l l also be inversely p r o p o r t i o n a l to the i n i t i a l v o l u m e raised to a power of 1/3 since the expected value of the g rowth rate is i tself inversely p r o p o r t i o n a l to i/(v\u00C2\u00B0y/\ Assumption 9 can be jus t i f ied by i n v o k i n g the C e n t r a l L i m i t T h e o r e m . T h e numerous variables in f luenc ing the g rowth rate unaccoun ted for i n the second-stage m o d e l are assumed to be add i t ive w i t h their average effect bel ieved to be negl igible . 3.2.3.4 Matrix Notation M o d e l (3.27) and its assumptions can be wr i t t en us ing m a t r i x no ta t ion as: 6 = X0 + e (3.29) w i t h the fo l lowing matrices: \" Ol\" X11 X21 0 = X = e = 9N X\N where Xn = Sli/'(V\u00C2\u00AE)1?3 and X2i = \u00E2\u0080\u0094Sli. F r o m the assumpt ions about the error t e rm, the d i s t r i bu t i on of 9 is easily found: 9 ~ N(X0,(f>2V) (3.30) where V is a d iagona l m a t r i x whose i t h d i agona l element is l / ( V i \u00C2\u00B0 ) 1 / ' 3 . T h i s concludes the second-stage mode l . T h i s stage provides a b io logica l ly- reasonable m o d e l to predict the growth rate. T h e est imate g iven by the second-stage m o d e l can be seen as the predic ted g rowth rate of a given p lo t . Chapter 3. Statistical Methods 55 3.2.4 Combined Model The second-stage model (3.27) can be substituted into the first-stage model (3.15). 3.2.4.1 Assumption A new assumption is needed. Assumption 10. 22[v{p|ej] = -E[rp|e;] == 0 Assumption 10 implies that the error terms from the first-stage model are uncorre-lated with the error terms from the second-stage model. This means that the individual-dependent variables unaccounted for in the model are assumed to be uncorrelated with the unaccounted time-dependent variables or the residual error term. It should be noted that the error terms from the first-stage model are also uncorrelated with the independent variables since the independent variables are assumed to be fixed. 3.2.4.2 Matrix Notation The combined model can be written with matrix notation as: Yip = Qx - 82 Sli + e; + r p + v;. (3.31) Y = FX/3 + Fe + Gr + v (3.32) Based on the assumptions about the error terms, the distribution of Y is: Y ~ N(FX3,V). (3.33) where \Er = $ >2FVF' + n. Chapter 3. Statistical Methods 56 3.2.5 Summary Growth rate and interest rate are very similar in nature, and the mathematics relevant to interest rates can be adapted to permanent sample plot growth rates. On recent permanent sample plots, the growth rate can be approximated by a constant. When the time period effect is removed, an estimate of the constant growth rate can be obtained. This can be called the observed instantaneous growth rate. The instantaneous growth rate can also be predicted from independent variables using a relationship derived from Von Bertalanffy's growth curve. This relationship presents the growth rate as a function of site index and volume at the beginning of the growing season. The second-stage estimate of growth rate can be called the predicted instantaneous growth rate. Both models can be joined to express the growth rate observations as a function of initial volume, site index, an individual effect, a time effect and an unexplained error term. For predictions of future volumes on a permanent sample plot, the parameter of interest is the random variable and its prediction is discussed in the next section. 3.3 Parameter Estimation Traditionally, the parameter of interest has been 0. However, this parameter does not recognize the individual character of each plot, represented by e;. The parameter 6i does recognize this individual character. The group of permanent sample plots provides the information to estimate the global parameter 0. Repeated observations on a given plot yield the information to predict the individual term e; of this plot. The estimate of 0 and the predictor e; can be joined to predict the actual parameter of interest. It is a common technique in parameter estimation to assume that the covariance ma-trix of the model random terms are known in order to derive the parameter estimators. Chapter 3. Statistical Methods 57 This simplification often allows a better understanding of the estimators. The estima-tors of the parameters composing 6 are first derived for the case where the covariance matrix of the model random terms are assumed to be known. Second, estimators for the unknown covariance matrices are defined. Finally, the estimated covariance matrices are substituted in the parameter estimators to yield estimates of the best linear unbiased predictor of 0 and of its variance. 3.3.1 Estimation When Covariance Matrices Are Known When the covariance matrix of the various random terms in the combined model (3.32) are assumed known, only the parameters 8 and e need to be estimated. These parameter estimates are random variables, and their variances can also be estimated. The parameter estimates are discussed first followed by a discussion of the estimate of the variance of the predictor of 0;. Deriving and simplifying the estimators involved a lot of matrix algebra. Most of the algebraic steps were carried out to provide a proof of the estimators. While complex at first glance, these steps are easy to follow using the matrix inversion rules given in section 2.4. Of course, the first and the final lines of any proof are sufficient to understand the text. The intermediate lines are complimentary, for the pleasure of the mathematically-inclined reader. 3.3.1.1 Parameter Estimators Estimators for 8 and e are first derived from a theorem given by Harville (1976), and combined to form a predictor of 6. The simplification and interpretation of this predictor using matrix algebra leads to an easily interpretable form. The suggested predictor is shown to be a weighted average between the observed instantaneous growth rate and Chapter 3. Statistical Methods 58 the predicted instantaneous growth rate. Interpretation of 9 is discussed following the mathematical derivation of 3 and e. Harville (1976) gave an extension of the Gauss-Markov theorem for prediction of random terms. Using theorem 1 from Harville (1976), the best Hnear unbiased predictor of e can be found. This leads to 9, the individual growth rate. The random parameter 9 is the sum of a fixed parameter, X8 and a random term, e (see (3.29)). The B L U P of 9 will be the sum of the B L U E of X3 and the B L U P of e. From the full model (3.32), the B L U E of 3 is the GLS estimate $ = [ X ' F ' ^ F X J ^ X ' F ' t f ^ Y (3.34) which leads to X/3, the B L U E of X8. FoUowing Harville (1976), the B L U P of e is e = 2VF'*_1(Y - FX/3) (3.35) which can be developed further using the definition of ^ and matrix inversion rules e = ^ V F ' ^ Y - 0 2 V F ' \u00C2\u00A5 - 1 F X / ? = ^ v F ' j n - 1 - n-^KF'n^F) + ^V-^ F^'O-^ Y - ^ 2 V F ' { n _ 1 - n ^ F K F ' n ^ F ) + ^ v ^ ^ F ' n - ^ F x / ? = 0 2 vF'n _ 1Y - ^ v F ' n ^ F K F ' n ^ F ) + ^-'v-^^F'n^Y - ^ V F ' n ^ F X / S + ^ V F ' n - ^ K F ' n ^ F ) + t-'v-^F'n^Fxp = < ?6 2V{I-(F'n _ 1F)[(F'fl\" 1F) + ^- 2 V- 1 ] - 1 }F , n - 1Y -4> 2v{i - (F 'n^FjKF'n^F) + ^ v - ^ - ^ F ' n ^ F X / ? = ^ v K F ' n ^ F j - ^ ^ v i - ^ F ' n \" ^ ) - 1 ^ - ^ Chapter 3. Statistical Methods 59 - ( /^[(F'ft - 'F)- 1 + 0 2V]- 1X/3 (3.36) This predictor is unbiased in the sense that E[e] \u00E2\u0080\u0094 E[e] = 0 where 0 is a column-vector of zeros. It is best in the sense that E[{e - e)'(e - e)] < E[{e - e)'(e - e)] (3.37) where e is any other predictor of e. The BLUP of 6 can now be defined as 9 = X/3 + e (3.38) = X/^ + ^ V ^ F ' O - ' F ) - 1 + ^ 2 V ] - 1 ( F ' f t - 1 F ) - 1 F ' \u00C2\u00A3 r 1 Y - ^ [ ( F ' f i - ' F ) - 1 + 2V]-lXp = ^ V K F ' O - ' F ) - 1 + ^ v j - ^ F ' n ^ F ^ F ' n ^ Y +{I - ^ [ ( F ' O ^ ' F ) - 1 + ^ vj-^x/S = ^ V K F ' J T ' F ) - 1 + 2\]-1(F'n-1F)-1F'ii-1Y +(F'n _ 1 F)- 1 [(F'fi\" 1 F)- 1 + ^ 2 V ] _ 1 X / 3 (3.39) This predictor of 9 can further be simplified; however, the terms it includes first must be recognized. The term ( F ' f i ^ F ^ F ' f i ^ Y is the GLS estimate of 9, 9, from the first-stage model, model (3.18). Its distribution is conditional on 9 since 9, a random variable, was assumed fixed in the first stage. This distribution is 9\9 ~ N(9, [F'fi^F]- 1). (3.40) Therefore, the term [F'fi^F]- 1 is the variance of the variable 9 \ 9. From the distribution of 8, shown in (3.30), the variance of 9 is 2V, and its expected value is X/3. Combining Chapter 3. Statistical Methods 60 this information with the distribution shown in (3.40), the unconditional distribution of 9 is: 0~iV(X/3 ,$) (3.41) where # = [F'ft^F]- 1 + cfV. The GLS estimate of 3 from the full model, defined in (3.34), is equivalent to the GLS estimate of 3 from the unconditional distribution of 9, shown in (3.41) (Amemiya 1971). This can be shown using the matrix inversion rules. J3 = [ X ' F ' S ^ P X ^ X ' P ' * \" ^ (3.42) = [ X ' F ' J T X F X - X ' F / f i \" 1 F ( F / n \" 1 F + ( ^ - 2 v - 1 ) - 1 F ' f i ~ 1 F X ] - 1 \u00E2\u0080\u00A2 [ X ' F ' f t _ 1 Y - x ' F ' n ~ 1 F ( F ' n - 1 F + ^ - ' v - ^ ^ F ' n ^ Y = [X'{I- ( F ' n _ 1 F ) ( F ' f i _ 1 F + 0 - 2 v - 1 ) - 1 } ( F ' n _ 1 F ) X ] - 1 \u00E2\u0080\u00A2[X'{I - ( F ' f i - 1 F ) ( F ' n _ 1 F + ^ V - ^ - ^ F ' I T ' Y ] = [X'^F ' fi - 'F)- 1 +rV]-\F'fl-1'E)-1F'Ql-1Y = [ X ' S ^ X j ^ X * - ^ (3.43) This demonstrates that X/3 in 9 (see model (3.39)) is equivalent to the GLS estimate of the expected value of the unconditional distribution of 9. The BLUP of 9, equation (3.39), can therefore be rewritten: 9 = var(0)[var(0 \6) + vai(0)] _ 10 + var(0 \ 9)[vav(9 \ 9) + var(0)]_1\u00C2\u00A3[0] (3.44) Chapter 3. Statistical Methods 61 It is almost identical to the expected value of the conditional variable 9 given 9. Using Bayes' rule, the expected value of 9 \ 9 is shown to be: E[9\6) = var(0)[var(0 19) + v a r ^ ) ] \" ^ + var(0 | 0)[vai(01 9) + v a r ^ ) ] \" 1 ^ ] (3.45) The only difference is that the unknown expected value of the unconditional variable 9 is replaced by an estimate in the BLUP. The B L U P can be interpreted as a weighted average between the growth rate observed on a plot and the growth rate predicted from some plot attributes. This can easily be understood. The individual growth rates can be predicted from some independent variables with some error, and it can also be observed with some error. The best predictor will be an average between the prediction and the observation with the weights proportional to the confidence in each estimate. This confidence is inversely proportional to the respective variance of each estimate. The B L U P is this weighted average. As mentioned before, the B L U P is a random variable. The precision of the B L U P can be obtained by estimating the variance of the predictor. This variance estimator follows. 3.3.1.2 Estimation of the Variance of the Predictor The B L U P is a random variable, and estimation of its variance will be needed to build a confidence interval around this predictor. The individual B L U P of 9{ can be defined by: \u00C2\u00A7i = Xi/3 + e; (3.46) where X; = [Sli/(V\u00C2\u00AE)1/3 \u00E2\u0080\u0094 Sli]. Because 9i is a predictor of 9it the variance of 6{ can be expressed as: EiiBi-Bif] = E[(ynl3 + ei - Xi3 - erf] = EidxiP-Xi^ + ei-ei)2] Chapter 3. Statistical Methods 62 = E[(xi$ - X i/3) 2] + E[el\ + E[4\ + 2E[(xi/3 --2E[(xiJ3 - Xif3)ei] - 2E[eiei] = var(x;/3) + var(ej) + var(ei) + 2cov(xi/3, hi) \u00E2\u0080\u0094 2cov(0;, e^ ) (3.47) The various components can be estimated separately, using the mathematics of linear combinations. First, it can be shown that the covariance between 3 and e is zero which means that cov(x;/3, \u00C2\u00A7;) is zero. cov(/3,e) = cov ( [X'$- 1 X]- 1 X'$- 1 ^, 0 2 V * _ 1 [ 0 - X/3]) = cov([x ,$- 1x]- 1x'$- 1^, { ^ v s - ^ i - x p r ^ x ^ x ' s - 1 ] } ^ = [X'#- 1X]- 1X ,#- 1var(fl){^ 2V*- 1 [I - X p C ' ^ X l ^ X ' * - 1 ] } ' = [ X ' S ^ X I ^ X ' S - ^ I J I - ^ - ^ [ x ' * - 1 ^ - ^ ' ] * - 1 ^ } = [ X ' s - ' x ] - 1 ^ * - V 2 v - s ^ x p c ' s ^ x i ^ x ' s - V ' v ] = [ x ' * - 1 ^ - 1 ^ * - 1 ^ - [X'#- 1x]- 1x , *-V 2v = 0 (3.48) The variance of the B L U P of e; can be derived from the definition of e in (3.36). Let {02V*I>-1}; be the i t h row of the matrix defined by ^> 2V$ - 1. Then the variance of e; is: var(8i) = var({^ 2 V*- 1 } i [f l -X / S]) = {(^ 2 V#- 1 } i var(^-X /3){^ 2 V$- 1K = .{0 2V*- 1} ivar([I-X[X ,*--- 1X]- 1X ,*- 1]fl){^ 2V*- 1K = {^V*\" 1 }^ - X[X'* - 1 X]- 1 X ,#- 1]var(tf) .[i - x t x ' * - 1 ^ - 1 ^ * - 1 ] ^ ^ * - 1 } ; = {^v* - 1 }^ - x[x , #- ' 1 x]- 1 x'$- 1 ]$ Chapter 3. Statistical Methods 63 \u00E2\u0080\u00A2[i - $- 1 x[x'#- 1 x]- 1 x , ]{^ 2 v#- 1 } i }; = { ^ V S - ^ i ^ - X p C ' ^ ^ X j ^ X ' l ^ V * - 1 } ; (3.49) To define the covariance between 0; and e;, it will be useful to express in terms of Y . From equation (3.38), 9{ can be written as: = [ ^ - { ^ V ^ - ^ X ^ + ^ V * - 1 } ^ = [xj - { ^ V ^ - ^ i X K X ' * - 1 ^ - 1 ^ * - 1 ^ + {^V*-1}^ = ([X i - {^ 2 v#- 1 } i x][x '$- 1 x]- 1 x'*- 1 + { ^ V * \" 1 } , ) ( F ' n ~ 1 F ) - 1 F , n _ 1 Y = T i Y (3.50) where Ti = ([x; - { ^ V * - 1 } ^ ] ^ ' * - 1 ^ - ^ ' * \" 1 + { ^ V * \" 1 } ; ) ( F ' n ^ F J ^ F ' n - 1 . Therefore, the covariance of and e; is cov(0;,e;) = cov(riY,ei) = cov(r i(FX/5 + Fe + G r + v),ei) = cov(riFe,ej) = cov(riFjej, e,-) = r i F i var(e i ) (3.51) The variance of e{ is given by var(ei) = 0 2 / ( ^ \u00C2\u00B0 ) 1 / 3 (3.52) The term var(Xj/3) is easily obtained from the variance of 0: var(x;/3) = x;var(/3)x-= X i f X ' S ^ X ] - 1 ^ (3.53) Chapter 3. Statistical Methods 64 Using these different estimators, the variance of 8{ can then be defined as E[{\u00C2\u00A7{ - Q i f ] = x j x ' s ^ x ] - 1 ^ + o ' v * - 1 } ^ * - x p c ' s ^ x i ^ x ' K ^ v * - 1 } ; . +(1 - 2TiFi)2V. The best strategy for estimating the unknown variances is to use the two-stage struc-ture of the model; each stage can be dealt with independently. In the first stage, estimates of 8, 82 and cr2,... based on the observed growth rates are obtained. In the second stage, the variable 8 is unobservable, but the unconditional distribution of 8, an estimate of 8 from the first stage, can be derived. From this distribution, estimates of 8 and 2 are computed. 3.3.2.1 First-Stage Estimators Different estimators for error component models have been suggested in the literature. However, none of these estimators were found adequate for the assumptions made in the Chapter 3. Statistical Methods 65 first stage. Instead, theoretical estimators had to be derived based on moment estimators of the unobserved residuals. It is a common practice to replace the unobserved residuals by their OLS estimates in theoretical estimators (Wallace and Hussain 1969, Zellner 1962, Swamy 1971, among others). For small sample sizes, when replacing the unobserved residuals by their OLS estimates in an estimator, it is probably worthwhile to adjust the degrees of freedom. A) Estimator of o~2 Following Wallace and Hussain (1969), the estimator for a2 can be derived using the squared residuals from period p with the following reasoning: E \u00C2\u00A3 n i p Sr ip E E \u00C2\u00A3 i E (5 p r p + E i v i p ) V i p Sr ip E f Vi\u00C2\u00BB ~ V' ^ p v pJ 2 = E = E ^ ^ p - v (3.55) Therefore, the theoretical estimator for a2 can be 1 Jp i Hp E i u\u00C2\u00BBp Sr> (3.56) which, after correcting for the degrees of freedom and replacing the unknown residuals by their OLS estimates ri ; p , becomes 1 E u. E t ^tp tp 'V J (3.57) B) Estimator of 82 Chapter 3. Statistical Methods 66 U s i n g a s imi la r s t rategy as for a2, an est imator for 8 can be der ived us ing the squared residuals , the estimates of p = ^ Y Y sPSQX p p>q (3.58) (3.59) (3.60) F r o m the s u m of these three est imators , the es t imator of 82 is easi ly ob ta ined 1 v CY Y Y Y uiruj* - Y sPaD (3.61) i P J 9 where u = W + E p <&=fi1 + E P Eq>P SpSqXm. Afte r correct ing for the degrees of freedom a n d replacing the u n k n o w n residuals by their O L S residuals , the es t imator becomes 1 (qnP)#o (3.62) C ) E s t i m a t o r of 9 T h e estimates of a2 a n d 82 are replaced i n il to get the est imate il. T h e F G L S est imate of 9 becomes: 9 = [ F ' n ^ F J ^ F ' n ^ Y (3.63) A g a i n , the double hat above 9 impl ies that i t was ob ta ined after first es t imat ing il. T h e F G L S est imate is a s y m p t o t i c a l l y n o r m a l l y d i s t r ibu ted . W h e n W is large, the n o r m a l Chapter 3. Statistical Methods 67 distribution should be a relatively good approximation. The approximate distribution of \u00C2\u00A7\6,ia: l\8 ~ N(8, [F'fWF]-1) (3-64) The inversion of matrix fi, a W X W matrix, can be simplified using matrix inversion rules. The inversion of fi is equivalent to: f T 1 = [^GAG' + S] - 1 . = I T 1 - S ^ G ' G'S 1G + l/S2A~1 G'S 1 (3.65) G'S G + 1/S2A' which requires the inversion of S, a diagonal matrix, and A and two M x M matrices. 3.3.2.2 Second-Stage Estimators In the second stage, an estimate of the unconditional distribution of 6, the GLS estimate of 6 is needed. Unfortunately, the GLS estimate of 8 cannot be obtained because the covariance matrix of the error terms is unknown. This difficulty can be overcome using the unconditional distribution of 8, a FGLS estimate of 9. From the conditional distribution of 8 given 8, shown in (3.64), and the distribution of 8, shown in (3.30), the unconditional distribution of 8 is easily found: 8 ~ N(X0, [F'fV'F]-1 + 2V) (3.66) In the second-stage, estimates of 8 and (f)2 are needed. GLS estimates are not possible because they would depend on the unknown value of each other. An iterative FGLS estimate is then required. Theoretically, the GLS estimate of 0 would be J3 = [X'[2V + (F'n~ 1 F)- 1 ] - 1 X]- 1 XV 2 V + (F'fV'F)-1]-1^ (3.67) Chapter 3. Statistical Methods 68 which is also the maximum likelihood estimate since the error terms are assumed to be normally distributed. Similarly, a possible estimate of c62 is its maximum likelihood estimate. However, equating the differential of the pdf of 6 with respect to 2 with zero does not yield a simple solution. No satisfactory estimate was found in the literature. An alternative estimate had to be derived. Based on the following reasoning: E E[(6 -X0)'(6 -X0)] = tr[cA2V + ( F ' f i - ' F ) - 1 ] = c62tr[V] + t r K F ' n ^ F ) \" 1 ] (9 - X0)'(e - X0) - t rKF'SWF)\" 1 ]\" tr[V] (3.68) the following estimate was found adequate 12 N [6 - X0)'{6 - X/3) - trKF'iWF)-- 1 ] ^ = (iV^2) iriVJ ( 3 ' 6 9 ) Because the unknown 0 is replaced by an estimate, 2 degrees of freedom are lost for the two coefficients 0\ and 02-Since both estimates, <^ 2 and fi, depend on the unknown value of the other parameter, FGLS estimates are needed. To find these FGLS estimates, an iterative solution must be used. The algorithm stops when both $ and 2 is used to get the FGLS estimate of 0: $ = [x'[02v + ( F ' n 1F)-1]-1x]-1x,[^2v + [F'n 1F]~1]-1e (3.70) Chapter 3. Statistical Methods 69 3.3.3 The Best Linear Unbiased Predictor and Its Variance The estimates of o~2, S2, (f>2, and 0 can now be replaced in the BLUP, defined in equation (3.39), to obtain what could be called an \"estimated best linear unbiased predictor\". Equation (3.39) becomes: \u00C2\u00A7 = ^vKF ' f t^F)- 1 + 4>2v\-H + ( F ' n ^ F ^ K F ' n ^ F ) - 1 + ^ V ^ X / f l (3.71) Similarly, an estimate of the variance of the BLUP of di is obtained by replacing the unknown parameters in (3.54) by their estimates. This yields: E[(\u00C2\u00A7i - 6i)2] = X i p c ' * - 1 ^ - 1 ^ + {^v*-1}^* - x tx ' ^x j^x 'K^v*- 1 } ; +(1 - 2 f i F i ) ^ 2 / ( ^ \u00C2\u00B0 ) 1 / 3 (3-72) where $ = (F'fi 2V ti = ({^-{^V^-^iX^X'^X]-^'^-1+ {4>2V^'1}i)(Fn 1 F ) - 1 F ' n 1 The estimate of the BLUP of 9i can now be used to predict future growth rates and compound the current volume of a permanent sample plot. Confidence intervals around the prediction can be built with the variance estimate. 3.4 Model Predictions The annual growth rate on a given plot can be predicted from model (3.46) 6i=Xii3 + -ei (3.73) Once the unknown parameters 0 and e; have been estimated, the instantaneous growth rate can be predicted from two plot attributes: the site index and the initial Chapter 3. Statistical Methods 70 volume. The growth rate can be assumed to be constant in a given year and to be chang-ing from year to year. The volume can thus be updated annually. Each updated volume becomes the new initial volume, [i.e. the volume at the beginning of the next growing season). The growth rate that will occur in this particular year is then predicted leading to a new updated volume. This process is repeated, say, K times to get the volume K years from now. This recursive system offers more flexibility to predict future volumes because the growth rate is updated annually. Assuming a constant instantaneous growth rate over a period of one year is close to reality. The annual volume update involves computing the expected value of the antilog of the predictor 9{. The technique for estimating this expected value is discussed in section 3.4.1. If a confidence interval is desired around the predicted volume, the confidence interval around each intermediate volume will be needed too since the volume at the end of the growing season is a function of the predicted volume at the beginning of the season. In section 3.4.2, the problem of variance estimation is analyzed, and a solution is suggested. 3.4.1 Volume Predictor Various techniques exist to estimate the expected value of the antilog of an estimate (Ung and Vegiard 1988). However, no technique was found to estimate the expected value of the antilog of a predictor. The subtle difference can be explained by recalling that 0i is a predictor of 9{, a random variable, while the expected value of an estimate is a constant. The technique used by Baskerville (1972) was simple enough to be adapted for estimating the expected value of the antilog of a predictor and was, consequently, the chosen technique. Future volume prediction is driven by the initial volume in the first iteration. The last volume measurement observed on the sample plot can serve as the initial volume since Chapter 3. Statistical Methods 71 this volume represents the best knowledge of initial volume to predict future growth rates. A possible estimator of the updated volume could be the volume at the beginning of the growing season multiplied by the antilog of the predicted growth rate: V/+l = V? x exp$] (3.74) where V? is the predicted volume on plot i at iteration j and Q\ is the predicted growth rate. However, this estimator is biased. Based on the moment generating function, it can be shown that E[X] = exp which can be approximated by E[X] = exp E[\nX] + var(ln X) E[\nX] + v a r i (lnX) (3.75) (3.76) Goldberger (1968) pointed out that this result was asymptotically unbiased.2 The true model expressed in a logarithmic scale is lnV?'+ 1 =lnV? +61 + uf (3.77) where u\ is the error term on plot i at iteration j. The first term on the right-hand side of the equation can also be rewritten in a similar form. Hence, the true model really is lny / + 1 =lnVf-fE^ + Eu* fc=0 k=0 (3.78) Given that the volume is estimated on plot i and that annual observations are used, if the first volume is assumed to be measured without error, then i E[\nVi0+1}=\xvV? + Y2 was 0. The second-stage was the poorest Chapter 5. Discussion 97 estimator of observed rates. To make comparison possible with the first-stage model, the second-stage model was used to predict not the individual observations. Better results could have probably been achieved by using the observed initial volumes and by trying to predict the observed rates. Another reason for the weak results in the second stage is that it could not predict extreme growth rates. This poor performance was also expected since the number of plots was not large and the estimators of 3 and, especially, (j)2 were based on large sample theory. Surprisingly, the second-stage model could predict future growth rates with respectable precision but a rather large bias. Overall, the second-stage model underestimated the volume, suggesting that the initial volume used to fit the model was possibly underestimated. Writing the code for the prediction of future volumes was simple using a matrix al-gebra language. When no variance estimates are computed, the predictions are rapidly obtained. However, computing the confidence intervals of the predictions can be cum-bersome. The overall results of the weighted predictor were always between the overall results from the first-stage and the second-stage models. The weighted predictor also did poorly in estimating the observed mean growth rates. This is not surprising since it is derived from the best estimate for observed rates, the first-stage estimate, to which a certain quantity was added. Combining the first- and the second-stage was of no use for fitting the observed growth rates. If fitting past observations had been the only purpose, the first-stage estimator would have been sufficient. When future volumes were considered, the weighted predictor was more accurate than either of the individual estimates. The variation of the prediction errors was considerably better than the first-stage estimate and slightly smaller than the second-stage model. The weighted predictor was by far the best option for predicting future growth rates. Analyzing the detailed results, it is seen that the relative bias was more important for predictions on plots less than 20 years of age initially. For projection time less than Chapter 5. Discussion 98 10 years, on plots with initial ages less than 20 years, the bias was close to the negative value of the mean absolute deviation, suggesting that the weighted predictor consistently overestimated the actual volumes. However, because these plots have less volume, they carry only a small weight when the overall bias is computed. For plots with initial age between 31 and 40 years, when only the projection time was analyzed, there was a strong bias. This probably occurred by chance, no rationale could be found to explain such a behaviour. Curiously, the length of the total lag did not seem to have an effect on the accuracy of the predictions. The projection time had a stronger influence on accuracy. The precision of the estimates was dependent on both the projection time and the total lag. There was a major difference in precision between plots with initial age less than 40 years and plots with initial age over 40 years. This probably indicated where the curve of growth rate over time for Douglas-fir in the data set began to level off. For plots with initial age over 40, the variability was less than 10% most of the time. For these plots, the first-stage estimate yielded precise estimates which explains the good precision of the weighted predictor. The mean absolute deviation is the arithmetic average of the error terms in absolute value. For plots over 40 years of initial age, it was usually less than 5%, even for projection times greater than 11 years. It was also relatively larger on young plots. This lack of precision was mainly caused by the first-stage estimates, as confirmed by looking at the results in detail. The precision of the second-stage estimator was relatively constant. Chapter 6 Conclusions Statistical models based on a collection of remeasured entities are common in a few scientific fields. However, permanent sample plots in forestry have particular complexities that models developed in other fields do not address. Adequate techniques must be developed to predict growth in forestry permanent sample plots. These techniques must be adaptable to plots with only a few measurements since this is the situation most often encountered in Canada. Above all, these techniques must rely on a clear understanding of the purposes of permanent sample plots, and their possibilities and limitations as sources of information. The framework presented in this dissertation reflects this understanding. Growth rates were explained by site index and initial volume in a model inspired from Von Bertalanffy's growth curve. The individual relationship characterizing a sample plot was expressed as a random variation around a general model. Observations on a plot were assumed to be random outcomes of the plot's individual model. The current volume of a plot was updated using the growth rate predicted from the individual model. This technique was shown to be more precise than using the general model common to all plots. The weighted predictor also emphasizes the importance of the error structure in forestry statistical models, as pointed out by Gregoire (1987). Forestry information is not equally valuable; this fact must be recognized in a model. Separating the error term into components facilitates assigning a specific weight to each piece of information. 99 Chapter 6. Conclusions 100 A greater weight must be given to reliable information when computing an estimator. Precision of an estimator is as important as its bias. However, too often foresters have not paid attention to precision. Accurate mathematical models are sought to the detri-ment of statistical models. The power of statistics lies in its recognition of the difference between a variable and a constant. Statistical models will always be more precise than mathematical models when the error structure is correctly specified. Because models are used to make decisions, precise models are needed to make well-informed decisions. The results of the numeric example illustrated the advantages of the weighted pre-dictor. Despite the approximations made in both stages, the weighted predictor was reasonably accurate and precise. It is probable that using a larger number of plots would have resulted in more precision with the second-stage model and consequently to the weighted predictor as well. Permanent sample plots with a longer history would also have provided more precise results. Fortunately, the weighted predictor can be adapted to data sets where all sample plots have more than three measurements. The weighted predictor is certainly not perfect. A major problem with the weighted predictor is the absence of statistical theory dealing with random coefficient modeling in finite samples. The estimator is approximately precise, but how approximate is un-known. The first-stage approximation can be justified for a short period of time only. No stochastic term is included for predicting catastrophic losses. Despite these weaknesses, the weighted predictor can provide useful information for solving problems facing many forest agencies in Canada. The proposed framework, which led to the weighted predictor, opens new horizons for modeling using permanent sample plots. Developments in the theory related to ran-dom coefficient models in finite samples can be expected in the near future and should answer many important questions. New computer technology will facilitate large simu-lation problems. Research on techniques to include probabilities of catastrophic losses in Chapter 6. Conclusions 101 prediction could improve the decision-making process based on growth models. Foresters wil l benefit from these improvements if they adopt the right framework. Perfect models might still be far ahead, but we are definitely getting closer. References Cited Aitchison, J. and Brown, J. A . C. 1957. The lognormal distribution. Cambridge Univ. Press, Cambridge, England, 176 p. Aitken, A. C. 1935. On least squares and linear combinations of observations. Proceed-ings of the Royal Society of Edinburgh 55:42-48. Amemiya, T. 1978. A note on random coefficient model. Inter. Econom. Rev. 19:793\u00E2\u0080\u0094 796. 1971. The estimation of the variances in a variance-component model. In-ternat. Econom. Rev. 12:1-13. Arner, S. L. and Seegrist D. W. 1980. Missing: a computer program for the maximum likelihood estimates of the parameters of the multivariate hnear model with incomplete measurements. USDA For. Serv. Gen. Tech. Pap. NE-56, 19 p. 1979. A computer program for the maximum likelihood estimator of the general multivariate Hnear model with correlated errors. USDA For. Serv. Gen. Tech. Pap. NE-51, 10 p. Arora, S. S. 1973. Error components regression models and their apphcations. Ann. Econom. Social Measurement 2:451-461. 102 References Cited 103 Axelsson, E. and Axelsson, B. 1986. Changes in carbon allocation patterns in spruce and pine trees following irrigation and fertilization. Tree Physiology 2:189-204. Baltagi, B. H. , 1981. Pooling: an experimental study of alternative testing and estimation procedures in a two-way error component model. J . Econometrics 17:21-49. Baskerville, G. L. 1972. Use of logarithmic regression in the estimation of plant biomass. Can. J . For. Res. 2:49-53. Bayes, T. R. 1763. An essay towards solving a problem in the doctrine of chances. Philo. Trans. Roy. Soc. 53:370-418 (reprinted in Biometrika 45:296-315). Berkey, C. S. 1982. Bayesian approach for a nonlinear growth model. Biometrics 38:953-961. Biging, G. S. 1985. Improved estimates of site index curves using a varying-parameter model. For. Sci. 31:248-259. Buckman, R. R. 1962. Growth and yield of red pine in Minnesota. USDA For. Serv. Tech. Bull. 1272, 50 p. Carter, R. L. and Yang, M . C. K. 1986. Large sample inference in random coefficient regression models. Comm. Statist. A:Theory-Methods 15:2507-2525. Chapman, D. G. 1961. Statistical problems in population dynamics, In Proc. Fourth Berkeley Symp. Math. Stat, and Prob., Univ. Calif. Press, Berkeley, p. 153-168. References Cited 104 Clutter, J . L. 1963. Compatible growth and yield models for loblolly pine. For. Sci. 9:354-371. Curtis, R. 0. 1967. A method of estimation of gross yield of Douglas-fir For. Sci. Monogr. 13, 24 p. Davis, L. S. and Johnson, K. N . 1987. Forest management. Third edition. MacGraw-Hill Book Company, New York, 790 p. Davis, A. W. and West, P. W. 1981. Remarks on \"Generalized least squares estimation of yield functions\" by I. S. Ferguson and J . W. Leech. For. Sci. 27:233-239. Dempster, A. P., Laird, N . M . , and Rubin, D. B. 1977. Maximum likelihood from incomplete data via the E M algorithm.(with discussion) Jour. Roy. Stat. Soc. series B 39:1-38. Dhrymes, P. J . 1971. Equivalence of iterative Aitken and maximum likelihood estimators for a system of regression equations. Austral. Econom. Papers 10:20-24. Dielman, T. E. 1989. Pooled cross-sectional and time-series data. Marcel Dekker Inc., New York, 249 p. _ 1983. Pooled cross-sectional and time-series data: a survey of current statis-tical methodology. Am. Stat. 37:111-122. Fearn, T. 1975. A Bayesian approach to growth curves. Biometrika 62:89-100. Ferguson, I. S. and Leech, J . W. 1978. Generalized least squares estimation of yield functions. For. Sci. 24:27-42. References Cited 105 1981. Reply to remarks by A. W. Davis and P. W. West on \"Generalized least squares of yield functions.\" For. Sci. 27:589-591. Fuller, W. A. and Battese, G. E. 1974. Estimation of linear models with crossed-error structure. J . Econometrics. 2:67-78. Furnival, G. M . and Wilson, R. W. Jr. 1971. Systems of equations for predicting forest growth and yield. In Statistical Ecol. 3:43-57. G. P. Patil, E. C. Pielou, W. E. Waters, eds. Penn State Univ. Press, University Park. Garcia, 0. 1983. A stochastic differential equation model for the height growth of forest stands. Biometrics 39:1059-1072. Geisser, S. 1970. Bayesian analysis of growth curves. Sankhya A 32:53-64. Gerhardt, E. 1930. Ertragstafeln fur reine und gleichartige hochwaldbestande von eiche, buche, tanne, fichte, kiefer, gruner Douglasie und larche. Second edition. Springer-Verlag, Berlin. Goldberger, A. S. 1968. The interpretation and estimation of Cobb-Douglas functions. Econometrica 36:464-472. Gregoire, T. G. 1987. Generalized error structure for forestry yield models. For. Sci. 33:423-444. Grizzle, J . E . and Allen, D. M . 1969. Analysis of growth and dose response curves. Biometrics 25:357-381. References Cited 106 Gumpertz M . and Pantula, S. G. 1989. A simple approach to inference in random coefficient models. Am. Stat. 43:211-215. Harville, D. 1976. Extension of the Gauss-Markov theorem to include the estimation of random effects. Ann. Stat. 4:384-395. Hui, S. L. and Berger, J. 0. 1983. Empirical Bayes estimation of rates in longitudinal studies. J . Am. Stat. Ass. 78:753-760. Judge, G. G., Griffiths, W.E. , Hill , R.C., Liitkepohl, H. , and Lee, T. C. 1985. The theory and practice of econometrics. Second edition. John Wiley and Sons, New York, 1019 p. Keyes, M . R. and Grier, C. C. 1981. Above- and below-ground net production in 40-year old Douglas-fir stands on low and high productivity sites. Can. J . For. Res. 11:599-605. King, J . E. 1966. Site index curves for Douglas-fir in the Pacific Northwest. Weyerhaeuser Co. For. Res. Center, For. Paper 8. 50 p. Kmenta, J . and Gilbert R. F. 1970. Estimation of seemingly unrelated regressions with autoregressive disturbances. J . Am. Stat. Ass. 65:186-197. Kurz, W. A. and Kimmins, J . P. 1987. The influence of site quality on tree resource allocation to fine roots and its effect on harvestable productivity of coastal Douglas-fir stands. Can. For. Service, Victoria, B .C. F R D A report 034, 102 p. (+ appendices) References Cited 107 Laird, N . M . and Ware, J . H . 1982. Random-effects models for longitudinal data. Bio-metrics 38:963-974. Leak W. 1966. Analysis of multiple systematic remeasurements. For. Sci. 12:69-73. Lee, L. F. and Griffiths, W. E. 1979. The prior likelihood and best linear unbiased prediction in stochastic coefficient linear models. Univ. of New England Working Papers in Econometrics and Applied Statistics No. 1, Armidale, Australia. LeMay, V . 1988. Comparison of fitting techniques for systems of forestry equations. Unpublished Ph.D. dissertation, Univ. of B. C , 172 p. Lindley, D. V . and Smith, A. F. M . 1972. Bayes estimates for the linear model. J . Roy. Stat. Soc. B 34:1-41. Maddala, G. S. 1977. Econometrics. McGraw-Hill, New York, 516 p. 1971. The use of variance components models in pooling cross-section and time-series data. Econometrica 39:341-358. Maddala, G. S. and Mount, T. D. 1973. A comparative study of alternative estimators for variance component models used in econometric applications. J . Am. Stat. Ass. 68:324-328. Marshall, P. L. and Jahraus, K. V . 1987. Growth and yield availability and use in British Columbia. Unpublished report, 96 p. References Cited 108 Mundlak, Y . 1978. Models with variable coefficients: integration and extension. Ann. de 1'INSEE 30-31:483-509. _ 1961. Empirical production functions free of management bias. J . Farm Econom. 43:44-56. Parks, R. W. 1967. Efficient estimation of a system of regression equations when dis-turbances are both serially and contemporaneously correlated. J. Am. Stat. Ass. 62:500-509. Phillips, P. C. B. 1985. The exact distribution of the SUR estimator. Econometrica 53:745-756. Pienaar, L. V . and Turnbull, K . J . 1973. The Chapman-Richards generalization of Van Bertalanffy's growth model for basal area growth and yield in even-aged stands. For. Sci. 19:2-22. Pothoff, R. F. and Roy, S. N . 1964. A generalized multivariate analysis of variance model useful especially for growth curve problems. Biometrika 51:313-326. Rao, C. R. 1987. Prediction of future observations in growth curve models. Stat. Sci. 2:434-471. . 1975. Simultaneous estimation of parameters in different hnear models and applications to biometric problems. Biometrics 31:545-554. 1972. Estimation of variance and covariance components in hnear models. J . Am. Stat. Ass. 67:112-115. References Cited 109 _ 1965. The theory of least squares when the parameters are stochastic and its application to the analysis of growth curves. Biometrika 52:447-458. Rao, U . L. G. 1982. A note on the unbiasedness of Swamy's estimator for the random coefficient regression model. J . Econometrics 18:395-401. Richards, F. J . 1959. A flexible growth function for empirical uses. J . Exp. Bot. 10:290\u00E2\u0080\u0094 300. SAS Institute Inc. 1982. SAS user's guide: statistics, 1982 edition. SAS Institute Inc., Cary, NC, 584 p. Schmidt, P. 1977. Estimation of seemingly unrelated regressions with unequal number of observations. J . Econometrics 5:365-377. Schumacher, F. 1939. A new growth curve and its applications to timber-yield studies. J . For. 37:819-820. Seegrist, D. W. and Amer, S. L. 1978. Statistical analysis of linear growth and yield models with correlated observations from permanent plots remeasured at fixed intervals, p 209-233 In Growth models for long-term forecasting of timber yields, (eds. J . Fries , H . Burkhart, T. A. Max), Sch. For. Wild. Res. VPISU, Publ. FWS 1-78. Smith, A. F. M . 1973. A general Bayesian linear model. J . Roy. Stat. Soc. B 35:67-75. Srivastava, V . K. 1973. The efficiency of an improved method of estimating seemingly unrelated regression equations. J . Econometrics 1:341-350. References Cited 110 Sullivan A. D. and Clutter, J . L. 1972. A simultaneous growth and yield model for loblolly pine. For. Sci. 18:76-86. Sullivan A. D. and Reynolds M . R. Jr. 1976. Regression problems from repeated mea-surements. For. Sci. 22:382-385. Swamy, P. A. V . B. 1974. Linear models with random coefficients, In Frontiers in econo-metrics (ed. P. Zarembka), New York, Academic Press, p. 143-168. . 1971. Statistical inference in random coefficient regression models. Springer-Verlag, New York, 209 p. . 1970. Efficient inference in a random coefficient regression model. Econo-metrica 38:311-323. Swamy, P. A. V . B. and Arora, S. S. 1972. The exact finite sample properties of the estimators of coefficients in the error components regression models. Econo-metrica 40:261-275. Swamy P. A . V . B. and Metha J . S. 1979. Estimation of common coefficients in two regression equations. J . Econometrics 10:1-14. _ 1975. On Bayesian estimation of seemingly unrelated regressions when some observations are missing. J . Econometrics 3:157-169. 1973. Bayesian analysis of error components regression. J . Am. Stat. Ass. 68:648-658. References Cited 111 Swindel, B. F. 1968. On the bias of some least squares estimators of variance in a general Hnear model. Biometrika 55:313-316. Tait, D. E. N . , Cieszewski, C. J . , and Bella, I. E. 1988. The stand dynamics of lodgepole pine. Can. J. For. Res. 18:1225-1260. Taylor, W. E. 1980. Small sample considerations in estimation from panel data. J . Econometrics 13:203-223. Ung, C. H. and Vegiard, S. 1988. Problemes d'inference statistique rehes a la transfor-mation logarithmique en regression. Can. J . For. Res. 18:733-738. Wallace, T. D. and Hussain, A. 1969. The use of error component models in combining cross-section with time-series data. Econometrica 37:55-72. ZeUner, A. 1987. An introduction to Bayesian inference in econometrics. R. E. Krieger Publishing Co., Malabar, FL , 431 p. . 1962. An efficient method for estimating seemingly unrelated regressions and tests for aggregation bias. J . Am. Stat. Ass. 57:38-368. Zellner, A. and Vandaele, W. 1975. Bayes-Stein estimators for k-means, regressions, and simultaneous equation models In Studies in Bayesian Econometrics and Statistics (ed. S. E. Fienberg and A. Zellner), Amsterdam, North HoUand, p. 62 7-653. "@en . "Thesis/Dissertation"@en . "10.14288/1.0075366"@en . "eng"@en . "Forestry"@en . "Vancouver : University of British Columbia Library"@en . "University of British Columbia"@en . "For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use."@en . "Graduate"@en . "Growth prediction of recent permanent sample plots for forest inventory projection"@en . "Text"@en . "http://hdl.handle.net/2429/31159"@en .