Abstract

Forest ecosystem dynamics are driven by a complex array of simultaneous cause-and-effect relationships. Understanding this complex web requires specialized analytical techniques such as Structural Equation Modeling (SEM). The SEM framework and implementation steps are outlined in this study, and we then demonstrate the technique by application to overstory-understory relationships in mature Douglas-fir forests in the northwestern USA. A SEM model was formulated with (1) a path model representing the effects of successively higher layers of vegetation on late-seral herbs through processes such as light attenuation and (2) a measurement model accounting for measurement errors. The fitted SEM model suggested a direct negative effect of light attenuation on late-seral herbs cover but a direct positive effect of northern aspect. Moreover, many processes have indirect effects mediated through midstory vegetation. SEM is recommended as a forest management tool for designing silvicultural treatments and systems for attaining complex arrays of management objectives.

1. Introduction

Ecosystem processes involve complex interactions of many cause-and-effect relationships [1]. With sufficient understanding of these processes and relationships, forest management activities can be designed to produce specific outcomes and desired results, within the limits of natural ecosystem dynamics. The relationship between forest structure and function assumes particular importance because silvicultural knowledge exists to produce a variety of structural outcomes in many forest types [2]. Structural diversity is widely accepted as an important driver of biodiversity despite the fact that many mechanisms linking the two are still poorly understood. Classical silvicultural systems such as shelterwood with reserves [3] and innovations like variable-retention harvesting [4] have been proposed as part of a strategy to meet diverse forest management objectives in late-seral Douglas-fir (Pseudotsuga menziesii (Mirb.)) forests on public lands in the Pacific Northwest USA. The heterogeneous structures retained in the proposed systems have been hypothesized to maintain taxa and ecological processes characteristic of late-seral forests in the short-term and accelerate their recovery in the long-term [5]. The success of these systems depends on complex responses to prescribed management activities, on our understanding of the cause-and-effect pathways of these responses, and on our ability to utilize this understanding to design efficient systems for meeting defined objectives. However, for silvicultural strategies that retain late-seral or old-growth attributes, routine statistical methods may not be sufficient for understanding the mechanisms linking stand structure to biodiversity and for testing the interactions of various taxa and ecological processes [6, 7]. In the context of silvicultural management of biodiversity, it is crucial that research on ecosystem responses employ statistical methods that imply causal mechanisms among structural components and processes that operate simultaneously through complex and often indirect pathways.

Manipulative experiments based on the principle of randomization are the historical foundation for inferring causality [8, 9]. Analysis of Variance (ANOVA) models ascribe treatment effects to orthogonal factors in these designed experiments [10]. The effects represent the net influence of factors on a response variable, but do not distinguish between the direct and indirect pathways that link cause and effect. ANOVA has the potential for making causal inference, but lacks the capability to uncover much of the information about underlying processes or response mechanisms [10].

Most data from forest ecosystems are observational and multivariate. Data reduction by categorization (classification) or by creation of synthetic continuous variables (ordination) is very efficient for identifying general patterns in the data from an exploratory and descriptive perspective [11]. Their capability for understanding functional links is limited because they cannot explicitly incorporate complex interrelationships between variables based on a priori knowledge [12]. However, the methods are useful for discovering underlying trends and for supplementing the process of model building and hypothesis formulation. Lastly, multiple regression can test whether an estimated coefficient is significantly different from zero in a statistical sense, but cannot draw causal conclusions from the test [13]. Additionally, regression models ignore the possibility that a predictor may indirectly influence the response through other predictors [13]. In short, multiple regression allows us to predict but not to explain [13].

Structural Equation Modeling (SEM) is an alternative method for testing our understanding of complex ecological processes. SEM is a collection of procedures that tests hypothesized relationships among observed variables [14, 15]. Complex interactions are first translated into a network of directional paths linking variables and are then evaluated against multivariate data [16]. These paths postulate direct and indirect effects among ecosystems components, as well as spurious associations between variables that may be attributable to common causes. A direct effect describes direct regulation of a response variable (effect) by a causal variable, while an indirect effect implies that the regulation is mediated through other variables. Hence, SEM is often related to causal modeling [17]. It is philosophically a confirmatory data analysis, but its application extends to testing alternative a priori models or to model building [18], and can therefore be regarded as blending confirmatory and exploratory analyses [19]. The key to successful SEM rests on the competence of a researcher to posit initial cause-and-effect models drawing from accumulated knowledge, prior experience, and published results.

Comprehensive assessments by Shipley [13], Pugesek et al. [20], and Grace [12] have brought SEM into the context of natural systems. Laughlin and Abella [21] and Laughlin et al. [22] applied SEM to an observational study of abiotic and biotic factors influencing plant community composition and species richness in a ponderosa pine (Pinus ponderosa P. and C. Lawson) forest ecosystem. Youngblood et al. [23] studied the effects of experimental thinning and prescribed burning treatments on mortality of ponderosa pine with ANOVA and SEM. They found that SEM was able to represent and detect cascading effects of several fire-related factors on mortality, a process that the ANOVA model could not represent. These studies have demonstrated the flexibility of SEM for accommodating a variety of study designs. With increasing recognition of the benefits possible from embedding existing knowledge in data analysis and of the limitation of common statistical methods for large-scale ecological experiments [24], SEM may offer a promising alternative for investigating causal relationships.

