Landsat near-infrared (NIR) band and ELM-FATES sensitivity to forest disturbances and regrowth in the Central Amazon

Forest disturbance and regrowth are key processes in forest dynamics, but detailed information on these processes is difficult to obtain in remote forests such as the Amazon. We used chronosequences of Landsat satellite imagery (Landsat 5 Thematic Mapper and Landsat 7 Enhanced Thematic Mapper Plus) to determine the sensitivity of surface reflectance from all spectral bands to windthrow, clear-cut, and clear-cut and burned (cut+ burn) and their successional pathways of forest regrowth in the Central Amazon. We also assessed whether the forest demography model Functionally Assembled Terrestrial Ecosystem Simulator (FATES) implemented in the Energy Exascale Earth System Model (E3SM) Land Model (ELM), ELM-FATES, accurately represents the changes for windthrow and clearcut. The results show that all spectral bands from the Landsat satellites were sensitive to the disturbances but after 3 to 6 years only the near-infrared (NIR) band had significant changes associated with the successional pathways of forest regrowth for all the disturbances considered. In general, the NIR values decreased immediately after disturbance, increased to maximum values with the establishment of pioneers and early successional tree species, and then decreased slowly and almost linearly to pre-disturbance conditions with the dynamics of forest succession. Statistical methods predict that NIR values will return to pre-disturbance values in about 39, 36, and 56 years for windthrow, clear-cut, and cut+ burn disturbances, respectively. The NIR band captured the observed, and different, successional pathways of forest regrowth after windthrow, clear-cut, and cut+ burn. Consistent with inferences from the NIR observations, ELMFATES predicted higher peaks of biomass and stem density after clear-cuts than after windthrows. ELM-FATES also predicted recovery of forest structure and canopy coverage back to pre-disturbance conditions in 38 years after windthrows and 41 years after clear-cut. The similarity of ELM-FATES predictions of regrowth patterns after windthrow and clearcut to those of the NIR results suggests the NIR band can be used to benchmark forest regrowth in ecosystem models. Our results show the potential of Landsat imagery data for mapping forest regrowth from different types of disturbances, benchmarking, and the improvement of forest regrowth models.

