Improving simulated Amazon forest biomass and productivity by including spatial variation in biophysical parameters

. Dynamic vegetation models forced with spatially homogeneous biophysical parameters are capable of produc-ing average productivity and biomass values for the Amazon basin forest biome that are close to the observed estimates, but these models are unable to reproduce observed spatial variability. Recent observational studies have shown substantial regional spatial variability of above-ground productivity and biomass across the Amazon basin, which is be-lieved to be primarily driven by a combination of soil physical and chemical properties. In this study, spatial heterogeneity of vegetation properties is added to the Integrated Biosphere Simulator (IBIS) land surface model, and the simulated productivity and biomass of the Amazon basin are compared to observations from undisturbed forest. The maximum RuBiCo carboxylation capacity ( V cmax ) and the woody biomass residence time ( τ w ) were found to be the most important properties determining the modeled spatial variation of above-ground woody net primary productivity and biomass, respectively. Spatial heterogeneity of these properties may lead to simulated spatial variability of 1.8 times in the woody net primary productivity (NPP w ) and 2.8 times in the woody above-ground biomass (AGB w ). The coefﬁcient of correlation between the modeled and observed woody productivity improved from 0.10 with homogeneous parameters to 0.73 with spatially heterogeneous parameters, while the coefﬁcient of correlation between the and observed woody above-ground biomass improved our with IBIS dynamic key in the biome severely Clearer understanding of the mechanisms that drive of and cmax is to achieve further to properties is essential for useful prediction of the future carbon cycle of Amazon forests. In this work we investigate the importance of representing spatial heterogeneity in vegetation properties in a commonly used DGVM.


Introduction
Tropical forests play an important role in the global carbon cycle, accounting for about one-third of the global net primary productivity and 55 % of total global forest carbon (Pan et al., 2011;Malhi, 2010). The Amazon contains about 50 % of the world's tropical forests (Pan et al., 2011). It is highly diverse in terms of climate, soil physical and chemical properties, and species composition (Davidson et al., 2012;Quesada et al., 2011;Fyllas et al., 2009;Phillips et al., 2004). However, most global vegetation models represent Amazon rainforests through a single set of parameters, which do not vary in space and thus fail to represent its complex spatial diversity. Dynamic Global Vegetation Models (DGVMs) are powerful tools for understanding past and potential future carbon fluxes and stocks. An accurate representation of spatial and temporal variability of the forest biophysical properties is essential for useful prediction of the future carbon cycle of Amazon forests. In this work we investigate the importance of representing spatial heterogeneity in vegetation properties in a commonly used DGVM.

A. D. A. Castanho et al.: Improving simulated Amazon forest biomass and productivity
Field observations from undisturbed old-growth Amazon forest plots have recently quantified the regional variation of many forest attributes such as geological history (Higgins et al., 2011), soil properties , foliar physiological properties , above-ground live biomass (Malhi et al., 2006;Baker et al., 2004b), above-ground wood productivity , net primary productivity Malhi et al., 2009), and residence time of plant components Galbraith et al., 2013). These analyses have shown a general east-to-west gradient of tree structure and dynamics in Amazon forests. Forests in the west tend to have higher above-ground productivity and lower aboveground biomass, while those in the east and central Amazon are slower growing with higher above-ground biomass Malhi et al., 2004;Baker et al., 2004a). The east-west productivity and biomass gradient appears to be related to soil fertility and soil physical properties rather than to climate variations . The soil properties appear to be related to geological history and the exposure of more fertile Miocene-age marine or lacustrine sediments in western Amazonia (Higgins et al., 2011).
Although the spatial variation in biomass and productivity has been measured and described and the relationship to edaphic properties noted, the underlying mechanisms are still not well understood (Malhi, 2012).
With respect to woody net primary productivity (NPP w ) a number of mechanisms have been discussed in the literature. It has been suggested that the variability in wood productivity could be directly related to a shift in the balance of carbon allocation between roots and wood and/or respiration rate . It is also possible that the NPP w spatial variability could be driven by variability in gross primary productivity (GPP) due to a potential limitation of photosynthesis by soil P availability Mercado et al., 2011;Lloyd et al., 2010). Data analyses by Malhi et al. (2004) also suggest that carbon use efficiency (CUE = NPP/GPP) is higher in the west. This suggests that, in addition to a potentially higher GPP in the west there are also potentially higher autotrophic respiration rates in the less fertile soils of central Amazonia, making the CUE less efficient there. The direct effect of climate factors such as temperature, shortwave radiation and precipitation have been argued to be too small to explain the observed variability in the productivity Senna et al., 2009;Quesada et al., 2012;Galbraith et al., 2013).
The spatial variability of above-ground biomass across the Amazon basin has also been discussed in the literature. Field data have indicated the importance of wood density, canopy height and plant allometry in biomass estimates (Baker et al., 2004b;Malhi et al., 2006;Feldpausch et al., 2011). Woody biomass residence time has been suggested to be an important factor for better representation of above-ground biomass in vegetation models (Delbart et al., 2010). Field data show that the fast-growing forests in the western Amazon have lower wood density, lower residence time and therefore lower above-ground biomass than the slow-growing forests in the central and eastern Amazon . The mechanisms that determine the spatial variability in plant residence time are still unresolved (Galbraith et al., 2013). It may be due to external factors such as soil physical properties, disturbance (e.g. landslides and erosion on steep slopes), or climate impacts on tree mortality rates, or intrinsic factors such as high growth rates intensifying light competition and associated tree mortality. Residence time is strongly correlated with soil physical properties, and forest growth rates are strongly related to available soil P and climate. However, the large-scale variation in biomass appears not to be explained by any of the edaphic properties alone .
Dynamic global vegetation models predominately characterize the Amazon forest with tropical broadleaf evergreen trees as the plant functional type (PFT), which is represented by a set of parameters that are invariant in space and time. There is a growing awareness that such an approach is unable to capture spatial variations in key biophysical properties (Senna et al., 2009;Delbart et al., 2010;Fyllas et al., 2012). In this paper we address this issue through data integration and model improvement. We derive the most important parameters for simulating the spatial variability of above-ground woody net primary productivity and biomass. We use the Integrated Biosphere Simulator (IBIS) DGVM with spatially varying observational estimates of key biophysical parameters (woody biomass residence time (τ w ), maximum ribulose-1,5-bisphosphate carboxylase/oxygenase (RuBisCO) carboxylation capacity (V cmax ), and NPP allocation to wood) to simulate AGB w and NPP w and evaluate model performance, in comparison to field data. We create basin-wide raster data sets of the key parameters by extrapolation of the site-specific heterogeneous parameters. Finally, we evaluate the impact of using spatially varying parameters on simulated AGB w and NPP w throughout the Amazon.

