Extrapolating leaf CO2 exchange to the canopy: a generalized model of forest photosynthesis compared with measurements by eddy correlation

Over the last 4 years, two data sets have emerged which allow increased accuracy and resolution in the definition and validation of a photosynthesis model for whole forest canopies. The first is a greatly expanded set of data on the nitrogen-photosynthesis relationship for temperate and tropical woody species. The second is a unique set of long-term (4 year) daily carbon balance measurements at the Harvard Forest, Petersham, Massachusetts, collected by the eddy-correlation technique. A model (PhET-Day) is presented which is derived directly from, and validated against, these data sets. The PnET-Day model uses foliar nitrogen concentration to calculate maximum instantaneous rates of gross and net photosynthesis which are then reduced for suboptimal temperature, photosynthetically active radiation (PAR), and vapor pressure deficit (VPD). Predicted daily gross photosynthesis is closely related to gross carbon exchange at the Harvard Forest as determined by eddy-correlation measurements. Predictions made by the full canopy model were significantly better than those produced by a multiple linear regression model. Sensitivity analyses for this model for a deciduous broad-leaved forest showed results to be much more sensitive to parameters related to maximum leaf-level photosynthetic rate (A max) than to those related to light, temperature, VPD or total foliar mass. Aggregation analyses suggest that using monthly mean climatic data to drive the canopy model will give results similar to those achieved by averaging daily eddy correlation measurements of gross carbon exchange (GCE).


