Numerical Terradynamic Simulation Group 4-2014 A satellite data driven biophysical modeling approach for estimating northern peatland and tundra CO 2 and CH 4 fluxes

The northern terrestrial net ecosystem carbon bal­ ance (NECB) is contingent on inpnts from vegetation gross primary prodnctivity (GPP) to offset the ecosystem respi­ ration (Reco) of carbon dioxide (CO2 ) and methane (CH4 ) emissions, bnt an effective framework to monitor the re­ gional Arctic NECB is lacking. We modified a terrestrial car­ bon llnx (TCF) model developed for satellite remote sens­ ing applications to evalnate wetland CO2 and CH4 llnxes over pan-Arctic eddy covariance (EC) llnx tower sites. The TCF model estimates GPP, CO2 and CH4 emissions nsing in sitn or remote sensing and reanalysis-based climate data as inpnts. The TCF model simnlations nsing in sitn data explained > 7 0 % of the r^ variability in the 8 day cnmnlative EC measmed llnxes. Model simnlations nsing coarser satellite (MODIS) and reanalysis (MERRA) records acconnted for approximately 69 % and 75 % of the respective r^ variability in the tower CO2 and CH4 records, with cor­ responding RMSE nncertainties of < 1 .3 g C m -2 d i (CO2 ) and 18.2mgCm“ ̂d“ ̂ (CH4 ). Althongh the estimated annnal CH4 emissions were small (< 18gC m “ ^yr“ )̂ relative to Reco (> 180gC m “^yr“ ^), they rednced the across-site NECB by 23 % and contribnted to a global warming potential of approximately 165 ± 128gC02eqm “ ^yr“ ̂when consid­ ered over a 100 year time span. This model evalnation indi­ cates a strong potential for nsing the TCF model approach to docnment landscape-scale variability in CO2 and CH4 llnxes, and to estimate the NECB for northern peatland and tnndra ecosystems.


Introduction
Northern peatland and tundra ecosystems are important components of the terrestrial carbon cycle and store over half of the global soil organic carbon reservoir in seasonally frozen and permafrost soils (Hugelius et al., 2013). However, these systems are becoming increasingly vulnerable to carbon losses as CO 2 and CH 4 emissions, resulting from climate warming and changes in the terrestrial water balance (Kane et al., 2012;Kim et al., 2012) that can increase soil carbon decomposition. Recent net CO 2 exchange in northern tundra and peatland ecosystems varies from a sink of 291 Tg C yr −1 to a source of 80 Tg C yr −1 , when considering the substantial uncertainty in regional estimates using scaled flux observations, atmospheric inversions, and ecosystem process models (McGuire et al., 2012). The magnitude of a carbon sink largely depends on the balance between carbon uptake by vegetation productivity and losses from soil mineralization and respiration processes. High-latitude warming can increase ecosystem carbon uptake by reducing cold-temperature constraints on plant carbon assimilation and growth (Hudson et al., 2011;Elmendorf et al., 2012). Soil warming also accelerates carbon losses due to the exponential effects of temperature on soil respiration, whereas wet and inundated conditions shift microbial activity towards anaerobic consumption pathways that are relatively slow but can result in substantial CH 4 production (Moosavi and Crill, 1997;Merbold et al., 2009). Regional wetting across the Arctic Zhang et al., 2012a) may increase CH 4 emissions, which have a radiative warming potential at least 25 times more potent than CO 2 per unit mass over a 100 year time horizon (Boucher et al., 2009). The northern latitudes already contain over 50 % of global wetlands and recent increases in atmospheric CH 4 concentrations have been attributed to heightened gas emissions in these areas during periods of warming (Dlugokencky et al., 2009;. Northern peatland and tundra (≥ 50 • N) reportedly contribute between 8 and 79 Tg C in CH 4 emissions each year, but these fluxes have been difficult to constrain due to uncertainty in the parameterization of biogeochemical models, the regional characterization of wetland extent and water table depth, and a scarcity of ecosystem-scale CH 4 emission observations (Petrescu et al., 2010;Riley et al., 2011;Spahni et al., 2011;McGuire et al., 2012;Meng et al., 2012).
Ecosystem studies using chamber and tower eddy covariance (EC) methods continue to provide direct measurements of CO 2 and CH 4 fluxes and add valuable insight into the environmental constraints on these processes. However, extrapolating localized carbon fluxes to regional scales has proven difficult and is severely constrained by the limited number of in situ observations and the large spatial extent and heterogeneity of peatland and tundra ecosystems. Recent approaches have used satellite-based land cover classifications, photosynthetic leaf area maps, or wetness indices to "upscale" CO 2 (Forbrich et al., 2011;Marushchak et al., 2013) and CH 4 (Tagesson et al., 2013;Sturtevant and Oechel, 2013) flux measurements. Remote sensing inputs have also been used in conjunction with biophysical process modeling to estimate landscape-level changes in plant carbon assimilation and soil CO 2 emissions (Yuan et al., 2011;Tagesson et al., 2012a;Yi et al., 2013). Previous analyses of regional CH 4 contributions have ranged from the relatively simple modification of CH 4 emission rate estimates for wetland fractions according to temperature and carbon substrate constraints (Potter et al., 2006;Clark et al., 2011) to the use of more complex multi-layer wetland CH 4 models with integrated hydrological components (McGuire et al., 2012;Wania et al., 2013). Yet, most investigations have not examined the potential for simultaneously assessing CO 2 and CH 4 fluxes, and the corresponding net ecosystem carbon balance (Sitch et al., 2007;McGuire et al., 2012;Olefeldt et al., 2012) for peatland and tundra using a satellite remote sensing-based model approach.
It is well recognized that sub-surface conditions influence the land-atmosphere exchange of CO 2 and CH 4 production. However, near-surface soil temperature, moisture and carbon substrate availability play a crucial role in regulating ecosystem carbon emissions. Strong associations between surface soil temperature (≤ 10 cm depth) and CO 2 respiration have been observed in Arctic peatland and tundra permafrost systems (Kutzbach et al., 2007). Significant relationships between CH 4 emissions and temperature have also been reported (Hargreaves et al., 2001;Zona et al., 2009;. Although warming generally increases the decomposition of organic carbon, the magnitude of CO 2 production is constrained by wet soil conditions (Olivas et al., 2010) that instead favor CH 4 emissions and decrease methantrophy in soil and litter layers (Turetsky et al., 2008;Olefeldt et al., 2012). Oxidation by methanotrophic communities in surface soils can reduce CH 4 emissions by over 90 % when gas transport occurs through diffusion (Preuss et al., 2013), but this constraint is often minimized when pore water content rises above 55-65 % (von Fischer and Hedin, 2007;Sjögersten and Wookey, 2009). Despite increases in the availability of organic carbon and accelerated CO 2 release due to soil warming and thickening of the active layer in permafrost soils (Dorrepaal et al., 2009), anaerobic communities have shown a preference for light-carbon fractions (e.g., amines, carbonic acids) that are more abundant in the upper soil horizons (Wagner et al., 2009). Similarly, labile carbon substrates from recent photosynthates and root exudates have been observed to increase CH 4 production relative to heavier organic carbon fractions (Ström et al., 2003;Dijkstra et al., 2012;Olefeldt et al., 2013) that require longer decomposition pathways to break down complex molecules into the simple compounds (i.e., acetate, H 2 + CO 2 ) used in methanogenesis (Le Mer and Roger, 2001).
The objective of this study was to evaluate the feasibility of using a satellite remote sensing data driven modeling approach to assess the daily and seasonal variability in CO 2 and CH 4 fluxes from northern peatland and tundra ecosystems, according to near-surface environmental controls including soil temperature, moisture and available soil organic carbon. In this paper we incorporate a newly developed CH 4 emissions algorithm within an existing terrestrial carbon flux (TCF) CO 2 model framework Yi et al., 2013). The CH 4 emissions algorithm simulates gas production using near-surface temperature, anaerobic soil fractions and labile organic carbon as inputs. Plant CH 4 transport is determined by vegetation growth characteristics derived from gross primary production (GPP), plant functional traits and canopy/surface turbulence. Methane diffusion is determined based on temperature and moisture constraints to gas movement through the soil column, and oxidation potential. Ebullition of CH 4 is assessed using a simple gradient method (van Huissteden et al., 2006). The integrated TCF model allows for satellite remote sensing information to be used as primary inputs, requires minimal parameterization relative to more complex ecosystem process models, and provides a framework to monitor the terrestrial net ecosystem carbon balance (NECB). Although the NECB also encompasses other mechanisms of carbon transport, including dissolved and volatile organic carbon emissions and fire-based particulates, the NECB is limited in this study to CO 2 and CH 4 fluxes, which often are primary contributors in high-latitude tundra and peatland ecosystems (McGuire et al., 2010).
To evaluate the combined CO 2 and CH 4 algorithm approach, we compared TCF model simulations to tower EC records from six northern peatland and tundra sites within North America and Eurasia. For this study, baseline simulations driven with tower EC-based GPP and in situ meteorology data were first used to assess the capability of the TCF model approach to quantify temporal changes in landscapescale carbon (CH 4 and CO 2 ) fluxes. Secondly, CO 2 and CH 4 simulations using internal TCF model GPP estimates (Yi et al., 2013) and inputs from satellite and global model reanalysis records were used to evaluate the relative uncertainty introduced when using coarser scale information in place of in situ data. These satellite and reanalysis driven simulations were then used to determine the annual CO 2 and CH 4 fluxes at the six tower sites, and the relative impact of CH 4 emissions on the NECB.

TCF model description
The combined TCF model CO 2 and CH 4 framework regulates carbon gas exchange using soil surface temperature, moisture and soil organic carbon availability as inputs, and has the flexibility to run simulations at local and regional scales. TCF model estimates of ecosystem respiration (R eco ) and net ecosystem CO 2 exchange (NEE) have been evaluated against tower EC data sets from boreal and tundra systems using GPP, surface (≤ 10 cm depth) soil temperature (T s ) and volumetric moisture content (θ) inputs available from global model reanalysis and satellite remote sensing records (Kimball et al., 2009;McGuire et al., 2012). A recent adjustment to the TCF model Yi et al., 2013) incorporates a light-use efficiency (LUE) algorithm that provides internally derived GPP calculations to determine R eco and NEE fluxes at a daily time step. The adjusted TCF CO 2 model also allows for better user control over parameter settings and surface meteorological inputs . The CO 2 and newly added CH 4 flux model components are described in the following sections. A summary of the TCF model inputs, parameters, and the associated parameter values used in this study are provided in the Supplement (tables S1 and S2; Fig. S1).