Material and methods
The study area is the Amazon region (Fig. 1). The spatial analyses applied are at 1 • × 1 • horizontal resolution. The 1 • × 1 • spatial resolution has been chosen as a compromise between the spatial resolution of the model drivers (e.g. climate and soil properties) and computer run-time. In this section we describe (1) the IBIS dynamic vegetation model used, (2) the field data used in calibration and validation, and (3) the sequence of model simulation exercises.

IBIS 2.6. Integrated Biosphere Simulator
IBIS is a comprehensive model of terrestrial biospheric processes (Foley et al., 1996;Kucharik et al., 2000). The model uses an integrated framework based on land surface biophysics (canopy and soil physics, plant physiology), A. D. A. Castanho et al.: Improving simulated Amazon forest biomass and productivity 2257 Fig. 1. Locations of the main field observation data: woody above-ground biomass (Malhi et al., 2006); woody net primary productivity ; maximum carboxylation capacity of RuBisCO and specific leaf area index ; total soil Pposphorus ; woody carbon allocation ; woody residence time (Galbraith et al., 2013). Shaded areas include the Amazonian sensu stricto (Amazon basin below 700 m a.s.l., light gray) with an estimated area of ∼ 5.65 million km 2 , Amazon River basin (light gray including southeast Planalto, western Andes in dark gray) and tropical forest areas in the north (Guiana, dark gray) (Eva et al., 2005). Each field site that provided data is marked by a hexagon, which is divided into 6 wedges. Each wedge corresponds to a particular variable (see key in lower right). If a particular field variable is available at a site, the wedge corresponding to it is black. For example if woody above-ground biomass (AGB w ) was collected at a particular site, the lower left wedge is black. vegetation phenology, vegetation dynamics and competition, terrestrial carbon and nutrient cycling. IBIS has been previously validated and applied to the Amazon (Senna et al., 2009;Delire and Foley, 1999;Foley et al., 2002;Coe et al., 2007). In those studies the model adequately simulated the carbon, energy and water budgets of the basin. However, the authors in those studies pointed out the need for better spatial representation of parameters to improve model performance in comparison with observations across the Amazon. Detailed descriptions of the model can be found in those publications.
The IBIS land surface module is derived from the land surface transfer model (LSX) (Thompson and Pollard 1995a, b). Land surface processes are represented by two vegetation layers (woody and herbaceous plants), and six soil layers (to simulate the diurnal and seasonal variations of heat and moisture in the total soil depth). The dynamics of soil volumetric water content are simulated for each layer. The soil water infiltration rate is based on the Green-Ampt formulation (Green and Ampt, 1911;Li et al., 2005Li et al., , 2006. The model has 12 PFTs that compete for light and water, using different ecological strategies. Nutrient competition is not currently included. The model allows for one or more PFT per grid cell that combined define a vegetation type. In IBIS the Amazon basin is predominantly represented by the tropical broadleaf evergreen tree PFT. The goal of this study is to clarify the importance for carbon cycle simulation of using more detailed spatially-varying parameters within this single plant functional type.

Model review
In this section we perform a brief review of the main processes that involve NPP allocation, τ w and V cmax in 2258 A. D. A. Castanho et al.: Improving simulated Amazon forest biomass and productivity simulated woody above-ground productivity and aboveground biomass in the IBIS numerical model. NPP allocation refers to the partitioning of new growth into different plant tissues including wood, leaf and fine root. Allocation is very important for simulating the carbon cycle as it directly influences long-term carbon storage . Furthermore, the amount of carbon allocated to leaves influences the total canopy photosynthesis, and the amount of carbon allocated to roots influences the amount of water uptake and nutrient acquisition, among other processes. Carbon partitioning varies between numerical models. Some models use a dynamic carbon allocation, while others are based on a predefined ratio between main plant compartments fixed by each PFT . The original configuration of IBIS used a fixed partitioning of net C of 50 % to wood (α w ), 30 % to leaves (α l ) and 20 % to roots (α r ) for the tropical broadleaf evergreen trees PFT (Eq. 1).
The biomass residence time (τ ) defines the lifetime of a unit of biomass in the plant. Many global vegetation models assume a predefined and constant value of τ for each PFT and for each plant compartment (wood (τ w ), leaf (τ l ) and fine root (τ r )). For tropical broadleaf evergreen trees in IBIS, τ w is set to 25 yr, while (τ l ) and (τ r ) are set to 1 yr. Other global vegetation models assume a constant τ w for tropical forests, ranging from 20 to 200 yr. The woody biomass residence time is a key parameter for accurately simulating biomass stocks in an ecosystem. The change in the biomass (M) of an individual plant compartment (i: wood, leaf or fine root) over a period of time is described in Eq. (2): α represents the fraction of net primary productivity (NPP) allocated to biomass pool i, and τ is the residence time of that pool, expressed in years (Foley et al., 1996). V cmax refers to the photosynthetic capacity of the plant. It is the carboxylation capacity of the enzyme RuBisCO, which catalyzes the CO 2 reaction during its assimilation process in leaves (RuBisCO is the CO 2 receptor molecule in the Calvin cycle). It is directly related to the GPP of the plant, and in IBIS it is defined initially for tropical broadleaf evergreen trees as 65 (µmol CO 2 m −2 s −1 ).

Field observation database and basin-wide extrapolations
We have assembled a wide range of published data from field observations at several sites across the Amazon basin (Fig. 1). The sites are all in undisturbed old-growth forest, with most of them being part of the RAINFOR network ("Rede Amazónica de Inventarios Forestales", Amazon Forest Inventory Network; http://www.rainfor.org/). The RAIN-FOR project is an international effort to monitor structure, composition and dynamics of the Amazonian forest in order to better understand their relationship to soil and climate (Malhi et al., 2002;Peacock et al., 2007). The RAINFOR field data are in generally based on one-hectare plots (see references for more detailed information). In this study, plot data are aggregated to the 1-degree horizontal resolution used by IBIS.

Carbon allocation
There are few plot measurements of carbon allocation to stems, roots, and leaves reported in the literature.  compiled a carbon allocation database for tropical forests worldwide. They report the partition of carbon between wood, fine roots and leaves for 10 plots in the Amazonian basin that represent 6 sites at the 1 • ×1 • grid cell resolution of the model ( Table 1). The authors showed that in general there is nearly equal allocation of new carbon between wood, leaves and fine roots. Aragão et al. (2009) suggested that the C allocation partition appears related to soil texture rather than soil fertility. The authors identified that carbon allocation to roots decreases with increasing soil clay content. They argue that this happens because, in well-drained sandy soils, roots grow with less resistance from the soil and have faster water absorption. Following the hypothesis of Aragão et al. (2009) we tested the correlation between percent sand content with both fine-root and leaf carbon. We obtained a similar correlation to their study (Fig. 2, Eq. 3 in Table 1). The carbon allocation between wood, leaves and fine roots for the whole basin is estimated based on the regressions (Eq. 3, Table 1) applied to the soil texture map . The correlations between carbon allocation and soil texture are available for a small number of sites, and may be limited by other factors that are either not well known or are not well represented by this limited database. Applying Eq. 3 to the entire basin, the estimated woody carbon allocation for the region varies between 30 and 40 % (Fig. 5a, background map). This estimate does not reproduce the amplitude of the site-specific measurement of carbon allocation (25-50 %) (Fig. 5a, bullets), which suggests that the assumption of this hypothesis does not critically affect the final results (Supplement Table B).

Woody residence time (τ ω )
Woody biomass residence time (τ ω ) has been estimated from field measurements, and a strong spatial variability has been reported (25-100 yr) . Spatially varying τ ω is included as an input parameterization in the model. In this work we use the compiled data on τ ω from Galbraith et al. (2013), which is in terms of carbon residence time, more appropriate for this study. The data set includes analyses of 129 plots across Amazonia for 5-25 yr time series between 1971 and 2011 (Table 1) (Galbraith et al., 2013). There are 34 1 • × 1 • grid cells associated with these sites (Fig. 1). There are strong indications that τ ω is correlated to soil physical properties ; however the mechanisms that would explain the spatial variability of τ ω are not completely understood , Galbraith et al., 2013. For this reason we scale up τ w to the entire basin, using simple kriging interpolation of the field data points (Fig. 5b).

Total soil phosphorus (Ptot)
Phosphorus (P) is known to be a limiting factor for productivity of mature tropical forests (Vitousek, 1984;Lloyd et al., 2010;Mercado et al., 2011); therefore it is used in this study to represent the soil fertility limitation in our model. The total available P in this work is used to estimate the maximum carboxylation capacity of RuBisCO. Quesada et al. (2010) performed extensive collection and analyses of soil data at 71 sites with varying soil properties throughout Amazonia. The 71 sites are grouped into 26 1 • × 1 • grid cells (Table 1, Fig. 1). Based on field data  we defined a relationship between total soil P measured (average to depth from 0 to 0.3 m) at the site level and the respective soil class (Fig. 3). The relationship between total P and the labile P pools is not linear and could depend on external factors such as soil texture for example. The P map estimate we present here is focused on the large-scale variability of P in soil where total P varies from 50 -500 mg Kg −1 . The total soil P correlates significantly (r 2 = 0.65, p < 0.005) with P in leaves. Quesada et al. (2011) presented a map of basin-wide distribution of soil coverage for each reference soil group. Based on the relationship presented in Fig. 3 and the soil class map, we created a spatial map of total soil P content (the average by depth from 0 to 0.3 m) of each 1 • ×1 • grid cell in the Amazon (Fig. 4). Due to the large variability of soil types within the grid cell, we expect to see discrepancies between the site level measurement and the P content derived for the grid cell. The derived total soil P map qualitatively reproduces the east (100 mg kg −1 ) to west (450 mg kg −1 ) increase in fertility as observed in the independent site level measurements. Lowest total soil P occurs in northern Brazil and southern Venezuela, which coincides with the highest values of soil sand content ( Fig. 4). Total derived P values are estimated to exceed 300 mg kg −1 in a portion of central Amazonia, but there are no observations to corroborate these estimates. Lacking further field measurements to validate our map, we use it cautiously in this study as a means of understanding the sensitivity of simulated biomass to fertility variation.

Maximum carboxylation capacity of RuBisCO (V cmax )
Maximum carboxylation capacity of RuBisCO (V cmax ) and specific leaf area index (SLA) are important properties for simulating photosynthesis. We collected the existing data on these to explore their spatial distribution in the Amazon Basin and to use as input parameterizations in the model. Fyllas et al. (2009) analyzed leaf properties at 62 RAIN-FOR plots across the Amazon Basin. These data, when averaged to our grid cell of 1 • × 1 • , represent 22 data points ( Fig. 1, Table 1). The authors present data for leaf mass per unit area (the inverse of SLA), and leaf concentration of the main growth-limiting nutrients such as N and P. Their analyses showed that soil fertility is one of the most important predictors for observed higher nutrient concentration in leaves. Mercado et al. (2009) noted a correlation between V cmax observed from the field and the concentrations of P in leaves.
We developed a similar regression equation to that of Mercado et al. (2009Mercado et al. ( , 2011 but between V cmax and total P concentration in soil instead of the P concentration in leaves ( Fig. 3b, Eq. 4 in Table 1). The advantage of this empirical regression with respect to soil P is the ability to estimate V cmax for the whole basin based on Eq. (4) ( Table 1) and the map of total soil P concentration (Fig. 4). The V cmax spatial distribution shows the same spatial structure as the P map, with a gradient from west to east (Fig. 5c). There is in general a good agreement between field estimates of V cmax  and the estimated map in this work (Fig. 5c). The San Carlos do Rio Negro (SCR, Venezuela) site represents a significant outlier, as the observed V cmax at this site (ranging around 65 µmol m −2 s −1 ; Mercado et al., 2011) is considerably higher than the estimated value based on soil P content (ranging around 35 µmol m −2 s −1 , Fig. 5c). The reason for this difference is not clear but may be due to the large differences between foliar P ) and soil P for this specific site . The limitation of the linear regression between V cmax and total soil P is that it does not reproduce the saturation in V cmax due to high levels of P content. One example of this is the high V cmax value estimated in Cuzco Amazonico (CUZ, Peru) ( Fig. 5c) due to the elevated total soil P in this site (Fig. 4). Mercado et al. (2011) suggested the use of a modified photosynthesis model that includes both P and N limitation of the main photosynthetic The dots represent averaged field plot measurements in one-degree grid cells . The total soil phosphorus map (background) is derived based on soil class map and the relationship between site level total soil P content and soil class (Fig. 3).
parameters, as described in Domingues et al. (2005Domingues et al. ( , 2010. More detailed physiological analyses are important for a better definition of the relation between V cmax and the P limitation. There are no clear relationships between SLA and other biophysical properties; therefore we interpolated the site values to the entire basin using the kriging interpolation method (Fig. 5d).

Above-ground woody net primary productivity (NPP w )
The above-ground wood net primary productivity (NPP w ) field database is used in this work for comparison to the simulated NPP w . IBIS, like many other ecosystem models, simulates a generic woody biomass pool that includes all aboveground wood and coarse roots. Therefore, to facilitate comparison with the field data, which is above-ground woody productivity only, the simulated woody net primary productivity was divided by 1.21 to remove the below-ground coarse-root fraction of the simulated NPP wood as suggested in Houghton et al. (2001). Malhi et al. (2004) present a large data set of above-ground coarse wood productivity in 104 neotropical forest plots. Their sites are all located at an elevation lower than 1000 m, in mixed-age old-growth humid forest, and with no human disturbance. For comparison with our simulated results these data are aggregated to the 1 • × 1 • grid cell resolution of this study, representing 25 grid-points across Amazonia (Fig. 1, Table 1). There is high productivity in west and low productivity in central and eastern Amazonia, varying in space from 0.15 up to 0.55 kg-C m 2 yr −1 , with an overall variability of 260 % (or 0.55/0.15 = 3.6 times) (Fig. 6a).   Malhi et al. (2006) present a synthesis of data on woody above-ground live biomass of old-growth lowland tropical forest for 227 plots across South America. This data was rescaled to the one-degree grid resolution, resulting in 69 sites for comparison with our simulation results (Fig. 1, Table 1). The spatial distribution of biomass shows high biomass in the slow-growing central Amazonia forest and Guyana, with low biomass in the western faster-growing forests and the dryer southern and eastern margins. The observed woody biomass ranges from 9 up to 20 kg-C m −2 , with a spatial variability of 120 % (or 20/9 = 2.2 times) among forest sites (Fig. 6b).

Model configuration
In order to quantify the response of IBIS to spatially varying parameters, we performed a suite of simulations in 4 different categories ( Table 2). The first category is a simulation over the entire Amazon basin that uses the original configuration of the IBIS model. This simulation serves as a reference to the other experiments and is referred to as the control simulation (CA). The second simulation assumes the original configuration CA but alters the allocation of NPP to wood,  Malhi et al. (2004) and (2006), respectively. Dots represent the average from measurement plots averaged to 1 • × 1 • grid cell.
foliage and roots so that it is more consistent with observations, allocating one-third for each component (34 %, 33 % and 33 %, respectively) (SA3a) . In these tests, constant parameter values are assigned and fixed in space for the entire Amazon basin (homogeneous parameterization). The third simulation category referred to as the site level simulation (SS) tests the importance of using the spatially heterogeneous field data to represent the fundamental parameters. In these simulations the model is run only at the site level with the parameters from our database (Fig. 5, dots and Table 1). Comparison of these simulations (SS) with CA and SA3a and with observations provides an understanding of the ability of the model to simulate the productivity and biomass at individual locations. The fourth simulation category uses the basin-wide spatially varying parameter values (Fig. 5, background map and Table 1) derived in this study. The results of this exercise (referred to as RS), when compared to CA, SA3a and the observations, quantify the value of using best estimates of basin-wide parameters to derive productivity and biomass values across the Amazon. We considered spatial variation in carbon allocation to wood, leaves and fine roots, woody biomass residence time, maximum carboxylation capacity of RuBisCO (V cmax ) and SLA ( Table 2). The specific parameter values and simulation runs are summarized in Table 2.
The model was forced with prescribed climate based on the Sheffield et al. (2006) database, which is a combination of global observation-based data sets and reanalysis data from the National Center for Environmental Prediction -National Center for Atmospheric Research (NCEP-NCAR). The data set is available from 1970 to 2008 (39 yr), has one-degree spatial resolution and 3 h time resolution, which was linearly interpolated to one hour.
The model simulations were run for a total of 680 yr (1329-2008). The long simulation was required to allow the slow carbon pools to come to equilibrium. There was an initial spin-up of 386 yr (from 1329-1715) under constant preindustrial atmospheric CO 2 values (278 ppm). The spin-up simulation started from near bare ground until soil carbon, vegetation structure and biomass achieved an equilibrium state. The runs were continued from 1715 up to 2008 with increasing prescribed atmospheric CO 2 concentrations (from 278 to 386 ppm). During the entire 680 yr run the prescribed climate was applied cyclically.
Soil texture data was based on the IGBP-DIS global soil and Quesada et al. (2010) data set. The control simulation (CA) and regional simulation (RS) used the regional map of texture, while the site level simulations used site level soil texture information from Quesada et al. (2010). The soil depth is considered homogeneous with 10 m in all simulations. There are 6 soil layers with thicknesses from the top layer to the bottom of 0.25, 0.375, 0.625, 1.25, 2.5, 5 m depth.
No land use changes or other disturbances (e.g. fire) were incorporated in the simulations. Therefore, the results are for potential vegetation conditions (the vegetation in equilibrium with the prescribed soil and climate). Potential vegetation simulations were chosen because they should be most comparable to the field data, which were collected in undisturbed old-growth tropical forest plots

Statistical Analyses
The simulated variables are averaged for the last 10 yr of simulation (1999-2008) and compared to the field data within a grid cell, which represent an average of the period of sampling. For statistical comparison we use the correlation coefficient, linear regression and the index of agreement (Willmott, 1982). The index of agreement provides information on how correlated and how distant the simulated data points are from the reference (observation), while the correlation coefficient might have a high value just because the data are well correlated even if not necessarily close in absolute values.
The index of agreement varies from 0 to 1 where 0 means a very poor agreement and 1 the maximum agreement.

Comparison of simulations and field observations
We performed a series of sensitivity analyses with the model to identify the factors affecting the spatial variability of IBISsimulated woody productivity and biomass (Supplement A and B). The sensitivity analyses indicated that spatially varying values of maximum RuBisCO were important for accurate simulation of NPP w (Supplement Table B). Woody residence time was identified as the most relevant factor for reproducing the spatial variability in above-ground biomass (Supplement Table B). In this section we quantify how the site-specific observed parameter information (SS, Table 2) affected simulated productivity and biomass compared to the homogeneous parameter assumption (CA and SA3a, Table 2). These heterogeneous parameterization simulation results (SS) were compared to observations available for wood productivity and above-ground biomass (Fig. 7).

Homogeneous parameterization
The simulated NPP w in the control simulation (CA, where carbon allocation to wood is defined as 50 %,) is systematically overestimated compared to the observations (Fig. 7a, dark square). The simulated NPP w in (SA3a, where carbon allocation to wood is defined as 34 %) is in better agreement with the average observed values (Fig. 7a, dark triangles). The spatial variation of the observed NPP w was not reproduced in these simulations.
Similar analyses were performed for the above-ground biomass comparison (Fig. 7b). The simulated biomass in the control (CA) does not reproduce the spatial variation, but average values are similar to the observed (Fig. 7b, dark square). If only one of the parameters is corrected, for example carbon allocation (SA3a, dark triangle), the estimated biomass deviates strongly from the observations (Fig. 7b), which is consistent with results with the DVM ORCHIDEE (Delbart et al., 2009).

Heterogeneous parameterization
The results from simulations that include the combined effects of all site level heterogeneous parameterizations are presented in this section (Fig. 7, red and black).
The NPP w analyses are made for two series of data. The first considers only the data points where field estimates of V cmax are available ( Fig. 7a and c, in red), hereinafter Series A. The second data series includes data points where the V cmax has been extrapolated to all other locations where there are NPP w measurements (Table 1), hereinafter Series B (Fig. 7c, in black). The full set of NPP w that include sites with known V cmax from field and extrapolated V cmax in this work is called Series A + B (Fig. 7c, in red and black).
The NPP w from simulation SS with heterogeneous parameterization (Fig. 7a, red circles) shows better agreement with the observations when compared to the previous homogeneous assumption SA3a (Fig. 7a, triangles). The coefficient of correlation improves from 0.10 to 0.71, the regression slope (error) improves from −0.04(0.03) to 0.6(0.2), and intercept (error) from 0.35(0.05) to 0.16(0.07) ( Table 3, Series A). The slope coefficient improves but is still low, indicating that NPP w is overestimated by the model where observed values are low (Fig. 7a). The NPP w simulated by IBIS with the new spatial parameterization agrees better with the observed spatial variability. The index of agreement (Fig. 7a, Table 3) increases from 0 (SA3a, Series A) to 0.7 (SS, Series A).
The comparison of results of Series B (extrapolated parameterization) in simulation RS with SA3a (Table 3, Fig. 7c in black only) shows that some locations are closer to the observations (index of agreement 0.4 compared to 0.3) but the correlation is not improved (coefficient of correlation 0.66 compared to 0.79). This indicates, as expected, that the regionalized parameter data are not as good as site-specific field estimates. However, the conclusions for the regionalization of NPP w are still limited to a validation against only 10 data points that are not representative of the entire basin.
Similar to NPP w , AGB w comparisons were made for two series: one including only the data points where woody residence time estimates are known (Fig. 7b, in red), referred to as Series C. The second data series includes locations where the woody residence time is extrapolated to all other AGB w data collection points, called Series D (Fig. 7d, in black). The full set of simulated AGB w values derived from both direct measurements of residence time and the extrapolations is called Series C + D (Fig. 7d, in red and black). The heterogeneous woody residence time data has the strongest influence on the simulated spatial variability of AGB w (Fig. 7b, red circles). The data series SS captures the spatial variability of the AGB w much better than the SA3a simulation with homogeneous parameterizations (dark triangles). The   Table 3. Statistical summary of the comparison of woody net primary productivity between IBIS-simulated results and field estimates. The table presents the number of data points within the studied area, mean and standard deviation, regression slope, intercept and correlation coefficient, and index of agreement (d). The statistical analyses were made for all sites excluding four outliers (JEN, CAQ, SCR, CUZ, discussed in Supplement D). The statistical analyses were divided in groups of data point as follows: Series A and B represent the series of all data points that have available NPP w field information; Series A represents the series of data where the V cmax was estimated based on field information; and Series B represents all other data points where NPP w field was known and V cmax was extrapolated based on the methodology described in this work. coefficient of correlation of simulated woody above-ground biomass improves from 0.22 with homogeneous parameterization (SA3a, Series C) to 0.80 with spatial varying parameters (SS, Series C). The regression analyses show a significant improvement of the slope(error) from 0.05(0.05) to 1.06(0.18), and intercept(error) from 6.5(0.8) to −2.8(2.8) ( Table 4, Series C). The index of agreement in Series C improved from 0 with homogeneous parameterization (Table 4, Series C, SA3a) to 0.8 with heterogeneous parameterization (Table 4, Series C, SS). Some outliers were identifyed and are discussed in detail in Supplement D.
The statistical analyses of AGB w from Series D (with the extrapolated woody residence time) provide a measure of the value of the extrapolation adopted in this work (Table 4, Series D). The results show improvement of all statistical parameters comparing the regional simulation (RS, Series D) to the homogeneous assumption (SA3a, Series D). The coefficient of correlation improved from −0.006 to 0.52 and the index of agreement from 0 to 0.6. As expected, the results Table 4. Statistical summary of the comparison of wood above-ground biomass between IBIS-simulated results and field estimates. The table presents the number of data points within the studied area, mean and standard deviation, regression slope, intercept and correlation coefficient, and index of agreement (d). The statistical analyses were made for all sites excluding three (CUZ, CHN and AMB, in Supplement D). The statistical analyses were divided in groups of data point as follows: Series C and D represent the series of all data points that have available AGB w field information; Series C represents the series of data where the woody residence time was estimated based on field information; and Series D represents all other data points where AGB w field was known and woody residence time was extrapolated based on the methodology described in this work. derived from extrapolated parameters (RS, Series D) are in poorer agreement with the results derived from SS Series C where residence time is site measured, with a lower slope (0.44) and intercept of 6.9 and larger variance of the distribution differences (6.5 compared to 3.3 in Series C).

Regional Simulation Analyses (RS)
In this section we present the basin-wide simulated woody above-ground productivity and biomass (simulation RS), based on our maps of spatially varying parameters (Table 1, Fig. 5). Quantitative validation of simulated NPP w and AGB w was discussed in Sect. 3.1 (Table 3, Series B; Table 4, Series D) at the specific locations in the basin. Qualitatively, there is much greater spatial variation of the NPP w and AGB w across the basin in the RS simulation ( Fig. 8c and d) compared to CA ( Fig. 8a and b). The simulated NPP w spatial variability (Fig. 8c) follows the soil total phosphorus map patterns (Fig. 4). There is higher productivity in the west where the fertility is higher and also in central west Amazonia where P content is higher and the soil is silt. The productivity decreases in central and east Amazonia and increases again in the northeast. There is a region of low NPP w in northern Brazil and southern Venezuela due to the low estimated soil fertility in that area. A qualitative comparison of this simulation with a published satellite-derived map of NPP in Amazonia (Nunes et al., 2012) suggests similarities in spatial variability. A detailed analysis of the comparison of these simulations including land use change and other disturbances and comparison with satellite products is in progress and will be addressed in a different study.
The biomass map shows a west-to-east trend of higher biomass in central and northeastern Amazonia (RS, Fig. 8d) that is qualitatively consistent with the observed field data. The simulated spatial variability is linked to the woody residence time map. The biomass map (Fig. 8d) can be qualita-tively compared to the estimated biomass map from Malhi et al. (2006). Both maps suggest higher biomass occurs in the central and northeastern Amazon and lower biomass in the west and south. In the central Amazon the highest biomass values in the Malhi et al. (2006) data set are clustered around the sites of measurement, which is most likely an effect of extrapolation from the few data points available (Malhi et al., 2006). The gradient of biomass in the transition to the savanna environment in the southeast diverges in the absolute values from the Malhi et al. (2006) estimates. There is little field data for this region, making validation difficult.

Discussion
From this modeling exercise we have identified six major points regarding the controls on and importance of spatial variation in above-ground productivity and biomass.

Homogeneous parameterizations and climate alone do not capture the spatial variability of woody above-ground productivity (NPP w ) and biomass (AGB w ) present in the field data
The results suggest that models may be able to simulate average NPP w and AGB w in reasonable agreement with observations but sometimes because of compensation of errors. For example, in simulation CA (the woody carbon allocation parameter is 50 % and the woody residence time is 25 yr) the resulting average biomass is in relatively good agreement with observations (Fig. 7b, dark square). As discussed earlier, those parameter values, although commonly used, have been shown by observations to be incorrect. Using a high woody allocation rate (50 % vs 34 % observed) and low residence time (25 yr vs. 25-100 yr observed) causes a large carbon allocation each year but low turnover rate and coincidently total biomass in good agreement with observations. In addition, the simulated results using homogeneous parameters do not capture the observed spatial variability of NPP w (260 %, the highest value is 3.6 times higher than the lowest value; Malhi et al., 2004) and AGB w (120 %, highest value is 2.2 times higher than the lowest value; Malhi et al., 2006) across the Amazon basin ( Fig. 7a and b dark square). There is relatively large spatial variation in the Amazon climate. However, field data analyses have shown poor correlation between climate variables and large-scale variation of productivity and biomass patterns in the basin Quesada et al., 2012;Galbraith et al., 2013). An analysis of the spatial patterns of productivity  found no obvious relationship between the spatial distribution of wood productivity and precipitation, dry season length, or radiation, but the authors did find some decline in woody productivity with increasing temperature. However, the lower temperature in western Amazonia is indirectly correlated to higher soil fertility in that region, making it difficult to directly correlate temperature and productivity given the strong correlation between productivity and soil fertility. In analyses of basal area spatial patterns, Malhi et al. (2006) found some correlation with dry season length and precipitation. The decline in basal area, however, was evident only with extreme water stress and long dry season length (more than 4 months).
The spatial variation of climate in the Amazon basin imparts simulated productivity and biomass variation of 35 % and 45 % from the basin average minimum (Table B). IBIS and most numerical models underestimate or do not explicitly consider mortality rates due to short-or long-term disturbances such as temperature extremes, drought or flooding Galbraith et al., 2013). However, because the observed data were collected at sites specifically chosen to be free of recent disturbance, lack of climate-stress mortality is not a factor in the simulated low sensitivity to climate variation. Therefore, the low variability as a function of climate in the IBIS simulations appears to be a robust indication of the scale of the actual variability imparted by climate.

Maximum carboxylation capacity of RuBisCO activity (V cmax ) as a function of soil fertility is the primary variable controlling the simulated variation of woody above-ground productivity across the Amazon
Observations suggest that soil fertility plays an important role in creating spatial variation of productivity across the Amazon basin, but little is known about the specific pathways Aragão et al., 2009. Some of the possible factors that have been explored are differences in gross primary productivity, respiration rates, and carbon allocation between carbon pools.

2268
A. D. A. Castanho et al.: Improving simulated Amazon forest biomass and productivity Malhi et al. (2004) suggest that carbon allocation may be related to the spatial variability of NPP w , due to shifts in allocation to roots in less fertile soils. Alternatively the author's data analyses also suggest that carbon use efficiency (CUE = NPP/GPP) is higher in the west, which suggests that in addition to a potentially higher GPP in the higher fertility soils of the western Amazon, there are also potentially higher autotrophic respiration rates in the lower fertility soils of central Amazonia. We found that alteration of GPP and autotrophic respiration achieved through a change in leaf photosynthetic capacity (related to maximum carboxylase capacity of RuBisCO) is the strongest candidate to explain the spatial variability of productivity. The sensitivity tests show that spatial variation of RuBisCO (from 75-40 µmol CO 2 m −2 s −1 , SA5, Table A) leads to a simulated 60-80 % change in NPP w (Table B). A shift in woody carbon allocation from 50 % to 25 % imparts a 60 % change in NPP w (Table B). Our analyses suggest that V cmax driven by soil fertility plays a stronger role than carbon allocation in the spatial variability of NPP w .

Heterogeneous parameterizations of key biophysical properties based on site-specific data improves the simulated woody above-ground productivity compared to homogeneous parameterizations
We estimated the maximum carboxylation capacity of Ru-BisCO as a function of site-observed total soil P from Quesada et al. (2010). This estimate was based on the established correlation between V cmax and leaf P content (Mercado et al., , 2011, and the understanding that leaf nutrients are directly related to soil nutrients , The spatial heterogeneity of other parameters of minor effect in the NPP w calculation, such as carbon allocation, woody residence time, and specific leaf area, were also derived based on published site field data. The IBIS-simulated NPP w using heterogeneous parameterization was compared to NPP w field data . As presented in the Results (Sect. 3) the use of heterogeneous parameterizations in IBIS based on field data significantly improved the simulation of NPP w (correlation coefficient −0.1 to 0.71, respectively, for the homogeneous and heterogeneous parameterization, Table 3). The simulated NPP w values at the lower end of the observations were not well captured by the model as indicated by a slope of 0.58 ± 0.20 and an offset of 0.16 ± 0.07 (Table 3). This could be due to an overestimation of the V cmax for low-fertility sites or by other factors that are not captured by the model. For example, lowland areas may have higher respiration rates than what the model predicts, due to the higher temperatures and/or the higher respiratory costs due to slower plant growth in less fertile soils . These effects may not be fully represented in IBIS and may contribute to the overestimation of the low end of NPP w , in central and east Amazo-nian sites. A detailed analysis of field information on the respiration rates, available biomass in decomposition state and how the model reproduces these processes needs to be addressed in the future. A factor that may be contributing to a general overestimation of NPP w is the simulated Leaf Area Index (LAI). The LAI in the model is systematically higher than the observations, which would cause an increase in simulated NPP. Each year, because of increasing atmospheric CO 2 there is more carbon to allocate, and the allocation in the model is constant in time as is the residence time of the leaves. Therefore, the LAI in the model increases with time.
The assumption of constant allocation and residence time is probably too simplistic, but the right balance of changes in carbon allocation between the plant components in space and time is poorly known.
In our analyses we observed that the differences of simulated NPP w from the observations do not appear to be related to misrepresentation of carbon allocation. For example, if we define woody carbon allocation to a value that minimizes the NPP w error, then the carbon allocation to wood would have to be unrealistically low (0.15-0.25 compared to the observed values 0.25-0.5). There are most likely other unknown factors contributing to the lack of good agreement, including shorter timescale variability of parameter response to drought, fertility and disturbance for which we do not yet have data. Further improvements of these parameters are clearly required and are discussed in suggestion for future work below.

Woody residence time is the most important mechanism affecting the magnitude and spatial distribution of simulated AGB w
Our results of the IBIS sensitivity analyses suggest that woody residence time is the most important mechanism affecting the magnitude and spatial distribution of simulated AGB w . This result is in agreement with the analyses of Senna et al. (2009) and Delbart et al. (2010). Spatial variation of woody residence time within the range of 25 to 100 yr changes simulated AGB w by 180 % from the basin minimum (Table B).

Heterogeneous parameterizations of key biophysical properties based on site-specific data improves the simulated woody above-ground biomass when compared to homogeneous parameterizations
We used field data available for 21 sites across the Amazon basin to represent the spatial heterogeneity of woody residence time Galbraith et al., 2013). The spatial variability of other parameters of minor effect in AGB w such as carbon allocation, V cmax , and specific leaf area index were also included in the simulation of AGB. The comparison of the simulated AGB w with observations showed an improvement in the coefficient of correlation from 0.22 in the homogeneous experiments (SA3a) to 0.80 in the heterogeneous ones (SS) ( Table 4). The regression fit is much closer to the 1 : 1 relation. The good agreement comes from two factors: one is that the AGB simulated by the model was based on residence times estimated from field data; the other is that most of the sites for comparisons of biomass come from direct measurement of individual trees (Baker et al., 2004b). As a result, both model and field measurement methodologies used the most accurate information available. 4.6 The regional maps of key parameter values developed in this study significantly improve simulated AGB w compared to simulations using homogeneous parameterization assumptions

Biogeosciences
The quality of the regional maps of the physiological properties depends on three factors: (1) the number of site level data points available; (2) how representative these sites are of the larger scale; and (3) how well we understand what drives the spatial variability of these properties. We developed a set of extrapolations of the observed field data to the entire Amazon basin in order to introduce greater heterogeneity into the simulations of NPP w and AGB w . Unfortunately there are relatively few field data sites (Table 1) and the processes that govern the spatial variability are not completely understood, which makes our regional maps somewhat speculative. The V cmax was extrapolated to the entire basin using a soil map of P and the regression equation between field-based V cmax and P (Table 1, Fig. 5c). The comparison of simulated NPP w (from the extrapolated V cmax ) with observations was represented by a few points (n = 9). Although the simulated results for some of the individual data points were improved, it was not significant (p > 0.5) considering the small size of the sample.
We estimated woody biomass residence time for the entire basin based on a simple kriging interpolation (Fig. 5b). The large-scale spatial variability showed a trend of lower τ w in west than in central and east Amazonia, which is in qualitative agreement with the observations. The comparison of the basin-wide simulations suggests that this first attempt at spatially varying biophysical parameters improved the simulations. For example, the coefficient of correlation of simulated AGB w with observations at 46 locations increased (0.52, RS) compared to the more standard homogeneous parameters (∼ 0, CA) ( Table 4). The results showed that the proposed map of woody residence time, despite its uncertainties, was a better option than a homogeneous assumption for the entire basin.
Given the importance of the woody residence time in simulating the AGB w , it is fundamental to understand the factors that govern its spatial distribution. However there is no single mechanism that is known to control its geographic pattern. Many authors have explored this topic including  who have argued that woody residence time is correlated with soil physical properties such as soil depth, soil structure and topography. However, residence time has also been found to be well correlated with soil fertility within the Amazon Galbraith et al., 2013). There are also other soil conditions, such as rooting depth limitations, low drainage capacity, poor soil structure, and topographic position that might affect woody residence time . Delbart et al. (2010) presented an alternative solution of dynamic estimation of woody residence time as a function of NPP, whereby the function was defined based on the empirical correlation between them. This NPP-based assumption is, however, only valid for a forest in equilibrium, and the NPP needs to be well estimated by the model. This combination of poorly known governing factors makes it difficult to mechanistically determine residence time across the Amazon. We believe that our prescribed approach using a simplified interpolation is a good starting point, although we cannot assume residence times will remain unchanged under future scenarios of climate change.

Conclusions
This work has presented efforts at incremental improvement of numerical modeling of tropical broadleaf evergreen forests in the Amazon. It has identified some of the most important and relevant parameters for the simulation of C fluxes and stocks by a DGVM and shown the importance of spatial and temporal representation of functional diversity of tropical ecosystem modeling. This project has also helped identify several areas of research in data collection and model development that could lead to further improvements in model representation of tropical forest environments. Accuracy of regional simulations is still limited by data scarcity of field data. Data needs include (1) expansion of the network of field data monitoring to a better characterization of the basin, (2) examination of the physiological processes that govern the spatial variability of the main parameterizations, and (3) development of new satellite and airborne measurements of key biophysical properties and mechanisms.
Accuracy is also limited by poor numerical representation of several important biophysical properties. Model improvements include (1) dynamic representation of disturbances such as land use change and fire, (2) improvement of model biophysical processes such as photosynthetic response to nutrient limitation and tree mortality as a function of drought stress, and (3) development of more-dynamic representations of PFTs to better characterize the functional diversity of tropical forests.
Continued stepwise improvement in models and data collection such as those described here will yield greater understanding of the spatial and temporal variations of forest carbon stocks throughout the tropics.