SEM has been under development for decades in human science [25], but its application in the field of natural science has not been as widespread. SEM applications in natural science have provided only a brief account of the methodology (e.g., [26, 27]). In one exception, Mitchell [28] illustrated path analysis, one component of SEM, using data from a hummingbird pollination study, with few remarks on other important properties and capabilities. Grace [14] summarized general aspects of applying SEM in observational ecology studies with two examples. On an advanced level, Grace and Bollen [29] provided details about the SEM framework with emphasis on a specific issue related to composite variables. Grace et al. [25] extended the basic SEM framework to Structural Equation Meta-Models (SEMMs) to represent multiple orders of latent and composite variables. Lastly, Shipley [9] generalized SEM to accommodate multilevel data. To the best of our knowledge a succinct summary of SEM components and properties is lacking in the silvicultural literature. This type of summary could benefit scientists who wish to apply SEM to silvicultural research but are insufficiently familiar with the basic concepts. The goal of this paper was to cover major aspects of SEM that are essential for an introductory understanding of the methodology. Three specific objectives were to (1) systematically present the technical framework of SEM in a silvicultural context, (2) illustrate the application of SEM to understanding the mechanisms of overstory-understory relationships in mature Douglas-fir forests in the Pacific Northwest region of the United States, and (3) discuss the advantages of SEM in the context of this specific example. Topics were presented in a manner that facilitated easy access of the information to readers with different specific interests in the technique.

2. Structural Equation Modeling (SEM)

The two central components of SEM are the path model and the measurement model. The path model or path analysis quantifies specific cause-and-effect relationships between observed variables [30, 31]. The measurement model quantifies linkages between (1) hypothetical constructs that might be known but unobservable components of forest ecosystems and (2) observed variables that represent a specific hypothetical construct in the form of a linear combination. LISREL (LInear Structural RELations) was developed as a unifying and flexible mathematical framework to specify these linkages [32]. Hayduk [33] and Kelloway [34] described the framework in detail. The following summary is primarily drawn from Grace [12], Kline [19], and Schumacker and Lomax [15], and the presented mathematical framework is a summary of the LISREL system.

2.1. Path Model

Specification of a path model involves hypothesized cause-and-effect relationships between observed variables. These relationships are usually based on theoretical considerations or evidence from prior studies. However, certain conditions must be met for a variable to be designated as cause versus effect. Kenny [17] proposed three conditions: time precedence, functional relationship, and nonspuriousness. First, for variable to cause variable , has to precede in time, so time precedence implies an asymmetric relationship between the two. Second, and should be functionally related because there is no causal relationship if they are independent. Third, if the relationship between and is spurious due to a common cause, it will disappear once the common cause is identified and represented in the model. Pearl [35] reviewed recent advances in causal inference, and stressed the paradigmatic shift needed to move from traditional statistical analysis to causal analysis.

Consider a hypothetical model for representing the response of small mammal abundance at the forest floor to stand structural characteristics (Figure 1). In this scenario, overstory tree density and coarse woody debris are exogenous observed variables whereas understory tree cover, herbaceous cover, and abundance of small mammal species are endogenous variables. An exogenous variable is always considered only a cause; the causes of the exogenous variable itself are generally unknown or of little interest; therefore, it is not represented in the model. An endogenous variable is an effect, but it may also be a cause to other endogenous variables. A path model may look similar to multiple regression [13], where an exogenous variable is analogous to a predictor variable and an endogenous variable is the response; however, the difference is that the endogenous variable can be both a predictor and a response in a system of equations. An important assumption is that exogenous variables are measured without error. Any measurement errors, including those attributable to data entry, field recording, or other causes [33], could lead to bias in estimated path coefficients [16]. On the other hand, each endogenous variable is assumed to have error. The error term reflects all unmeasured causes of variation in the associated endogenous variable as well as true measurement error [19].

A unidirectional arrow is a path that represents a direct effect of the causal variable on a response variable, whereas a bidirectional arrow represents an unanalyzed association. Unanalyzed associations are customarily specified between exogenous variables because no hypothesis is included to explain why they covary. In general not all causal processes acting on a system are specified because the objective of the model has limited scope, so the model represents a simplification to address the specific objective [22]. A basic assumption underlying the path model is that a path represents a linear relationship between variables. Nonlinear relationship must be accommodated by transformations or higher-order terms.

Once the model is fitted to data, the path coefficient (direct effect) of each path is estimated and interpreted similarly to a regression coefficient. An indirect effect of any causal variable is estimated as the product of the chain of direct effects, and a total effect is the sum of all direct and indirect effects. In the example, overstory tree density is hypothesized to have two indirect effects on abundance, mediated by understory trees and herb cover. Thus, the latter indirect effect is the product of two path coefficients, one representing the effect of tree density on herb cover and the other the effect of herb cover on small mammal abundance. The total effect on abundance is then the sum of the two indirect effects.

Direct and indirect effects should be interpreted with caution [37]. As with multiple regression, an effect is interpreted as the change induced by fixing other variables in a model and changing only the subject variable. A direct effect would occur if all other variables in a model remained constant [13]. In estimating an indirect effect, all other variables in the model are controlled except for the mediating variables in the path representing the indirect effect of interest [13]. In other words, direct and indirect effects are interpreted as if they are orthogonal, but in reality variables from systems as complex as managed or unmanaged forests are multicollinear and require the same precautions as interpreting multiple regression coefficients under differing levels of multicollinearity.

Path coefficients can be expressed as standardized or unstandardized coefficients. Unstandardized coefficients are more intuitive because they represent the change in the response variable per unit change in the causal variables; that is, both response and causal variables remain on a scale consistent with original units of measurements. Furthermore, unstandardized results are important because significance tests on the coefficients are based on standard errors of the unstandardized solutions [12]. Standardized coefficients are typically computed either by first standardizing all variables (subtracting their mean and dividing by their standard deviation) and then computing the path coefficients, or by the ratios of the standard deviations of the variables. These coefficients are therefore expressed in units of standard deviations for the corresponding variables. Standardized coefficients allow direct comparison of the magnitude of effects of two causal variables measured on different scales [12].

2.2. Measurement Model

A common method for evaluating a measurement model is confirmatory factor analysis (CFA). CFA is a process of specifying the number and types of observed variables associated with one or more hypothetical constructs and analyzing how well the observed variables measure the constructs. A hypothetical construct is a conceptual variable which cannot be directly measured. Conversely, an observed variable can be measured and is used to infer the construct. Two ways of representing constructs are with latent and composite variables. A latent variable is a cause of its corresponding observed variables, whereas a composite variable is a collective effect of the variables [29]. For brevity, only latent variables are further discussed. Grace and Bollen [29] provided more details about the theory and application of composite variables.

