Process-Based Modeling to Assess the Effects of Recent Climatic Variation on Site Productivity and Forest Function across Western North America

A process-based forest growth model, 3-PG (Physiological Principles Predicting Growth), parameterized with values of soil properties constrained by satellite-derived estimates of maximum leaf area index (LAImax), was run for Douglas-fir (Pseudotsuga menziesii) to contrast the extent to which site growth potential might vary across western North America between a cool, wet period (1950–1975) and a more recent, generally warmer and drier one (2000–2009). LAImax represents a surrogate for overall site growth potential, as demonstrated from a strong correlation between the two variables, with the latter based on the culmination of mean annual increment estimates made at 3356 ground-based U.S. Forest Service survey plots across the states of Oregon and Washington. Results indicate that since 2000, predicted LAImax has decreased more than 20% in portions of the Southwest USA and for much of the forested area in western Alberta. Similar percentage increases in LAImax were predicted for parts of British Columbia, Idaho and Montana. The modeling analysis included an assessment of changes in seasonal constraints on gross primary production (GPP). A general reduction in limitations caused by spring frost occurred across the entire study area. This has led to a longer growing season, along with notable increases in summer evaporative demand and soil drought for much of the study area away from the maritime influence of the Pacific Ocean. OPEN ACCESS


Introduction
Since the turn of the century, the public's perception that significant changes in climate may be underway has increased following unparalleled outbreaks of fires, insects and pathogens in western North America [1][2][3].Satellite-borne instruments document the extent of these disturbances and quantify the rates of change since the acquisition of Landsat imagery first began in 1972 [4].With the launch of the National Aeronautics and Space Administration's Moderate Resolution Imaging Spectroradiometers (MODIS), daily changes in the Earth's terrestrial surface are currently being monitored.Satellite observations allow global monitoring of the state of vegetation, both as gradual changes occur over time and following disturbances [5].
The question arises, how much disturbance is natural and how much can be attributed to changes in climate?Without periodic disturbances, many forests become overly dense and predisposed to fire, insects and disease [6,7].A lengthening in the frequency of a disturbance, such as the interval between wildfires, by itself, may allow longer-lived species to replace more short-lived, fire-adapted ones, altering forest composition for centuries [8].Alternatively, too short of a fire-return interval may prevent fire-adapted species from reaching maturity and producing seed.Ideally, it would be desirable to know the composition, density and age class distribution of trees in every forest [9].
More realistically, we might select a property of forests that is known to respond quickly to changing environmental conditions, such as maximum leaf area index, LAI max [10].Forest ecologists recognize that there is a close relationship between LAI max and the growth potential of vegetation in a particular region [10][11][12].The natural range in both variables is large in western North America.Across a 250-km transect in Oregon, LAI max varies from <1.0 in a juniper woodland to 12.0 m 2 m −2 in a coastal Sitka spruce-western hemlock forest [13,14] with a similar range and spatial distribution of the maximum mean annual increment (MAI max ) from 1.5 to 15 m 3 ha −1 yr −1 [15].
A number of process-based models have been developed that quantify the link between LAI max and productivity by simulating gross photosynthesis, plant respiration and the allocation of growth to leaf area, as well as to stems and roots [16,17].Most process-based models take into account constraints imposed by seasonal variation in incident radiation, temperature and the evaporative demand of the atmosphere.Some also include a soil water balance to account for observed differences in growth and mortality during periods of extended drought [18,19].Few models, however, account for spatial variation in soil fertility [20], although the above-ground growth responses to rising atmospheric CO 2 are recognized to be severely constrained on infertile soils [21], as well as on soils with extreme water deficits [22].
Nightingale et al. [23] investigated the implications of variation in both of these soil properties by performing sensitivity analyses on modeled gross photosynthesis (GPP) among nine broadly defined ecoregions of the contiguous USA.Most forested areas in the western part of the United States were shown to be moderately to highly sensitive to variation in estimates of water availability and fertility.
Even within a small area of southwestern Oregon, accurate information of soils proved necessary to ensure that model predictions were in close agreement with the measured growth [12].The above cited authors concluded that conventional soil maps were insufficient to parameterize forest growth models across large areas.As an alternative to intensive field surveys, Coops et al. [24] applied model inversion and optimization techniques to derive and map soil properties across western North America at a spatial resolution of 1 km.Their mapping of soil properties was dependent on convergence with satellite-derived estimates of LAI max , which saturated at 7.0 m 2 m −2 .We performed sensitivity tests, whereby baseline values for the available water storage capacity (ASWC) of soils defined by Coops et al. [24] were varied by ±50% to assess the effect on modeled predictions of the current distribution of 20 western tree species.In three cases, altering ASWC values shifted the predicted ranges of species by >20% from their known distributions on surveyed field plots [25].
In this paper, we accept the baseline values of soil properties delineated for forested regions by Coops et al. [24] and expand on the approach to assess the extent that soil properties might buffer or magnify the effects of climatic variation on LAI max and potential forest productivity in western North America.Specifically, we contrast conditions in a cool, wet period  with a notably generally warmer and drier one (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009).If climatic variation were sufficient to affect GPP and growth potential between the two periods, seasonal changes in the frequency of frost, the intensity of drought and evaporative demand should be evident, along with changes in LAI max .With the exception of the presence or absence of a heavy snowpack, which favor some species over others, climatic effects should be most important during the growing season when solar irradiance is high and canopy development peaks.We recognize, however, that in places where the temperature remains above freezing, considerable photosynthesis is also possible during early spring and late autumn [26].To account for photosynthesis that might occur during the dormant season, we parameterize a forest growth model, 3-PG (Physiological Principles Predicting Growth) developed by Landsberg and Waring [27] for the most widely distributed evergreen conifer in the study area, Douglas-fir (Pseudotsuga menziesii) [28].