CO 2 flux component
The internal TCF model GPP algorithm estimates daily fluxes according to a biome-dependent vegetation maximum LUE coefficient (ε max ; mg C MJ −1 ) that represents the optimal conversion of absorbed solar energy and CO 2 to plant organic carbon through photosynthesis . To account for daily minimum air temperature (T min ) and atmospheric vapor pressure deficit (VPD) constraints on photosynthesis (Running et al., 2004), ε max is reduced (ε) using dimensionless linear rate scalars ranging from 0 (total inhibition) to 1 (no inhibition) that are described elsewhere (i.e., Kimball et al., 2012;Yi et al., 2013). In this study we also account for the sensitivity of shallow rooted vegetation and bryophytes, which lack vascular tissues for water transport, to changes in surface volumetric soil water (Wu et al., 2013), where θ min and θ max are the specified minimum and maximum parameter values: (1) where SW rad (W m −2 ) is incoming shortwave radiation and FPAR is the fraction of daily photosynthetically active solar radiation (PAR; MJ m −2 ) absorbed by plants during photosynthesis. For this approach, PAR is assumed to be 45 % of SW rad (Zhao et al., 2005). Remotely sensed normalized difference vegetation index (NDVI) records have been used to estimate vegetation productivity (Schubert et al., 2010a;Parmentier et al., 2013) and changes in growing season length (Beck and Goetz, 2011) across northern peatland and tundra environments. Daily FPAR is derived using the approach of Badawy et al. (2013) to mitigate potential biases in low biomass landscapes (Peng et al., 2012): This approach uses NDVI or simple ratio (SR; i.e., (1 + NDVI)/(1 − NDVI)) indices as input index values. The results are then averaged to obtain FPAR. Index range corresponds to the difference between the 2nd and 98th percentiles in the NDVI and SR distributions (Badawy et al., 2013). Biome-specific autotrophic respiration (R a ) is estimated using a carbon use efficiency (CUE) approach that considers the ratio of net primary production (NPP) to GPP (Choudhury, 2000). Carbon loss from heterotrophic respiration (R h ) is determined using a three-pool soil litter decomposition scheme consisting of metabolic (C met ), structural (C str ) and recalcitrant (C rec ) organic carbon pools with variable decomposition rates. The C met pool represents easily decomposable plant residue and root exudates including amino acids, sugars and simple polysaccharides, whereas the C str pool consists of litter residues such as hemi-cellulose and lignin (Ise et al., 2008;Porter et al., 2010). The C rec pool includes physically and chemically stabilized carbon derived from the C met and C str pools and also corresponds to humified peat. A fraction of daily NPP (F met ) is first allocated as readily decomposable litterfall to C met and the remaining portion (1 − F met ) is transferred to C str (Ise and Moorcroft, 2006;Kimball et al., 2009). To account for reduced mineralization in tundra and peatland environments, approximately 70 % of C str (F str ) is reallocated to C rec (Ise and Moorcroft, 2006;Ise et al., 2008): Daily CO 2 loss from the C met pool (i.e., R h,met ) is determined as the product of C met and an optimal decomposition rate parameter (K p ). The realized decomposition rate (K met ) results from the attenuation of K p by dimensionless T s and θ multipliers (T mult and W mult , respectively), that vary between 0 (fully constrained) and 1 (no constraint): The temperature constraints are imposed using an Arrheniustype function (Lloyd and Taylor, 1994;Kimball et al., 2009) where decomposition is no longer limited when average daily T s exceeds a user-specified reference temperature (T ref ; in K) that can vary with carbon substrate complexity, physical protection, oxygen availability and water stress (Davidson and Janssens, 2006). The W mult modifier accounts for the inhibitory effect of dry and near-saturated soil moisture conditions on heterotrophic decomposition (Oberbauer et al., 1996). For this study, θ opt is set to 80 % of pore saturation to account for ecosystem adaptations to wet soil conditions (Ise et al., 2008;Zona et al., 2012) and near-surface oxygen availability provided by plant root transport (Elberling et al., 2011). Decomposition rates for C str and C rec (K str , K rec ) are determined as 40 % and 1 % of K met , respectively (Kimball et al., 2009), and R h is the total CO 2 loss from the three soil organic carbon pools: Finally, the TCF model estimates NEE (g C m −2 d −1 ) as the residual difference between R eco , which includes R a and R h respiration components, and GPP. Negative (−) and positive (+) NEE fluxes denote respective terrestrial CO 2 sink and source activity:

CH 4 flux component
A CH 4 emissions algorithm was incorporated within the TCF model to estimate CH 4 fluxes for peatland and tundra landscapes. The model estimates CH 4 production according to T s , θ, and labile carbon availability. Plant CH 4 transport is modified by vegetation growth and production, plant functional traits, and canopy aerodynamic conductance that takes into account the influence of wind turbulence on moisture/gas flux between vegetation and the atmosphere. The CH 4 module is similar to other process models (e.g., Walter and Heimann, 2000;van Huissteden et al., 2006), but reduces to a one-dimensional near-surface soil profile following Tian et al. (2010) to simplify model parameterization amenable to remote sensing applications. For the purposes of this study, the soil profile is defined for surface (≤ 10 cm depth) soil layers as most temperature and moisture retrievals from satellite remote sensing do not characterize deeper soil conditions. Although this approach may not account for variability in carbon fluxes associated with deeper soil constraints, field studies from high-latitude ecosystems have reported strong associations between CH 4 emissions and near-surface conditions including T s and soil moisture (Hargreaves et al., 2001;Sachs et al., 2010;von Fischer et al., 2010;Sturtevant et al., 2012;Tagesson et al., 2012b).

CH 4 production
Soil moisture in the upper rhizosphere is a fundamental control on CH 4 production and emissions into the atmosphere. Methanogenesis (R CH 4 ) within the saturated soil pore volume (ϕ s ; m −3 ; the aerated pore volume is denoted as ϕ a ) is determined according to an optimal CH 4 production rate (R o ; µM CH 4 d −1 ) and labile photosynthates: For this study, CH 4 production was driven using the soil C met pool to reflect contributions by lower weight carbon substrates (Reiche et al., 2010;Corbett et al., 2013) in labile organic carbon-rich environments. Carbon from the C str pathway may also be allocated for CH 4 production in ecosystems with lower labile organic carbon inputs and higher contributions by hydrogenotrophic methanogenesis (Alstad and Whiticar, 2011). The Q 10p temperature modifier is used as an approximation to the Arrhenius equation and describes the temperature dependence of biological processes (Gedney and Cox, 2003;van Huissteden et al., 2006). The reference temperature (T p ) typically reflects mean annual or non-frozen season climatology. Both Q 10p and T p can be adjusted, in addition to R o , to accommodate varying temperature sensitivities in response to ecosystem differences in substrate quality and other environmental conditions (van Hulzen et al., 1999;Inglett et al., 2012). Methane additions from R CH 4 are first allocated to a temporary soil storage pool (C CH 4 ) prior to determining the CH 4 emissions for each 24 h time step; C met is also updated to account for carbon losses due to CH 4 production.

CH 4 emission
The magnitude of daily CH 4 emissions (F CH 4 ) from the soil profile is determined through plant transport (F plant ), soil diffusion (F diff ) and ebullition (F ebull ) pathways: Vegetation plays an important role in terrestrial CH 4 emissions by allowing for gas transport through the plant structure, avoiding slower diffusion through the soil column and often reducing the degree of CH 4 oxidation (Joabsson et al., 1999). Daily F plant is determined using a rate constant (C p ) modified by vegetation growth and production (f grow ), and an aerodynamic term (λ) and a rate scalar (P trans ) that account for differences in CH 4 transport ability according to plant functional type: A fraction of F plant is oxidized (P ox ) prior to reaching the atmosphere and can be modified according to plant functional characteristics (Frenzel and Rudolph, 1998;Ström et al., 2005;Kip et al., 2010). Plant transport is further reduced under frozen surface conditions to account for pathway obstruction by ice and snow or bending of the plant stem following senescence (Hargreaves et al., 2001;Sun et al., 2012).
The magnitude of f grow is determined as the ratio of daily GPP to its annual maximum and is used to account for seasonal differences in root and above-ground biomass (Chanton, 2005). Aerodynamic conductance (g a ) represents the influence of near-surface turbulence on energy/moisture fluxes between vegetation and the atmosphere (Roberts, 2000;Yan et al., 2012) and gas transport within the plant body (Grosse et al., 1996;Sachs et al., 2008;Wegner et al., 2010;Sturtevant et al., 2012): Values for z m and d are the respective anemometer and zero plane displacement heights (m); z om and z ov are the corresponding roughness lengths (m) for momentum, heat and vapor transfer. The von Karman constant (k; 0.40) is a dimensionless constant in the logarithmic wind velocity profile (Högström, 1988), µ m is average daily wind velocity (m s −1 ), d is calculated as 2/3 of the vegetation canopy height, z om is roughly 1/8 of canopy height (Yang and Friedl, 2002), and z ov is 0.1 z om (Yan et al., 2012). The estimated g a is then scaled between 0 and 1 to obtain λ using a linear function for sites with a lower observed sensitivity to surface turbulence; for environments with a higher sensitivity to surface turbulence, a quadratic approach is used when µ m exceeds 4 m s −1 : Although this approach focuses on the influence of wind turbulence on plant gas transport within vegetated wetlands, it is also applicable for inundated microsites where increases in surface water mixing can stimulate CH 4 degassing . In addition, Eq. (15) reflects near-neutral atmospheric stability and adjustments may be necessary to accommodate unstable or stable atmospheric conditions (Raupach, 1998). The upward diffusion of CH 4 within the soil profile is determined using a one-layer approach similar to Tian et al. (2010). The rate of CH 4 transport (D e ; m −2 d −1 ) is considered for both saturated (D water ; 1.73 × 10 −4 µM CH 4 d −1 ) and aerated (D air ; 1.73 µM CH 4 d −1 ) soil fractions: Potential daily transport through diffusion (P diff ) is estimated as the product of D e and the gradient between C CH 4 and the concentration of CH 4 in the atmosphere (Air CH 4 ). This is further modified by soil tortuosity (τ ; 0.66), which increases exponentially for T s < 274 K to account for slower gas movement at colder temperatures and barriers to diffusion resulting from near-surface ice formation (Walter and Heimann, 2000;Zhuang et al., 2004), and pathway constraints within the saturated pore fraction (1 − θ ): A portion of diffused CH 4 is oxidized (R ox ) before reaching the soil surface, using a Michaelis-Menten kinetics approach that is scaled by φ a : where V max is the maximum reaction rate and K m is the substrate concentration at 0.5 V max (van Huissteden et al., 2006). Oxidation during soil diffusion is modified by soil temperature Q 10 constraints (Q 10d ); T d is the reference temperature and can be defined using site-specific mean annual T s (Le Mer and Roger, 2001). Total daily CH 4 emission (F diff ) from the soil diffusion pathway is determined by substracting R ox from P diff . The CH 4 algorithm uses a gradient-based approach to account for slow or "steady-rate" ebullition from inundated micro-sites in the landscape (Rosenberry et al., 2006;Wania et al., 2010), whereas episodic events originating deeper within the soil require more complex modeling techniques and input data requirements (Kettridge et al., 2011) that are beyond the scope of this study. Emission contributions due to ebullition occur when C CH 4 exceeds a threshold value (v e ) of 500 µM (van Huissteden et al., 2006). The magnitude of gas release is determined by steady-rate bubbling (C e ) applied within the saturated soil pore space (ϕ s ): (20)