Choice of the observed variables for a measurement model must consider the validity and reliability of the observed variables. Validity refers to the accuracy of an observed variable for representing the effect of a latent variable. Because a latent variable can be multifaceted, the observed variables selected should measure these different facets [38]. Reliability refers to the consistency in measurement of an observed variable or the amount of random measurement error. The idea is similar to estimating the precision of a measuring device by repeatedly measuring the same observation under similar conditions.

Consider a hypothetical measurement model for species diversity (Figure 2). Conventionally, a circle represents a latent variable and a rectangle represents an observed variable. Following Krebs [39], the model postulates that species diversity is a multifaceted concept with three observed variables measuring its different facets. The total number of species (observed and unobserved), which could be estimated from a Jackknife estimator [39, 40], measures the total species richness. Shannon diversity index measures the diversity and the evenness of a community. The Simpson evenness index measures only the evenness of a community, which is hypothesized to be measuring a different facet of diversity than the Shannon diversity index. The single arrowheads pointing from the latent variable to the observed variables are factor loadings. They represent the direct effects and are interpreted similarly to regression coefficients. The arrowheads to left of the observed variables depict the random measurement errors.

A general recommendation is to have three or more observed variables per latent variable to adequately account for measurement error and to meet the requirement for model identification, which ensures model convergence and proper solutions [19]. However, it may not always be possible to achieve the ideal of three or more observed variables; in this case, one or two observed variables per latent variable are still feasible with some constraints on model specification (see [19]). Observed variables should be both valid and reliable with respect to measurement of the latent variables. Hayduk [33] commented that observed variables with greater than 40% measurement error were prone to estimation problems. Bollen [16] provided further details on the use of latent variables.

2.3. Structural Regression Model

A structural regression (SR) model is a path model with latent variables, thus combining principles of path and measurement models. The goal is to take measurement errors of observed variables into account when evaluating a path model. This is the most general kind of core model that is widely applied in SEM [19]. A fully latent SR model has only latent variables in the path model whereas a partially latent SR model is a mix of observed and latent variables. Similar to path and measurement models, an important phase of analysis in an SR model is model identification. Model identification is a property that determines whether the model allows for unique parameter estimates. The two basic conditions for identification are (1) model degrees of freedom equal to or greater than zero () and (2) a known scale for every latent variable.

The total degrees of freedom for an SR model is , where is the number of observed variables. This total corresponds to the number of variances and unique covariances in a variance-covariance matrix for variables. The model degrees of freedom, , is the total degrees of freedom minus the number of estimated parameters. A just-identified model () will have unique parameter estimates, but an overidentified model () is desirable for model testing and assessment. An underidentified model () will not have unique solutions for all parameters. Empirical underidentification occurs when the effective degrees of freedom are reduced due to two or more observed variables that have a very high correlation (e.g., correlation ≥ 0.90), in turn leading to problems with parameter estimation [17, 19]. Kline [19] proposed two basic ways to deal with extreme collinearity: removing one of the variables or combining redundant ones into a composite variable. Scale is a property of each latent variable. Because a latent variable is not measurable, it must take on the scale or units of measure from one of its observed variables. A way to assign scale is by imposing a unit loading identification (ULI) constraint by fixing a factor loading of one observed variable at a value of 1.0.

An SR model that meets these conditions on scale and degrees of freedom does not guarantee identification. Bollen [16] suggested a two-step rule for checking identification of an SR model. The first step is viewing an SR model as a CFA model. The resulting CFA model is identified if it meets the following sufficient requirements and assumptions: (1) at least two observed variables per latent variable, (2) independence between measurement errors and latent variables, and (3) independence between measurement errors. The second step is then to check for recursiveness of the path model part of the SR model, ignoring any observed variables used to measure latent variables. A path model is identified if it meets the following requirements for recursiveness: (1) errors associated with endogenous latent variables are uncorrelated, and (2) all causal effects are unidirectional. If both the CFA model and path model are identified by these respective sets of requirements, then the whole SR model is identified, and model fitting can proceed.

The two-step rule from Bollen [16] can be relaxed, and some SR models that fail the rule can still be identified. One example is correlated measurement errors in a CFA model. Another common example is a nonrecursive path model, for example, a model with feedback loops or correlated errors, which requires special rules for identification. A special case is an SR model that has one observed variable per latent variable and therefore requires a priori assignment of measurement errors.

2.4. Model Estimation

The LISREL framework can be summarized into three matrix equations, two for the measurement model component and one for the path model component [12]. For the measurement model component, where is a vector of observed exogenous variables, and it is a linear function of a vector of exogenous latent variables and a vector of measurement error . is a matrix of factor loadings relating to . Similarly, is a vector of observed endogenous variables, is a vector of endogenous latent variables, is a vector of measurement error for the endogenous variables, and is a matrix of factor loadings relating to . Associated with (1) and (2), respectively, are two variance-covariance matrices, and . The matrix is a matrix of variances and covariances among measurement errors , and is a matrix of variances and covariances among measurement errors .

For flexibility, LISREL describes the path model component as relationships among latent variables, where is a matrix of path coefficients describing the relationships among endogenous latent variables, is a matrix of path coefficients describing the linear effects of exogenous variables on endogenous variables, and is a vector of errors of endogenous variables. Associated with (3) are two variance-covariance matrices: is a variance-covariance matrix of latent exogenous variables, and is a matrix of covariances among errors of endogenous variables. With only these three equations, LISREL is a flexible mathematical framework that can accommodate any specification of a SEM model.

