Manish Verma and Mark Friedl ,

1. Within and Across-Biome Seasonal Variation in GPP: An Assessment of Remote Sensing Methods Using the La Thuile Dataset Manish Verma and Mark Friedl, Department of Earth and Environment, Boston University, 675 Commonwealth Avenue, Boston, MA 02215, USA. Andrew D. Richardson Department of Organismic and Evolutionary Biology, Harvard University, HUH, 22 Divinity Avenue, Cambridge, MA 02138, USA.

2. Summary • Objective: • To compare the relative suitability of 9 different remote-sensing-based models for capturing intra- and inter-annual variations in daily GPP in different biomes. • Methods: • Used 134 sites (intra-annual analysis) and 94 site years (inter-annual analysis) of daily GPP data from the La Thuile dataset to examine the agreement between • (i) modeled vs tower GPP at seasonal scale (intra-annual), and • (ii) Interannual anomalies in modeled vs tower GPP • Key findings: • Models were significantly more powerful capturing site-specific intra-annual variation than inter-annual variation in every biome in GPP. • More complex models were better able to capture intra-annual variation in daily GPP than simple models in all biomes, except croplands, grasslands and mixed forests; But, model complexity did not improve model results at inter-annual time scale. • Light use efficiency models performed just as well as a non-parametric, empirical neural network model in both intra- and inter-annual analysis, suggesting no real advantage for the latter approach.

3. Background Terrestrial gross primary productivity (GPP) is the largest component of the global carbon cycle (Beer et al., 2010). Regular monitoring of terrestrial GPP is required to understand the global carbon cycle and project future climate (Bunn et al., 2006). Different models that use satellite data have been proposed to estimate seasonal variation in GPP both within and across biomes.

4. Background Average intra-annual variation in daily GPP based on 6 years of data Inter-annual variation (anomalies) in daily GPP for 6 different years. We compared the ability of 9 different remote-sensing-based models for capturing site-specific intra- and inter-annual variation in daily GPP. The figures below illustrate the two types of variability at a site (FR-Hes) in daily (8-day mean) GPP .

5. Background & Justification Intra- and inter-annual variations in daily GPP have different magnitudes. Typically, intra-annual variation is significantly greater than inter-annual variation in daily GPP for most biomes. The biotic and abiotic drivers of the two types of variations are also different (e.g. in deciduous biomes, leaf area is the dominant control on intra-annual variations, but not on inter-annual variations). Here we examine the ability of models to capture these two types of variability.

6. Background & Justification The 9 models use different input variables and parameters, have different functional forms, and cover a range of model complexity.

7. Background & Justification They are also based on different hypotheses about controls on GPP Increasing complexity

8. Motivating Questions While photosynthesis at leaf and canopy level is well understood, a great deal of uncertainty remains about the drivers and controls of variations in ecosystem level productivity (Beer et al., 2010). It would therefore be useful to understand how different satellite-based models (which represent different hypotheses about controls on GPP) vary in their ability to capture site specific intra- and inter-annual variations in seasonal GPP.

9. Data & Methods: Selection of Sites from the La Thuile Dataset We first identified sites where for each site-year more than 95% of the days had daily GPP data, and the mean daily quality flag was more than 0.75 (Richardson et al., 2010). We then removed sites with excessive local heterogeneity in land cover surrounding the flux site. To do this, we extracted MODIS land cover type data (Friedl et al., 2010) at 500 m spatial resolution for 7 by 7 windows centered on each site, and excluded sites where fewer than 7 pixels were classified as belonging to the same land cover type as the tower site. The final data set included 144 sites with 422 site-years of data. This data was further screened for uncertainty in intra- and inter-annual analysis (described on slides 13 to 15).

10. Data & Methods: MODIS Data • We used the following MODIS data and associated quality flags • Normalized difference vegetation index (NDVI; 500m), • Enhanced vegetation index (EVI; 500m), • Land surface water index (LSWI; 500m; Xiao et al, 2004), • Growing period length (GPL; 500m; Ganguly et al., 2010) • Land cover (500 and 1000m; Friedl et al., 2010) • Fraction of absorbed photosynthetically active radiation (FPAR; 1000m; Myneni et al 2002), • Day and night land surface temperature (LST; 1000m; Wan et al., 2002), • MOD17 GPP (1000m; Running et al., 2004). • We downloaded 500 and 1000 meter products for 7 by 7 and 3 by 3 pixel windows (respectively) centered on each site (Heinsch et al., 2006; Schaefer et al., 2012). • Using the MODIS land cover type product, we then selected (i) the center pixel, and (ii) other pixels in the window with land cover labels equivalent to the land cover type at each flux tower. MODIS data at each site were then averaged over the selected pixels to produce a single value for each MODIS product at each time step. • Spatial averaging over surrounding pixels with similar land cover is likely to minimize random noise, maximize correspondence between remote sensing and tower measurements, and reduce the effects of gridding artifacts (Tan et al., 2006).

11. Data & Methods: Model Calibration • We used daily GPP, PAR, temperature and precipitation data from the La Thuile dataset. • We downloaded data for the MOD17 product from NASA’s Land Processing Distributed Active Archive Center (LP DAAC) and followed Sims et al. (2008) to estimate TG model results. • The remaining seven models (see Slide 6) were calibrated to tower GPP at daily time step using linear or non-linear least squares. For each model in every biome, we used a leave-one-site-out cross-validation method and the used calibrated parameters to predict daily GPP at each test site.

