SITE DESCRIPTION

The Abisko tundra site (hereafter AT) is 6.2 km SSE of Abisko, Sweden, in mixed tundra vegetation. An eddy covariance tower was located on a moraine plateau at UTM coordinates 411180.486 E, 7577776.226 N, (zone 34W) at an elevation of 752.3 m a.s.l. A meteorological tower was located nearby at 411191.840 E, 7577785.175 N, 751.9 m a.s.l. Sensor location was measured by real-time kinetic GPS using a Trimble 5800 rover with Trimble 5700 base station (Trimble Navigation Ltd., Sunnyvale, California, U.S.A.). Six types of tundra vegetation were previously classified at AT: fell field, open Empetrum heath, closed Empetrum. heath, poor fen, shrub tundra, and grey willow scrub (mostly Salix glauca, noting frequent hybridization) (Shaver et al., 2007; Street et al., 2007; Fox et al., 2008). The same assemblages existed during the study period in 2007.

EDDY COVARIANCE

NEE was measured using an open path eddy covariance system comprised of a LI-COR 7500 infrared gas analyzer (Li-Cor, Lincoln, Nebraska, U.S.A.) and a Gill R3 Sonic Anemometer (Gill Instruments Limited, Lymington, Hampshire, U.K.). Sensors were mounted at 3.0 m on a telescopic mast. Half-hourly fluxes were calculated using FluxView (Centre for Ecology and Hydrology, Wallingford, U.K.). Processing involved correction for the angle of attack (Gash and Dolman, 2003), adjustment for lag times between sonic anemometer and infrared gas analyzer, coordinate rotations to align the horizontal wind vector along the half-hourly mean wind direction and to force the vertical component to zero, correcting observed heat flux for humidity effects (Schotanus et al., 1983), correction for spectral flux losses, the Webb-Pearman-Leuning correction for the effects of air density fluctuations on open path eddy covariance measurements (Webb et al., 1980), and the Burba cor- rection to account for heating of the open-path infrared gas analyzer (Burba et al., 2008).

METEOROLOGICAL MEASUREMENTS

Photosynthetic photon flux density (PPFD) was measured using a SKP 215 Quantum Sensor (Skye Instruments, Llandrindod Wells, U.K.) at a height of 2 m on the meteorological tower. Missing measurements were replaced using the linear relationship with incident shortwave radiation measured by a co-located CM11 Pyranometer (Kipp & Zonen, Delft, Netherlands). Air temperature was measured using a HMP45 sensor (Vaisala, Helsinki, Finland), also at a height of 2 m.

TOWER AND AIRCRAFT NDVI MEASUREMENTS

NDVI was measured from the micrometeorological tower using 20-min averaged output of SKR1800 two-channel reflectance sensors (Skye Instruments) pointing both upward and downward and mounted at 1.5 m on the meteorological tower. The SKR1800 as configured has a field of view of 60° and a footprint of 1.7 m in diameter. NDVI was converted to the formulation used by van Wijk and Williams (2005) and Williams et al. (2008) as discussed in Stoy et al. (2009b). L was estimated using the NDVI conversion for sensors at 1.5 m height following Williams et al. (2008):

Representative daily NDVI values were chosen by selecting periods with solar zenith angles between 55° and 75° to ensure that different days of observation were comparable given the low solar zenith angles experienced in Abisko before mid-June and after early September.

Aircraft-based NDVI measurements were used for the spatially explicit flux footprint analysis. NDVI was calculated for 4 m pixels in a 2 km × 2 km area including the AT tower footprint using observations from an Azimuth Systems AZ-16 Airborne Thematic Mapper (ATM) flown on 17 July 2005 (day of year [DOY] 198) as described in Stoy et al. (2009b).

The tower-based and aircraft-based NDVI measurements were compared against surface reflectance measurements made every 20 m in 200 m transects in the four cardinal directions from the meteorological tower during the peak of the growing season (DOY 185; 2007) using a FieldSpec2 spectroradiometer (ASD Inc., Boulder, Colorado, U.S.A.), henceforth called an “ASD.” The present analysis concentrates on the 200 × 200 m area centered around the AT tower that corresponds to the flux source area during periods with adequate convective and mechanical production of turbulent kinetic energy for accepting the half-hourly eddy covariance flux values (Aubinet et al., 2000). A further check on NDVI surrounding AT was provided by ASD measurements from a University of Edinburgh G-GEOS aircraft flyover on DOY 198, 2008.