SEM has been typically implemented through covariance structure modeling where the variance-covariance matrix is the basic statistic for modeling. Model fitting is based on a fitting function that minimizes the difference between the model-implied variance-covariance matrix and the observed variance-covariance matrix , where is estimated from observed data, is predicted from the causal and noncausal associations specified in the model, and is a generic function of the difference between and based on an estimation method that follows. As Shipley [13] concisely stated, causation implies correlation; that is, if there is a causal relationship between two variables, there must exist a systematic relationship between them. Hence, by specifying a set of theoretical causal paths, one can reconstruct the model-implied variance-covariance matrix from total effects and unanalyzed associations. Hayduk [33] outlined a step-by-step formulation under the LISREL mathematical framework, specifying the following mathematical equation for : where . Note that in (5) the derivation of does not involve the observed and latent exogenous and endogenous variables (i.e., , , , and ).

A common method in SEM for estimating parameters in is maximum likelihood (ML). In ML estimation, the algorithm iteratively searches for a set of parameter values that minimizes the deviations between elements of and [12]. This minimization is accomplished by deriving a fitting function (4) based on the logarithm of a likelihood ratio, where the ratio is the likelihood of a given fitted model to the likelihood of a perfectly fitting model. The maximum likelihood procedure requires the endogenous variables to follow a multivariate normal (MVN) distribution, and to follow a Wishart distribution. Hayduk [33] described the steps in the derivation and expressed the fitting function as where refers to the trace of a matrix and and are defined as above. Proper application of (6) also requires that observations are independently and identically distributed and that matrices and are positive definite [33]. After minimizing (6) through an iterative process of parameter estimation, the final results are the estimated variance-covariance matrices and path coefficients for the specified model.

2.5. Model Assessment

Kline [19] and Schumacker and Lomax [15] provided a comprehensive listing of indices and criteria to assess model fit, but four basic fit statistics are summarized here. The goal of model assessment is to test the causal implications of a model [13]. The first is the overall model chi-square test based on a test statistic that is a function of the mentioned fitting function (6) as follows: where is sample size and follows a chi-square distribution with degree of freedom as defined above. Subsequently, a value is estimated and evaluated against a significance level.

The overall model chi-square test is only applicable for an overidentified model, that is, when . A just-identified model (), for example, a path model representation of a multiple regression, does not have the required degrees of freedom for model testing [13]. The null (causal) hypothesis associated with the test is that there is no difference between model estimates and the data, and the alternative hypothesis is otherwise. Therefore, failure to reject the null hypothesis is the ultimate objective of the modeling process. Although it may seem to be contrary to the intent of common hypothesis testing in ANOVA, this approach is consistent with the accept-support context where the null hypothesis represents a researcher’s belief [41]. Nonetheless, as with common hypothesis testing, failure to reject the fitted model does not prove the specified causal relationships in the model. One should be particularly aware of existing equivalent models, that is, models that have different hypothesized causal relationships but fit the data equally well.

The second fit statistic to consider is the Root Mean Square Error of Approximation (RMSEA), which is parsimony-adjusted index that accounts for model complexity. The index approximates a noncentral chi-square distribution with the estimated noncentrality parameter as where is computed from (7) and is defined above. The magnitude of reflects the degree of misspecification of the fitted model. The RMSEA is then defined as Thus, RMSEA measures the degree of misspecification per model degree of freedom, adjusted for sample size. RMSEA also reflects the view that the fitted model is an approximation of reality, so that RMSEA measures the error of approximation [42]. Browne and Cudeck [43] suggested that RMSEA ≤ 0.05 indicates a close approximation or fit, a value between 0.05 and 0.08 indicates a reasonable approximation, and a value ≥ 0.1 suggests a poor fit.

The third index is the standardized root mean square residual (SRMR), which is relatively easy to compute. Both and are transformed into correlation matrices, and the residual matrix is the difference between the two. Hence the mean square of the elements in the residual matrix is the SRMR. In general, SRMR <0.10 is considered a good fit of as an approximation to .

Lastly, the Jöreskog-Sörbom Goodness of Fit Index (GFI) is a measure of relative amount of variances and covariances jointly accounted for by the model, and it is defined as [44] where is identity matrix. GFI ranged from 0 to 1.0 with 1.0 indicating the best fit.

In general, statistical tests for the overall model fit and values of parameter estimates are less important in SEM than in univariate regression models. One reason is that all parameters are simultaneously estimated in SEM, so the significance of a parameter estimate should be viewed in the context of the whole model. Second, the confirmatory aspect of the model is weakened if model modification is based on the significance of estimates rather than the theory behind the model structure. Finally, SEM is still a large-sample technique, and hypothesis testing is generally affected by sample size, especially the chi-square test and to a lesser extent RMSEA, SRMR and GFI.

3. Overstory-Understory Relationships

An appropriate topic for SEM application is modeling the relationship between overstory structure and understory plants. The range of implications includes forage production [45], conservation of functional groups such as late-seral herbaceous species under variable-retention harvesting [46], and controls on vertical structure under mature forest conditions [47]. Overstory tree cover [48], species composition [49], and tree density [50], in addition to shrub cover [51] and soil litter [52], are among the variables found strongly associated with richness and abundance of the understory community. These variables serve as surrogates for ecological factors such as light attenuation, throughfall precipitation, and soil water and nutrient availability. Many researchers have underscored the complex interactions and multivariate nature of these factors [53]. A SEM approach to the overstory-understory relationship could apply existing knowledge to the process of proposing functional links and could thereby test and expand our understanding of relationships. Toward this end, the following example illustrates the application of SEM in quantifying ecological mechanisms influencing the abundance of late-seral herb species in mature Douglas-fir forests of the Pacific Northwest.

3.1. Hypotheses and Processes

Application of SEM requires a set of well-defined hypotheses generated from theoretical considerations, previous knowledge, and personal observation. In this study, the processes of interest are light attenuation through the overstory canopy, competition for belowground resources, effects of forest floor litter, and topographic effects.