Materials and Methods
This section contains brief descriptions of the study area, the 3-PG model and the sources of data required to parameterize, to run it and to validate a general relationship between LAI max and field-derived estimates of maximum mean annual increments.More detailed descriptions are available in the publications cited below.

Study Area
Across western North America, all forested areas were delineated using MODIS-derived land cover information [29].In temperature and precipitation, the mountainous terrain enhances local variation in both climate and soils.As a result, the study area includes portions of 51, Level III, US Environmental Protection Agency's defined ecoregions [30], within which there are four major subregions.The latter are described below based on information extracted from several sources [13,14,23,24].
The Marine West Coast Forest zone, which extends from Alaska in a progressively narrowing band to the San Francisco Bay area, is the most productive of the four subregions, with the highest percentage of forested land (82%).Western hemlock (Tsuga heterophylla) and Sitka spruce (Picea sitchensis) are present throughout most of this maritime-affected zone, with Douglas-fir (Pseudotsuga menziesii), coastal redwood (Sequoia sempervirens), Alaska yellow cedar (Chamaecyparis nootkatensis), western red cedar (Thuja plicata) and grand fir (Abies grandis) well represented in some areas.Maximum values of LAI generally range from 8 to 12 m 2 m −2 .
The Northwest Forested Mountains is the second most productive subregion in the study area and is 55% forested.Douglas-fir and western hemlock are abundant in mixtures with Pacific silver fir (Abies amabilis), noble fir (Abies procera) and western larch (Larix occidentalis).At higher elevations and latitudes, lodgepole pine (Pinus contorta), whitebark pine (Pinus albicaulis), mountain hemlock (Tsuga mertensiana), subalpine fir (Abies lasiocarpa) and Engelmann spruce (Picea engelmannii) are characteristic species.Maximum LAI varies from 3.0 at some high elevation subalpine forests to 10.0 m 2 m −2 for the most productive Douglas-fir and hemlock stands [13,14].
On drier sites, ponderosa pine (Pinus ponderosa) and incense cedar (Calocedrus decurrens) are usually present.Ponderosa pine extends its range southward through the Temperate Sierra Ecoregion (20% forested) into the fringe of the North American Desert Ecoregion, with further reductions in productivity.In Oregon and Northern California, in the rain shadow of the Cascade and Sierra Mountains, open pine forests grade into western juniper (Juniperus occidentalis) woodlands and, finally, to sagebrush steppe [13,31].Some species, such as Jeffrey pine (Pinus jeffreyi), grow predominantly on soils derived from ultramafic parent materials in the northern extent of its range [32].Maximum LAI ranges between 0.5 and 2.5 m 2 m −2 .
In the North American Desert of the Southwestern U.S., tree species make up only 2% of the flora, of which pinyon pines (mainly Pinus edulis and P. monophylla) and junipers (Juniperus spp.) are the most common [33], with LAI max generally less than 1.5 m 2 m −2 [19].

