Advances in Space Research 39 (2007) 125–130 www.elsevier.com/locate/asr
Estimating net primary production for Scandinavian forests using data from Terra/MODIS Pontus Olofsson a
a,*
, Lars Eklundh a, Fredrik Lagergren a, Per Jo¨nsson
b,c
, Anders Lindroth
a
Department of Physical Geography and Ecosystems Analysis, Geobiosphere Science Centre, Lund University, So¨lvegatan 12, SE-223 62 Lund, Sweden b Department of Physics, Lund University, So¨lvegatan 14, SE-221 00 Lund, Sweden c ¨ stra Varvsgatan 11 H, SE-205 06 Malmo¨, Sweden Division of Mathematics, Natural Sciences and Language, Malmo¨ University, O Received 1 November 2004; received in revised form 28 May 2005; accepted 15 February 2006
Abstract A model for estimating net primary production (NPP) across Scandinavia has benn developed. The model is based on the light-use efficiency (LUE) concept, where NPP is calculated as a product of the amount of absorbed photosynthetically active radiation (APAR) and a LUE-factor (e) controlling the efficiency by which vegetation transforms photosynthetically active radiation (PAR) into biomass. The fractional APAR (FAPAR) is obtained by a linear transformation of 250 m normalized difference vegetation index (NDVI) data from the Moderate Resolution Imaging Spectroradiometer (MODIS). Prior to the transformation, the NDVI time series were seasonally adjusted by fitting local asymmetric Gaussian curves to the data thereby minimizing cloud contamination and other noise factors such as BRDF-related variation. The APAR is then calculated as a product of FAPAR and incident PAR, where the latter is obtained from the ‘‘Swedish Regional Climate Modelling Programme’’ (SWECLIM). The LUE-factor is modeled as a daily function of temperature, latitude, and day of the year (DOY). The FAPAR dataset was validated against measurements of fractional intercepted PAR (FIPAR) that were carried out at Norunda experimental site in central Sweden (605 0 N, 1729 0 E) between August and October, 2001. The calculated FAPAR time series are in good agreement with the measurements. The modeled NPP is evaluated against flux measurements carried out in Norunda 1997–1999. The determination coefficients obtained when comparing modeled with measured data are R2 = 0.82 for 2000 (RMSE = 2.71 g C m2 d1), and R2 = 0.68 for 2001 (RMSE = 3.57 g C m2 d1) (time series averaged every 10th day). 2007 Published by Elsevier Ltd on behalf of COSPAR. Keywords: NPP; FAPAR; LUE; MODIS; Forest; Scandinavia
1. Introduction It is necessary to strengthen our knowledge regarding the role of ecosystems as sinks for carbon dioxide to meet the needs of information with the fulfilment of the Kyoto protocol of the UN convention of climate change (UNFCCC).
*
Corresponding author. Tel.: +46 46 2229683; fax: +46 46 2224011. E-mail addresses:
[email protected] (P. Olofsson), lars.
[email protected] (L. Eklundh),
[email protected] (F. Lagergren),
[email protected] (P. Jo¨nsson). 0273-1177/$30 2007 Published by Elsevier Ltd on behalf of COSPAR. doi:10.1016/j.asr.2006.02.031
Due to thier geographical location, the high-latitude ecosystems of Scandinavia plays an important role in the carbon cycle. Soils in these areas store high amounts of carbon, and northern-latitude ecosystems play a significant role in the global carbon balance (Schlesinger, 1991). A common and appealing way to determine the carbon dioxide exchange from remotely sensed data is to estimate the photosynthetic rate through the light absorption of the vegetation canopy. This method is based on the principle that net primary production (NPP) can be estimated from the product of absorbed photosynthetically active radiation (APAR) and a light-use efficiency coefficient (e) (Monteith, 1977; Ruimy et al., 1994; Seaquist, 2001):
126
P n ¼ eEa ;
P. Olofsson et al. / Advances in Space Research 39 (2007) 125–130
ð1Þ
where Pn is the NPP and Ea the APAR. The fraction of photosynthetically active radiation (FAPAR) absorbed by vegetation is usually estimated with spectral vegetation indices such as the Normalized Difference Vegetation Index (NDVI) (Asrar et al., 1984; Tucker and Sellers, 1986; Goward and Huemmrich, 1992). The efficiency of the conversion of radiation into carbon, e, may vary substantially in time and space. The most common approach to determining epsilon is through its relation with various environmental and vegetation-type specific factors that describe the limits to the carbon uptake (Prince, 1991; Lagergren et al., 2005). In our current model PAR data and epsilon values have been estimated with meteorological model data from the SWECLIM project (Rummukainen et al., 2000). 2. Material and methods 2.1. Measurements 2.1.1. Measurement sites The Norunda area, located in central Sweden about 100 km northwest of Stockholm (6050 0 N, 1729 0 E), mainly consists of old and middle aged (50–100 yr) Scots pine and Norway Spruce with a leaf area index (LAI) estimated to range from 3 to 7 using the Li-Cor LAI-2000 (Lundin et al., 1999). The mean air temperature is 5.5 (1961–1990) and the soil is a sandy glacial till. More information about the Norunda area can be found in Lundin et al. (1999). Skyttorp is located approximately 20 km east of Norunda. The two plots where the measurements were carried out consist of a 30-year-old Scots pine stand and a 60-year-old mixed Scots pine and Norway spruce stand. The latter is located approximately 4 km south of the former. 2.1.2. PAR Incident Photosynthetically Active Radiation (PAR), corresponding to wavelengths between 0.4 and 0.7 lm, was measured below canopy at two stands in Norunda and at two in Skyttorp: five loggers (‘‘MiniCube VV’’ from Environmental Measuring Systems, Brno, Czech Republic) labelled L1–L5 each logged data from twelve PAR sensors (‘‘JYP 1000’’ from SDEC France, Reignac Sur Indre, France), labelled L1:1, L1:2, . . . , L5:12, were connected. To logger L1 only 10 sensors were connected. Loggers L2 and L4 were placed in the same satellite pixels as L3 (both at 250 m and 1 km resolution), and the sensors connected to L4 were repositioned 17 September, 09:35. The loggers were placed in the center of the Norunda area and in Skyttorp: L1. Skyttorp, Scots Pine, 30-year-old stand (WGS 84: 6007.698 0 N, 1750.135 0 E); MODIS tile h18v02, R 4738, C 4263 for 250 m resolution (R 1184, C 1065 for 1 km resolution)
L3. Norunda, Norway Spruce (WGS 84: 6005.033 0 N, 1728.500 0 E); MODIS tile h18v02, R 4759, C 4183 for 250 m resolution (R 1189, C 1045 for 1 km resolution) L5. Skyttorp, mixed Scots pine and Norway spruce, 60year-old stand (WGS 84: 6005.730 0 N, 1749.593 0 E); MODIS tile h18v02, R 4754, C 4266 for 250 m resolution (R 1188, C 1066 for 1 km resolution) The incident PAR above canopy, required to calculate the fraction of intercepted PAR (FIPAR), was measured by a sensor placed on the Norunda tower at a height of 102 m (Lundin et al., 1999). The correlation between sensor on the ground and sensor on the mast was R2 = 0.984. The observations recorded by the sensors, both above and below the canopy, were averaged every fifth minute. 2.1.3. FIPAR and FAPAR The fraction of photosynthetically active radiation intercepted by the canopy is the sum of the absorbed and reflected fractions and is often estimated as 1 t, where t is the canopy transmittance (Russell et al., 1989). Using the terminology in Gower et al. (1999), FIPAR (fi) is estimated using the following equation: Ei ; ð2Þ E where Ei is the incident PAR below canopy (transmitted PAR; in our case measured by L1:1, L1:2, . . . , L5:12), and E the incident PAR above canopy (in our case measured by the sensor mounted on the Norunda tower). Accordingly, FIPAR can be calculated for the four stands from 14 August to 25 October 2001. In order to acquire FAPAR (fa) from FIPAR the reflected fraction of the intercepted radiation must be subtracted from the FIPAR-values. Hence, by setting the reflected fraction to 3% (Russell et al., 1989), FAPAR is calculated according to: fi ¼ 1
fa ¼ fi 0:03:
ð3Þ
2.1.4. NPP In order to validate the modeled NPP, daily flux measurements of net ecosystem exchange (NEE) for 1997– 1999 collected at Norunda were recalculated into NPP by adding the heterotrophic respiration. Heterotrophic respiration was calculated from forestfloor respiration and the root-share of the forest-floor respiration (Lagergren et al., 2005; Wide´n and Majdi, 2001). The forest-floor respiration was calculated from soil temperature at 5 cm from a relationship established in Norunda (Wide´n, 2001). 2.2. Model Input The model requires three input datasets: FAPAR (from NDVI), PAR and e.
P. Olofsson et al. / Advances in Space Research 39 (2007) 125–130
2.2.1. NDVI The raw NDVI data was acquired from the MODIS VI product, MOD13Q1. NDVI is calculated from MODIS surface reflectances that have been corrected for molecular scattering, ozone absorption, and aerosols (MODIS VI User Guide). The spatial resolution of the NDVI dataset is 250 m and data is composed over periods of 16 days, with the final version available from 24 February 2000. In order to cover Scandinavia, two MODIS tiles are required: h18v02 and h18v03, where h and v denote the horizontal and vertical tile number, respectively. Prior to the scaling of NDVI to FAPAR (see Section 2.2.2), the data was seasonally adjusted by nonlinear least-squares fitting of local asymmetric Gaussian model functions to the time-series, using the computer program TIMESAT (Jo¨nsson and Eklundh, 2002, 2004). The data points in the time-series were classified into six different quality levels, assigning each point a weight, w, between 1 and 0. The levels were based on the NDVI quality dataset included in the MODIS VI product. Using the weights, a measurement uncertainty was calculated and used when obtaining the parameters for the local functions. An example of the function fitting is displayed in Fig. 1 (since NDVI data for the second part of 2002 were lacking at the time of the analysis, corresponding data for 2000 were added to the 2002 dataset). 2.2.2. FAPAR Because of a linear relationship between NDVI and FAPAR (see e.g. Asrar et al., 1984; Goward and Huemmrich, 1992; Myneni and Williams, 1994), the FAPAR dataset was obtained from a linear scaling of the NDVI. Since we lack data regarding which percentiles to use for the scaling, the figures used were taken from Sellers et al. (1994): the 5th percentile for the yearly NDVI is assumed to represent no vegetation cover with a FAPAR value of approximately 0, and the 98th percentile assumed to represent a maximum vegetation cover with a FAPAR value of 0.95. The percentiles were calculated for five different land cover classes: (1) needleleaf
250
NDVI
200
tree cover, (2) broadleaf tree cover, (3) mixed tree cover, (4) open fields (pastures, etc.), (5) lichens and mosses. For the needleleaf FAPAR time series the percentiles used were calculated to 0.0165 and 0.8871, respectively. The calculations were made for 2001. Accordingly, a relationship between FAPAR and NDVI was obtained for each year and land cover class. The land cover was obtained from the Global Land Cover 2000 database (GLC2000) based on data from VEGETATION on SPOT-4 (Bartalev et al., 2003). 2.2.3. PAR A dataset of downward surface solar radiation for 2000 and 2001 was obtained from the SWECLIM project (Bergstro¨m et al., 2001). The data is generated every day at 37 · 35 points across northern Europe. In order to obtain PAR from solar radiation, the latter values were multiplied by 0.48 (Guyot, 1998). Fig. 2 shows the geographical coverage of the original point data draped on a map of the Nordic countries. 2.2.4. Light-use efficiency factor, e The LUE factor, e, was modeled using the equations developed by Lagergren et al. (2005) as a function of temperature and day of the year (DOY). The model deploys two sets of equations, one for spring and one for autumn. In both cases the model calculates a seasonal factor for each day describing the seasonal variation of the LUE. The seasonal factor was calibrated against measured values of e in order to obtain the e used in the NPP modeling. Accordingly, LUE is modeled with daily resolution. The date deciding when to deploy the spring and autumn equations, respectively, is here treated as a function of latitude so that spring equations are deployed earlier farther south (approximately 2 days per degree, with DOY = 222 at latitude 60). Accordingly, the model deployed here is also a function of latitude besides temperature and DOY. The temperature data used in the model was obtained from the SWECLIM project. In order to obtain sufficient spatial coverage, the epsilon values were interpolated in the same manner as the PAR data: kriging interpolation with a new set of interpolation parameters approximately every fiftieth day.
150
100
50
0 10
20
30
40
50
NDVI Gauss weight = 0.0 weight = 0.4 weight = 0.9 start and end 60
16 day periods, 2000–2002 Fig. 1. Raw and TIMESAT-adjusted NDVI at Flakaliden (64.07N, 19.27E).
127
Fig. 2. Geographical extent of the SWECLIM dataset.
P. Olofsson et al. / Advances in Space Research 39 (2007) 125–130
3. Results 3.1. FAPAR Figs. 3–5 show the time series of NDVI-derived FAPAR plotted together with measured FAPAR, and FAPAR from the MODIS product MOD15A2 between DOY 200 (19 July) and DOY 300 (27 October) 2001. Measurements were carried out in a 50 · 50 m area, the NDVI used for obtaining FAPAR has a spatial resolution of 250 m, and MOD15A2 a spatial resolution of 1 km; a comparison of three introduces a scaling error because of subpixel heterogeneity, making the measurement plot not fully representative for the satellite pixels of 250 and 1000 m. Still, the 250 m NDVI pixels are classified as boreal forest (by GLC2000), and looking at the plots it is obvious that the NDVI-derived FAPAR is stable throughout the measured FAPAR time series; and more importantly, is able to correctly estimate FAPAR (the agreement is not as good for the plot named L5 as for the other two, with a slight overestimation throughout the measurement period).
FAPAR: Mixed con., Skyttorp
128
100
80
60
40
20
0 200
Measured FAPAR Calculated FAPAR MOD15A2 220
240
260
280
300
DOY 2001 Fig. 5. Measured, modeled and MODIS-generated FAPAR at plot L5, Skyttorp, mixed Scots pine and Norway spruce, 60 y old stand.
FAPAR: Pine, Skyttorp
100
80
60
40
20
0 200
Fig. 6. Ten day averaged NPP at Norunda 2000 plotted with measured NPP, averaged 1997–1999. Measured FAPAR Calculated FAPAR MOD15A2 220
240
260
280
300
DOY 2001 Fig. 3. Measured, modeled and MODIS-generated FAPAR at plot L1, Skyttorp, Scots Pine, 30-year-old stand.
FAPAR: Spruce, Norunda
100
3.2. Net primary production
80
60
40
20
0 200
The MODIS-generated FAPAR, MOD15A2, performs well for the most part during the growing season, but fails certain days (e.g. during the period DOY 245–290 in plot L3). By consulting the quality control datasets attached in the product it is possible to derive the cause of failure. For example, the value for day 249 in L3 (Fig. 4) with an FAPAR of 0.09 suffers from cloud presence, while the value for day 265 in the same figure with an FAPAR of 0.08 is flagged as being derived by the backup algorithm.
Measured FAPAR Calculated FAPAR MOD15A2 220
240
260
280
300
DOY 2001 Fig. 4. Measured, modeled and MODIS-generated FAPAR at plot L3, Norunda, Norway Spruce.
NPP is calculated each day as a product of e and APAR. A daily resolution of the APAR datasets was obtained by letting each value in the FAPAR time series represent the following 16 days. Figs. 6 and 7 show modeled NPP at Norunda for 2000 and 2001, respectively, plotted with the mean of measured NPP for 1997-1999 (all averaged every 10th day). For daily data, the determination coefficients are R2 = 0.50 for both 2000 and 2001. If averaged every 10th day the determination coefficients increase to R2 = 0.81 for 2000, and R2 = 0.67 for 2001. The root mean squared error (RMSE)
P. Olofsson et al. / Advances in Space Research 39 (2007) 125–130
Fig. 7. Ten day averaged NPP at Norunda 2001 plotted with measured NPP, averaged 1997–1999.
between modeled and observed NPP for the averaged time series was calculated to 2.72 g C m2 d1 for 2000 and 3.57 g C m2 d1 for 2001. The total NPP at Norunda is 893 g C m2 for 2000, and 664 g C m2 for 2001. The figure for 2000 is slightly higher than the one reported by Lagergren et al. (2005) who estimated the total NPP at Norunda to, on average for 1997– 1999, 810 g C m2 y1. For 2001 an unrealistic dip for modeled NPP during spring 2001 lowers the modeled yearly production. The dip occurs around DOY 80-90 (late March), which is followed by a steep increase in production during the beginning of April. The reason for this is that APAR increases during spring while the e-values, due to low temperature, remain negative until the beginning of April (see Lagergren et al. (2005) for a discussion regarding the e-model). When parameterizing the model using data measured at Norunda, this dip is not present. A comparison between the temperature time series from SWECLIM, used in the e-model, and the measured time series reveals that the former underestimates the daily air temperature by about 1.7 C in average during the first 100 days in 2001. This phenomena is present in the data for 2000 but less pronounced (0.7 C underestimation). The PAR time series on the other hand exhibits a slight overestimation throughout both years, which in combination with the low temperatures is responsible for the dip in modeled NPP during spring. Since the modeled time series is evaluated against averaged data from preceding years the exact behavior of NPP during spring 2001 is unknown. However, a production decrease as seen in the modeled time series is unrealistic, although a slight decrease in production could be expected as a result of lower temperature. More validation is needed to make further conclusions about the performance of the model. 4. Discussion In summary, the presented NPP model with a day to day LUE variation, although rather simple and requiring just a
129
few input datasets, performs well when compared with measured data. The model is able to correctly trace the production decrease, especially during the second part of the year. As can be seen in Fig. 1, the NDVI time series are rather noisy due to the presence of clouds. Accordingly, a seasonal adjustment of the NDVI data by the use of TIMESAT, eliminating the cloud disturbance, is critical if the data is to be used for extracting other biophysical parameters such as FAPAR. Since the transformation of NDVI into FAPAR is linear, the cloud-induced distortions present in the NDVI data will be present also in the resulting FAPAR. This adjustment of the NDVI data is probably the primary reason for the good agreement between the obtained FAPAR and the measured FAPAR. Accordingly, it is suggested that a seasonal adjustment of the NDVI is performed prior to the transformation into FAPAR in areas suffering from cloud presence for long periods of time. In this study the obtained FAPAR is validated against data from only three plots, all located within the same area and with the same land cover. In order to conclude that the method is robust and fully functional, more measurements are needed in areas representing different forest types and throughout the year. Using single pixel values for central Sweden of MOD15A2 is hazardous and it is recommended to spatially average the data; or consulting the QC datasets. It might also be possible to perform the same kind of season adjustment for FAPAR as was done with NDVI. Acknowledgements Thanks to Ranga Myneni and Wenze Yang for assisting in the quality assessment of MOD15A2. The Terra/ MODIS data used in this study were provided by Earth Observing System Data and Information System (EOSDIS), Distributed Active Archive Center at Goddard Space Flight Center which archives, manages and distributes the data. The Global Land Cover 2000 database was provided by the European Commision Joint Research Centre. Incoming radiation and temperature data were obtained from the Swedish Regional Climate Modelling Programme. This study was financed by the Swedish National Space Board. References Asrar, G., Fuchs, M., Kanemasu, E.T., Hatfield, J.L. Estimating absorbed photosynthetic radiation and Leaf Area Index from spectral reflectance in wheat. Agron. J. 76, 300–306, 1984. Bartalev, S., Erchov, D., Isaev, A., Belward, A. A new SPOT4-Vegetation derived land cover map of Northern Eurasia. Int. J. Remote Sens. 24, 1977–1982, 2003. Bergstro¨m, S., Carlsson, B., Gardelin, M., Lindstro¨m, G., Pettersson, A., Rummukainen, M. Climate change impacts on runoff in Sweden – assessments by global climate models, dynamical downscaling and hydrological modelling. Clin. Res. 16, 101–112, 2001. Goward, S.N., Huemmrich, K.F. Vegetation canopy PAR absorptance and the normalized difference vegetation index: an assessment using the SAIL model. Remote Sens. Environ. 39, 119–140, 1992.
130
P. Olofsson et al. / Advances in Space Research 39 (2007) 125–130
Gower, S.T., Kucharik, C.J., Norman, J.M. Direct and indirect estimation of leaf area index, fAPAR, and net primary production of terrestrial ecosystems. Remote Sens. Environ. 70, 29–51, 1999. Guyot, G. Physics of the environment and climate. Wiley & Sons Ltd and Praxis Publishing Ltd (Wiley-Praxis series in atmospheric physics and climatology), Chichester, 1998. Jo¨nsson, P., Eklundh, L. Seasonally extraction by function fitting to timeseries of satellite sensor data. IEEE Trans. Geosci. Remote Sens. 40, 1832–1834, 2002. Jo¨nsson, P., Eklundh, L. TIMESAT – a program for analyzing time-series of satellite sensor data. Comput. Geosci. 30, 833–845, 2004. Lagergren, F., Eklundh, L., Grelle, A., Lundblad, M., Mo¨lder, M., Lankreijer, H., Lindroth, A. Net primary production and light use efficiency in a mixed coniferous forest in Sweden. Plant Cell Environ. 28, 412–423, 2005. Lundin, L.-C., Halldin, S., Lindroth, A., Cienciala, E., Grelle, A., Hjelm, P., Kellner, E., Lundberg, A., Mo¨lder, M., More´n, A.-S., Nord, T., Seibert, J., Sta¨hli, M. Continuous long-term measurements of soilplant-atmosphere variables at a forest site. Agric. For. Meteorol. 98– 99, 53–73, 1999. MODIS VI User Guide (last visited 25th May 2005) http://tbrs.arizona.edu/project/MODIS/UserGuide_doc.php. Monteith, J.L. Climate and the efficiency of crop production in Britain. Philos. Trans. R. Soc. Lond. Ser. B: Biol. Sci. 281, 277–294, 1977. Myneni, R.B., Williams, D.L. On the relationship between FAPAR and NDVI. Remote Sens. Environ. 49, 200–211, 1994. Prince, S.D. A model of regional primary production for use with coarse resolution satellite data. Int. J. Remote Sens. 12, 1313–1330, 1991.
Ruimy, A., Saugier, B., Dedieu, G. Methodology for the estimation of terrestrial net primary production from remotely sensed data. J. Geophys. Res. 99, 5263–5283, 1994. Rummukainen, M., Bergstro¨m, S., Ka¨llen, E., Moen, L., Rodhe, J., Tjernstro¨m, M. SWECLIM – The first three years. SMHI Reports Meteorology and Climatology, No. 94, December 2000. Russell, G., Jarvis, P.G., Monteith, J.L. Absorption of radiation by canopies and stand growth, in: Russell, G., Marshall, B., Jarvis, P.G. (Eds.), Plant Canopies: Their Growth, Form and Function. Cambridge University Press, Cambridge, 1989. Schlesinger, W.H. Biogeochemistry An Analysis of Global Change, second ed Academic Press, San Diego, 1991. Seaquist, J. Mapping primary production for the West African Sahel with satellite data. Meddelanden fra·n Lunds Universitets Geografiska Institutioner. Doctoral Thesis 140, Lund University, Lund, 2001. Sellers, P.J., Los, S.O., Tucker, C.J., Justice, C.O., Dazlich, D.A., Collatz, G.J., Randall, D.A. A global 1 by 1 NDVI data set for climatic studies, pt. 2: the adjustment of the NDVI and generation of global fields of terrestrial biophysical parameters. Int. J. Remote Sens. 17, 3519–3546, 1994. Tucker, C.J., Sellers, P.J. Satellite remote sensing of primary production. Int. J. Remote Sens. 7, 1395–1416, 1986. Wide´n, B. CO2 exchange within a Swedish coniferous forest, spatial and temporal variation. Swedish University of Agricultural Sciences, Uppsala, 2001. Wide´n, B., Majdi, H. Soil CO2 efflux and root respiration at three sites in a mixed pine and spruce forest: seasonal and diurnal variation. Can. J. For. Res. 31, 786–796, 2001.