Hypothesis 1. Light interception by the overstory tree canopy directly affects late-seral herb cover, with an indirect effect through its influence on midstory shrub and tree cover.

Overstory effects on understory vegetation are most commonly assumed to be mediated by shading. Light levels under the canopy can be determined by the structure and distribution of overstory trees that absorb incident light [54]. Overstory canopy cover in this analysis is applied as a measurable surrogate for the amount of light intercepted by overstory trees. Likewise, understory plants occupying a given vertical layer have been shown to influence plants of even lower stature, for example, by controlling germination, growth, or survival of tree seedlings and thereby functioning as ecological filters [55, 56]. This filtering effect of understory vegetation can be at least partly attributed to shading effects at the forest floor [57]. In this analysis, understory shrub and tree cover served as a proxy for light attenuation by this middle layer of forest vegetation.

Hypothesis 2. Belowground competition from overstory trees, as represented by overstory stem density, directly affects late-seral herb cover, with an indirect effect through its influence on midstory shrub and tree cover.

Root trenching experiments at the forest floor have demonstrated positive responses of diversity, abundance, and biomass of the shrub and herb community, most likely due to increased availability of soil water and nutrients after eliminating roots from overstory trees [58, 59]. Overstory stem density was chosen to represent belowground competition as it has been found to positively correlate with belowground biomass (e.g., [60]). Overstory stem density is probably correlated to some degree with overstory canopy cover. Conditional on a given level of overstory canopy cover, it was hypothesized that overstory stem density would serve as a surrogate for belowground resource competition. Overstory competition for belowground resources was assumed to have a direct effect on late-seral herbs and an indirect effect through midstory shrubs and trees.

Hypothesis 3. Cover of coarse and fine woody litter directly affects late-seral herb cover, with the former having an indirect effect through its influence on midstory tree and shrub cover.

These model components represented the mechanical interference of forest floor litter on understory plant establishment [61]. Forest floor litter can directly impede the growth and germination of shrubs and herbs by reducing light and temperature and preventing root contact with mineral soil [62]. Two different processes of coarse and fine woody litter accumulation were considered. Coarse woody litter is produced by tree mortality in response to a number of interrelated factors such as stand density, site quality, age class, and disease intensity [63]. Overstory stem density and stand age were hypothesized as the main driving factors for tree mortality in mature Douglas-fir ecosystems. Overstory and midstory canopy cover contribute to fine woody litter, and the litter accumulation was assumed proportional to rate of crown recession of these two tree layers. Because higher stand density generally hastens crown recession [64], overstory stem density and average diameter were proposed as the main factors influencing fine woody litter accumulation. Lastly, we assumed that there was a myriad of unknown factors contributing to both coarse and fine woody litter production.

Hypothesis 4. Topographic aspect directly affects late-seral herb cover, with an indirect effect through fine woody litter and midstory shrub and tree cover.

Topographic aspect affects the amount of incident solar radiation, which strongly influences the microclimatic conditions such as air and soil temperature [65]. In the northern atmosphere, southwest aspects are often the most severe sites for vegetation establishment and growth [66]. Sariyildiz et al. [67] found that, in the northern hemisphere, fine litter deposited in stands on northern aspects had a higher decomposition rate than litter deposited on southern aspects, likely due to drier and hotter microclimatic conditions at the latter. Hence, the difference in decomposition rate due to aspect might influence late-seral herb cover.

Additional factors could be proposed that might influence late-seral herb cover. In their reviews, Barbier et al. [49] and Wagner et al. [68] suggested that overstory tree species, soil water availability, and phytotoxic compounds could be expected to affect the forest herb community. Their effects specifically on late-seral herb species are not quite clear, however, and would require observable surrogates that were not available in this study.

3.2. Data Collection

The data for modeling overstory-understory relationships came from the Demonstration of Ecosystem Management Options (DEMO) study [69], a large-scale operational research experiment implemented in western Oregon and Washington, USA. The DEMO study looked at the effects of variable-retention harvesting on various aspects of biodiversity, microclimate, and human perceptions. The analysis presented here was based on preharvest data collected in 1994/1995 from 36 mature Douglas-fir stands ranging in age from 65 to 170 years (based on stand examination records from the USDA-Forest Service and Washington Department of Natural Resources). Maguire et al. [2] provided details on the experimental and treatment designs and initial stand structures.

A permanent 8 × 8 or 7 × 9 grid of 40 m spacing was installed in each 13-ha stand. Overstory and understory vegetation were observed on approximately half the grid points of a permanent grid, that is, 32–37 sample points depending on the stand. Thus, the total number of sample points was 1181 for all 36 stands. Detailed sampling protocols were described by Halpern et al. [46]. Percent cover of herbaceous species (typically <1 m tall at maturity) was recorded for each of 24 microplots (0.2 × 0.5 m) clustered at each sample point. Percent cover of coarse and fine woody litter cover was also recorded for the same microplots. Percent cover of tall shrub species (typically >1 m tall at maturity) and of understory coniferous and hardwood trees (<5.0 cm dbh) was measured by the planar intercept method on four 6 m long transects radiating from each sample point. At the 0 m and 6 m marks of each transect, a Moosehorn densiometer was used to estimate percent overstory tree cover. Overstory trees were sampled with a set of nested circular plots: 0.01 ha plot for trees with dbh ≥ 5 and <15 cm and 0.04 ha plot for trees with dbh ≥ 15 cm [2]. At each sample point, aspect was recorded as azimuth in the downhill direction.

In the DEMO study, 48 herb species were classified as late-seral herb species, that is, species that reached maximum abundance in old-growth forest conditions and that were sensitive to canopy removal or disturbance [46]. The nine observed variables for modeling late-seral herb cover with SEM were (1) mean percent late-seral herb cover (LSHERB, %), (2) mean percent tall shrub and understory coniferous and hardwood tree cover (UNDER, %), (3) mean percent coarse woody litter cover (CLITTER, %), (4) mean percent fine woody litter cover (FLITTER, %), (5) mean percent overstory tree cover (TREE, %), (6) tree stem density expressed as number of trees per hectare (TPH), (7) quadratic mean diameter at breast height (dbh) of overstory trees (QMD, cm), (8) cosine-transformed aspect (ASPECT), and (9) stand age (AGE, years).