Introduction
Predicting the effects of global change on forest ecosystem function requires the development of simple, generalizable, well-validated, data-based models that can be run for large regions using only simple driving variables. Such models should be derived directly from existing physiological data on component processes, such as photosynthesis and respiration, and should accurately predict the measured rates of function of intact systems. Rigorous parameterization and validation increase the likelihood of accurate predictions for modified climate regimes. In a previous paper (Aber and F?d?rer 1992), a simple model of whole forest carbon and water balance (PnET) was presented and validated against annual net primary productivity data for ten forests throughout North America, and for monthly or annual water balance data for three forested watersheds.
The photosynthesis routine in this model was based on data available at that time relating maximum rates of net leaf-area-based photosynthesis (Amax) to weight-based foliar ? concentration. Limitations in this model included a mixing of units in the prediction of Amax and the limited number of observations available for parameterizing this relationship. Validations of seasonal carbon flux were not possible because of a lack of data on monthly carbon exchange by forest systems.
Over the last four years, two unique data sets have emerged which overcome these limitations. The first is a greatly expanded set of data on the nitrogen-Amax relationship for temperate and tropical woody species (Reich et al. 1991b(Reich et al. , 1992 The purposes of this paper are: (1) to present a simple model (PnET-Day) of seasonal changes in whole forest canopy photosynthesis driven by daily climatic data, (2) to validate this model against the 4-year, eddy-correlation data set available for the Harvard Forest; (3) to examine the effects of four different methods of aggregating daily climate data and model output to estimate monthly carbon balances; and (4) to test the sensitivity of the model to changes in both driving variables (e.g., temperature) and input parameters (e.g., foliar ? concentration).

General
The structure of the model presented here (PnET-Day, Fig. 1) is similar to the photosynthetic routines in the Forest-BGC (Running and Coughlan 1988; Running and Gower 1991) and TCX (Bonan 1993) models, but differs from those in the explicit use of foliar ? levels to determine Amax. Response functions for radiation intensity, temperature and vapor pressure deficit are used with daily mean climate drivers for these variables to calculate realized Amax for leaves at the top of the canopy. A layered canopy is then simulated, with both radiation intensity and specific leaf weight (SLW) declining with canopy depth. Leaf respiration is a function of Amax and temperature, and is calculated separately for daytime and nighttime temperatures. Maximum (summer) and minimum (winter) leaf mass are input parameters. The onset of canopy development in spring is driven by a growing degree day sum algorithm, and canopy senescence results from negative carbon balances in autumn. Both of these are responsive to weather patterns unique to specific years. The model assumes no significant water stress but is intended to replace the photosynthesis routine in the original PnET which performs full water balance and water stress calculations. Amax which is then separated into potential gross photosynthesis and dark respiration. Potential gross photosynthesis is reduced for suboptimal conditions of light, temperature and vapour pressure deficit (VPD) to give realized gross photosynthesis. Light levels in the simulated layered canopy are determined by ambient photosynthetically active radiation (PAR), cumulative leaf area index (LAI) and the light attenuation constant. Respiration is modified by temperature using a ?10 function Model algorithms and parameters (Table 1) Model parameters can be divided between those that are generalized and should apply to any species within the broad-leaved deciduous and needle-leaved evergreen groups, and those that need to be specified for an individual site or canopy.

Generalized parameters
Instantaneous Amax as a function of foliar ? (AmaxA, AmaxB). In a broad context, all wild C3 species demonstrate a common linear relationship between foliar ? concentration (g N.g1 leaf) and Amax (nmol C02 g"1 leaf s1, Field and Mooney 1986; Reich et al. 1991aReich et al. , 1992, making the use of photosynthesis-N relationships in modeling canopy gas exchange rates (e.g., Reich et al. 1990; Aber and Federer 1992) a powerful, generalized approach. However, this relationship varies between species, or functional groups of species, and species-specific or group-specific relationships predict observed patterns better than a single, general one (Reich et al. 1994(Reich et al. , 1995. Given this, we could attempt to aggregate species specific curves (if available) for each tree species in each simulated system for use at the canopy level. However, in keeping with the goal of model generality, data used here are summarized for two groups, broad-leaved deciduous and needle-leaved evergreen, that show different relationships ( Daily averaged Amax (AmaxErac). Maximum instantaneous rates of net photosynthesis are not generally maintained throughout an entire day. In order to run PnET-Day on daily to monthly timesteps an average daily Amax is needed. Reductions in gas exchange rates correlated with increasing evaporative demand are well known, and addressed in the model through the DVPD variable (see below, "Effects of vapor pressure deficit"). However, even on days of apparently non-limiting VPD, maximum early morning rates are not maintained throughout the day. In addition to periods of less than saturating irradiance, several poorly understood factors, which may include end-product inhibition, more negative xylem water potentials and inherent circadian rhythms, combine to yield a daily averaged Amax which is below the maximum, early morning instantaneous rate. For example, Ellsworth and Reich (1992) observed on days of low (and ostensibly non-limiting) evaporative demand that achieved daily leaf level carbon gain for sugar maple was 77% of that possible based solely on light limitation to photosynthesis, while on days of high evaporative demand achieved carbon gain was 57% of that possible. Other studies show similar values for this ratio under non-VPD limiting conditions (Table 2). In the PnET model AmaxFrac is set to 0.76, the mean of measured values for 11 eastern deciduous species.     1990; Ellsworth and Reich 1992; Kruger and Reich 1993) and found a consistent decline in net photosynthetic rate and leaf conductance with increasing VPD, especially above 1 kPa. To effectively but simply capture this non-linear pattern, and to aggregate the effects of partial days above 1 kPa, we use a power function (DVPD1 ? VPDDVPD2) with the values given below. This formula-tion also allows for a linear response from 0 kPa as described for western conifers (e.g.

Site-specific parameters
The following variables are site-specific and ideally should be determined for the area of the forest sampled by the eddy flux system. However, the dimensions of the average tower footprint are not known with precision, and so we rely on parameters which were measured either spatially across the Harvard Forest ( Timing of leaf-out and senescence (GDDFolStart, GDDFolEnd, SenescStart). These variables determine the timing of leaf out and the earliest time at which foliar senescence can occur. The first two are expressed as total accumulated growing degree days calculated as all mean temperatures above 0?C. There is very little quantitative data on the timing of initiation and completion of foliar expansion as a function of climatic variables although the potential for obtaining these relationships from simultaneous climate monitoring networks and satellite remote sensing (e.g., AVHRR data) is high. The values used here (Table 1)  , extending 6-10 m above a mixed oakmaple canopy, was instrumented with an array of meteorological sensors. Supplemental instrumentation and data acquisition equipment were installed in a climate-controlled hut 20 m from the tower base. Throughout the investigation, air temperature at 30 m (aspirated thermistor), incident PAR (silicon quantum sensor), and soil surface temperature (potted thermistors), were logged at 0.5 Hz. Half-hour means were calculated from these data and further aggregated to daily maximum and minimum temperature, and total daily incident PAR.
For input to PnET-Day from either the NOAA or tower data sets, mean daily temperature is taken as the average of the maximum and minimum, and day and night temperatures are calculated as the average of the mean and the maximum and minimum temperatures, respectively.

Hours of daylight per day (hr)
This variable is calculated using an algorithm for daylength as a function of day of year and latitude drawn from Smith (1974). Total daily PAR (moles nr2 d1) is divided by hr*.0036 to give mean daily instantaneous PAR (umoles nr2 s1)?

Vapor pressure deficit (VPD)
This variable is calculated assuming that the atmosphere is saturated at the daily minimum temperature. VPD then the difference between the saturated vapour pressure at the daytime (not maximum) temperature and the minimum temperature (as in PnET, Aber and F?d?rer 1992).

Validation
Two approaches are taken to determining how well the PnET-Day model predicts carbon balances at the Harvard Forest. First, predicted gross carbon exchange (GCE, equivalent to gross photosynthesis) is compared with the tower GCE, calculated as the sum of measured net carbon exchange (NCE) and estimated total forest respiration, derived as a function of temperature using nighttime tower flux measurements. As a test of the degree to which a full canopy model improves upon a simple statistical model, a statistical model is derived and also tested for goodness of fit with the tower data.

The Harvard Forest eddy-correlation data set
The net exchange of C02 by the forest surrounding the meteorological tower was measured from January 1991 through December 1994 using the eddy correlation method (

min1 of air down a 50-m tube and through a closed path infrared gas analyzer (IRGA), a process that introduced a lag of several seconds. Following an adjustment for this lag, we calculated the net C02 exchange as the 30-min covariance of vertical wind and linearly detrended C02 concentration. We compensated for errors associated with sonic alignment and local topography by rotating the flux to the plane with zero mean vertical velocity (McMillen 1988). A series of simulations, laboratory tests, and spectral analyses indicated a small underestimation of flux due to the loss of high-frequency C02 fluctuations. (The 90% response determined by C02 addition on the tower was faster than Is.) We corrected for this error by increasing the measured C02 flux in proportion to the underestimation of sensible heat flux associated with a simulated reduction in the high frequency response of the temperature detector (Leuning and King 1992; Goulden et al. 1995). This correction was generally small (< 10%). The eddy-correlation method measures the net exchange of C02 through a plane at 30 m height. This flux may differ from that into and out of organisms if the quantity of C02 stored between the forest floor and the plane at 30 m changes. In order to more directly assess the flux into plants and soils, we measured the quantity of C02 stored below 30 m by frequently sampling the mixing ratio at eight heights through the canopy (Wofsy et al. 1993). This change in C02 storage was then added to the eddy flux to calculate net carbon exchange (NCE; positive values represent net movement of carbon into the ecosystem). This flux is similar to the total ecosystem respiration during well-mixed nocturnal periods, and the sum of photosynthesis and respiration during the day. Respiration from soil heterotrophs and plant maintenance is strongly controlled by temperature (Jarvis and Leverenz 1983). Seasonal relationships were derived between soil temperature and ecosystem respiration using an overall Q{0 of 2.2 (Goulden et al. 1995). Hourly estimated respiration rates were summed for each day and added to net ecosystem flux measurements to estimate GCE. GCE is equal to the combined rate of RUBISCO carboxylation and oxyg?nation (Goulden et al. 1995).
Predictions for total daily respiration from the entire ecosystem obtained by this method can be compared with similar estimates for soil-only respiration using the equations of Kicklighter et al. (1994) and mean daily air temperature (Fig. 4). The tower-based estimates average 34% higher than the generalized equation from Kicklighter et al. (1994). This could result from site-specific conditions surrounding the tower area, and from the inclusion of above-ground plant respiration in the tower data. The Kicklighter equation predicts respiration from soils only.
The tower measurements were frequently interrupted by rainfall, calibration, maintenance and data collection, and occasionally for extended periods by equipment failure. In the present analysis we use 538 days with uninterrupted observations for the years 1991-1994. This represents the longest and most continuous eddycorrelation data set available for any forest ecosystem. For the pur-

One of the purposes of both the tower measurements and the PnET-Day model is to predict total carbon balances over long time periods (months to years). A critical step in arriving at final values
for these periods is solving the aggregation problem: how to extrapolate even this relatively complete data set to cover days for which direct measurements are not available. Four methods of aggregating to a monthly total were compared. The first is to simply take the average of all daily tower GCE measurements within each month. The second also uses the partial data, but applies the model to the daily climate drivers measured at the tower and averages all model predictions within a month. The third averages all of the tower climate data for a month and uses that average data to run one day of the model, which is then applied to the entire month. The fourth uses the NOAA record of mean monthly climatic to run the model for the average day in each month. Monthly data from this source are derived from a more complete set of daily observations, but have no direct link with tower data.

Sensitivity analyses
Sensitivity analysis allows a clearer understanding of the relative importance of the different parameters and algorithms in the model in controlling daily GCE. It provides insight into which factors are most important in controlling model predictions, and also shows the degree to which errors in parameters or input data will result in errors in prediction. It will also show the extent of shortterm change in ecosystem function expected from different changes in climate.
To test the sensitivity of the model to both driving variables (climate) and input parameters (Table 1), the effect of a 10% increase in each variable on predicted GCE over the entire simulation period was recorded. A preliminary test showed that sensitivity responses were symmetrical, that either an increase or decrease resulted in the same relative change in GCE, but with opposite sign. Therefore, results are presented for the 10% increase only.
At a larger scale, the effects of changing species group, and of predicted changes in climate were also tested. For the former, the parameter set was altered to represent a Harvard Forest pine stand and the model was rerun with the monthly averaged NOAA climate data. For the latter, temperature was increased by either 3?C or 5?C (day and night temperature, all months).

PnET-Day versus eddy correlation data
All parameter values used in the validation exercise are as derived from field data (Table 1)

Aggregation
The four aggregation methods (Fig. 7) produced very similar results, with differences occurring mainly in the early spring leaf-out period of 1991. Average daily GCE for the 4-year period ranged only from 3.92 to 4.27 g C irr2 day1 for the four treatments. These results suggest that accurate predictions of monthly GCE can be derived using only monthly mean data from nearby weather stations. This type of climatic information is available over wide areas and can be extrapolated to regional data planes for use with geographic information systems (GIS; e.g., Ollinger et al. 1994). Combining PnET-Day with such a GIS would yield spatially explicit estimates of GCE over the region involved.

Sensitivity analyses
GCE predictions were more sensitive to changes in parameters related to Amax than to any others (Table 3). Increases of 10% either the foliar ? concentration or the AMaxB parameter yielded increases of 14% in GCE. The ratio between realized daily Amax and the instantaneous rate (AMaxFrac) was the next most sensitive parameter. These were followed by parameters related to light interception per unit leaf mass and the shape of the photosynthetic response curve (k, HalfSat and SLWmax). Parameters related to foliar respiration, temperature effects, and the timing of leaf out and senescence all showed responses of 3% or less to a 10% change in value.
Of particular interest in a global remote sensing context are the relative sensitivities of GCE to total leaf biomass and foliar ? concentration.
Many current models of biosphere-atmosphere exchange (see discussion by Running and Hunt 1993) are based on NDVI (normalized difference vegetation index), which is often assumed to represent spatial differences in green biomass. These results suggest that, at least in closed canopy, broad-leaved deciduous forests, gross photosynthesis is 7 times more sensitive to the foliar ? concentration in the foliage dis- Altering all of the vegetation parameters to run a pine canopy rather than a deciduous canopy (Table 1) reduces predicted gross carbon exchange by about 20% in mid summer (Fig. 8). This results from the lower Amax values which are only partially offset by higher SLW. A slight increase in early-spring carbon gain for the pines relative to the hardwood simulation results from the presence of over-wintered foliage. It should be noted that the much lower slope of the Amax-N relationship for needle-leaved evergreen species will result in a much lower sensitivity to changes in foliar ? concentration. Predicted GCE for the Harvard Forest is also less sensitive to changes in climate variables than to those determining Amax (Table 3). Increases of 10% in PAR and temperature result in 5% increases in GCE. However, if changes in temperature are expressed as absolute increases on the order of those predicted for the next century (3-5?C), increases of 16-21% in GCE are predicted. It must be noted that this assumes no dilution in foliar ? content by increased carbon availability and biomass production, no water stress due to either a longer growing season or reduced rainfall, no acclimation of Amax or respiration to altered temperature, and no photoperiodic controls on carbon acquisition and allocation. More complete models which include a full water balance, and carbon and nitrogen allocation routines are required to assess the importance of these system-level feedbacks. Not only did the PnET-Day and statistical models differ in terms of accuracy of validation, they also gave different results for responses to warming (Fig. 9). While both models predicted increased GCE, the statistical model failed to capture interactions between temperature and other variables, yielding a flat, linear increase across all months when iave is above 0?C. In contrast, the full model shows both sharp increases in GCE by the extension of the growing season into both late winter and late fall (greatest increases in March-April and October-November), with possible depressions in mid-summer due to temperatures above the optimum for gross photosynthesis. Predicted increases in average daily GCE were 0.81 and 1.76 g C nr2 day1 for the PnET-Day and statistical models, respectively.

Conclusions
The results presented here demonstrate that a simple, daily time step model based on physiological measurements at the leaf level can accurately predict seasonal changes in gross carbon exchange by a forest canopy. Predictions made by the full canopy model (PnET-Day) were significantly better than those produced by a multiple linear regression model. Sensitivity analyses predicted that a deciduous stand at the Harvard Forest should be more sensitive to parameters related to maximum photosynthetic rate (Amax) than to those related to light, temperature, VPD or total foliar mass. Aggregation analyses suggest that using monthly mean climatic data to drive the canopy model will give results similar to those achieved by averaging daily measurements of GCE within a month.