NET ECOSYSTEM EXCHANGE MODEL

NEE was modeled using the observed relationship between photosynthesis, L, and irradiation (i.e. PPFD) and between respiration, air temperature, and L (“PLIRTLE”; Shaver et al., 2007). The gross primary productivity (GPP) submodel follows the aggregated canopy photosynthesis model of Rastetter et al. (1992). GPP, in µmol CO₂ m^-2 s^-1, is assumed to follow a saturating response to PPFD, in µmol photons m^-2 s^-1, and is integrated through the canopy using the Beer-Lambert law for light attenuation:

where k is the Beer's law light extinction coefficient, assumed here to be 0.5 after Shaver et al. (2007), E_o?? (µmol CO₂ µmol photons^-1) is the initial slope of the light sensitivity response of photosynthesis, and P_max is the light saturated rate of photosynthesis per unit leaf area (µmol m^-2 s^-1).

Shaver et al. (2007) found that ER models for tundra that include two sources, one sensitive to L and T_a and one not (the latter presumably from deeper soil horizons), fare better than models that assume a single substrate pool. We use model “ER₂” (Shaver et al., 2007):

where R₀ (µmol m^-2 s^-1) is base respiration at 0 °C, R_x (µmol m^-2 s^-1) is the component of ecosystem respiration that is insensitive to observed L and T_a, and β (°C^-1) is the temperature sensitivity of R₀. The global parameter set using observations from Sweden and Alaska as described in Shaver et al. (2007) was chosen for this analysis to investigate the proposed pan-Arctic applicability of PLIRTLE and the generality of the chamber-based model upscaling approach (Table 1). The net CO₂ flux from atmosphere to biosphere is denoted as negative following the micrometeorological convention.

FLUX FOOTPRINT MODEL

A two-dimensional semi-analytic footprint model after Hsieh et al. (2000), extended to two dimensions by Detto et al. (2006) (see also Novick et al., 2004; Stoy et al., 2006; Oishi et al., 2008) was used to sample NDVI pixels from the ATM-generated NDVI map (see Fig. 3 later herein). AU NDVI pixels that contributed less than 0.01% to the footprint were excluded for efficient computation, and the total remaining footprint was rescaled to equal one. These pixels were weighted by flux footprint source area to derive a probability distribution of flux-weighted NDVI values, which were converted to L following Equation (1). This derived distribution of L entered the PLIRTLE model prior to averaging following Stoy et al. (2009b). Modeled NEE based on L estimates from the tower-mounted NDVI sensor (NEE_mod,tower) and modeled NEE based on L estimates from the source-area weighted footprint analysis (NEE_{mod,footprint}) were compared against eddy covariance—measured NEE. We focused the model comparison on the 83-day period between DOY 174 when continuous eddy covariance measurements began and DOY 257 when tower-mounted NDVI measurements became unreliable. We additionally highlight a period between DOY 186 and DOY 200 when PPFD and T_a exhibited characteristic mid-summer variability but tower-mounted NDVI measurements were to a first order constant.

TABLE 1

The pan-Arctic parameter set for the PLIRTLE model with ecosystem respiration model (ER) 2 and a fixed light extinction coefficient (k) of 0.5 following Shaver et al. (2007). Units for E_o are defined in the text.

FIGURE 1.

(A) Air temperature (T_a), and (B) photosynthetically active photon flux density (PPFD) between day of year (DOY) 120 and 250 of 2007, with a focus on the model comparison period between DOY 186 and DOY 200 (C: T_a; D: PPFD).

We note that chamber observations and subsequent modeling work were made at the same site during previous years (Shaver et al., 2007; Street et al., 2007), and that the ATM flyover also took place during previous years (Stoy et al., 2009b). If this previous measurement and modeling work can be used as a basis for upscaling NEE during subsequent years, it may be argued that the ability to upscale models based on chamber measurements to the eddy covariance flux footprint is further simplified.

MICROCLIMATE

The summer of 2007 was characterized by brief warm spells extending from multiple days to over one week (Fig. 1). The highlighted period (DOY 186–200) was characterized by variable PPFD and T_a, which was generally >8 °C.

TOWER AND AIRCRAFT-MEASURED NDVI