The values of all variables except AGE were estimated at the sample point level and were used for fitting the models described below. LSHERB, CLITTER, and FLITTER were the average of 24 microplots, and UNDER and TREE were the average of four transects at each sample point. As recommended by Beers et al. [66], aspect was cosine-transformed with a predetermined phase shift of 45° resulting in a variable that would range from one at a northeast aspect and negative one at a southwest aspect [70].

3.3. Structural Equation Modeling

The four hypotheses described above were translated into an SR model depicted in Figure 3. The SR model failed the two-step rules for model identification by Bollen [16] due to one observed variable per latent variable and the correlated errors of the fine and coarse litter. The errors of endogenous latent variables reflected omitted causes [19]. Hence, the correlated errors between the latent variables fine and coarse litter, depicted as the bidirectional arrow (Figure 3), reflected the assumption in Hypothesis 3 that a myriad of unknown factors contributing to both coarse and fine woody litter production. Hayduk [33] suggested that the solution for the former identification problem was to assign a priori measurement errors and the ULI constraint. Measurement error expressed as the percent of variance of each observed variable was estimated as the average of three expert opinions (Table 1). The amount of variance attributed to measurement error was the product of this percentage and the observed variance and was then entered into the SR model (Figure 3). Subsequently, the ULI constraint assigned a value of 1.0 for each factor loading (Figure 3). Although latent variables seemed identical to observed variables in the SR model (Figure 3), particularly in cases having the same labels (e.g., aspect), they were conceptually different because latent variables were treated as true measurements after accounting for measurement errors. The SR model had a bow-free pattern; that is, errors of fine and coarse woody litter were correlated but there was no direct effect (unidirectional arrow) between the two latent variables [19]. In practice, a path model with bow-free pattern is considered recursive and thus is identified [19]. Consequently, the whole SR model was identified.

Exploratory data analysis showed nonlinear bivariate relationships, ill-scaled observed variance-covariance matrix, and nonnormality in the data (see [71]). The bivariate relationship between LSHERB and TPH was nonlinear; thus, TPH was transformed by the natural logarithm to satisfy the linearity assumption. An ill-scaled observed variance-covariance matrix, that is, extremely large difference between the largest and smallest variances, may cause a problem in the iterative model estimation process [19]. Hence, ASPECT and ln(TPH) were each multiplied by 10 to improve the properties of the observed variance-covariance matrix (Figure 3). The correlation matrix and standard deviations of the nine observed variables with the transformed TPH and ASPECT are presented in Table 2. The correlation between ln(TPH) × 10 and QMD was relatively high but was below the extreme value suggested by Kline [19] (−0.804, Table 2). Hence, we assumed that this level of correlation would not severely affect parameter estimation. Most observed variables were not univariate normally distributed due to moderate or extreme skewness and/or kurtosis. SEM is sensitive to violation of the multivariate normality assumption, particularly kurtosis. Therefore, following the recommendations of Jöreskog et al. [72], both observed and asymptotic variance-covariance matrices were fitted by Robust Maximum Likelihood (RML) to obtain standard error estimates that accounted for the nonnormality in the data.

Because sample points were nested within each stand as mentioned, observations within the stand were not independent. If the nesting structure of the data was ignored, the parameter estimates would be unbiased but the associated standard errors might be underestimated [10]. Following methods discussed by Asparouhov and Muthén [73] and described in LISREL documentation [74], standard errors were adjusted and the model chi-square statistic (7) was scaled to reflect the nesting structure of the data. This final formulation of the SR model (Figure 3) was fitted with LISREL 8.8 [75].

The model converged to an admissible solution, but the resulting fit was poor; the scaled- was 33.75 with and value < 0.001, and the RMSEA was 0.049 with 95% confidence interval between 0.032 and 0.067. Sets of modification indices were provided by LISREL with suggested paths and unanalyzed associations to improve model fit. One suggested path, Stand Age Late-Seral Herb Cover, was plausible because stand structure and microclimatic conditions gradually become more favorable to the development of late-seral herb community through the course of stand development. By including the additional path, the final model converged with indices indicating relatively good overall model fit; that is, the scaled- was 6.30 with and a value = 0.61 (sample size = 1181), and the RMSEA was < 0.001 with 95% confidence interval between 0.0 and 0.029.

3.4. Results

Estimated unstandardized path coefficients depicted the effects among variables in the SR model (Figure 4). To facilitate comparison between relative magnitudes of effects, the width of each arrow in the diagram corresponds to the value of the respective standardized solution; that is, wider arrows indicated larger standardized coefficients and thus stronger relative effects.

Overstory cover affected late-seral herb cover both directly and indirectly, but the indirect effect through Understory cover was weaker. The indirect effect was estimated as the product of the direct effects of overstory cover on understory cover and understory cover on late-seral herb cover or approximately −0.09. In contrast, the direct effect of overstory cover on late-seral herb cover was −0.19 (Figure 4). In other words, a 1% increase in overstory cover directly reduced late-seral herb cover by 0.19% ( value = 0.011, Table 3) and indirectly by 0.09% (−0.09 = −0.51 × 0.18, Figure 4). Hence, the total effect of overstory cover, which included the direct effect and all possible indirect effects, was significantly negative (−0.24, value < 0.001, Table 5). All possible indirect effects included non-significant indirect paths such as through fine litter.