12. Data & Methods: Uncertainty in Tower GPP To explore agreement between tower and modeled GPP we used 8-day mean daily GPP values. We assume that at the scale of 8-day time steps, uncertainty in tower GPP is 10% (± 1 standard deviation) of GPP. Based on this assumption we selected sites (for intra-annual analysis) and site years (for inter-annual analysis) where variance in tower data was large relative to uncertainty (see next 3 slides).

13. Data & Methods: Site Selection for Intra-Annual Analysis We excluded sites where the standard deviation in daily GPP was less than 3 times the mean uncertainty (except in EBF, see the next slide). Hainich (Germany): Seasonal variance is many times larger than uncertainty Santarem (Brasil): Seasonal variance is low relative to uncertainty.

14. Data & Methods: Site Selection for Intra-Annual Analysis In EBF mean daily GPP is large but seasonal variability is small. Thus, the ratio of the standard deviation in daily GPP to mean uncertainty in daily GPP is relatively small. EBF however are a special case where GPP seasonality is not synchronous with the seasonality of leaf area, and it is important to know if remote sensing models succeed in this case. Thus, in EBF we relaxed the criteria and excluded sites where the standard deviation in daily GPP was < 2 times the mean uncertainty. Relationship between mean and standard deviation in daily tower GPP in photosynthetically active period (daily GPP > 0.33 gC/m2/day)

15. Data & Methods: Site-Year Selection for Inter-Annual Analysis We selected sites based on 2 criteria: (1) sites were required to have 3 or more years of data; (2) only site-years with “large anomalies”, defined based on biome-specific thresholds, were included. “Large anomalies” have higher signal to noise ratio and are therefore robust relative to uncertainty. Seasonal anomalies from mean fluxes at Harvard Forest. We selected only the year corresponding to blue circles. Horizontal lines show ± 0.75 gC/m2 Seasonal variation in fluxes for three years at Harvard Forest

16. Data & Methods: Summary of Sites and Site Years from the La Thuile Dataset Used in the Study CSH: Shrubland; CRO: Croplands; DBF: Deciduous Broadleaf ; EBF.: Evergreen Broadleaf; ENF: Evergreen Needleleaf; GRA: Grassland; MF: Mixed Forest; SAV: Savannas and Woody Savannas

17. Data & Methods: Statistical Analysis Intra-Annual Variation For each model 'm' and site 's', we calculated the R2 , relative-RMSE and slope of regression (Beta) between daily (8-day mean) modeled and tower GPP. Relative-RMSE is the RMSE normalized by the mean daily tower GPP. Using the site specific R2, relative-RMSE and Beta values we calculated biome level means and 95% confidence interval for each of the three metrics for each model. Inter-Annual Variation Similarly, for each model ‘m’ and site-year ‘sy’ we calculated the the R2 , relative-RMSE and Beta between tower and modeled anomalies to calculate biome specific means and bootstrap confidence intervals for each metric. EVI-Linear is the simplest of the nine model and functions as the “null model” in our analysis.

18. Results: Intra-annual Analysis CSH: Shrubland; CRO: Croplands; DBF: Deciduous Broadleaf ; EBF.: Evergreen Broadleaf; ENF: Evergreen Needleleaf; GRA: Grassland; MF: Mixed Forest; SAV: Savannas and Woody Savannas

19. Results: Intra-annual Analysis R2: Models capture a large proportion of intra-annual variability in tower GPP. More complex models (the VPRM, MOD17-Tower and neural network model) perform better (greater R2 and/or smaller variations around mean) than the “null model” in all biomes except CRO, GRA and MF. RMSE: Relative RMSE for more complex models (the VPRM, MOD17-Tower and neural network model) tended to be around 0.5 in all biomes except CRO, where it was around 0.7. The VPRM, MOD17-Tower and neural network models, recorded lower relative RMSE (smaller mean and/or variations around mean) than the “null model” in all biomes except CRO, GRA and MF. Beta: In CRO, DBF and MF there was almost no difference in slopes between the other models and the “null” model. In the remaining biomes, the slope between predicted and tower GPP was closer to 1 for the VPRM, MOD17-Tower and neural network model relative to the null model.

20. Results: Inter-annual Analysis CRO: Croplands; DBF: Deciduous Broadleaf ; EBF.: Evergreen Broadleaf; ENF: Evergreen Needleleaf; GRA: Grassland; SAV: Savannas and Woody Savannas

21. Results: Inter-annual Analysis R2: Models were generally unsuccessful in capturing interannual variations in daily tower GPP. R2 was close to or greater than 0.5 Only in GRA. More complex models (the VPRM, MOD17-Tower and neural network model) did not perform better than the “null model”, except in GRA. RMSE: Relative RMSE was generally close to or more than 1.5. The VPRM, MOD17-Tower and neural network models, did not record lower relative RMSE than the “null model” except in GRA. Beta: There was little difference in across the different models.

22. Conclusions Models were relatively more successful in capturing intra-annual variations in GPP. In capturing intra-annual variation in tower GPP, more complex models performed better than the “null model” in shrubland, deciduous broadleaf forest, evergreen needleleaf forest, evergreen broadleaf forest, and savannas but not in cropland, grassland and mixed forest. More complex models did not perform better than the “null model” in capturing inter-annual variation. There was little difference between the performance of the two light use efficiency models (MOD17-Tower and VPRM), and the neural network model in both intra- and inter-annual analysis.