Mean tower-measured daily NDVI was ∼0.55 preceding the period of rapid vegetation growth, which began around DOY 150. Mean daily NDVI reached its peak of ∼0.75 around DOY 180 and began a slow decline around DOY 210 (Fig. 2, part A). The mean daily L estimates that resulted from the NDVI observations following Equation (1) ranged from <0.2 m² m^-2 to ∼1 m² m^-2 (Fig. 2, part B).

Mean NDVI from the ASD-sampled transects within the flux footprint (0.683) and NDVI measured by the tower-mounted sensor (0.675) were not significantly different (two-sided t-test, p < 0.05) on DOY 185, 2007, when the ASD measurements were made. Mean tower-measured NDVI during DOY 186–200 was 0.73 and the mean weighted NDVI captured by the dynamic flux footprint (e.g. Fig. 3) over the same period was 0.71. A representative distribution of NDVI and L within the flux footprint for a half-hourly period is shown in Figure 4. Mean ATM-measured NDVI from the Airborne Remote Sensing Facility flyover on 17 July (DOY 198) 2005 in the vicinity of the tower representing the flux footprint was 0.72. Mean ASD-measured NDVI from the G-GEOS flyover in the vicinity of the tower was 0.73 on DOY 198, 2008.

MODEL-MEASUREMENT COMPARISON

Across the observation period, NEE_mod,tower explained 61% of the variance of NEE measurements with a modeling efficiency (Loague and Green, 1991; Meyer and Butler, 1993) of 0.58 and root mean squared error (RMSE) of 15.7 µmol CO₂ m^-2 s^-1. Mean NEE_mod,tower was -1.74 µmol CO₂ m^-2 s^-1, which is less than the mean measured NEE of -1.27 µmol CO₂ m^-2 s^-1, indicating that the model overestimated ecosystem C uptake by 37%. NEE_{mod,footprint} explained only 53% of the variability of NEE and had a lower modeling efficiency (0.53), but also had a far lower RMSE (1.3 µmol CO₂ m^-2 s^-1). Mean NEE_{mod,footprint} differed from observed NEE by less than 3%, far within the range of uncer-tainty on the order of 10% or more commonly attributed to eddy covariance observations (Moncrieff et al., 1996).

The differences in goodness-of-fit between NEE_{mod,footprint} and NEE_mod,tower extended beyond single-day examples; the RMSE of NEE_mod,tower was generally lower during the earlier part of the measurement period before DOY 210, and the RMSE of NEE_mod, footprint was generally lower during the latter part of the measurement period (Fig. 6). Both approaches replicated observed variability poorly during a period centered around DOY 203. This period was marked by relatively low retrieval of eddy covariance flux data that averaged less than 50% of available measurements per day owing in part to large rain events. Approximately 8.5 mm of rain was measured at the Abisko Research Station during the week beginning DOY 200. The RMSE of NEE_{mod,footprint} tended to be lower during the latter part of the observation period after DOY 210. This period corresponded to larger footprint volumes; the peak of the footprint source weight function averaged 41m before DOY 210 and was over 60% longer (66 m) after DOY 210, on average. The increase in footprint dimension corresponded to a decrease in the sensible heat flux toward the end of the measurement period (r = -0.76; p < 0.001); wind speed was not statistically related to DOY (p = 0.59). Both modeling approaches explained over 60% of the variability of measured NEE when excluding the period between DOY 200 and DOY 204 during the period of low NEE measurement retrieval, and after DOY 240 when L estimates were frequently lower than 0.4 m² m^-2 and vegetation was senescing (Fig. 2, part B).

FIGURE 2.

(A) Daily measurements of the normalized difference vegetation index (NDVI) using a tower-mounted sensor, and (B) subsequent estimates of the leaf area index (L) at the Abisko Tundra flux tower during the 2007 growing season. Circles represent the daily mean of acceptable NDVI measurements and subsequent L estimates. Error bars represent 1 standard deviation about the daily mean.

When focusing on the highlighted period between DOY 186 and DOY 200, NEE_{mod,footprint} matched the diurnal patterns of observed NEE better during certain days (e.g. DOY 186 and DOY 196; Fig. 5, part A) and tended to better replicate periods with large flux magnitude (Fig. 5, part C). NEE_mod,tower better matched the diurnal pattern or observed NEE on other days (e.g. DOY 195 and DOY 197). Both modeling approaches frequently missed large positive and negative excursions of observed NEE (e.g. DOY 190). The local Abisko-area PLIRTLE parameterizations (Shaver et al., 2007) did not improve model fit versus the pan-Arctic parameterization (data not shown).