Tree density effects on late-seral herb cover were also primarily indirect and mediated by understory cover (Figure 4). This was supported by an insignificant direct effect (−0.08, value = 0.62, Table 3). The combined indirect effect, which included all possible indirect paths, suggested that a twofold increase in tree density reduced late-seral herb cover by 1.92% (−1.92 = −0.277 × 10 × ln(2); value = 0.003; Table 4). Of all the indirect effects, the path Tree Density → Understory Cover → Late-Seral Herb Cover contributed the most to the combined indirect effect, amounting to a reduction of 1.24% cover for a doubling of tree density (; Figure 4).

Coarse litter and fine litter had insignificant direct, indirect, and total effects on late-seral herb cover (all values ≥ 0.13, Tables 3, 4, and 5). As hypothesized, both tree density and stand age were significant predictors for coarse litter ( values = 0.002 and 0.007, resp., Table 3), and both contributed positively to its accumulation. Similarly for fine cover, all effects from overstory cover, understory cover, tree density and tree size were significant (all values ≤ 0.007, Table 3) The direct effects of both overstory cover and understory cover were positive on fine litter accumulation, but the effects of the other two were negative. However, the fitted model could only explain a small portion of the observed variation in coarse litter and fine litter ( = 0.04 and 0.12, resp., Figure 4). Finally, the errors of coarse litter and fine litter were negatively correlated with large error covariance of −42.28 (Figure 4) and an error correlation of −0.707.

Aspect had positive direct and indirect effects on late-seral herb cover ( values = 0.009 and 0.017, resp., Tables 3 and 4), with a total effect amounting to 0.73 ( value < 0.001, Table 5). As aspect shifted from southwest towards northeast, the late-seral herb cover increased. The significant indirect effects were dominated by the path through understory cover. In contrast, the path Aspect Fine Litter Late-Seral Herb Cover had an estimated indirect effect of 0.0 (Figure 4). As expected, stand age had a strong direct effect on late-seral herb cover, with the late-seral herb cover increasing by 0.16% for every one-year increase in Stand Age ( value < 0.001, Table 3).

There were many other pathways and cascading indirect effects that one could explore in the SR model (Figure 4), but we presented only those most pertinent to the formulated hypotheses. With all the predicted direct and indirect effects, the fitted model explained approximately half of the observed variance in late-seral herb cover (, Figure 4).

3.5. Interpretation

The results supported the first hypothesis that greater light attenuation negatively affected late-seral herb cover. The significance of both direct and indirect effects supported the common perception that light was the most limiting parameter for ground vegetation coverage and abundance [68]. The majority of sample points in the DEMO study had dense overstory canopy cover ranging from 65% to 100%. Late-seral herbs were unlikely to thrive directly under an extremely dense overstory canopy because of limiting light availability [76]. The model also concurrently accounted for a cascading effect of overstory canopy density on small tree and shrub cover and the latter on late-seral herb cover. The low shade cast by midstory shrubs and trees [77] probably eliminated more light demanding species and favored the relatively shade-tolerant late-seral species, but increasing overstory canopy was unlikely to benefit both. In short, late-seral herbs occupying the forest floor responded to the cumulative effects from successively higher layers of vegetation.

The second hypothesis on belowground resource competition was supported by the results that decreasing tree density had a positive net effect on late-seral herb cover. Lindh and Muir [78] found that thinning dense 20-year-old Douglas-fir stands had a positive impact on the frequency of late-seral herbs 20 years later. On the other hand, Ares et al. [79] found that richness of early seral herbs and shrubs increased but late-seral herb response was minor six years after thinning in 40- to 60-year-old Douglas-fir stand. However, their results might not readily discriminate between release from belowground resource competition and effects from increasing light availability to forest floor. The SEM model in fact suggested that the indirect effect through understory shrub and tree cover was the dominant effect. Possible differences in the general depth of root systems between layers of forest vegetation might partially explain the observed response. Smith and Huston [80] suggested that shade-tolerant plants would invest resources into aboveground development at the expense of a limited root system. Hence, late-seral herbs might have smaller and shallower root system than understory trees and shrubs (see [81]). In areas with lower overstory stem density, understory shrubs and small trees may be able to exploit the deeper belowground resources such as water otherwise used by overstory trees (see [68]), thereby minimizing the effects of belowground competition on late-seral herbs that would be utilizing resources closer to the surface. Lastly, the bivariate relationship was observed to be nonlinear, which suggested that late-seral herbs were more responsive to release from belowground resource competition in a highly dense forest stand, perhaps due to immediate availability of belowground resources from a decrease in density. Nevertheless, the effect of root competition from successively higher layers of vegetation on late-seral herbs is poorly known because most work has focused on the composite response of all understory species (e.g., [82]).

Contrary to the common perception that forest floor litter might negatively affect late-seral herbs by physically obstructing growth and establishment [62], the third hypothesis was not supported by the fitted model. Although the suggested influence of overstory and understory crown dynamics on rates of litterfall was supported, rates of decomposition differed between the two litter types. Hence, the net and cumulative effects on late-seral herbs could be more complex than depicted by the model. It was also apparent that a myriad of unknown processes controlled litter accumulation as indicated by the low and large negative correlation between errors of coarse and fine litter. Future research could measure other dimensions of forest floor litter accumulation such as litter depth and thereby improve our understanding of mechanisms that have been proposed in the literature.

The fitted model indicated that conditions in Douglas-fir stands on northeast aspects were more favorable to late-seral herb establishment and maintenance than those on southwest aspects. This was supported by Lindh [83], who found that presence of old-growth associated herbs was positively associated with north-facing aspect in second-growth Douglas-fir stand. The benefit of northeast aspects to late-seral herb species can probably be attributed to the lower solar angle and greater light attenuation for a given overstory canopy, resulting in more mesic conditions at the forest floor.

Given sufficient time, natural stand dynamics gradually led to stand conditions that were favorable for late-seral herb community. Many aspects of the biotic and abiotic environment change throughout the course of stand development, including changes in the chemical and structural properties of the soil [49]. In some forest types, age might be the strongest integrator of the composite effects on understory community structure [84]. In the fitted models, average age of overstory trees emerged as a strong predictor of late-seral herb cover, but this variable alone offered limited insight into mechanisms associated with increasing age.