Study sites and in situ data records
Tower EC records from six pan-Arctic peatland and tundra sites in Finland, Sweden, Russia, Greenland and Alaska were used to assess the integrated TCF model CO 2 and CH 4 simulations ( Fig. 1; Table 1). The Scandinavian tower sites include Siikaneva (SK) in southern Finland and Stordalen Mire (SM) in northern Sweden near the Abisko Scientific Research Station. The Lena River (LR) delta site is located on Samoylov Island in northern Siberia and EC measurements from the Kytalyk (KY) flux tower were collected near Chokurdakh in northeastern Siberia. The Zackenberg (ZK) flux tower is located within the Northeast Greenland National Park, and tower data records for Alaska were obtained from a water table manipulation experiment (Zona et al., 2009;2012;Sturtevant et al., 2012) approximately 6 km east of Barrow (BA). With the exception of Siikaneva, the EC tower footprints represent wet permafrost ecosystems with complex, heterogeneous terrain that includes moist depressions, drier, elevated hummocks and inundated microsites. Vegetation within the tower footprints Riutta et al., 2007;Sachs et al., 2008;Jackowicz-Korczyński et al., 2010;Parmentier et al., 2011a;Zona et al., 2011;Tagesson et al., 2012b) consists of Carex and other sedges, dwarf shrubs (e.g., Dryas and Salix), grasses (e.g., Arctagrostis) and Sphagnum moss (with the exception of Zackenberg). Mean daily T s and θ site measurements corresponding to near-surface (≤ 10 cm) soil depths were selected when possible (Table 1), to better coincide with the soil penetration depths anticipated for upcoming satellite-based microwave remote sensing missions . For Siikaneva, reanalysis θ was used in place of in situ measurements as only water table depth information was available to describe soil wetness . At the Lena River site T s and θ (≤ 12 cm) observations were obtained from the nearby Samoylov meteorological station and represent tundra polygon wet center, dry rim and slope conditions Sachs et al., 2008). Although θ was also measured during the summer of 2006, the in situ records are limited to the wet polygon center location (Boike, personal communication, 2012) and were not used in this study due to the potential for overestimating saturated site conditions. For Zackenberg, site T s measurements were obtained at a 2 cm depth (Tagesson et al., 2012a, b) within the tower footprint, while near-surface θ (< 20 cm) and ≥ 5 cm T s measurements were collected adjacent to the site (Sigsgaard et al., 2011). At Stordalen, site θ measurements were not avail- able at the time of this study (Jackowicz-Korczyński et al., 2010). Barrow (Zona et al., 2009;Sturtevant et al., 2012) includes southern (S), central (C) and northern (N) tower locations; in 2007 only CO 2 and CH 4 EC measurements from the northern tower were used in the analysis, due to minimal EC data availability for the other tower sites following data processing (Zona et al., 2009). Many of the Barrow CO 2 measurements were also rejected for the 2009 period; as a result NEE was not partitioned into R eco and GPP (Sturtevant et al., 2012).