FIGURE 3.

A normalized difference vegetation index (NDVI) map at 4 m resolution of the area around the Abisko Tundra eddy covariance tower (AT), overlaid with a flux footprint for a half-hourly flux measurement at 10:00 (local time) on DOY 186, 2007. The outermost ellipsoid represents the 80% footprint, and the innermost ellipsoid represents the 20% footprint.

FIGURE 4.

(A) The probability distribution function of tower-measured normalized difference vegetation index (NDVI) at the Abisko, Sweden, tundra research site with a kernel density estimate of footprintweighted NDVI for the time period shown in Figure 3 (10:00 on DOY 186,2007). (B) The leaf area index (L) distributions that result from converting the NDVI distributions in part A using the transfer functions described in Williams et al. (2008) (see Equation 1).

FIGURE 5.

(A) The measured net ecosystem exchange of CO₂ (NEE) at a tundra ecosystem near Abisko, Sweden, with PLIRTLE (Shaver et al., 2007) model outputs that employ leaf area index (L) estimated from tower-based normalized difference vegetation index (NDVI) measurements (NEE_mod,tower, B) and NDVI sampled from a dynamic flux footprint (NEE_mod,_footprint, C).

FIGURE 6.

Goodness-of-fit statistics for the PLIRTLE model with leaf area index estimated by a tower-mounted normalized difference vegetation index (NDVI) sensor (“Tower”) and with leaf area index estimated using NDVI from a flux footprint analysis following Hsieh et al. (2000) and Detto et al. (2006) using a 4 m NDVI map (“Footprint”, see Fig. 3). (A) EF: modeling efficiency; (B) r ²: coefficient of determination; (C) RMSE: root mean squared error.

UPSCALING CO₂FLUX IN TUNDRA

A model (PLIRTLE) based on microclimate and L (Shaver et al., 2007) that was parameterized using pan-Arctic chamber-based flux measurements (Shaver et al., 2007; Street et al., 2007) explained up to 75% of the variance of NEE for select two-week periods and explained over 60% of the variance of NEE across most of the growing season. The 4% difference in magnitude between NEE and NEE_{mod,footprint} is well within NEE error from eddy covariance observations, which is usually considered to be on the order of 10–15% (Goulden et al., 1996; Moncrieff et al., 1996; Stoy et al., 2006) depending on flux magnitude (Hollinger and Richardson, 2005; Richardson et al., 2006). The maximum percent of variance explained by the models is similar to but slightly less than the explanatory power of the PLIRTLE model based on chamber measurements (80% ± 5%) in pan-Arctic tundra ecosystems (Shaver et al., 2007).

Our results contrast with those of Fox et al. (2008), who found that upscaling the chamber measurements of Shaver et al. (2007) and Street et al. (2007) based on a map of dominant species gave NEE values 60% greater than eddy covariance-based measurements, and suggested that the chamber sampling design was not systematic. This interpretation is correct; the sampling design of Shaver et al. (2007) and Street et al. (2007) was meant to capture representative examples of dominant tundra ecosystem types for the purposes of a pan-Arctic comparison, not to capture characteristic rates of ecosystem CO₂ exchange at the eddy covariance footprint scale. In other words, assuming that vegetation type controlled surface fluxes, rather than vegetation characteristics (here via L), resulted in erroneous upscaled fluxes. We recommend upscaling biosphere-atmosphere flux observations through validated models of ecosystem function (Rastetter et al., 1992) rather than species composition.

FOOTPRINT MODELING

Adding a dynamic flux footprint via NEE_{mod,footprint} decreased the explanatory power of the model for many periods. NEE_{mod,foot print} represented an improvement over NEE_modtower in certain instances, including large CO₂ uptake events (Fig. 5). On occasion, NEE_{mod,footprint} estimated greater CO₂ sequestration than measurements supported. We suggest that this overestimate is likely due to the frequent presence of high NDVI outliers (Figs. 3 and 4) when scaling using the flux footprint. These high NDVI patches correspond largely to Salix shrubs, which were also frequent outliers when parameterizing the PLIRTLE model (Shaver et al., 2007). Further modeling work with chamber observations suggested that the complex light environment near the edges of Salix vegetation patches caused GPP to deviate from a simple, linear relationship with L (Fletcher et al., 2012). NEE_{mod,footprint} often included NDVI pixels that corresponded to the Salix patches to which PLIRTLE does not fit as well. Future research should focus on modeling NEE in Salix, especially if shrub growth and encroachment across the circumpolar tundra biome continue (Forbes et al., 2010; Macias-Fauria et al., 2012).