4. Discussion

The overstory-understory model illustrated an application of SEM to understanding the implications of vertical forest structure on late-seral herb species in a forest ecosystem. The overstory-understory relationships involved cascading effects of higher layers of forest vegetation on successively lower layers. The interrelationships could be complex and difficult to manipulate and test experimentally. Development of the SEM drew on published univariate models, exploratory analyses, theory and expert knowledge to establish a set of working hypotheses on the interactions of different layers of forest vegetation. The hypotheses were translated into a mathematical model, and the latter was fitted to a large dataset from mature Douglas-fir forests. The process was confirmatory at this stage of development. However, the initial structure was altered based on modification indices and a new tentative structure was imposed to improve the model fit. Although supported by additional theory, causal inference on the overstory-understory relationship was theoretically weakened by changes to the a priori model structure. Nonetheless, the blend of confirmatory and exploratory SEM analyses provided insights into the complexity imposed by cascading effects of different layers of forest vegetation.

The overstory-understory model was fitted by maximum likelihood. Because maximum likelihood was a full-information method, all processes and interactions represented in the model were considered simultaneously during model estimation [19]. This aspect of SEM was particularly important as the mechanisms for understory vegetation responses to overstory characteristics were interactive and rarely dominated by one specific ecological factor [68]. As a result, light attenuation, belowground resource competition, microclimatic effects associated with aspect, and stand development over time were identified simultaneously as important drivers of the observed patterns. Physical obstruction from coarse and fine woody litter, however, had a relatively small influence on late-seral herbs in the presence of other apparently more dominating processes. The power of SEM derived from the fact that processes and components operating in these complex forest ecosystems were simultaneously accounted for in a way that tentatively identified causal mechanisms and multiple pathways that contributed to the responses of interest.

Given that there was only one indicator per latent variable, the overstory-understory SR model might alternatively be specified as a path model. This would simplify the estimation process but would assume that observed variables were measured without error. A sensitivity analysis on the effects of measurement errors on estimated path coefficients was implemented for the overstory-understory SR model and suggested that, when measurement errors were large, some path coefficients increased in strength and absolute value while others decreased in both [71]. Field observation was often subjected to certain degree of measurement error that included both observer and instrument error. Fiala et al. [85] found that the Moosehorn densiometer generally provided more conservative estimates of canopy cover compared to other measurement methods. Given the implied error associated with estimating canopy cover, an SR model that accounted for measurement error would be more appropriate than a path model under the circumstances.

The final overstory-understory SR model could be revised by removing predicted paths that are statistically insignificant. Unlike multiple regression, however, model parsimony was not the ultimate goal of SEM. Empirical respecification of a model by dropping paths that were not significant increased the chance of overfitting a model to a specific dataset; for example, the statistically insignificant path might be due to random variation or unique attributes of the analyzed dataset [19]. As mentioned above, a path was evaluated simultaneously in the context of other paths; thus, its significance was meaningful only in the presence of other paths represented in the model. Empirical respecification was therefore not recommended, and a researcher should not feel compelled to drop insignificant paths from a model until replication of the results with other datasets or sufficient accumulation of evidence from related research indicated otherwise [19, 86].

Although the final overstory-understory model fitted the DEMO data reasonably well, the results did not prove the causal relationships hypothesized in the model. Causal inference can be supported by replicating the model with independent datasets, considering equivalent models before model fitting [19], or by very specific manipulative experiments. Finally, it was important to note that ecosystem effects and responses were explained in terms of mechanisms that we did not measure directly; for example, light attenuation was represented by overstory canopy cover rather than a direct measure of light intensity at the forest floor. To the extent that more direct mechanisms and ultimate variables are becoming more directly observable, for example, the light filtering effects quantified spatially and temporally in terms of PAR (photosynthetically active radiation) irradiance [68], the dominating processes and mechanisms can be inferred with greater confidence.

5. Conclusions

The aim of this paper was to provide a conceptual framework for applying SEM to understand mechanisms influencing one component of forest structure. The primary motivation was the assumption that a successful SEM should facilitate forest management by providing insight into causal mechanisms and thereby assisting in the design of silviculture treatments and systems that best met specified management objectives. The larger goal was to establish a basic understanding of principles, limitations, and assumptions of SEM and to promote potential applications. Data from natural systems typically create analytical challenges, such as nonnormality in endogenous variables, nonlinear relationships, multicollinearity, ill-scaled variance-covariance matrices, and nested study designs. SEM continues to evolve and address special issues such as modeling binary, ordinal and categorical data [87, 88], and accommodating interactions and curvilinear relationships between latent variables [89, 90]. Recent years have seen the development of different model estimation procedures (see, [15]), including a Bayesian approach to SEM [91]. The latter has been an attractive alternative for modeling hierarchical SEMs. With these advancements, SEM will become an important and complementary analytical tool in understanding the processes that drive the dynamics of complex natural ecosystems.

Acknowledgments

The authors gratefully acknowledge the Edmund Hayes Endowment for Silviculture Alternatives for providing the funding to complete this analysis. They thank Paul Anderson and Doug Mainwaring for providing their expertise and knowledge in quantifying the measurement errors used in the analysis. They thank the many contributors to the Demonstration of Ecosystem Management Options (DEMO) study for allowing them to use DEMO data in this study. Finally, the authors thank the reviewers for excellent reviews that have hugely improved early versions of the paper. The DEMO study is a joint effort of the USDA Forest Service Region 6 and Pacific Northwest Research Station. Research partners include the University of Washington, Oregon State University, University of Oregon, Gifford Pinchot and Umpqua National Forests, and the Washington State Department of Natural Resources. Funds were provided by the USDA Forest Service, PNW Research Station (PNW-93–0455, PNW-97-9021-1-CA, and PNW-01-CA-11261993-091).