Forest Inventory and Analysis (FIA) Data
To demonstrate that a close relation exists between simulated LAI max and site growth potential, we used data collected from 3356 U.S. Forest Service FIA (Forest Inventory and Analysis) plots in the states of Oregon and Washington, where Douglas-fir, ponderosa pine and/or western hemlock were present and site productivity was measured [15].The FIA plots, each consisting of four subplots totaling 0.07 ha in size, are part of the national inventory of forests for the United States [34,35].A tessellation of hexagons, each approximately 2400 ha in size, is superimposed across the nation's forests, with one field plot randomly located within each hexagon.In the western states, approximately one tenth of the plots are measured annually.The data includes the ages and heights acquired over the period from 2000-2006 for the three dominant tree species in the region.Hanson et al. [36] describes site tree selection procedures, which, since 2001, requires at least 3, and sometimes 5 or 10 trees be sampled, depending on the size and diversity of species present.They also present site index equations for each species, as well as the method by which yield tables [37][38][39][40] were employed to estimate the culmination of the mean annual increment (MAI max ) for different species, site productivity classes and regions.
We note the potential of obtaining more direct measurements of periodic annual increment (PAI) following a complete decadal resurvey of all plots (next, in 2020), but such an analysis would need to be restricted to fully stocked stands with trees of ages between 20-30 years when PAI peaks in all of the above-cited yield tables.We assessed the relationship between modeled LAI max (described below) and MAI max , using regression analysis.

Process-Based Growth Model
There are a wide variety of physiologically-based process models, but only a few have been designed to scale projections of photosynthesis, structural growth and mortality across landscapes (see the review by Mä kelä et al. [17]).Among the most widely used is the 3-PG model [27].The 3-PG model differs from others primarily in a number of simplifying assumptions: (1) that monthly mean climatic data are adequate to capture major trends; (2) that autotrophic respiration (R a ) and net primary production (NPP) are approximately equal fractions of gross photosynthesis (GPP); and (3) that the proportion of NPP allocated to roots increases from 25% to 60% as nutrients (particularly nitrogen) become progressively less available.3-PG also differs from most other models by including subroutines for thinning and defoliation, which are important attributes as outbreaks of insects and disease become more frequent.
The model calculates gross photosynthesis, canopy evaporation and transpiration, growth allocation and litter production at monthly intervals (∆t).It reduces potential photosynthesis and transpiration by imposing restrictions on stomatal conductance through modifiers, taking values between 0 = complete restriction and 1.0 = no restriction, that account for the effects of frost, high vapor pressure deficits and deficits in soil water content.The soil water modifier is determined from the ratio of the amount of water available in the root zone of the trees (ASW) to the maximum value (ASW max ).ASW max is the water holding capacity of the soil, which is the difference between the water content in the root zone at field capacity and at the wilting point.For any given month, ASW is calculated from: (1) where ASW(t) is the value at the beginning of the month and P, E and T denote the monthly values of precipitation, evaporation and transpiration, respectively.The model includes a term to account for rainfall interception by the forest canopy; this water is assumed to be lost by evaporation, giving E. Transpiration is calculated from the Penman-Monteith equation, which incorporates a canopy conductance term derived from stomatal conductance and LAI.If the value of ASW on the left-hand-side of Equation (1) exceeds ASW max , i.e., the whole root zone is at field capacity, the excess is assumed to be lost as runoff or drainage.
At monthly time steps, the model is unable to compute a snow water balance accurately, although one may assume that precipitation in months with average temperatures below freezing is largely in the form of snow.At annual time steps, the model sums monthly changes in tree number, mean diameter, stand basal area, above-ground volume and biomass and updates the changes in LAI.
To account for seasonal adjustments in temperature optima [41] and to incorporate the large genetic variation among populations of Douglas-fir, we broadened the range for which photosynthesis could remain above 50% of maximum to lie between 0 °C and 35 °C by setting minimum, optimum and maximum temperatures at −7 °C, 18 °C and 40 °C, respectively.The photosynthetic response at temperatures below −2 °C was truncated to zero, because, below that threshold, stomata are closed [42].
The fertility-dependent growth modifier in the 3-PG model is a function of the soil fertility rating, FR, which ranges between zero, for the poorest soils, to unity for most fertile soils [27].Landsberg and Sands [43] in their book "Physiological Ecology of Forest Production" define soil fertility as a site property and, therefore, the same for all species.Differences among species in their response are incorporated into the model by altering how photosynthetic efficiency is changed and the way that growth is allocated above and below ground.As previously mentioned, we generated estimates of soil fertility (FR) and available soil water storage capacity (ASWC) at 1-km resolution at all forested sites by inverting the 3-PG model to achieve close agreement between modeled LAI max and MODIS-derived observations over the range from 0.5 to 7.0 m 2 m −2 , the latter value representing current sensor saturation [24].