forest degradation (Silverio et al., 2018;Schwartz et al., 2017). Under natural conditions, disturbed forests recover to their pre-disturbance conditions through complex interactions that vary across spatial and temporal scales (Chazdon, 2014). In general, it is known that forest pathways of regrowth (i.e., pattern of regrowth) are initiated with fastgrowing and shade-intolerant species (pioneers) that are established from seeds and dominate a few years after disturbance, followed by the recruitment and establishment of shade-tolerant species, and finally the establishment of a closed-canopy old-growth forest (Chazdon, 2014;Denslow, 1980;Mesquita et al., 2001;Swaine and Whitmore, 1988). Understanding of the dynamics of forest regrowth following natural and anthropogenic disturbances in the Amazon, however, has so far been limited by lack of long-term observational data showing different stages of forest regrowth.
Remote sensing data can be used to assess forest regrowth via changes in spectral characteristics (Frolking et al., 2009;Roberts et al., 2004;Schroeder et al., 2011;DeVries et al., 2015;Lucas et al., 2002;McDowell et al., 2015). Landsat satellite imagery is appropriate for examining land surface changes due to its long-term record availability and spatial resolution of 30 m (Loveland and Dwyer, 2012;Wulder et al., 2012;Alcantara et al., 2011;Woodcock et al., 2008;Cohen and Goward, 2004;Hansen et al., 2013). Landsat imagery has been used to detect forest disturbance and pathways of regrowth in temperate and boreal forests in the United States and Canada (Kennedy et al., 2007(Kennedy et al., , 2012Pickell et al., 2016;Schroeder et al., 2011;Dolan et al., 2009Dolan et al., , 2017 and for the detection of forest disturbance and regrowth of biomass in the Amazon (Vieira et al., 2003;De-Vries et al., 2015;Lucas et al., 2002;Powell et al., 2010;Lu and Batistella, 2005;Steininger, 2000;Shimabukuro et al., 2019). These studies suggest that Landsat may be sensitive to different types of disturbances and their subsequent pathways of forest regrowth in the Amazon, but these applications have not yet been assessed.
The ability to forecast future trajectories of forests depends upon the fidelity with which disturbance and regrowth processes are represented within terrestrial biosphere models. These models capture processes operating between the leaf and landscape scales and can represent regrowth changes over large regions (Fisk, 2015), long time periods (Holm et al., 2017;Putz et al., 2014), a range of disturbance intensities (Powell et al., 2013), and interactions between multiple disturbance types and disturbance histories (Hurtt et al., 2006). But how well these models simulate and capture the diverse array of successional pathways of forest regrowth after anthropogenic or natural disturbances needs to be more thoroughly evaluated, given observed increases in disturbance rates (Lewis et al., 2015). The few modeling studies analyzing tropical disturbances have focused on the effects of fragmented edges or the regrowth of specific tree species (Dantas de Paula et al., 2015;Kammesheidt et al., 2002).
Cohort-based dynamic vegetation demographic models (VDMs) are particularly suitable tools for expanding upon these studies . In contrast to traditional land surface models, VDMs include ecological demographic processes, such as discretized vegetation height, with different plant types competing for light within the same vertical profile, and heterogeneity in light availability along disturbance and recovery trajectories, all of which facilitate direct simulation of regrowth dynamics during succession. This structured demography in VDMs allows for simulation of canopy gap formation, competitive exclusion, and coexistence of vegetation, thus producing variability in forest stand age and composition (Fisher et al., 2010;Moorcroft et al., 2001;Longo et al., 2019). VDMs are designed for vegetation to dynamically respond to variation in traits (Fyllas et al., 2014) leading to differences in plant mortality, growth, and recruitment rates (Shugart and West, 1980). These attributes influence the ecosystem fluxes of carbon, energy, and water (Bonan, 2008). Despite their potential for simulating regrowth processes, there has been limited VDM testing of regrowth following tropical forest disturbances. Importantly, projections of future climate using Earth system models (ESMs) are strongly influenced by the terrestrial carbon cycle in the tropics (Arora et al., 2013;Friedlingstein et al., 2014), which is strongly regulated by disturbance and regrowth Trumbore et al., 2015;Magnabosco Marra et al., 2018).
Observational studies have shown that Amazon forests follow a range of successional regrowth pathways after clearcutting and burning (Mesquita et al., , 2015. Thus, the type of disturbance and the pre-disturbance ecosystem state are important determinants of the successional pathways of forest regrowth. Nonetheless, this information is difficult to obtain in remote forests of the Amazon. In this study, we addressed this issue in the context of windthrow, clear-cut, and clear-cut and burned (cut + burn) disturbances to analyze (i) the sensitivity of Landsat to detect and distinguish these relevant disturbances and their pathways of forest regrowth and (ii) the time span of forest regrowth. This understanding of forest regrowth was used to (iii) test the modeled forest regrowth of the Functionally Assembled Terrestrial Ecosystem Simulator (FATES) model (Fisher et al., 2015) implemented in the Energy Exascale Earth System Model (E3SM) Land Model (ELM) Zhu et al., 2019), ELM-FATES. This study provides insights into the use of remote sensing to identify drivers of forest disturbance in the Amazon, and a better understanding of the pathways of forest regrowth provides insights into the resilience of these forests to repeated disturbances and can help improve land models. For the spectral characteristics before and after disturbances we used cells of 3 × 3 pixels (blue squares) over disturbed and undisturbed areas. For the pathways of forest regrowth after clear-cutting and burning sites, we analyzed areas with different distances from the disturbance edge (A 1 , A 2 , and A 3 in yellow). The background image in (a) is from © Google Earth Pro. The background images in (b), (c), and (d) are from Landsat 5 on 12 July 1987, 1 June 1984, and 12 July 1987 and were composed in RGB color using bands 5, band 4, and band 3, respectively.
2 Study area and methods

Study area and sites
Forests in the Central Amazon affected by windthrow, clearcut, and cut + burn were addressed in this study. Windthrows (Mitchell, 2013) in the Amazon are caused by strong descending winds that uproot or break trees (Garstang et al., 1998;Negrón-Juárez et al., 2018;. In clear-cut areas, forests are cut and cleared, and in cut + burn areas forests are cleared and burned (Mesquita et al., , 2015Lovejoy and Bierregaard, 1990). The windthrow, clearcut, and cut + burn sites used in this study were selected based on the following conditions: (a) prior to disturbance they were upland (no flooding) old-growth forest and located in the same region, with similar climatic, edaphic, and floristic composition; (b) long-term records of satellite imagery and corresponding field data before and after disturbance are available; and (c) no subsequent disturbance has occurred.
The three windthrow, clear-cut, and cut + burn sites analyzed in this study are located near the city of Manaus, Central Amazon (Fig. 1a). The windthrow (centered at 3 • S, 60.75 • W; Fig. 1b) was located near the village of Tumbira, about 80 km southwest of Manaus, occurred in 1987, and covered an area of ∼ 75 ha. At this site, data on forest regrowth including forest structure and species composition for trees ≥ 10 cm DBH (diameter at breast height of 1.3 m) have been collected since 2011 covering disturbed and undisturbed areas and found that the genus Cecropia is one of the dominant species in the most disturbed areas (Magnabosco . The clear-cut and cut + burn sites were experimentally created within the Biological Dynamics of Forest Fragments Project (BDFFP), which encompasses an area of ∼ 1000 km 2 (centered at 2.5 • S, 60 • W) located 80 km north of the city of Manaus, Brazil. The BDFFP was established and managed in the early 1980s by Brazil's National Institute of Amazonian Research (INPA) and the Smithsonian Institution and is the longest-running experiment of forest fragmentation in the tropics (Bierregaard et al., 1992;Lovejoy et al., 1986;Laurance et al., 2011Laurance et al., , 2018Tollefson, 2013). Further details of the BDFFP are in Bierregaard et al. (2001). The selected clear-cut site is located on the Porto Alegre farm (centered at 2.35 • S, 59.94 • W; Fig. 1c). This site was clear-cut in 1982 without subsequent use and was dominated by the pioneer tree genus Cecropia 6-10 years after abandonment (Mesquita et al., 1999. The cut + burn site is located on the Dimona farm (centered at 2.33 • S, 60.11 • W; Fig. 1d), which was clearcut and burned in September 1984, maintained as pasture for 2-3 years, and then abandoned. By 1993 this site was 6 years old and dominated by the pioneer tree genus Vismia (Mesquita et al., 1999. In the Manaus region the mean annual temperature is 27 • C (with higher temperatures from August to November and a peak in October) and the mean annual rainfall is 2365 mm with the dry season (rainfall < 100 mm per month; Sombroek, 2001) from July to September . The topography is relatively flat with landforms ranging from 50-105 m a.s.l. (Laurance et al., 2007(Laurance et al., , 2011Renno et al., 2008), and the mean canopy height is ∼ 30 m, with emergent trees reaching 55 m (Laurance et al., 2011;Lima et al., 2007;Da Silva, 2007). The soil in this region comprises Ferrosols (Quesada et al., 2011;Bierregaard et al., 2001;Ferraz et al., 1998) following the Food and Agriculture Organization (FAO) classification and with similar floristic composition (Bierregaard et al., 2001;Carneiro et al., 2005;Vieira et al., 2004;Higuchi et al., 2004). In the BDFFP and for old-growth forest trees with DBH ≥ 10 cm, there are 261 ± 18 species ha −1 , the stem density is 608 ± 52 stems ha −1 , and the basal area is 28 m 2 ha −1 (Laurance et al., 2010), which is representative of the region (da Silva et al., 2002;Vieira et al., 2004;Carneiro et al., 2005;Magnabosco Marra et al., 2014. In this region 93 % of stems are between 10 and 40 cm in DBH (Higuchi et al., Landsat NIR band and ELM-FATES sensitivity 2012) and the annual tree mortality is 8.7 trees ha −1 for trees ≥ 10 cm in DBH (Higuchi et al., 1997).

Landsat satellite data and procedures
The Landsat Ecosystem Disturbance Adaptive Processing System (LEDAPS) (Schmidt et al., 2013;Masek et al., 2006Masek et al., , 2008Masek et al., , 2013 surface reflectance (SR) from the Landsat 5 Thematic Mapper (TM) was used in this study to characterize the type of disturbance and the subsequent pathways of forest regrowth over our study areas. LEDAPS was developed to ensure that spectral changes in Landsat are associated with regrowth dynamics Schmidt et al., 2013) and to facilitate robust studies of land surface changes at different temporal and spatial scales in tropical forests (Kim et al., 2014;Valencia et al., 2016;Alonzo et al., 2016). LEDAPS SR Landsat 5 TM (L5 hereinafter) is generated by the United States Geological Survey (USGS) using the Second Simulation of the Satellite Signal in the Solar Spectrum (6S) that corrects for the influences of, among others, water vapor, ozone, aerosol optical thickness, and digital elevation on spectral bands (USGS, 2017;Vermote et al., 1997). L5 bands are derived using per-pixel solar illumination angles and generated at a 30 m spatial resolution on a Universal Transverse Mercator (UTM) mapping grid (USGS, 2017). LEDAPS in the Landsat 7 Enhanced Thematic Mapper Plus (ETM+) sensor (L7 hereinafter) was also used to corroborate our predictions (described below). Though L5 and L7 use the same wavelength bands, they are different sensors and differences in surface reflectance may exist, especially in tropical forests, due to large atmospheric effects (Claverie et al., 2015). Landsat 8 was not used since comparison between Landsat 8 and both L5 and L7 is not straightforward due to differences in the spectral bandwidth of these sensors. We used LEDAPS since a long time series of data is available with high spectral performance (Claverie et al., 2015) and it is suitable for ecological studies in the Amazon (van Doninck and Tuomisto, 2018;Valencia et al., 2016). L5 and L7 are available in Google Earth Engine (Gorelick et al., 2017), which we used to retrieve and analyze these data.
The L5 and L7 spectral bands used in this study were blue (0.45-0.52 µm), green (0.52-0.62 µm), red (0.63-0.69 µm), near infrared (NIR) (0.76-0.90 µm), shortwave infrared 1 (SWIR1) (1.55-1.75 µm), and shortwave infrared 2 (SWIR2) (2.08-2.35 µm). L5 and L7 measurements provide the fraction of energy reflected by the surface and range from 0 (0 %) to 10 000 (100 %). Only scenes from June, July, and August were used since these dry-season months present less cloud cover over our study area . This procedure also reduces effects associated with illumination or phenology since images correspond to the same period each year. Only images that were cloud free, cloud shadow free, and haze free over our disturbed areas were used to eliminate errors associated with these elements. For this procedure, visual inspection of visible bands and quality infor-mation from L5 and L7 were used. No further corrections were applied due to the robustness of L5 imagery over the Amazon (Valencia et al., 2016). All of the disturbances are in the Landsat scene Path 231, Row 062. The dates of L5 images used (Landsat 5 operational imaging ended in 2011) were 1 June 1984 (except for the windthrow), 6 July 1985, 12 July 1987, 2 August 1989, 20 July 1990, 8 August 1991, 31 July 1994, 21 June 1997, 26 July 1998, 13 July 1999, 24 July 2003, 4 August 2007, 6 August 2008  The spectral characteristics of old-growth forests and their changes after disturbances were investigated using 19 cells of 3 × 3 pixels (Fig. 1b, c, d). The average of each cell was used in our analysis. Spectral characteristics for oldgrowth forests for each site were determined from cells located in the same position as the disturbance but previous to the disturbance and/or from adjacent areas. Five cells of old-growth forests were located from 1 to 2 km away from the windthrow site. Though closer distances may also represent old-growth forests, we were conservative since Landsat is not sensitive to clusters of downed trees comprising fewer than eight trees (Negrón-Juárez et al., 2011). For the clear-cut and cut + burn sites the spectral characteristics of their respective old-growth forests (control) were studied from three cells per site located 500 to 800 m away from the edge of the disturbance to minimize edge effects that are relevant in the first 100 m (Lovejoy et al., 1986;Laurance et al., 2007;Mesquita et al., 1999). The spectral characteristics for the windthrow were acquired from two cells containing the highest level of SWIR1 values in 1987 which is associated with the maximum disturbance (Negrón-Juárez et al., 2011;Magnabosco Marra et al., 2018;. For the clear-cut site three cells were located 400-500 m from the edge, and for the cut + burn site, three other cells were located at a distance of 100-300 m with respect to the edge. For the clear-cut site we also selected four areas: A 1 , A 2 , A 3 , and A T (A T = A 1 +A 2 +A 3 ), shown in Fig. 1c. For the cut + burn site, we selected three areas: A 1 , A 2 , and A T (A T = A 1 +A 2 ), shown in Fig. 1d.
L5 data for the windthrow, clear-cut, and cut + burn sites encompass a period of 28 years with 13 years of missing data due to cloud cover or lack of image. In order to assess the forest regrowth to spectral levels similar to old-growth forests (control), we applied a gap-filling method (Gerber, 2018) for time series to obtain estimates for missing years using the R package "zoo" (Zeileis et al., 2018). The gap-filled datasets were analyzed using the smoothing-spline technique (R package, 2017). To determine whether L5 bands were sensitive to regrowth, we analyzed changes in the slope (β) of the bands across our chronosequence. A t test on the slope coefficient was used to test the null hypothesis that β is zero (H 0 , β = 0) against the alternative hypothesis (H 1 , β = 0) at a 5 % significance level (α = 0.05). If the computed test statistic (t stat) was inside the critical values, then H 0 was not rejected. The critical values (±t 1−α/2 , n − 2; n is the number of data points) were obtained from statistical tables (Neter et al., 1988). Forests in the Manaus region affected by windthrows were dominated by tree species from the genera Cecropia and Pourouma in about 3-5 years (Magnabosco Nelson and Amaral, 1994), and the clearcut and cut + burn sites were dominated by Cecropia and Vismia about 6 years after the disturbances (Mesquita et al., 1999. The slopes of the time series were determined after these periods, i.e., 1991, 1987, and 1990 for windthrow, clear-cut, and cut + burn sites, respectively. A comparison of successional pathways of forest regrowth among studied disturbances was conducted which was feasible due to the similar conditions of climate, soils, and structure and composition of the old-growth forests. Time series of L5 bands were analyzed using the statistical nonparametric function (univariate fit), with the smoothing spline and the Gaussian regression ANOVA (analysis of variance) model. Calculations were conducted on the R 3.5.2 software platform (R package, 2017) using the package "gss" (General Smoothing Splines) (Gu, 2018). We calculated the smooth spline (using the cubic fit algorithm) of observed data and the associated standard errors, from which we calculated Bayesian 95 % confidence intervals. Predictions of the time after disturbance needed to reach old-growth forests values are based on these data using the function "ssanova" (Fitting Smoothing Spline ANOVA Models) of the R package gss (General Smoothing Splines), version 2.1-9. The predictions were compared with published field observations (Sect. 2.1) where data were available, and L7 images were used to assess the reliability of our predictions.

Forest regrowth simulation in ELM-FATES
Time series of L5 bands sensitive to disturbances and the pathways of forest regrowth were compared with modeling results from ELM-FATES (Fisher et al., 2010(Fisher et al., , 2015Holm et al., 2020). The underlying model structure and concepts in FATES are based on the Ecosystem Demography (ED) concept (Moorcroft et al., 2001) and are described in detail at https://github.com/NGEET/fates, last access: 1 July 2019. A major development is the modularization of the model structure in FATES so that boundary conditions and vegetation can be coupled with ESM land models. FATES is integrated into the E3SM Land Model (ELM) Zhu et al., 2019) and the Community Land Model (CLM) Lawrence et al., 2019) coupled to the Community Earth System Model (Hurrell et al., 2013). In this study we used ELM-FATES. ELM-FATES simulates vegetation that varies in terms of successional age and size; plant competition; and dynamic rates of plant mortality, growth, and recruitment, all on landscapes partitioned by areas of disturbance. The main updates and modifications in ELM-FATES compared with ED include changes to carbon alloca-tion and allometry and the introduction of the perfect plasticity approximation (PPA) (Purves et al., 2008;Fisher et al., 2010) used for the accounting of crown spatial arrangements throughout the canopy and organizing cohorts into discrete canopy layers. Photosynthesis and gas exchange physiology in ELM-FATES follows the physics within the CLM v4.5 (Bonan et al., 2011) and, unlike ED, uses the original Arrhenius equation from Farquhar et al. (1980). ELM-FATES tropical forest simulations conducted here were based on parameter and demography sensitivity analysis at a site 40 km from the BDFFP (Holm et al., 2020), at the ZF2 research station (Magnabosco Marra et al., 2014). Holm et al. (2017) found that with the improved parameterization, ELM-FATES closely matched observed values of basal area, leaf area index (LAI), and mortality rates but underestimated stem density for a Central Amazon old-growth forest near the BDFFP.
Model simulations were driven by climate-forcing data derived from measurements collected between the years 2000 and 2008 at the K34 flux tower located at 2.6 • S, 60.2 • W ( de Araujo et al., 2002) about 40 km from the BDFPP, at the ZF2 research station. ELM-FATES (using the git commit "4a5d626" and the version corresponding to tag "sci.1.0.0_api.1.0.0") was run and spun up for 400 years until a stable biomass equilibrium was reached within the modeled forest. We then simulated a one-time logging treatment of a near-complete mortality of all trees (98 % "clear-cut") with a remaining 2 % consisting of only small trees > 5 cm DBH to aid in recruitment (Mesquita et al., 2015). Second, we simulated a one-time windthrow disturbance that killed 70 % of trees ("windthrow") as was reported in a recent observational study on windthrows in the same region (Magnabosco . All dead trees were "removed" and therefore did not enter modeled soil pools. The fire module in ELM-FATES is currently under final development and testing, and therefore burned simulations are not included in this study. The old-growth forest simulated by ELM-FATES, used as a pre-disturbance metric, was based on previously validated tropical parameterization and sensitivity testing in the same region (Holm et al., 2020); see supplementary material of Fisher et al. (2015) for a description of plant-functionaltype specific carbon allocation and allometry schemes and updates from the ED model framework. Simulations of disturbance and subsequent vegetation regrowth were initiated from this old-growth forest state. The model design used here only allows for simulating intact forests with natural disturbances (e.g., gap dynamics or windthrows) or harvested forests and not both at the same time or in adjacent patches. Accounting for distance to intact forests was excluded due to the current limited understanding of seed dispersal mechanisms (i.e., spatial variability, dispersal limitation, etc.) in tropical forests . We use a more general form of seed production, such that the individual cohorts in ELM-FATES use a targeted fraction of net primary production (NPP) during the carbon allocation process (after accounting for tissue turnover and storage demands), which Figure 2. Imposed trait variation used in the parameterization of ELM-FATES tropical evergreen plant functional types (PFTs) for the 20-ensemble simulations and the resulting average growth rate average over the 50-year simulation period. Each simulation consisted of a single PFT varying by all four traits at once: wood density, V cmax , the canopy area allometric coefficient, and leafclumping index for leaf self-occlusion. Dark green points represent fast-growing evergreen pioneer PFTs, while light green points represent slow-growing late successional PFTs. adds to the site-level seed pool for recruitment of new cohorts. Field data were not used to simulate or calibrate the modeled forest regrowth post-disturbance.
To account for uncertainty in the representation of plant physiology within tropical evergreen forests, we analyzed an ensemble of 20 simulations varying in targeted plant functional traits. We prescribed each ensemble with a single tropical evergreen plant functional type (PFT) that varied in wood density (0.44 to 1.06 g cm −3 ) and maximum rate of carboxylation (V cmax ; 42 to 55 µmol m 2 s −1 ) (Table 1), via random sampling. To evaluate changes in canopy coverage of the forest stand, each PFT additionally varied by an allometric coefficient (1.35 to 1.65), determining the crown-area-todiameter ratio, and a leaf-clumping index (0.59 to 1.0 out of a 0-1 fraction) that determines how much leaf self-occlusion occurs and decreases light interception, and the direct and diffuse extinction coefficients in the canopy radiation calculations. The default values for these parameters are based on, or derived from, references given in Table 1. Each ensemble member represents a single PFT across the spectrum of fast-growing "pioneer" PFTs and slow-growing "late successional" PFTs to provide a reasonable spread across the trait uncertainty when assessing regrowth from disturbance. We characterized pioneer plants in our simulations as having low wood density (Baker et al., 2004) and high V cmax based on the inverse relationship between these two plant traits, as well as a low crown area coefficient and low leaf-clumping factor; i.e., a monolayer planophile distribution (Lucas et al., 2002). These correlated relationships were applied to the ensemble-selected traits (Fig. 2). The opposite relationship was applied for slow-growing, late successional PFTs (i.e., high wood density, low V cmax , high crown area coefficient, and high leaf-clumping factor).
In order to evaluate ELM-FATES performance during forest regrowth we compared the NIR band, the most sensitive band to regrowth (see results), with ELM-FATES outputs of aboveground biomass (AGB, Mg ha −1 ), total stem density of trees ≥ 10 cm DBH (stems ha −1 ), leaf area index (LAI; onesided green leaf area per unit ground surface area; m 2 m −2 ), and total live crown area (m 2 m −2 ) since these variables directly influence the surface reflectance (Ganguly et al., 2012;Lu, 2005;Masek et al., 2006;Powell et al., 2010;Ruiz et al., 2005). We suggest that testing an array of modeled forest variables (e.g., biomass structure, density coverage of vegetation, and proportion of the tree crown that has live foliage) provides a robust comparison to NIR values, due to multiple forests characteristics contributing to and affecting NIR reflectance (Ollinger, 2011), and reduces model unknowns and biases that can arise when using only one model variable. The usage of different stand structure and canopy processes can be helpful when evaluating ELM-FATES during different phases of forest regrowth. In addition, we averaged modeled outputs of the crown area, stem density, and LAI, since each of these variables influences the reflectance of forests, and defined this average as the modeled "canopy coverage". Measurements of forest canopy cover have been used to analyze plant growth and survival, and it is an important ecological parameter related to many vegetation patterns (Ganey and Block, 1994;Jennings et al., 1999;Paletto and Tosi, 2009). Modeled diameter growth rates (cm yr −1 ) for trees with DBH ≥ 10 cm are also shown to provide information on the successional dynamics within ELM-FATES.

L5 bands and disturbances
All L5 bands showed an increase in surface reflectance immediately after windthrow, clear-cut, and cut + burn sites except the NIR values which decreased (with higher decrease after burning) (Fig. 3a, b, and c). This decrease in NIR values was due to exposed woody material and dry leaves, typical after windthrow (Negrón-Juárez et al., 2010a;Negrón-Juárez et al., 2011) andclear-cutting (Sect. 7 in Adams andGillespie, 2006) or the dark surface following burning (Pereira et al., 1997). For windthrows, these effects last about 1 year, after which vegetation regrowth covers the ground surface (Negrón-Juárez et al., 2010a, 2011. The spectral characteristics of old growth and disturbances are shown in Fig. 3d-f with the error bands representing the standard deviation of all pixels from respective cells. About 1 year after the disturbance the bands that experienced increases in surface reflectance showed a decrease in surface reflectance (the opposite to the NIR band) due to the increases in vegetation cover. A similar response is expected for the clear-cut that occurred in 1982 and therefore before the beginning of our available data (L5 imagery is available from 1984; Fig. 3e). The simi- Table 1. The range (minimum to maximum) of four key model input parameters used in the 20-ensemble ELM-FATES simulations for both windthrow and clear-cut simulations, to account for uncertainty in the representation of plant traits, along with the default value used in the ELM-FATES model. Wood density value from Moorcroft et al. (2001); V cmax based on Oleson et al. (2013) and Walker et al. (2014); crown area : DBH derived from Farrior et al. (2016) and adjusted based on site-specific sensitivity tests; and the leaf-clumping index based on radiation transfer theory of Norman (1979

Pathways of forest regrowth
About 6 years after the disturbance, NIR values reached a maximum and then decreased slowly with time showing a significant negative trend (Table 2). SWIR1 values also showed a significant negative trend with time but only for the clear-cut site (Table 2). In general, the green, blue, red, SWIR1, and SWIR2 band values returned to pre-disturbance values (control) about 6 years after the disturbance (Fig. 3d, e, f and Table 2). Therefore, we used the NIR band (which remained higher than pre-disturbance values throughout the time series and is potentially sensitive to ecosystem properties of regrowing forest) to investigate the regrowth dynamics in comparison to our control forests. We used the relationships presented in Figs. 4, 5, and 6 to determine the time taken for NIR values from the disturbance sites to become similar to control NIR values. The average control NIR value was 28 ± 1 %. For the windthrow site the NIR values became similar to control levels after about 39 years (range 32 to 57 years). For the clear-cut and cut + burn sites, this period was estimated to be 36 years (range 31 to 42 years) and 56 years (range 42 to 93 years), respectively. From Figs. 4-6 it is evident that the type of disturbance has a clear effect on the pathways of NIR recovery. L7 data, in general, are within the 95 % CI of predictions.
During the first 12 years following the windthrow, the spline curve fitted to the NIR data decreased by ∼ 0.13 % yr −1 after which the rate of decrease doubled (0.26 % yr −1 ; Fig. 4). For clear-cutting, NIR values decreased faster, i.e., ∼ 0.4 % yr −1 . The decrease in NIR values for the clear-cut site appears to be independent of the distance from the edge of the disturbance since the changes in NIR values of all selected areas (A 1 , A 2 , A 3 , and A T ) are similar (Fig. 5). For the cut + burn site, the rate of change in NIR values to reach values similar to the control forests was the slowest among all disturbances considered (∼ 0.15 % yr −1 ) (Fig. 6). The cut + burn site showed differences with respect to the border of the disturbance (areas A 1 , A 2 , and A T ), which may be related to the spatial heterogeneity of burnings and forest responses.

FATES model and regrowth from forest disturbance
To address our goal of improving the connection between remote sensing, model benchmarking, and the fidelity of future predictions of forest regrowth processes, we examine the representation of such processes within ELM-FATES. The average of the ELM-FATES 20-member ensemble predicted a continuous, and almost linear, regrowth of biomass ( Fig. 7a) after clear-cut and windthrows. The modeled recovering biomass returned to modeled old-growth forest values more quickly for windthrows (37 years, range 21 to 83 years) compared to clear-cuts (42 years, range 27 to 80 years). However, the annual rate of change in biomass regrowth over 50 years was faster in the clear-cut simulation (2.5 Mg ha −1 yr −1 ) than in the windthrow simulations (2.0 Mg ha −1 yr −1 ), which was due to the clear-cut site recovering from initial biomass of near zero and a proportionally greater contribution of fast-growing pioneer species. The model simulation of stem density, LAI, and crown area are shown in Fig. 7b-d, respectively. For stem density, ELM-FATES (black line) predicted an average of up to 8 years before any new stems reached ≥10 cm DBH (a standdeveloping period; Fig. 7b) for the clear-cut. The simulated stem density for old-growth forests ( Fig. 7b; green line) was ∼ 200 stems ha −1 (≥ 10 cm at 1.3 m), ∼ 408 trees lower than observations. The model ensembles with typical early successional traits predicted a forest with many fast-growing, small-diameter stems < 10 cm, with maximum early successional stem densities reaching 1560 and 1414 stems ha −1 for clear-cut and windthrows, respectively. Once the canopy closes and self-thinning dominates (average of 15 years after disturbance), there are declines in stem density as trees gain biomass, and canopy closure forces some trees into the understory, where they die at faster rates due to shad-  ing. The modeled forests returned to old-growth stem density conditions 39 and 41 years after windthrows and clearcut, respectively (Table 3). Though the two disturbance types have very similar times of return to pre-disturbance conditions, they differ in the speed of recovery. ELM-FATES predicts faster diameter growth increments (1.3 cm yr −1 ) and canopy closure for a forest composed of all pioneer type PFTs and slower (0.5 cm yr −1 ) diameter growth and more open canopies for the late successional forest type (Fig. 8).
Diameter growth is an emergent model feature of dynamic plant competition for light and stand structure and is consistent with observational studies of secondary forests growing through succession (Brown and Lugo, 1990;Winter and Lovelock, 1999;Chapin et al., 2003).
The LAI of the modeled old-growth forest (4.0 m −2 m −2 ), prior to disturbances, was close to the observed LAI (4.7 m −2 m −2 ) measured near our study sites . Due to disturbance, the initial modeled LAI (Fig. 7c) and total crown area (Fig. 7d) decreased, as expected. During regrowth from disturbance both the LAI and total crown area rapidly recovered, and the LAI even surpassed pre-disturbance values. This pattern resembles the initial NIR spike due to fast-growing PFTs. These two canopy coverage attributes reached maximum values after 3 to 6 years, depending on the disturbance and response of forest attributes (Table 3). To evaluate model results against remote sensing observations, we compared the initial period after the disturbance of the spikes in NIR values to the canopycoverage metric (combination of LAI, stem density, total  . Changes in NIR values after clear-cuts for areas A 1 , A 2 , A 3 , and A T = A 1 + A 2 + A 3 (shown in Fig. 1c) and prediction of NIR values reaching pre-disturbance values. The plots show the data (circles), the fit (solid line), and the 95 % confidence interval (CI, dashed lines). Grey bar represents the control (old-growth forests) NIR value of 28 ± 1 %, and the horizontal dashed black line is 28 %. crown area) over the same modeled period. ELM-FATES predicted that after a windthrow the forest took 5.7 years to reach maximum values of canopy coverage, which was sooner than the clear-cut simulation (7 years). While the modeled time span for this initial period was similar to that inferred from NIR values, there was disagreement between which disturbance recovery occurred fastest (windthrow in ELM-FATES vs. clear-cut in the NIR band), similar to disagreement in the recovery of AGB (Table 3).
ELM-FATES provided a prediction of the values in each forest variable when the stand reached its production limit and full canopy closure, at which point there was a shift to a declining trend and decreases in forest attributes that outpaced any gains (Table 3). At the maximum peak recovery Figure 6. Changes in NIR values for cut + burn site in areas A 1 , A 2 , and A T = A 1 +A 2 (shown in Fig. 1) and prediction of NIR values reaching pre-disturbance values. The linear fit (solid line) and the 95 % CI (dashed line) are shown. Grey bar represents the control (old-growth forests) NIR value of 28 ± 1 %, and the horizontal dashed black line is 28 %. and carrying capacity limit, the highest forest values occurred in the clear-cut simulation, matching the higher NIR values from clear-cut a few years after the disturbance (Figs. 4-6). Over the longer self-thinning period the modeled LAI decreased and returned to modeled old-growth values 26 years after clear-cut and gradually over 53 years for windthrows (Table 3). LAI was the only variable that had a noticeable faster recovery in the clear-cut simulations. After both disturbances the total crown area permanently remained high (0.99 m 2 m −2 ) and slightly higher than the crown area of the simulated old-growth forests (0.98 m 2 m −2 ), suggesting that disturbances can generate a denser canopy, as discussed below.

Discussion
Our results show that Landsat reflectance observations were sensitive to the initial changes in vegetation following windthrows, clear-cut, and cut + burn, three common disturbances in the Amazon. Specifically, a decrease in NIR values and an increase in SWIR1 values were the predominant spectral changes immediately (within a few years) following disturbances. The increase in SWIR1 values was different among the disturbances with the maximum increase observed in the cut + burn, followed by clear-cut and then the windthrow site. The highest increase in SWIR1 values in cut + burn sites may be related to the high thermal emission of burned vegetation (Riebeek, 2014). Likewise, the relatively higher moisture content of woody material in the windthrow site decreases the reflection of the SWIR2 band. On the other hand, in our control (old-growth) forests, we observed typically high NIR reflectance due to the cellular Figure 7. Simulated regrowth of the Central Amazon forest after a clear-cut (98 % tree mortality; black line) and windthrow (70 % tree mortality; orange line) event using 20 simulations of the demographic model ELM-FATES, compared to the modeled old-growth values prior to disturbance (dashed green line) and against field data (solid green line) from close sites, except for crown area data that were taken from lidar data in Acre, Brazil (Figueiredo et al., 2016). The shaded grey and orange areas represent the spread across the ensembles, showing minimum and maximum values of each forest attribute over its regrowth. (a) Regrowth of aboveground biomass (AGB; Mg C ha −1 ). (b) Regrowth of stem density (stems ha −1 ) of stems > 10 cm DBH and years of returning to modeled old-growth values. (c) Regrowth of leaf area index (m −2 m −2 ) and (d) regrowth of total crown area. Table 3. Summary of different time spans of regrowth (years) to old-growth forest status after two disturbance types: windthrows and clearcuts from ELM-FATES model results and remote sensing. Additionally, the time (years) it takes forest attributes to reach maximum values during regrowth and the corresponding value at this maximum peak. AGB (Mg C ha −1 ), stem density (stems ha −1 ), LAI, and crown area (m 2 m −2 ) refer to simulation results, as compared against NIR remote sensing. The average of AGB and stem density is characterized as modeled "forest structure". The averages of crown area, stem density, and LAI are characterized as modeled canopy coverage in this study and are additionally compared against NIR values.  Adams and Gillespie, 2006), absorption of red radiation by chlorophyll (Tucker, 1979), and absorption of the SWIR1 band by the water content in leaves (Sect. 7 in Adams and Gillespie, 2006). While the SWIR1 band is frequently used to identify exposed woody biomass immediately after disturbances (Negrón-Juárez et al., 2010b, we found that the NIR band was more sensitive to the successional pathways of regrowth for all the disturbances considered. The NIR band has also been associated with succession (Lu and Batistella, 2005) and regrowth (Roberts et al., 1998) in naturally and anthropogenically disturbed tropical forests (Laurance, 2002;Chazdon, 2014;Magnabosco Marra, 2016;Laurance et al., 2011). Previous studies have shown that changes in NIR values are related to leaf structure and surface characteristics (Roberts et al., 1998;Xiao et al., 2014) with the youngest leaves having higher NIR values with respect to fully formed and older leaves (Roberts et al., 1998). Maximum values of the NIR band were observed about 6 years after clear cut, which is the time pioneers form a closed canopy (Mesquita et al., 1999(Mesquita et al., , 2015, and characterized by a relatively uniform distribution of tree diameter and heights (Vieira et al., 2003). This maximum in the NIR band was higher in the clear-cut site dominated by species from the genera Cecropia and Pourouma (Mesquita et al., 2015;Massoca et al., 2012) than the site affected by cut + burn dominated by Vismia species (Mesquita et al., 2015;Laurance et al., 2018). The higher NIR values in Cecropia and Pourouma are due to the genera's monolayer planophile distribution of large leaves that produces high reflectance compared to Vismia, species of which have rougher and denser canopies that trap more NIR reflectance (Lucas et al., 2002). The high values of the NIR band might be related to the low leaf wax (Chavana-Bryant et al., 2017) from new trees and/or scattering related to leaf and canopy water (Asner, 2008). NIR values decrease with the dynamics of succession due to an increase in the canopy roughness (Hallik et al., 2019).
After the establishment of pioneers, the NIR values decrease with time but with different rates depending on the type of disturbance. In windthrown areas, tree mortality and subsequent recruitment may continue for several decades, promoting changes in functional composition and canopy architecture . Cecropia and Pourouma trees grow relatively quickly, and after closing the canopy they limit light penetration due to their large leaves creating a dark, cooler, and wetter understory Jakovac et al., 2014). As a result, light levels in the understory decline faster with time and thus allow the recruitment and establishment of shade-tolerant species. The cohort of Cecropia and Pourouma species has a relatively short lifespan, and, ∼ 20 years after disturbance, secondary and old-growth forest species start to be established (Mesquita et al., 2015). With the self-thinning of Cecropia and Pourouma, the growing understory traps more light, and consequently albedo decreases (Roberts et al., 2004). This pattern is con- Figure 8. Change in predicted diameter increment growth rate (cm yr −1 ) for one simulation, from the 20 ensembles, that represents a fast-growing pioneer forest stand and a slow-growing late successional forest stand from a clear-cut disturbance (black and grey) and a windthrow disturbance (orange). Variations in the diameter increment are a result of differences in the following traits: wood density, V cmax , crown area, and a leaf-clumping index. sistent with the decline in NIR values and observed changes in canopy architecture (Mesquita et al., 2015), photosynthesis, and the LAI (Saldarriaga and Luxmoore, 1991). In contrast, the architecture of Vismia species that dominate cut + burn areas allows higher light levels in the understory and subsequent recruitment of Vismia or other genera with similar light requirements. As a consequence, species turnover and structural changes are slower than in clear-cut areas (as found by Jakovac et al., 2014, in a study conducted a few kilometers from the BDFFP) and windthrows, consistent with changes in NIR values. In the course of succession, Vismia tends to be replaced by Bellucia, which is a genus with a similar leaf and canopy structure to Vismia (Mesquita et al., 2015). This pattern favors the penetration of light through the canopy  for several decades before a more shaded understory allows the germination and establishment of old-growth species .
For the windthrow we estimate that the NIR values should become similar to pre-disturbance conditions in about 39 years. This value agrees with the 40 years of biomass regrowth found using ground-based data in the Central Amazon . This result also corroborates previous studies that the NIR band operates in the best spectral region to distinguish vegetation biomass (Tucker, 1979(Tucker, , 1980 and photosynthesis (Badgley et al., 2017). For clear-cutting and cut + burn, the regrowth time was about 36 and 56 years, respectively, but no ground-based estimates were available for comparison. Still, NIR values showed that the pathways of regrowth from clear-cut and cut + burn are divergent with time (Figs. 5 and 6), which is consistent with observational studies (Mesquita et al., 2015). . Total stem density (stems ha −1 ) separated into six diameter (cm) size classes from Central Amazon field data located at the nearby ZF2 site (green bars) averaged from 1996-2011, and predicted by ELM-FATES (grey bars). The percentages represent the proportion of stems in each size class relative to the total stem density.
In general, we found that NIR values may be used as a proxy in modeling studies aimed at addressing forest regrowth after disturbances. Though NIR values are useful to distinguish successional stages up to decades after the disturbance, they may not represent all the successional processes. As soon as the forest canopy becomes structurally similar to that of the mature forest, the NIR band will no longer be sensitive to changes in vegetation attributes (Lucas et al., 2002). Though L5 NIR data may be complemented with current Landsat measurements (the L5 NIR band has a comparable performance to the Landsat 8 NIR Operational Land Imager (OLI) algorithm; Vermote et al., 2016), it is important to emphasize that our estimates of recovered reflectance and biomass in disturbed areas do not capture the full recovery of diversity in floristic attributes and species composition that can take centuries (Rozendaal et al., 2019). The predominance of Cecropia, after clear-cut, and Vismia, after cut + burn, have also been found in the western (Gorchov et al., 1993;Saldarriaga et al., 1986) and the southern (Rocha et al., 2016) Amazon suggesting that our findings are applicable to other regions. However, an Amazon-wide study is beyond the scope of our work, which is to explore the sensitivity of Landsat to different disturbance types.
Our analysis demonstrates that this version of ELM-FATES has the capacity to reproduce the initial response to disturbance and regrowth patterns after the clear-cut and windthrow that occurred over similar time ranges compared to NIR values. The strongest agreement occurred when ELM-FATES predicted higher peaks of the post-disturbance stem density and LAI in clear-cuts than in windthrows, which can be used for future benchmarking, consistent with the higher peak of NIR values from clear-cuts (Fig. 5 vs. Fig. 4). This effect may be due to ELM-FATES having more ho-mogeneous canopies after clear-cuts as well as more open disturbed area for fast-growing plants, which is also an observed trend. In addition, the average regrowth times to predisturbance values were close between ELM-FATES and NIR results (Table 3), showing that pathways of forest regrowth in ELM-FATES are comparable to observed patterns in tropical forests. ELM-FATES predicted a continuous, and almost linear, regrowth of biomass for the first 50 years of simulation after both clear-cut and windthrows (Fig. 7a), consistent with NIR results and observational studies (Mesquita et al., 2015;Saldarriaga et al., 1988;Jakovac et al., 2014;Magnabosco Marra et al., 2018). In addition, the changes in biomass rates predicted by ELM-FATES were similar to biomass observations recorded after clear-cut (2.3 vs. 2.6 Mg ha −1 yr −1 ) (Mazzei et al., 2010), as well as there being a faster rate of AGB accumulation after clear-cut compared to windthrow, similar to a study reporting higher regrowth rates in more highly disturbed sites (Magnabosco . Landsat showed a faster recovery of NIR values to pre-disturbance conditions in clear-cuts compared to windthrows. Faster growth is characteristic of anthropogenically driven secondary forests that reflect rapid colonization and monodominance of adapted species and genera in changed environmental conditions; e.g., high growth rates, low self-competition, high leaf area index, low herbivory rates Mesquita et al., 2015;Rozendaal and Chazdon, 2015). ELM-FATES predicted a faster recovery of structural AGB and canopy coverage to predisturbance conditions for windthrows (70 % tree mortality) compared to clear-cuts (98 % tree mortality). A major contributing factor to this pattern resulted from larger modeled diameter increments after windthrows (0.92 cm yr −1 ) compared to clear-cuts (0.82 cm yr −1 ) in the first 20 years after the disturbance, setting the trajectory for faster regrowth to pre-disturbance after windthrows. Only the LAI had a faster recovery to pre-disturbance values after clear-cuts, which is expected due to the newly developed forest having a simplified forest structure and the canopy being more homogeneous after a clear-cut (Rosenvald and Lohmus, 2008). ELM-FATES predicted the timing of peak canopy coverage was marginally sooner after windthrows compared to clear-cuts, opposite to the NIR pattern. This discrepancy may be related to more biomass loss and open canopy coverage, followed by a lack of rapid colonization in the modeled clear-cut. Due to the higher stand disturbance that naturally occurs from clearcuts and the diverse complexities in tropical forest composition, we emphasize that the dynamics of different competing PFTs in ELM-FATES requires further investigation. Emerging modeling studies that include plant trait trade-offs, for example, in the leaf and stem economic spectrum or fastgrowth vs. slow-growth strategies, may help to better capture the drivers of forest productivity and demography, enabling improved modeled responses to global change scenarios (Fauset et al., 2019;Sakschewski et al., 2015). Here we test the basic representation of biomass demographics, prior to the more challenging aspects of representing interacting functional diversity in recovering systems (Fisher et al., 2015;Powell et al., 2018).
ELM-FATES predicts the total stem density for a closedcanopy forest to be very low compared to observations (200 simulated vs. 600 stems ha −1 ), but modeled total AGB was close to that reported for the same region (110 simulated vs. ∼ 150 Mg C ha −1 ; Chambers et al., 2013). This discrepancy is due to ELM-FATES predicting a disproportionately high number of large trees ( Fig. 9; with 8 % of stems > 60 cm and 4.5 % of stems > 100 cm), resulting in a crowded canopy, which out-compete smaller understory trees. Higuchi et al. (2012) report 93 % of trees ≥ 10 cm DBH in the study site to be below 40 cm DBH, while ELM-FATES predicted noticeably fewer trees below 40 cm DBH (Fig. 9). Low stem density could be attributed to multiple model assumptions, such as high density-dependent mortality and self-thinning due to the marginal carbon economics of understory trees, low branch-fall turnover, a need for greater limitation of maximum crown area than currently modeled, and increases in mortality rates with tree size (Johnson et al., 2018). Our findings here will guide future ELM-FATES and ecosystem modeling development efforts towards improving the representation of forests comprising dense canopies and how they shift during regrowth.
Land surface models do not typically simulate spectral leaf reflectance, but there is potential to include such output within radiative transfer schemes as is currently done in the CLM. That addition would greatly assist our ability to compare with Earth observation datasets. In lieu of this development, we show that with successional aging, modeled forest structure returns to pre-disturbed values (through canopy closure) with a similar recovery time to that inferred from NIR data, occurring with the process of canopy closure, all which can be compared against remote sensing vegetation indices (see Supplement Fig. S1) and metrics. Which vegetation index (e.g., normalized difference vegetation index - Rouse et al., 1973;enhanced vegetation index -Huete et al., 2002) or metric properly represents the successional pathways following disturbances remains an important area of study.

Conclusions
We tested the sensitivity of Landsat surface reflectance to windthrow, clear-cut, and cut + burn forest in the Central Amazon. The NIR band was more responsive to the successional pathways of forest regrowth years after the disturbance. NIR values showed that pathways of forest regrowth were different among the disturbances, with cut + burn being the most different in terms of spatial heterogeneity and regrowth time to old-growth status, in agreement with observational studies. Our results indicate that after disturbances the NIR values will reach old-growth forest values in about 39 years following windthrows (in agreement with observed biomass regrowth), 36 years for clear-cuts, and 56 years for cut + burn. These results were then compared with simulations of regrowth after windthrows and clear-cut from ELM-FATES. The simulated forest structure and the remote sensing NIR values from the windthrow and clear-cut have similar return times to those of old-growth forest conditions. Future studies applying ELM-FATES should focus on improving stem density predictions, which were underestimated, and on enhancing the capacity to compare with remote sensing observations through representation of canopy-spectralreflectance characteristics.
Code and data availability. The Landsat data used in this study are freely available through the Google Earth Engine platform. ELM-FATES is available at https://github.com/NGEET/fates. Observational data used to compare remote sensing and modeling results have been previously published and references provided (e.g. Mesquita et al., 2015;Jakovac et al., 2014).
Author contributions. RINJ designed the study and performed the remote sensing analysis. JAH performed the model simulations and analysis. RNJ and JAH formulated the research goals and wrote the initial draft of the manuscript. BF performed the statistical analysis of the remote sensing data. DMM provided the field data for windthrows and contributed to the discussion of biomass recovery. RAF and JKS contributed with modeling improvements. ACdA contributed with the observational input data to initialize the model. WJR and JQC were responsible for project funding and administration. All authors contributed to the final writing and editing of the manuscript.
Competing interests. The authors declare that they have no conflict of interest. Financial support. This research has been supported as part of the Next Generation Ecosystem Experiments-Tropics, funded by the US Department of Energy, Office of Science, Office of Biological and Environmental Research under contract no. DE-AC02-05CH11231. This research was also supported through the Reducing Uncertainties in Biogeochemical Interactions through Synthesis and Computation Scientific Focus Area (RUBISCO SFA) under contract to LBNL, which is sponsored by the Regional and Global Climate Modeling (RGCM) Program in the Climate and Environmental Sciences Division (CESD) of the Office of Biological and Environmental Research (BER) in the US Department of Energy Office of Science.
Review statement. This paper was edited by Christopher Still and reviewed by two anonymous referees.