Our findings suggest that adding complexity via the footprint model, and thereby uncertainty, to an upscaling approach that matches observations reasonably well does not necessarily improve model fit. NEE_{mod,footprint} did not represent an unambiguous improvement over NEE_mod,tower (Fig. 6) except during the latter part of the measurement period when the flux footprint tended to be larger. It is unclear if disturbance near the tower or the larger population of NDVI pixels sampled by the larger footprint can be attributed to the improved fit of NEE_{mod,footprint} during the latter part of the growing season. Regardless, fortuitous tower-mounted NDVI sensor placement in an area with an average NDVI that closely matched the NDVI in the area of the flux footprint (Fig. 3) contributed to the success of NEE_mod,tower.

OBSERVATIONAL CHALLENGES FOR TUNDRA CO₂ FLUX UPSCALING

The nonlinear relationship between NDVI and L creates a nonlinear relationship between NDVI and NEE in PLIRTLE (van Wijk and Williams, 2005; Williams et al., 2008; Stoy et al., 2009b, 2009c). Most growing season NDVI measurements at AT are on the order of 0.7 (Fig. 2, part A), which is near the point where small errors in NDVI create large biases in PLIRTLE-based NEE estimates (see Fig. 7, part A, in Stoy et al., 2009c). This sensitivity, coupled with the sensitivity of NDVI to soil reflectance (Rocha and Shaver, 2009), suggests that including multiple remote-sensing products including the enhanced vegetation index (EVI) may improve multi-scale carbon flux estimates in tundra (Boelman et al., 2003; Rocha and Shaver, 2009). Future work should add multiple observational constraints to improve multi-scale understanding of tundra CO₂ flux (Quaife et al., 2008; Hill et al., 2011).

ADDITIONAL APPLICATIONS OF THE PLIRTLE MODEL

PLIRTLE and similar models have been used extensively in recent studies of CO₂ exchange in tundra at the chamber (Shaver et al., 2007; Street et al., 2007), tower (Rastetter et al., 2010; Loranty et al., 2011), landscape (Stoy et al., 2009a, 2009b), and regional scales (Loranty et al., 2011). The present analysis demonstrates that it can likewise be used as a strategy for upscaling chamber-scale observations to larger scales in space. It is important to note potential pitfalls of such an approach, including recent findings of lower GPP per unit ground area in tundra patch edges (Fletcher et al., 2012) and the potential role of patterned ground in explaining submeter spatial patterns of ecosystem respiration (Sommerkorn, 2008). We also note that our analysis is restricted to the growing season, and that shoulder and cold-season fluxes from tundra ecosystems may be non-trivial (Larsen et al., 2007a, 2007b; Street et al., 2012); in fact, both NEE_mod,tower and NEE_{mod, footprint} performed more poorly later in the growing season when vegetation began to senesce (Fig. 6). Despite potential limitations, we suggest that our upscaling approach may represent an improvement over previous studies that were unable to find convergent estimates of NEE by combining inference from chamber studies and from eddy covariance towers (Fox et al., 2008).

CONCLUSIONS AND FUTURE STUDIES

Our findings emphasize that canopy characteristics (via L) and microclimate (via PPFD and T_a), rather than community composition, can be used for upscaling CO₂ flux in tundra. Measurements of plant communities, including the extent of vegetation patches and the distribution of species, are not necessary for explaining or upscaling carbon flux during the growing season (Fox et al., 2008). We note, however, that plant response to regional climate change occurs at the species level. These alterations to plant community composition and the emergent ecosystem consequences, including subsequent changes to L (Jia et al., 2003), must be further studied (Chapin et al., 2005; McGuire et al., 2007; Spadavecchia et al., 2008). Research on present tundra carbon cycling should quantify the magnitude, distribution, and uncertainty of L and microclimate in tundra across scales in space (Loranty et al., 2011). Future research should focus on the response of species distribution and ecosystem attributes to a changing climate (Stow et al., 2004), including the response of climate to changes to vegetation distribution and land cover characteristics (McGuire et al., 2006, 2007; Forbes et al., 2010).