Climatic Data
Long-term weather observations across the region were interpolated spatially using Climate-WNA (Western North America) [44].Climate-WNA is based on the Parameter-elevation Regressions on Independent Slopes Model (PRISM), which accounts for variations in precipitation and temperature associated with mountainous terrain through interpolation of a digital terrain model described by Wang et al. [45].We obtained digital terrain data by expanding the 90-m digital elevation model (DEM), acquired from the Shuttle Radar Topography Mission (SRTM), to 1 km to provide consistent resolution for all of the datasets.
Mean monthly daytime vapor pressure deficits (VPDs) were estimated by assuming that the water vapor concentrations present on the day would be equivalent to those held at the mean minimum temperature [46].The maximum mean VPD was calculated each month as the difference between the saturated vapor pressure at the mean maximum and minimum temperatures.Mean daytime VPD was calculated at two thirds of the maximum value.The number of days per month with subfreezing temperatures (≤2 °C) was estimated from empirical equations with mean minimum temperature.In modeling, we did not attempt to account for extremes in temperature that might kill native species, because representatives of most conifers present in the region have survived for at least the last 50 years and, in many cases, a half a millennium or more [26].
Monthly estimates of total incoming short-wave radiation were obtained by combining the synoptic and zonal variation captured by the North American Regional Re-Analysis (NARR) with topographically-driven variation based on Fu and Rich [47], similar to the approach applied by Schroeder et al. [48].Spatial biases found through comparison with station networks were then removed.

Model Simulations
Two sets of simulations were performed with 3-PG.For both simulations, the 3-PG model was initialized with 500 seedlings per hectare and allowed to grow trees for 50 years under averaged monthly climatic conditions for the 1950-1975 and 2000-2009 periods.There is little self-thinning before 50 years when the initial stocking is set at such a low value, and no age-related decrease in growth was imposed to allow the maximum leaf area index to maintain a plateau from 20 to 50 years.
Otherwise, the LAI would decrease in response to the physiological requirement to maintain hydraulic homeostasis as trees increase in size [49].
We first generated estimates of LAI max for all forested areas for the two distinct periods (1950-1975 and 2000-2008) and then compared them.We binned the LAI max data by units of 1 in the latter period to reduce variations caused by the large 1-km cells before regressing with FIA-derived values of MAI max .
All model simulations were undertaken at 1-km spatial resolution across the study area using 3-PG, programmed in C++.Maps were first produced using ARCMAP software, then transferred to the Institute for Conservation Biology's website [50], where open access software is available to enhance detail and add topographical features and political boundaries.