Remote sensing and reanalysis inputs
Daily input meteorology was obtained from the Goddard Earth Observing System Data Assimilation Version 5 (GEOS-5) MERRA archive (Rienecker et al., 2011) with 1/2 × 2/3 • spatial resolution. The MERRA records were recently verified for terrestrial CO 2 applications in highlatitude systems (Yi et al., 2011(Yi et al., , 2013Yuan et al., 2011), and provide model enhanced T s and surface θ information similar to the products planned for the NASA Soil Moisture Active Passive (SMAP) mission . In addition to near-surface (≤ 10 cm) T s and θ information from the MERRA-Land reanalysis  required for the R eco and CH 4 simulations, daily MERRA SW rad , T min and VPD records were used to drive the internal GPP calculations. The MERRA near-surface (2 m) wind parameters were also used to obtain mean daily µ m for the CH 4 simulations. The MERRA-Land records for Greenland are spatially limited due to land cover/ice masking inherent in the reanalysis product, and MERRA T s and θ were not available  Sturtevant et al. (2012) for the Zackenberg tower site. As a proxy, T s was derived from reanalysis surface skin temperatures by applying a simple Crank-Nicholson heat diffusion scheme that accounts for energy attenuation with increasing soil depth (Wania et al., 2010); for θ, records from a nearby grid cell were used to represent moisture conditions at Zackenberg. For the daily LUE-based GPP simulations, quality screened cloud-filtered 16 day 250 m NDVI values from MODIS Terra (MOD13A1) and Aqua (MYD13Q1) data records (Solano et al., 2010) were used as model inputs. Differences between the MOD13A1 and MYD13Q1 retrievals were minimal at the tower locations, and the combination of Terra and Aqua MODIS records reduced the retrieval gaps to approximate 8 day intervals. The NDVI retrievals correspond to the center coordinate locations for each flux tower site, and temporal linear interpolation was used to scale the 8 day NDVI records to daily inputs. Coarser (500-1000 m resolution) NDVI records were not used in this study due to the close proximity of water bodies at the tower sites, which can substantially reduce associated FPAR retrievals. In addition, 250 m MODIS vegetation indices have been reported to better capture the overall seasonal variability in tower EC flux records .

TCF model parameterization
A summary of the site-specific TCF model parameters is provided in the Supplement (Table S2). Parameter values associated with grassland biomes were selected for the LUE model VPD and T min modifiers used to estimate GPP (Yi et al., 2013), as more specific values for tundra and mossdominated wetlands were not available. Parameter values for θ max were obtained using growing-season maximum θ measurements for each site and θ min was set to 0.15 for scaling purposes. Model ε max was specified as 0.82 mg C MJ −1 for the duration of the growing season, although actual LUE can vary throughout the summer due to differences in vegetation growth phenology and nutrient availability (Connolly et al., 2009;King et al., 2011). The tundra CUE ranged from 0.45 to 0.55 (Choudhury, 2000); a lower CUE value of 0.35 was used for the moss-dominated Siikaneva site due to a more moderate degree of carbon assimilation occurring in bryophytes that has been observed in other sub-Arctic communities (Street et al., 2012). For the TCF model F met parameter, the percentage of NPP allocated to C met varied between 70 % and 72 % for tower tundra sites (Kimball et al., 2009) compared to 50 % and 65 % for Siikaneva and Stordalen, where moss cover is more abundant. The TCF model R o parameter ranged from 4.5 and 22.4 µM CH 4 d −1 (Walter and Heimann, 2000;van Huissteden et al., 2006). Values for Q 10p varied between 3.5 and 4 due to an enhanced microbial response to temperature variability under colder climate conditions (Gedney and Cox, 2003;Inglett et al., 2012). A Q 10d of 2 was assigned for CH 4 oxidation (Zhuang et al., 2004;van Huissteden et al., 2006). Parameter values for P trans , which indicates relative plant transport ability, ranged from 7 to 9 (dimensionless); lower values were assigned to tower locations with a higher proportion of shrub and moss cover, whereas higher P trans corresponds to sites where sedges are more prevalent (Ström et al., 2005;. For λ, the scaled conductance for lower site wind sensitivity was used in the CH 4 model simulations, except for the Lena River, which showed higher sensitivity to surface turbulence. Values for P ox ranged from 0.7 in tundra to 0.8 in Sphagnum-dominated systems to account for higher CH 4 oxidation by peat mosses (Parmentier et al., 2011c). Due to a lack of detailed soil profile descriptions and heterogeneous tower footprints, soil porosity was assigned at 75 % for sites with more abundant fibrous surface layer peat (i.e., Siikaneva and Stordalen) and 70 % elsewhere to reflect more humified or mixed organic and mineral surface soils (Elberling et al., 2008;Verry et al., 2011).

TCF model simulations
The TCF model was first evaluated against tower EC records using simulations driven with in situ environmental data including EC-based GPP, T s , θ and µ m . This step allowed for baseline TCF model R eco and CH 4 flux estimates to be assessed without introducing additional uncertainties from input reanalysis meteorology and LUE model derived GPP calculations. Four additional TCF model simulations were conducted using reanalysis θ, T s , µ m (in the CH 4 module), or internal model GPP in place of the in situ data. A final TCF model run included only satellite and reanalysis-based data, and was used to establish annual GPP, R eco and CH 4 carbon budgets for each site. Baseline carbon pools were initialized by continuously cycling ("spinning-up") the model for the tower years of record (described in Table 1) to reach a dynamic steady-state between estimated NPP and surface soil organic carbon stocks (Kimball et al., 2009). In situ data records were used during the model spin-up to establish baseline organic carbon conditions for the first five TCF model simulations, although it was often necessary to supplement these data with reanalysis information to obtain a continuous annual time series. The final model simulation did not include in situ data in the spin-up process.
The temporal agreement between the tower EC records and TCF model simulations was assessed using Pearson correlation coefficients (r; ±one standard deviation) for the daily, 8 day, and total-period (EC length of record) cumulative carbon fluxes and corresponding tests of significance at a 0.05 probability level. The 8 day and total-period cumulative fluxes were evaluated, in addition to the daily fluxes, to account for differences between the model estimates and tower EC records stemming from temporal lags between changing environmental conditions and resulting carbon (CO 2 , CH 4 ) emissions (Lund et al., 2010;Levy et al., 2012). The mean residual error (MRE) between the tower EC records and TCF modeled CO 2 and CH 4 fluxes was used to identify potential positive (underestimation) and negative (overestimation) biases in the simulations; root-mean-square-error (RMSE) differences were used as a measure of model estimate uncertainty in relation to the tower EC records.

Surface organic carbon pools
The TCF model generated surface soil organic carbon pools represent steady-state conditions obtained through the continuous cycling of in situ or satellite and reanalysis environmental data for the years of record associated with each tower site (described in Table 1). Approximately 600 and 1000 years of model spin-up were required for C rec to reach dynamic steady state conditions. Over 95 % of the resulting total carbon pool was allocated to C rec by the TCF model, with 2-3 % stored as C met and the remainder partitioned to C str . The estimated carbon pools from the in situ (reanalysis-based) model spin-up ranged from approximately 3.3 kg C m −2 (2.3 kg C m −2 ) for Zackenberg and Stordalen to 1.3 kg C m −2 (2.1 kg C m −2 ) for the other tower sites.
Differences in carbon stocks, resulting from the use of satellite remote sensing and reanalysis information in the TCF model, reflect warm or cold biases in the input T s records relative to the in situ data that modified the rate of CO 2 loss during model initialization. The larger carbon stocks at Zackenberg, compared to the other tundra sites, resulted from higher tower EC-based GPP inputs that often exceeded 5 g C m −2 d −1 in mid-summer, and a short (< 50 day) peak growing season (Tagesson et al., 2012a) that minimized TCF modeled R h losses. Although it was necessary to use internal LUE-based GPP calculations for Stordalen in the absence of available CO 2 records, the resulting C met and C rec carbon stocks were similar in magnitude to surface litter measurements at this site (Olsrud and Christensen, 2011). The TCF model simulated carbon stock for Lena River was less than a 2.9 kg C m −2 average determined from in situ (≤ 10 cm depth) measurements of nearby river terrace soils (Zubrzycki et al., 2013), but this could have resulted from site spatial heterogeneity and the use of recent climate records in the model spin-up that may not reflect past conditions.

LUE-based GPP
The GPP simulations using reanalysis and satellite-based inputs captured the overall seasonality observed in the tower records ( Fig. 2; Table 2) and explained 76 % (r 2 ; p < 0.05, N = 7) of variability in the total EC period-ofrecord fluxes (Fig. 3). The across-site RMSE and MRE were 1.3 ± 0.51 and −0.1 ± 0.7 g C m −2 d −1 , respectively. Although the 8 day cumulative flux correspondence between the tower EC and TCF model GPP estimates was strong (r 2 = 75±16 %), the model-tower agreement decreased considerably for daily GPP (r 2 = 57 ± 22 %). These differences may reflect a delayed response in vegetation productivity following changes in atmospheric and soil conditions (Lund et al., 2010), and short term fluctuations in the reanalysis SW rad inputs. For Kytalyk, the large RMSE (2.2 g C m −2 d −1 ) observed for the TCF model GPP simulations resulted from warm spring air temperatures that reduced T min constraints on carbon assimilation, although a similar increase in GPP did not occur in the EC-based records. This lack of response likely resulted from a shallow (< 14 cm) early season thaw depth at this site, that limited bud break activity in deeper rooted shrubs (e.g., Betula nana and Salix pulchra). To address this, an additional simulation was conducted using a temperature driven phenology model described in Parmentier et al. (2011a) to better inform the start of growing season in the TCF model. This step reduced the corresponding RMSE Fig. 3. Correspondence between TCF model and tower EC records for cumulative (g C m −2 ) GPP, R eco , NEE, and CH 4 fluxes from six pan-Arctic tower locations. The TCF model simulations include those derived from in situ measurements (open circles) or MODIS remote sensing and MERRA reanalysis inputs (MDMR; in black). A 1 : 1 relationship is indicated by the dashed line. The r 2 agreement is significant at a 0.05 probability level, except for MDMRbased R eco and NEE (p = 0.16 and 0.27), and excludes NEE fluxes for KY (circled) due to large differences in the CO 2 response relative to the other sites. difference for Kytalyk by 56 % (to 1 g C m −2 d −1 ) with an associated r 2 of 67 %.
Although previous LUE models (e.g., Running et al., 2004;Yi et al., 2013) have relied solely on VPD to represent water related constraints to GPP, our approach also considers soil moisture to better account for the sensitivity of bryophytes and shallow rooted vegetation to surface drying (Wu et al., 2013). Including this additional moisture constraint reduced the overall TCF model and tower GPP RMSE and MRE differences by approximately 14 % and 92 %. However, the model simulations continued to overestimate GPP fluxes for Siikaneva, Lena River (2003), and Kytalyk (MRE = −0.6 ± 0.8 g C m −2 d −1 ). This residual GPP bias could be influenced by inconsistencies between the coarse-scale MERRA reanalysis inputs and local tower meteorology, as reported elsewhere (e.g., Yi et al., 2013), although systematic biases for the high-latitude regions have not been identified. For instance, periods of warmer (3 to 4 • C) reanalysis T min inputs relative to in situ measurements at Lena River in 2003 led to seasonally higher TCF modeled GPP fluxes. In contrast, the reanalysis T min at Barrow was 2 to 7 • C cooler in mid-summer than the local meteorology; this resulted in significantly lower (p < 0.05) TCF model GPP estimates relative to the tower EC records (Table 2). It is also possible that differences in the light response curve and respiration models, used when partitioning the site EC NEE fluxes into GPP and R eco (i.e., Aurela et al., 2007;Table 2. Tower EC CO 2 records and TCF modeled gross primary production (GPP), ecosystem respiration (R eco ) and net ecosystem exchange (NEE) derived using in situ information (in parentheses) or satellite remote sensing and reanalysis inputs. The Pearson correlation coefficients (r) are significant at a 0.05 probability level, excluding Kytalyk 2009 NEE (r ≤ 0.11, p ≥ 0.17) and Barrow 2007N GPP and NEE (r < 0.1, p ≥ 0.16).

Site
Year  Kutzbach et al., 2007;Parmentier et al., 2011a;Tagesson et al., 2012a;Zona et al., 2012), may have contributed to differences between the TCF model simulations and tower CO 2 records. However, further investigation is needed to determine the expected range of GPP and R eco that might result from variability in the flux partitioning routines.

R eco and NEE
The in situ TCF model R eco simulations accounted for 59 ± 28 % and 76 ± 24 % (r 2 ) of the observed variability in the respective daily and 8 day cumulative tower EC fluxes ( Fig. 4; Table 2). As with GPP, the r 2 agreement increased to 89 % (p < 0.05, N = 6) when considering the total-period cumulative fluxes (Fig. 3). The overall RMSE difference for the in situ based TCF model R eco and NEE simulations was 0.74 ± 0.45 g C m −2 d −1 when using 5 cm depth T s inputs. A corresponding across-site MRE of −2.1 ± 5.7 g C m −2 d −1 indicated that the TCF model simulations overestimated R eco relative to the tower records, and slightly underestimated NEE (MRE = 0.1 ± 0.4 g C m −2 d −1 ). We also conducted TCF model simulations using 8-10 cm depth in situ T s inputs, instead of those from ≤ 5 cm (as reported in Table 2), to investigate the influence of deeper soil thermal controls on site R eco response; this step reduced the overall RMSE by approximately 12 %. Incorporating the TCF internal LUE model GPP estimates increased the overall RMSE for R eco and NEE by 23 % relative to the in situ based simulations, compared to a respective 3 % and 14 % increase when using reanalysis θ or T s inputs (Fig. 5). The model-tower daily and 8 day cumulative correspondence was also lower (r 2 = 32 and 56 %, respectively) for CO 2 simulations driven using internally derived GPP, relative to those using reanalysis θ or T s inputs (r 2 = 57 and 72 %) in place of the in situ records. Without the in situ inputs, the respective RMSE and MRE difference between the reanalysis-based R eco (NEE) simulations and the tower EC records averaged 0.9 ± 0.4 and −0.2 ± 0.9 g C m −2 d −1 (1 ± 0.5 and 0.3 ± 0.05 g C m −2 d −1 ).
Correspondingly, the reanalysis and remote sensing-based TCF model R eco (NEE) simulations accounted for 51 ± 29 (45 ± 34) % and 71 ± 17 (62 ± 34) % of the observed r 2 variability in the respective daily and 8 day tower EC records. and Kytalyk, which did not show significant (r ≤ 0.20; p ≥ 0.16) agreement with the site EC records (Table 2). For Barrow, it is likely that the water table manipulations at this site led to local temperature and moisture variability that was not reflected in the coarse reanalysis and remote sensing inputs. The minimal agreement at Kytalyk is attributed to higher R h losses driven by warmer reanalysis T s inputs, and increased R a contributions due to the overestimation of GPP relative to the tower EC records.

CH 4 fluxes
The in situ TCF model CH 4 simulations explained 64 ± 11 % and 80 ± 12 % (r 2 ) of the respective daily and 8 day cumulative variability observed in the tower EC records ( Fig. 6; Table 3), when excluding Kytalyk (p = 0.1). The r 2 correspondence increased to 98 % when considering the total periodof-record emissions across the six sites ( Fig. 3; p < 0.05, N = 9). At Kytalyk, Parmentier et al. (2011b) reported large differences in measured half-hourly CH 4 fluxes following shifts in wind direction, and larger emissions from portions of the tower footprint containing Carex sp., E. angustifolium and inundated microsites. Although this may have contributed to the observed discrepancy between the TCF model estimates and tower EC record, attempts to systematically screen the CH 4 observations based on wind direction, or to use daily EC medians instead of mean values, did not substantially improve the model results.
On average, the in situ TCF model simulations overestimated CH 4 fluxes relative to the tower EC records (MRE = −2.2 mg C m −2 d −1 ), with RMSE differences varying from 6.7 to 42.5 mg C m −2 d −1 . Without including µ m in the TCF model, the resulting RMSE increased by > 10 % and the mean daily correspondence decreased to r 2 < 40 %. The most substantial difference was observed for Lena River, where excluding µ m reduced the daily and 8 day emission correspondence by over 60 %. Unlike the TCF model R eco results, deeper (10 cm depth) T s measurement inputs did not improve the RMSE values, except for Barrow (2007N) where the RMSE decreased by 35 %. This sensitivity to deeper T s conditions may reflect changes in active layer depth following water table manipulations at this site (Zona et al., 2009, and associated changes in carbon substrate availability. In contrast, the RMSE for Lena River was 15 % higher when using in situ 10 cm T s records in the TCF model simulations instead of 5 cm depth measurements. A 6 % decrease in the RMSE occurred for Zackenberg (2008) when using the warmer (3 to 5 • C) 2 cm depth T s records, relative to model simulations using 5 cm T s inputs. Contrary to expectations, the 2 cm depth T s inputs did not improve RMSE differences for Zackenberg in 2009 when site moisture conditions were drier (Tagesson et al., 2012a). The reanalysis driven TCF model CH 4 simulations ( Fig. 6; Table 3) accounted for 48 ± 16 % and 79 ± 8 % (r 2 ) of the respective daily and 8 day variability in the tower EC records when excluding the less favorable results for Kytalyk (r 2 = 8 and 44 %, respectively). Although slightly lower than the in situ TCF model CH 4 estimates, the coarser reanalysis and remote sensing driven simulations explained 96 % (r 2 ) of the total period-of-record emissions at these sites (Fig. 3). The corresponding model RMSE was 18.2 ± 13.6 mg C m −2 d −1 , with an associated MRE difference of 1.8 ± 7.3 mg C m −2 d −1 that indicated the slight model underestimation of daily CH 4 emissions. The model RMSE differences increased by approximately 15 % when using reanalysis µ m records or internal GPP estimates in place of the in situ inputs, and by 10 % when incorporating reanalysis T s and θ inputs (Fig. 7).

Estimates of annual carbon budgets
The reanalysis and remote sensing driven TCF model simulations indicated a net CO 2 sink (NEE = −34.5 ± 18.5 g C m −2 yr −1 ) for the tower sites, excluding Barrow in 2009 (NEE = 7.3 g C m −2 yr −1 ) where the estimated R eco emissions exceeded annual GPP (Fig. 8).
Other studies near Barrow have also reported NEE losses from wet tundra communities, resulting from drier microscale surface conditions and warming within the hummocky landscape (Huemmrich et al., 2010b;Sturtevant and Oechel, 2013) that can strongly influence R eco . The corresponding TCF model R eco estimates ranged from 133 (Zackenberg in 2009) to 494 g C m −2 yr −1 (Stordalen in 2006) with lower CO 2 emissions occurring in the colder, more northern tundra sites. The strongest NEE carbon sink indicated by the model simulations was observed for the peat-rich Siikaneva site (−70.3 g C m −2 yr −1 ) due to high annual GPP (462.5 g C m −2 yr −1 ) relative to the other tower locations. Although tower EC CO 2 records were not available for Stordalen to verify the TCF model NEE results (−50.8 and −65.8 g C m −2 yr −1 , respectively), the estimates are slightly smaller (∼ 30 g C m −2 d −1 ) than other NEE approximations over the same time period , but are similar to observations reported for other years at this site (Olefeldt et al., 2012;Marushchak et al., 2013). The annual TCF model CH 4 estimates determined using the reanalysis inputs averaged 6.9 (± 5.5) g C m −2 yr −1 for the six tower sites. The highest CH 4 emissions were observed for Stordalen and Siikaneva (≥ 11.8 g C m −2 yr −1 ) due to higher model-defined CH 4 production rates and summer reanalysis T s records that were often 5 • C warmer than the other sites. In contrast, model CH 4 emissions were lowest for Barrow (1.8 g C m −2 yr −1 ) due to smaller GPP estimates and colder summer reanalysis T s records that did not reflect the unusually warm site conditions in 2007 (Shiklomanov et al., 2010). The annual TCF model CH 4 emissions for Lena River were relatively small (2.3 g C m −2 yr −1 , on average), but are similar in magnitude to site CH 4 estimates determined using more complex coupled biogeochemical and permafrost models (i.e., Zhang et al., 2012b). Although the TCF modeled CH 4 fluxes contributed only 1-5 % of annual carbon emissions (R eco + CH 4 ) at the tower sites, which is similar to previous reports (Schneider von Deimling et al., 2012), these CH 4 emissions reduced the NECB (−23.3 ± 19.6 g C m −2 yr −1 ) by approximately 23 % relative to NEE. The annual model estimates indicated that the site CO 2 and CH 4 fluxes, excluding Bar- row and Lena River, contributed to a net global warming potential (GWP) of 188 ± 68 g CO 2 eq m −2 yr −1 over a 100 year time horizon (Boucher et al., 2009) with total GWP influences by CH 4 at approximately 9 to 44 % that of R eco . Similarly, the Lena River and Barrow sites mitigated GWP at a mean rate of −40 g CO 2 eq m −2 yr −1 in 2006 and 2007, but were net GWP contributors in 2003 and 2009 (25 and 160 g CO 2 eq m −2 yr −1 , respectively). Although site CO 2 contributions from methantrophy during plant transport and soil diffusion were estimated to range from 3.8 to 58.3 g C m −2 yr −1 , these contributions represented < 14 % of total TCF model derived R eco .

Discussion and conclusions
The level of complexity in biophysical process models has increased considerably in recent years but there remain large differences in carbon flux estimates for northern highlatitude ecosystems (McGuire et al., 2012;Wania et al., 2013). An integrated TCF model CO 2 and CH 4 framework was developed to improve carbon model compatibility with remote sensing retrievals that can be used to inform changes in surface conditions across northern peatland and tundra regions. Although the TCF model lacks the biophysical and hydrologic complexity found in more sophisticated process models (e.g., Zhuang et al., 2004;Wania et al., 2010), it avoids the need for extensive parameterization by instead employing generalized surface vegetation growth, temperature, and moisture constraints on ecosystem CO 2 and CH 4 fluxes. Despite the relatively simple model approach and landscape heterogeneity at the tower sites, the TCF model simulations derived from local tower inputs captured the overall seasonality and magnitude of R eco and CH 4 fluxes observed in the tower EC records. Overall the R eco , NEE and CH 4 emission simulations determined using local site inputs showed strong mean correspondence (8 day r > 0.80; p < 0.05) with tower EC records, but the strength of agreement varied considerably for the daily fluxes due to temporal lags between changing environmental conditions and carbon emissions (Zhang et al., 2012b), and larger EC measurement uncertainty at the daily time step (Baldocchi et al., 2008;Yi et al., 2013). The respective RMSE differences from the in situ TCF model CO 2 and CH 4 simulations averaged 0.7 ± 0.4 g C m −2 d −1 and 17.9 ± 11.5 mg C m −2 d −1 , which is comparable to other site-based model results (e.g., Marushchak et al., 2013;Sturtevant and Oechel, 2013).
In this study, we used near-surface T s records in the model simulations to better coincide with the soil depths represented by upcoming satellite remote sensing missions, but acknowledge that deeper T s controls are also important for regulating high-latitude carbon emissions. This was evident in TCF model R eco results where RMSE differences between the in situ based simulations and tower EC fluxes generally improved when using deeper 10 cm T s inputs instead of those from shallower (≤ 5 cm) soil depths. However, the TCF model CH 4 simulations were more favorable when using near-surface (2 to 5 cm) T s inputs. The observed CH 4 emission sensitivity to surface soil warming may be influenced by cold temperature constraints on CH 4 production in the carbon-rich root zone where organic acids are more abundant (Turetsky et al., 2008;Olefeldt et al., 2013). Lightweight carbon fractions have been shown to be more suscep-tible to mineralization following soil thaw and temperature changes than heavier, more recalcitrant soil organic carbon pools in high-latitude environments (Glanville et al., 2012). However, the depletion of older organic carbon stocks may also become more prevalent in permafrost soils subject to thawing and physiochemical destabilization (Schuur et al., 2009;Hicks Pries et al., 2013a) in the absence of wet, anoxic conditions (Hugelius et al., 2012;Hicks Pries et al., 2013b). Seasonal changes in T s constraints were also evident in this study, especially in the Zackenberg records where the TCF model underestimated tower R eco and CH 4 emissions in autumn by not accounting for warmer temperatures deeper in the active layer that can sustain microbial activity following surface freezing (Aurela et al., 2002). Allowing the TCF model vegetation CUE parameter to change over the growing season instead of allocating R a as a static fraction of GPP may also improve model and tower R eco agreement. In Arctic tundra, R a can contribute anywhere from 40 to 70 % of R eco , with higher maintenance and growth respiration occurring later in the growing season when root systems expand deeper into the soil active layer (Hicks Pries et al., 2013a). Representing R a as a fixed proportion of daily GPP in the TCF model, and not accounting for the use of stored plant carbon reserves, may also have contributed to the lower R eco estimates during spring and autumn transitional periods when photosynthesis is reduced.
Our estimates of peatland and tundra CO 2 fluxes using TCF model simulations driven by MERRA reanalysis and satellite (MODIS) remote sensing inputs showed favorable agreement relative to the tower EC observations, with relatively moderate RMSE uncertainties of 1.3 ± 0.5 (GPP), 0.9 ± 0.4 (R eco ) and 1 ± 0.5 (NEE) g C m −2 d −1 . These model accuracies are similar to those reported in a previous TCF model analysis for the northern regions (Yi et al., 2013), and other Arctic LUE-based GPP studies (Tagesson et al., 2012a;McCallum et al., 2013). The associated modeltower RMSE for CH 4 was 18.2 ± 13.6 mg C m −2 d −1 , and is comparable to results from previous remote sensing driven CH 4 analyses (Meng et al., 2012;Tagesson et al., 2013). The larger observed differences between TCF model and tower EC-based GPP results may reflect seasonal changes in nutrient availability (Lund et al., 2010), although one peatland study reported that nutrient limitations to plant productivity could be detected indirectly by MODIS NDVI retrievals (Schubert et al., 2010b). It is more likely that this reduced correspondence resulted from fluctuations in the reanalysis SW rad inputs (Yi et al., 2011) and uncertainty associated with satellite NDVI and resulting FPAR inputs stemming from residual snow cover and surface water effects on optical-IR reflectances (Delbart et al., 2005). High-latitude studies have reported difficulty in using satellite NDVI to determine the start of spring bud burst and seasonal variability in leaf development (Huemmrich et al., 2010a). Evaluating other portions of the visible spectrum, including blue and green reflectances, in addition to NDVI has helped to alleviate this problem in remote sensing applications (Marushchak et al., 2013) and should be considered in subsequent studies. Incorporating phenological constraints into the TCF LUE model may also better characterize early season GPP, especially for plant communities such as E. vaginatum that are sensitive to changes in active layer depth (Parmentier et al., 2011a;Natali et al., 2012). Considering T s as an additional constraint in the TCF LUE model may also better account for autumn GPP activity under frozen air temperatures if plant-available moisture is still available within the root zone (Christiansen et al., 2012). Yi et al. (2013) attempted to address this condition by incorporating satellite passive microwave-based freeze/thaw records (37 GHz) to constrain GPP according to frozen, transitional, or non-frozen surface moisture states but did not report a significant improvement, likely due to the coarse (25 km) resolution freeze/thaw retrievals.
The TCF model assessment of annual NECB for the six northern tower EC sites indicate that CH 4 emissions reduced the terrestrial net carbon sink by 23 % relative to NEE. Although GPP at the Lena River and Barrow sites mitigated GWP additions from R eco and CH 4 in two of the years examined, in most years the tower sites were GWP contributors by approximately 165 ± 128 g CO 2 eq m −2 yr −1 when considering the impact of CH 4 on atmospheric forcing over a 100 year time span. These results are consistent with other model-based analyses of Arctic carbon fluxes (McGuire et al., 2010) and emphasize the importance of evaluating CO 2 and CH 4 emissions simultaneously when quantifying the terrestrial carbon balance and GWP for northern peatland and tundra ecosystems Olefeldt et al., 2012). However, ongoing efforts are needed to better inform landscape-scale spatial/temporal variability in soil moisture, temperature and vegetation controls on CO 2 and CH 4 fluxes for future model assessments using a combined network of in situ soil measurements and strategically placed EC tower sites (Sturtevant and Oechel, 2013), and regional airborne surveys. The upcoming SMAP mission may also help to determine landscape soil moisture and thermal constraints on northern carbon fluxes through relatively fine-scale (3 km resolution) and lower frequency (≤ 1.4 GHz) microwave retrievals with enhanced soil sensitivity (Entekhabi et al., 2010;Kimball et al., 2012), complimented by recent improvements in Arctic-specific reanalysis data (Bromwich et al., 2010). These advances, in conjunction with a suitable model framework to quantify ecosystem NEE and CH 4 emissions, provide the means for regional carbon assessments and monitoring of the net ecosystem carbon budget and underlying environmental constraints.