Results
Figure 1 provides strong support for the basic assumption that a general relationship exists between LAI max and site growth potential (MAI max ).The relationship applies over the full range of values in western North America [13].Figure 2 contrasts modeled changes in GPP and LAI max for Douglas-fir between the two periods.As expected, the spatial patterns are similar, because soil fertility, which exerts the major control on the partitioning of carbon above or below ground in the model, remained the same for both periods.An increase in drought or VPD between the two periods could result in reduced photosynthesis, causing a decrease in LAI when monthly litterfall exceeds new foliage production.A reduction in the occurrence of frost or an increase in solar irradiance associated with a reduction in cloud cover would have the opposite effect.The smallest percentage change in LAI max occurred in the Marine West Coast Forest Zone.However, more than a 20% increase in LAI max was simulated from the cooler to warmer period on the east side of the Sierra Mountains in California and Nevada.The model predicts >20% decreases in LAI max since 2000 for parts of Alberta and the northern Rocky Mountains.Figure 3 provides insights into the underlying causes for simulated changes in GPP and LAI max displayed in Figure 2. The most consistent change across the study area is a general decrease in spring frost (Figure 3a).The coastal region and also larger parts of the interior of the Pacific Northwest of the U.S. and Southwestern Canada experienced a particularly strong decline in spring frost across much of the natural range of Douglas-fir.   .The units represent the proportional change in environmental constraints on GPP in reference to the first period.The range of Douglas-fir is from [28].
A reduction in frost permits trees to start growth earlier, but at the same time, warmer temperatures melt snow faster and increase transpiration, which, on shallow soils, leads to soil water deficits.Since 2000, summer evaporative demand (Figure 3b) has increased most notably in the Southwest USA [51].Although summer soil water deficits were also predicted to increase across parts of the Southwest, the largest deficits were predicted for the forests situated along the cordillera through western Alberta into Montana, a condition that continued in some places into the fall (Figure 3d).

General Points
In general, model predictions of LAI max seem reasonable, as they agree closely with the independent analysis performed by Law et al. [14] across Oregon, Washington and Northern California and parallel values of MAI max for much of the same area (Figure 1).Unfortunately, it is not possible to assess whether stand growth differs significantly on FIA survey plots between the two periods, because measurements of site index, from which MAI max are calculated, integrate growth over many years.The US Global Change Program documents that climatic conditions have changed across the nation, with the most rapid warming, reduction in spring frost and winter snowpack occurring in the western states [52].The noted reduction in spring frost throughout most of the study area (Figure 3a) could result in extending the growing season and increasing productivity, assuming that there was no commensurate increase in evaporative demand and soil water deficits during the summer [53][54][55].
It is evident from our analysis that the constraints on GPP and productivity are far from uniformly distributed, even within a major forest zone (Figure 3).This variation is a reflection in part of local differences associated with mountainous topography and soils.For the most part, predicted soil properties seem reasonable based on extensive sensitivity analyses by Mathys et al. [25] and are, in any case, more reliable than can be obtained from conventional sources, at least at these spatial scales, as documented by Coops et al. [24].There may be an exception in the Maritime West Coast Forest Zone, where an increase in soil drought is indicated in some areas in Figure 3c.To date, there are few signs of reduction in LAI max across the subregion and little evidence that GPP has decreased, as monitored at eddy-flux tower sites in Oregon [56], Washington [57] and British Columbia [58].In fact, with the extension of the growing season associated with more frost-free days, we might expect GPP to increase slightly, along with LAI max .

Future Improvements in the Approach
Although the model inversion and optimization technique used by Coops et al. [24] provides considerable improvements in the quantitative assessment of soil properties at large spatial scales in comparison with conventional mapping approaches, independent measurements of LAI max , ASWC and FR would be desired.To this end, a number of remote sensing techniques show promise, although some are not yet available from orbiting satellites.Three-dimensional laser scanning (Light Detection and Ranging, LiDAR) could provide useful validation over the full range in LAI max [59].LiDAR is similar to radar, but uses pulses of light rather than microwaves to distinguish properties of the Earth's surface.
There are many ways that the modeling approach could be improved.Recent refinements in the 3-PG model have been made by Wei et al. [60].At decadal intervals, however, such refinements yield values of LAI max within 5% of those reported here [61].In plantations, where climatic and soils data are gathered on site, 3-PG in its current configuration generally predicts interannual growth and canopy dynamics within the error of field measurements [62][63][64][65].

Implications on Forest Health
Our decision to compare decadal or longer intervals of climatic records may obscure extreme events that can trigger wide-spread tree mortality over a few years, but high rates of tree mortality have been reported continually since the 1990s across much of the study area [1,3,[66][67][68].In the Southwest, much of the increase in tree mortality may be associated with consecutive years of drought [19,[53][54][55].Our modeling analysis indicates that sustained increases in evaporative demand in the Southwest may be equally important as periodic drought (Figure 2b,c), which supports the conclusions by Weiss et al. [51] drawn from longer-term climatic records.
On the other hand, in British Columbia, where forests of lodgepole pine (Pinus contorta) and, to a lesser extent, interior spruce (Picea glauca and P. engelmannii) have sustained heavy mortality from bark beetles and various pathogens [68], summer soil water deficits appear to have decreased since 1950-1975 as a result of increased rainfall [10] and reduced evaporative demand (Figure 2c).It is likely that a warming trend with fewer frost days, rather than drought, is the cause of the expansive mortality recently recorded in the forests of British Columbia [68,69].
We recognize that species differ in their abilities to access water and to survive drought [70].Furthermore, genetic differences within a species are important, particularly at the regional scale [71][72][73].The extent that the distribution of tree species may shift or are shifting is the subject of many papers (e.g., [74][75][76]) and will be the topic of a subsequent paper based in part on the data presented here.An advantage of incorporating a process-based model to predict shifts in the distribution of species is that seasonal constraints on photosynthesis can be identified and their influences assessed with respect to seasonal differences in the amount of light absorbed by the forest canopy.For example, the widespread reduction in spring frost is particularly significant for many subalpine species, which can break bud and still withstand temperatures of −5 °C.Warmer temperatures, along with a reduction in snowpack, will favor the invasion of faster growing, less snow and frost-hardy species into the subalpine forest zone.The extent of recent disturbances in forest cover also could create a positive feedback to warming the climate by reducing the transfer of water vapor to the atmosphere and by increasing heat exchange [77].

Conclusions
The approach illustrated in this paper, by incorporating improved maps of soil properties and by assessing seasonal constraints on productivity and LAI max in a more functional manner, provides predictions that can be validated: (a) at sites where eddy-flux measurements of CO 2 and water vapor exchange are being measured [78]; (b) at decadal revisited Forest Inventory and Analysis survey plots; and (c) through independent measurements acquired by various sensors borne on aircraft and satellites orbiting the Earth [59,79].With even moderate refinement in the mapping of soil properties, we achieved closer agreement between predicted and measured growth than previously reported for parts of the US Pacific Northwest [80].
Finally, the analyses presented in this paper show that conditions across the entire study area have been generally warmer since 2000 compared to the period of 1950-1975, a trend particularly evident in the boreal forests of Canada [71,72].The U.S. Global Change Program substantiates that rapid warming and melting of winter snowpack have already occurred over the last few decades and predict the trend to continue and to accelerate.The full implications of such changes on western forests are just beginning to be appreciated, with expanding areas of forests affected by fire, insects and disease (i.e., [67,68]).Without limitations on the rise of greenhouse gases, the forests of the future may be in a continuous state of transition, bringing to scientists and to society unforeseen challenges to predict and with which to cope.

Figure 1 .
Figure 1.Simulated maximum leaf area index (LAI max ) for the period 2000-2009 is closely related to the maximum mean annual increment (MAI max ) estimated from height and age measurements acquired from 2000 to 2006 on >3356 FIA (Forest Inventory and Analysis) survey plots in Oregon and Washington [15].

Figure 3 .
Figure 3. Modeled average changes from 2000 to 2009 in (a) spring frost frequency, (b) summer evaporative demand, (c) summer soil water stress and (d) fall soil water stress compared to the cooler and wetter period.The units represent the proportional change in environmental constraints on GPP in reference to the first period.The range of Douglas-fir is from[28].