Creating a regional MODIS satellite-driven net primary production dataset for European forests

remote sensing Article Creating a Regional MODIS Satellite-Driven Net Primary Production Dataset for European Forests Mathias Neumann 1,*, Adam Moreno 1, Christopher Thurnher 1, Volker Mues 2, Sanna Härkönen 3,4, Matteo Mura 5,6, Olivier Bouriaud 7, Mait Lang 8, Giuseppe Cardellini 9, Alain Thivolle-Cazat 10, Karol Bronisz 11, Jan Merganic 12, Iciar Alberdi 13, Rasmus Astrup 14, Frits Mohren 15, Maosheng Zhao 16 and Hubert Hasenauer 1 1Institute of Silviculture, Department of Forest and Soil Sciences, University of Natural Resources and Life Sciences, Vienna 1190, Austria; [email protected] (A.M.); [email protected] (C.T.); [email protected] (H.H.) 2Centre for Wood Science, World Forestry, University of Hamburg, Hamburg 21031, Germany; volker[email protected] 3 Department of Forest Sciences, University of Helsinki, Helsinki 00014, Finland; [email protected] 4Finnish Forest Research Institute, Joensuu 80101, Finland 5 Department of Bioscience and Territory, University of Molise, 86090 Pesche (IS), Italy; [email protected] 6geoLAB—Laboratory of Forest Geomatics, Department of Agricultural, Food and Forestry Systems, Università degli Studi di Firenze, Firenze 50145, Italy 7Facuty of Forestry, Universitatea Stefan del Mare, Suceava 720229, Romania; [email protected] 8Tartu Observatory, Tõravere 61602, Estonia; [email protected] 9Division Forest, Nature and Landscape, Department of Earth and Environmental Sciences, KU Leuven—University of Leuven, Leuven 3001, Belgium; [email protected] 10 Technological Institute, Furniture, Environment, Economy, Primary Processing and Supply, Champs sur Marne 77420, France; Alain.THIVOLLECAZA[email protected] 11 Laboratory of Dendrometry and Forest Productivity, Faculty of Forestry, Warsaw University of Life Sciences, Warsaw 02-776, Poland; karol.br[email protected].pl 12 Faculty of Forestry and Wood Sciences, Czech University of Life Sciences, Prague 16521, Czech Republic; [email protected] 13 Departamento de Selvicultura y Gestión de los Sistemas Forestales, INIA-CIFOR, Madrid 28040, Spain; [email protected] 14 Norwegian Institute for Bioeconomy Research, Ås 1431, Norway; rasmus.astr[email protected] 15 Forest Ecology and Forest Management Group, Wageningen University, Wageningen 6700, The Netherlands; [email protected] 16 Department of Geographical Sciences, University of Maryland, Collage Park, MD 20742, USA; [email protected] *Correspondence: [email protected]; Tel.: +43-1-47654-4059 Academic Editors: Lars T. Waser and Prasad S. Thenkabail Received: 5 April 2016; Accepted: 25 June 2016; Published: 29 June 2016 Abstract: Net primary production (NPP) is an important ecological metric for studying forest ecosystems and their carbon sequestration, for assessing the potential supply of food or timber and quantifying the impacts of climate change on ecosystems. The global MODIS NPP dataset using the MOD17 algorithm provides valuable information for monitoring NPP at 1-km resolution. Since coarse-resolution global climate data are used, the global dataset may contain uncertainties for Europe. We used a 1-km daily gridded European climate data set with the MOD17 algorithm to create the regional NPP dataset MODIS EURO. For evaluation of this new dataset, we compare MODIS EURO with terrestrial driven NPP from analyzing and harmonizing forest inventory data (NFI) from 196,434 plots in 12 European countries as well as the global MODIS NPP dataset for the years 2000 to 2012. Comparing these three NPP datasets, we found that the global MODIS NPP dataset differs from NFI NPP by 26%, while MODIS EURO only differs by 7%. MODIS EURO also agrees with NFI NPP across scales (from continental, regional to country) and gradients (elevation, location, tree age, dominant species, etc.). The agreement is particularly good for elevation, dominant Remote Sens. 2016,8, 554; doi:10.3390/rs8070554 www.mdpi.com/journal/remotesensing Remote Sens. 2016,8, 554 2 of 18 species or tree height. This suggests that using improved climate data allows the MOD17 algorithm to provide realistic NPP estimates for Europe. Local discrepancies between MODIS EURO and NFI NPP can be related to differences in stand density due to forest management and the national carbon estimation methods. With this study, we provide a consistent, temporally continuous and spatially explicit productivity dataset for the years 2000 to 2012 on a 1-km resolution, which can be used to assess climate change impacts on ecosystems or the potential biomass supply of the European forests for an increasing bio-based economy. MODIS EURO data are made freely available at ftp://palantir.boku.ac.at/Public/MODIS_EURO. Keywords: NPP; bioeconomy; forest inventory; NFI; climate; carbon; biomass; downscaling; increment; MOD17 1. Introduction Net primary production (NPP), the difference between Gross Primary Production (GPP) and plant autotrophic respiration, is the net carbon or biomass fixed by vegetation through photosynthesis. NPP represents the allocation rate of photosynthetic products into plant biomass and can be used to measure the quantity of goods provided to society by ecosystems [ 1 – 3 ]. NPP of forest ecosystems is essential to estimate the potential supply of biomass for bioenergy, fiber and timber supply. NPP is also a key variable to assess environmental change impacts on ecosystems [ 4 ] since any variation in the growing conditions influences the carbon cycle due to changes in carbon uptake and/or respiration. As interest grows in utilizing forests for a “bio-based economy” [ 5 , 6 ], more accurate and realistic forest productivity estimates become increasingly important. In addition, competing forest ecosystem services, such as biodiversity or and nature conservation, need to be considered to ensure sustainable use of our forests and to avoid unsustainable over-exploitation of renewable resources. Within the EU-28 160.9 million ha or 37.9% of the total land area are covered with forests [ 7 ]. These forests provide resources for the timber industry, the energy sector (24.3% of the energy in the EU-28 is generated from renewable sources of which 64.2% consists of forest biomass and waste [ 8 ]), but also for non-timber ecosystem services such as clean air, water, biodiversity or protection against natural hazards. Accurate and consistent forest information is a precondition for assessing the production and harvesting potential of forest resources in Europe. There are conceptually different data sources and methods to assess forest productivity like: (i) The MODIS algorithm MOD17 uses remotely sensed satellite-data and climate data to predict spatially and temporally continuous NPP and GPP (Gross Primary Production or carbon assimilation) based on an ecophysiological modelling approach [ 2 ]. In addition to satellite reflectance data and climate data, it requires the biophysical properties of land cover types, which are stored in the Biome Property Look-Up Tables (BPLUT) [9]. (ii) National forest inventory data can be used to assess the timber volume stocks as well as volume increment and removal, if repeated observations are available [ 10 ]. This terrestrial bottom-up approach collects forest information by measuring sample plots arranged on a systematic grid design across larger areas. In combination with biomass expansion factors or biomass functions, volume or tree information can be converted into biomass or carbon estimates to account for differences in wood densities, the carbon fraction and different allocation into compartments [11,12]. (iii) Flux towers record the gas-exchange in plant-atmosphere interactions [ 13 ], which can be used to derive GPP from Net Ecosystem exchange (NEE). NEE is estimated using eddy covariance data, climate measurements and other ancillary data [14]. Remote Sens. 2016,8, 554 3 of 18 Net Primary Production (NPP) from (i) top-down satellite-driven MOD17 algorithm and (ii) bottom-up NPP estimates using terrestrial forest inventory data were compared in a pilot study for Austria on national scale [ 15 ]. Top-down and bottom-up refer to the level of scaling of the primary recorded information (for MOD17 1-km remote sensing products and for Terrestrial NPP single tree observations). Our definition for top-down differs from traditional carbon cycle modelling [ 16 ]. This study wants to extend and test this concept for Europe on a continental scale. For this purpose, we obtain two wall-to-wall spatially-explicit and consistent MODIS NPP datasets by acquiring the global dataset using global climate driver and by creating a regional dataset MODIS EURO using 1-km European climate data. We evaluate these two datasets by comparing with the NPP derived from forest inventory data from 12 European countries. We assess the reliability and potential discrepancies of the MODIS satellite-driven top-down versus the terrestrial bottom-up NPP estimates from continental to national scale and across different gradients like location, elevation or stand density. This will provide a better understanding of the reliability of remote sensing based NPP estimates, which could be used also for regions, where no terrestrial measurements are available. 2. Materials and Methods We used two conceptually different methods to estimate NPP, (i) the MODIS NPP algorithm MOD17 and (ii) terrestrial forest inventory data and tree carbon estimation methods. Both have their respective strengths and weaknesses. MODIS NPP has the advantage of providing spatially continuous estimates with a consistent methodology, which is important for any large-scale studies. It incorporates biogeochemical principles in mechanistic modelling environment and the vegetation feedback to climate conditions through changes in Leaf Area Index and absorbed radiation [ 17 ]. It does not distinguish between different vegetation apart from general Land Cover types, has a coarse spatial resolution and might not be able to represent specific local conditions due to its calibration to global conditions. In contrast, terrestrial forest inventory NPP assesses the actual carbon allocation by trees and captures local small-scale effects (e.g., site conditions, tree age or forest management) as well as regional differences in estimating tree carbon [ 12 , 18 ]. It covers only the increment of trees assessed by the inventory system and might not capture local specifics of litter fall and fine root turnover very well, since broad model assumptions have to be used. 2.1. MODIS NPP Since the year 2000, the MOD17 product provides spatially and temporally continuous NPP estimates across the globe [ 17 ]. The algorithm behind uses the reflectance data from the sensor MODIS (MODerate resolution Imaging Spectroradiometer) of the TERRA and AQUA satellites operated by National Aeronautics and Space Administration of the United States (NASA). MOD17 provides GPP and NPP estimates at a 1-km resolution [ 2 , 17 ] and incorporates basic biogeochemical principles adopted from Biome-BGC [ 19 ]. It integrates a light use efficiency logic using remotely sensed vegetation information to estimate GPP (Equation (1)) with a maintenance and growth respiration module to derive NPP (Equation (2)). GPP “LUEmax ˆfTmin ˆfvpd ˆ0.45ˆSWrad ˆFPAR (1) NPP “GPP´RM´RG(2) LUEmax is the maximum light use efficiency, which get adjusted by f Tmin and f vpd to address water stress due to low temperature (Tmin) and vapor pressure deficit (VPD). SWrad is short wave solar radiation load, of which 45% is photosynthetically active. FPAR is the fraction of absorbed photosynthetic active radiation. R M is the maintenance respiration and is estimated using LAI (Leaf Area Index), climate data and biome-specific parameters. R G is the growth respiration and is estimated to be approx. 25% of NPP. The complete algorithm is documented in [ 18 ] and more details are found in the cited literature therein. Remote Sens. 2016,8, 554 4 of 18 The MOD17 algorithm requires climate data, FPAR and LAI (leaf area index) data as well as land cover data, which is derived from MODIS reflectance data [ 20 ]. We obtained the global MODIS NPP product (MOD17A3 Version 055) provided by the Numerical Terradynamic Simulation Group (NTSG) at University of Montana available at ftp://ftp.ntsg.umt.edu/pub/MODIS/NTSG_Products/. This data set (hereafter called MODIS GLOB) covers the period of 2000 to 2012, which is the time period covered by our terrestrial data (see next chapter), and provides the annual NPP in gC¨m´2¨year´1. The source of FPAR and LAI input is MODIS15 LAI/FPAR Collection 5, which was temporally gap filled to close data gaps due to unfavorable atmospheric conditions such as cloudiness or heavy aerosol presence [ 9 ]. For Land cover, we used the land cover product MOD12Q1 Version 4 Type 2 [ 21 ] representing the conditions in year 2001. Climate data are important input into the MODIS NPP algorithm and climate data have a strong impact on the MODIS NPP results [ 15 , 22 ]. MODIS GLOB uses the global climate data set NCEP2 [ 23 ] described in the following Section 2.2. In Europe, we have high quality daily climate data, the E-OBS data set [24], which was recently downscaled to a 1-km resolution [25]. We next ran the MOD17 algorithm with the downscaled European climate data [ 25 ] and obtained an additional MODIS NPP estimate for the period 2000–2012 (hereafter called MODIS EURO), which differ from MODIS GLOB provided by NTSG only in the used daily climate input data. We used the same FPAR, LAI and Land cover input, as used for the global NPP product, MODIS GLOB. MODIS EURO covers our study region, the EU-28 including Norway, Switzerland and the Balkan states (see Figure 1) and is made available under ftp://palantir.boku.ac.at/Public/MODIS_EURO. Remote Sens. 2016, 8, 554 4 of 18 The source of FPAR and LAI input is MODIS15 LAI/FPAR Collection 5, which was temporally gap filled to close data gaps due to unfavorable atmospheric conditions such as cloudiness or heavy aerosol presence [9]. For Land cover, we used the land cover product MOD12Q1 Version 4 Type 2 [21] representing the conditions in year 2001. Climate data are important input into the MODIS NPP algorithm and climate data have a strong impact on the MODIS NPP results [15,22]. MODIS GLOB uses the global climate data set NCEP2 [23] described in the following Section 2.2. In Europe, we have high quality daily climate data, the E-OBS data set [24], which was recently downscaled to a 1-km resolution [25]. We next ran the MOD17 algorithm with the downscaled European climate data [25] and obtained an additional MODIS NPP estimate for the period 2000–2012 (hereafter called MODIS EURO), which differ from MODIS GLOB provided by NTSG only in the used daily climate input data. We used the same FPAR, LAI and Land cover input, as used for the global NPP product, MODIS GLOB. MODIS EURO covers our study region, the EU-28 including Norway, Switzerland and the Balkan states (see Figure 1) and is made available under ftp://palantir.boku.ac.at/Public/MODIS_EURO. Figure 1. Our study region separated into four regions, countries with forest inventory data for estimating terrestrial National Forest Inventory (NFI) Net Primary Production (NPP) are marked with dots. 2.2. Climate Data As outlined, the two MODIS NPP estimates, MODIS GLOB and MODIS EURO, differ only in the daily climate data input: MODIS GLOB employs the global NCEP2 climate data set [23] and MODIS EURO uses European downscaled climate data [25]. We provide here a brief overview of the two climate data sets. The NCEP2 data set (NCEP-DOE Reanalysis 2) is a reanalyzed global daily climate data set with a spatial resolution of 1.875° × 1.875°. This corresponds to approx. 220 km at the equator at latitude 0° (approx. 136 × 220 km at latitude 50°). To compensate the coarse spatial resolution, for MODIS GLOB the climate data for the 1 km MODIS pixels was deduced with an bilateral interpolation method based on the neighboring NCEP2 pixels [9]. The downscaled climate data used for MODIS EURO provide daily climate data on a 0.0083° × 0.0083° resolution (approx. 1 × 1 km at the equator and approx. 0.6 × 1 km at 50° latitude) [25]. This data Figure 1. Our study region separated into four regions, countries with forest inventory data for estimating terrestrial National Forest Inventory (NFI) Net Primary Production (NPP) are marked with dots. 2.2. Climate Data As outlined, the two MODIS NPP estimates, MODIS GLOB and MODIS EURO, differ only in the daily climate data input: MODIS GLOB employs the global NCEP2 climate data set [ 23 ] and MODIS Remote Sens. 2016,8, 554 5 of 18 EURO uses European downscaled climate data [ 25 ]. We provide here a brief overview of the two climate data sets. The NCEP2 data set (NCEP-DOE Reanalysis 2) is a reanalyzed global daily climate data set with a spatial resolution of 1.875 ˝ˆ 1.875 ˝ . This corresponds to approx. 220 km at the equator at latitude 0 ˝ (approx. 136 ˆ 220 km at latitude 50 ˝ ). To compensate the coarse spatial resolution, for MODIS GLOB the climate data for the 1 km MODIS pixels was deduced with an bilateral interpolation method based on the neighboring NCEP2 pixels [9]. The downscaled climate data used for MODIS EURO provide daily climate data on a 0.0083˝ˆ0.0083˝ resolution (approx. 1 ˆ 1 km at the equator and approx. 0.6 ˆ 1 km at 50˝latitude) [25]. This data set was developed out of the E-OBS gridded climate data set (0.25˝resolution , using data from 7852 climate stations) [ 24 ] in conjunction with the WorldClim data set [26]. 2.3. Terrestrial NFI NPP Terrestrial forest data such as national forest inventory (NFI) data assess accumulated carbon on a systematic grid using a permanent plot design. From repeated observations of diameter at breast height (DBH) and/or tree height (H) in combination with biomass functions or biomass expansion factors the carbon accumulation of trees is estimated. Since this method is based on single tree measurements and local biomass studies, NPP derived from forest inventory data incorporates local effects such as weather patterns, climate anomalies, stand age, differences in biomass allocation, site and soil effects and different forest densities due to forest management [15,27]. We obtained 196,434 forest inventory plots covering 12 European countries. In Europe, each country has its own National Forest Inventory (NFI) system, which all have different measurement periods, sampling designs and methodologies [ 10 ] (Table S1 in the Supplementary Material). Thus, we first had to develop a harmonized and consistent terrestrial dataset for estimating Terrestrial NPP. We calculated NPP using the forest inventory data according to Equation (3). NPP “CARBINC `FRTO `CLF (3) CARB INC is the carbon increment of trees (gC ¨ m ´2¨ year ´1 ). FR TO is the carbon used for fine root turnover [ 28 , 29 ]. Fine root turnover FR TO is assumed to be equal to the carbon flow into litter C LF [ 27 , 30 ]. Both processes are controlled by the same factors and the assumption of similarity between the aboveand belowground turnover of short-living plant organs is supported by recently collected European data on fine root turnover [ 29 ] and litter fall [ 31 ]. C LF is the flow of carbon into litter (gC ¨ m ´2¨ year ´1 ) estimated using a climate-sensitive and species-dependent model [ 31 ] and is calculated as: Broadleaf-dominated : CLF “CFexpp2.643 `0.726LnpT`10q ` 0.181LnpPqq (4) Coniferous-dominated : CLF “CFexpp2.708 `0.505LnpT`10q ` 0.240LnpPqq (5) CF is the carbon fraction of dry biomass which is set equal to 0.5 [ 11 ]. Tis the mean annual temperature from the year 2000 to 2012 ( ˝ C). Pis the mean annual precipitation 2000 to 2012 [mm]. For temperature and precipitation we use the European climate data [ 25 ] to capture important small-scale regional effects such as elevation or topography in a more realistic way. Equation (4) is applied for all plots where broadleaf species contribute most to total basal area and Equation (5) is used for coniferous-dominated plots (see Table S2 of the Supplementary Material). We used data from nine National Forest Inventories (Austria, Czech Republic, Germany, France, Finland, Norway, Poland, Romania, Spain), and three Regional Forest Inventories (Belgium, Estonia, Italy). We grouped our 12 countries in four geographic regions, North Europe, Central-West Europe, Central-East Europe and South Europe [ 7 ], to address the large environmental, elevational and climatic Remote Sens. 2016,8, 554 6 of 18 gradients in Europe. Countries within a region should have similar climatic and edaphic conditions as well as similar tree allometries and allocation patterns [32]. The original locations of the inventory plots were falsified to the nearest pixel of the MODIS grid to guarantee the locations of the plots remain unknown. Temporal consistency with the MODIS data (available since year 2000) was ensured by using only inventory data, which provide CARB INC (Equation (3)) for the time period 2000 to 2012. Figure 1shows our study region with the four geographic regions completely covered by MODIS EURO, and the 12 countries, where we have NFI NPP. Although all our terrestrial forest inventory data assess properties of trees, there are different sampling methods and increment calculation by country in place, which may strongly affect the resulting estimates [ 33 , 34 ]. Four different methods to estimate tree carbon increment CARB INC are used in our data: (1) repeated observations of fixed area plots (used in Norway, Poland, Belgium); (2) repeated angle count sampling (for Austria, Germany, Finland); (3) increment cores (France, Romania, Italy); as well as increment predictions from (4) tree growth models (Czech Republic, Estonia, Italy). Tree growth model predictions were used if no increment observations, either from repeated observations or from increment cores, were available. In the Supplementary Material, we provide all details for our 12 inventory data sets, the local sampling system, the available data and the used increment method (Table S1 in the Supplementary Material). The tree carbon results for determining carbon increment CARB INC (Equation (3)) were estimated using the carbon calculation method applied by the local forest inventory organization and compiled in [ 32 ]. Local biomass functions and biomass expansion factors were used to derive tree biomass and carbon fractions to convert biomass into carbon. In the Supplementary Material, we provide a detailed description on processing the NFI data, the tree carbon estimates and stand variables to describe the represented forests (e.g., mean age, basal area or stand density index). Using this methodology, we processed the forest inventory data from the 12 countries (Table S1) and derived harmonized carbon stocks for all inventory plots. The forest inventory data set consists of 196.434 plots, harmonized across 12 European countries. We applied the carbon increment method for each country and calculated NPP by inventory plot (hereafter called NFI NPP) using Equations (3)–(5). 2.4. Analysis of NPP Results We thus have three NPP sources: two using the MOD17 algorithm with different daily climate data: (i) MODIS GLOB produced by the Numerical Terradynamic Simulation Group (NTSG) at University of Montana and (ii) MODIS EURO by running the original MOD17 algorithm and the latest BPLUTs parametrized by [ 9 ] with downscaled daily climate data from Europe [ 25 ] as well as (iii) Terrestrial NFI NPP using forest inventory data from the 12 countries (Table S1) and local carbon estimation methods [32]. We compared the three NPP datasets across Europe, by our 4 regions (Figure 1) and the 12 countries to analyze our results across different spatial scaling. We extracted for each forest inventory plot at the corresponding MODIS cell the average NPP from MODIS GLOB and MODIS EURO for 2000 to 2012. We next computed for all plots the difference between the two MODIS NPP estimates and the Terrestrial NFI NPP ( ∆ NPP GLOB = MODIS GLOB minus NFI NPP and ∆ NPP EURO = MODIS EURO minus NFI NPP). We used each NFI plot separately and did not compute average values for MODIS pixels. This avoided smoothing effects due to different spacing between inventory grid points and the plot clusters used in some countries (Table S1). To analyze the effect of gradients on the NPP results, we collected potentially meaningful meta-information such as plot location (Longitude and Latitude in WGS1984), Elevation (EU-DEM 30 m resolution), MODIS Land Cover type or forest characteristics (dominant tree species, mean age, stand density, tree height, etc.) and analyzed patterns of ∆ NPP GLOB and ∆ NPP EURO across these gradients. Remote Sens. 2016,8, 554 7 of 18 Terrestrial and remote sensing NPP estimates exhibited discrepancies in previous research [ 15 , 18 ] and as explanation the authors suggested changes in stand density, which are commonly caused by forest management and disturbances [ 15 , 18 ]. Since major parts of the forests in Europe are managed [ 7 ] and affected by natural disturbances such as wind damage or forest fire [ 35 ], they should have experienced changes in stand density as compared to unmanaged forests. Stand density directly affects terrestrial NPP estimates by its impact on the development of DBH and H of the remaining trees after forest management operations until canopy closure is reached. On the other hand MODIS NPP is based on the “big leaf” concept and assumes a full coverage of forest area. We thus use Stand density index (SDI) [36] in the analysis of our NPP estimates. 3. Results NPP estimated using the MOD17 algorithm has the advantage of providing spatialand temporal-continuous NPP estimates across Europe on a 1-km resolution and Figure 2illustrates this by showing MODIS EURO for the years 2000 to 2012. Note that MODIS EURO also covers not-forest land cover types such as crops, shrubor grassland. Remote Sens. 2016, 8, 554 8 of 18 Remote Sens. 2016, 8, x; doi:10.3390/ www.mdpi.com/journal/remotesensing Figure 2. MODIS EURO NPP on 1-km resolution representing average NPP for the period 2000–2012 using European daily climate data (available under ftp://palantir.boku.ac.at/Public/MODIS_EURO). Our NFI dataset covers the full elevational and latitudinal range of forest conditions in Europe including different site conditions, tree species, development stages or management practices. For most countries we have more than 5000 inventory plots (exception: Belgium with 512 plots) and in most cases a plot spacing of at least 4 by 4 km (Table S1). This dataset also provides information on forest properties such as tree age, carbon stocks or stand density and Table 2 indicates that these characteristics vary across Europe. Table 2. NPP and ∆NPP (always using median) for the whole dataset (“All Countries”), for each country separately and for each region (MODIS NPP using global climate data—MODIS GLOB; MODIS NPP using local European climate data—MODIS EURO and NPP using forest inventory data—NFI NPP); ∆NPP and Rel. ∆NPP both for MODIS GLOB and MODIS EURO. Positive differences indicate that MODIS NPP overestimates NFI NPP and vice versa. NPP and ∆NPP (gC·m−2·year−1) MODIS MODIS ∆NPP Rel. ∆NPP [%] GLOB EURO NFI NPP GLOB EURO GLOB EURO All Countries 680 577 539 141 38 26% 7% North Europe Finland 471 399 414 57 −15 14% −4% Norway 484 406 409 75 −3 18% −1% Estonia 534 504 492 42 12 9% 3% all 519 479 461 58 18 13% 4% Central-West Europe Austria 739 612 634 105 −22 17% −4% Belgium 732 599 644 88 −45 14% −7% France 787 666 604 183 62 30% 10% Germany 692 602 716 −24 −114 −3% −16% all 759 645 615 144 30 23% 5% Central-East Europe Czech Republic 696 618 553 143 65 26% 12% Poland 641 571 659 −19 −88 −3% −13% Romania 713 562 565 148 −3 26% −1% all 677 592 595 82 −3 14% −1% Figure 2. MODIS EURO NPP on 1-km resolution representing average NPP for the period 2000–2012 using European daily climate data (available under ftp://palantir.boku.ac.at/Public/MODIS_EURO). Terrestrial NFI NPP is driven by forest information collected by field crews. Thus it provides NPP and the carbon accumulation by forest stands during a certain time period. Table 1gives a summary of the forest inventory results by country, by region and the whole dataset, with the terrestrial NFI NPP at the right side. Remote Sens. 2016,8, 554 8 of 18 Table 1. Summary of the forest inventory results: Number of plots with data, Time period covered by NFI NPP, Mean elevation (range Minimum–Maximum) in meter above sea level (EU-DEM 30 m resolution). For the following plot statistics we provide mean and standard deviation: Mean quadratic DBH (cm), Mean Tree height (m), Basal area at 1.3 m height (m 2¨ ha ´1 ), Stem number (ha ´1 ), Tree carbon per hectare (gC ¨ m ´2 ), Median age class, SDI Stand Density Index [ 36 ] (for details on this variables see Supplementary Material), NPP is the NFI Net primary production (gC ¨ m ´2¨ year ´1 ) according to Equation (3), For Czech Republic we only have country means. Empty cells (-) indicate that this variable is not available from the NFI data set. At the end of each section, statistics of the region are given and at the bottom of the table summary statistics for whole Europe. Region Country Number of Plots Time Period Mean Elevation (min–max) (m) Mean DBH (cm) Mean Tree Height (m) Basal Area (m2¨ha´1) Stem Number (ha´1) Tree Carbon (gC¨m´2) Median Age (Years) SDI NPP (gC¨m´2¨year´1) North Europe Estonia 19930 2000–2010 66 (2–275) 17 ˘8 17 ˘7 19 ˘8 1540 ˘2554 5240 ˘2929 40–60 449 ˘192 509 ˘163 Finland 6442 2000–2008 141 (1–400) 18 ˘7 14 ˘5 18 ˘8 3522 ˘13251 4859 ˘3020 40–60 400 ˘236 446 ˘173 Norway 9562 2000–2009 391 (0–1253) 15 ˘6 9 ˘3 15 ˘12 930 ˘682 4003 ˘3691 60–80 368 ˘265 442 ˘143 all 35379 2000–2010 161 (0–1253) 16 ˘7 14 ˘7 18 ˘9 1736 ˘5983 4856 ˘3199 40–60 419 ˘224 482 ˘162 Central-West Europe Austria 9562 2000–2009 912 (113–2299) 32 ˘14 21 ˘7 32 ˘19 987 ˘1070 10364 ˘6973 60–80 688 ˘396 681 ˘251 Belgium 512 2009–2013 39 (2–278) 29 ˘12 18 ˘6 30 ˘13 660 ˘446 11507 ˘6475 40–60 648 ˘279 671 ˘195 France 33152 2001–2011 444 (0–2707) 23 ˘11 15 ˘7 23 ˘15 778 ˘602 8083 ˘6457 60–80 512 ˘298 649 ˘254 Germany 5894 2000–2008 344 (´5–1879) 28 ˘12 22 ˘7 31 ˘14 833 ˘814 11811 ˘6371 60–80 628 ˘302 754 ˘185 all 49120 2000–2013 514 (´5–2707) 25 ˘12 17 ˘8 25 ˘17 824 ˘749 9034 ˘6698 60–80 564 ˘328 667 ˘253 Central-East Europe Czech Rep. 13929 2001–2004 541 (138–1503) 25 20 33 812 17340 ˘10858 60–80 809 ˘441 643 ˘266 Poland 17281 2005–2013 193 (´4–1459) 23 ˘9 18 ˘5 29 ˘14 883 ˘614 10656 ˘6623 40–60 612 ˘263 720 ˘288 Romania 5509 2003–2011 542 (´1–1968) 24 ˘11 - 28 ˘15 878 ˘723 10355 ˘7256 40–60 582 ˘289 571 ˘164 all 36719 2001–2013 443 (´4–1968) 23 ˘10 18 ˘5 28 ˘15 881 ˘673 12376 ˘8793 40–60 652 ˘345 649 ˘248 South Europe Italy 15183 2002–2009 860 (7–2891) 20 ˘8 12 ˘4 22 ˘13 839 ˘636 6315 ˘4897 20–40 497 ˘293 635 ˘179 Spain 60033 2000–2008 842 (1–2549) 23 ˘13 10 ˘4 13 ˘11 491 ˘516 4003 ˘3918 40–60 288 ˘246 606 ˘293 all 75216 2000–2009 831 (1–2891) 22 ˘12 10 ˘4 15 ˘12 561 ˘560 4469 ˘4237 40–60 330 ˘269 578 ˘275 All countries - 196434 – 548 (´5–2891) 22 ˘11 13 ˘7 20 ˘15 900 ˘2646 7298 ˘6916 40–60 469 ˘325 597 ˘252 Remote Sens. 2016,8, 554 9 of 18 Our NFI dataset covers the full elevational and latitudinal range of forest conditions in Europe including different site conditions, tree species, development stages or management practices. For most countries we have more than 5000 inventory plots (exception: Belgium with 512 plots) and in most cases a plot spacing of at least 4 by 4 km (Table S1). This dataset also provides information on forest properties such as tree age, carbon stocks or stand density and Table 2indicates that these characteristics vary across Europe. Table 2. NPP and ∆ NPP (always using median) for the whole dataset (“All Countries”), for each country separately and for each region (MODIS NPP using global climate data—MODIS GLOB; MODIS NPP using local European climate data—MODIS EURO and NPP using forest inventory data—NFI NPP); ∆ NPP and Rel. ∆ NPP both for MODIS GLOB and MODIS EURO. Positive differences indicate that MODIS NPP overestimates NFI NPP and vice versa. NPP and ∆NPP (gC¨m´2¨year´1)MODIS MODIS ∆NPP Rel. ∆NPP [%] GLOB EURO NFI NPP GLOB EURO GLOB EURO All Countries 680 577 539 141 38 26% 7% North Europe Finland 471 399 414 57 ´15 14% ´4% Norway 484 406 409 75 ´3 18% ´1% Estonia 534 504 492 42 12 9% 3% all 519 479 461 58 18 13% 4% Central-West Europe Austria 739 612 634 105 ´22 17% ´4% Belgium 732 599 644 88 ´45 14% ´7% France 787 666 604 183 62 30% 10% Germany 692 602 716 ´24 ´114 ´3% ´16% all 759 645 615 144 30 23% 5% Central-East Europe Czech Republic 696 618 553 143 65 26% 12% Poland 641 571 659 ´19 ´88 ´3% ´13% Romania 713 562 565 148 ´3 26% ´1% all 677 592 595 82 ´3 14% ´1% South Europe Italy 862 657 635 227 22 36% 4% Spain 632 555 503 129 52 26% 10% all 691 584 519 172 65 33% 13% 3.1. NPP Estimates across Different Scales Comparing all our three NPP estimates on a European scale allowed us to explore the general behaviour and evaluate the agreement of the two remote sensing driven NPP products, MODIS GLOB and MODIS EURO, with the terrestrial driven NFI NPP estimates (Figure 3). Re-running the MOD17 algorithm with local climate data reduced the remotely sensed MODIS NPP in terms of median, mean and variation as compared to the global climate driver (Figure 3). NFI NPP is close to MODIS EURO regarding median and mean, but show larger variation. In addition, Figure 3confirms that our data is clearly right-skewed (NFI NPP in particular). Zooming in and examining the different NPP estimates by ecoregion and country allowed us to analyze our results on a higher spatial resolution and to assess local effects such as different regional growing conditions, the impact of local biomass allometries or tree species composition [ 32 ] as well as the potential effect of different forest management practices in Europe [7]. We provide in Table 2the median NPP for the three NPP sources (MODIS GLOB, MODIS EURO and NFI NPP) and the differences between MODIS and NFI NPP ( ∆ NPP GLOB and ∆ NPP EURO ), both in absolute values in gC ¨ m ´2¨ year ´1 and normalized in relation to NFI NPP (Rel. ∆ NPPi in %). Results are given in Table 2for Europe, by country and for the four eco-regions [7]. At the European level, the MODIS GLOB gives an NPP of 680 gC ¨ m ´2¨ year ´1 , the MODIS EURO resulted in 577 gC ¨ m ´2¨ year ´1 , and the NPP from the NFI data exhibit a value of 539 gC ¨ m ´2¨ year ´1 . The differences in NPP ( ∆ NPP GLOB ) using the global dataset MODIS GLOB are larger than ∆ NPP EURO using the regional dataset MODIS EURO (+26% vs. +7%). The same pattern is evident across all four Remote Sens. 2016,8, 554 16 of 18 measures for European forests. Since the lifetime of the satellites carrying the MODIS sensor is unknown, we strongly suggest the implementation and testing of this concept in the upcoming European satellite technologies such as the Copernicus Programme to ensure consistent and realistic productivity estimates also in the future. MODIS EURO data are made freely available for 2000 until 2012 under ftp://palantir.boku.ac.at/ Public/MODIS_EURO. Supplementary Materials: The following are available online at www.mdpi.com/2072-4292/8/7/554/s1, Table S1: Summary of the properties of the different forest inventory datasets, Table S2: Tree species groups used in this study, description and selected tree species, Figure S1: Direct pixel-to-plot comparison of MODIS EURO and NFI NPP, Figure S2: For North Europe ∆ NPP grouped by Elevation, Latitude and Longitude, Figure S3: For Central-West Europe ∆ NPP grouped by Elevation, Latitude and Longitude, Figure S4: For Central-East Europe ∆ NPP grouped by Elevation, Latitude and Longitude, Figure S5: For South Europe ∆ NPP grouped by Elevation, Latitude and Longitude, Figure S6: Difference ∆ NPP grouped by age classes, Figure S7: Difference ∆ NPP grouped by tree height classes, Figure S8: Difference ∆ NPP grouped by MODIS Land cover types, Figure S9: Difference ∆ NPP grouped by dominant species, Figure S10: MODIS EURO and NFI NPP by Stand density Index (SDI) classes. Acknowledgments: This work was conducted as part of the collaborative project “FORest management strategies to enhance the MITigation potential of European forests” (FORMIT). The research leading to these results has received funding from the European Union Seventh Framework Programme under grant agreement n ˝ 311970. Special thanks to all the field crews collecting the sample data for the inventory plots. We are also grateful to the responsible people from the various forest inventory organizations for providing us with the data of their forest inventory systems and thereby making this work possible in the first place. We also want to acknowledge the open data policy of NASA and the work of the MODIS land product science team providing us with the FPAR and LAI products. We further want to thank Loretta Moreno for proof reading the manuscript. We want to thank in particular the editor and the anonymous reviewer on their helpful comments on an earlier draft of the manuscript. Author Contributions: M.N. conceived and designed the study, coordinated compiling the NFI NPP dataset, calculated the forest inventory results for Austria and wrote the first draft of the manuscript, A.M. and C.T. developed the code for computing MODIS EURO and maintain the ftp-server, V.M. calculated the forest inventory results for Germany, S.H. for Finland, M.M. for Italy, O.B. for Romania, M.L. for Estonia, G.C. for Belgium, A.T. for France, K.B. for Poland, J.M. for Czech Republic, I.A. for Spain, R.A. for Norway, M.Z. provided the original MOD17 code and helped in preparation of the input data, F.M. and H.H. coordinated and supervised the analysis and the manuscript writing, all authors contributed equally in writing and revising the manuscript. Conflicts of Interest: The authors declare no conflict of interest. References 1. Gower, S.T.; Kucharik, C.J.; Norman, J.M. Direct and Indirect Estimation of Leaf Area Index, f APAR, and Net Primary Production of Terrestrial Ecosystems. Remote Sens. Environ. 1999,70, 29–51. [CrossRef] 2. Running, S.W.; Nemani, R.R.; Heinsch, F.A.; Zhao, M.; Reeves, M.; Hashimoto, H. A Continuous Satellite-Derived Measure of Global Terrestrial Primary Production. Bioscience 2004 ,54, 547–560. [CrossRef] 3. Imhoff, M.L.; Bounoua, L.; DeFries, R.; Lawrence, W.T.; Stutzer, D.; Tucker, C.J.; Ricketts, T. The consequences of urban land transformation on net primary productivity in the United States. Remote Sens. Environ. 2004 , 89, 434–443. [CrossRef] 4. Friedlingstein, P. A steep road to climate stabilization. Nature 2008,451, 297–298. [CrossRef] [PubMed] 5. Enriquez, J. GENOMICS: Genomics and the World’s Economy. Science 1998 ,281, 925–926. [CrossRef] [PubMed] 6. McCormick, K.; Kautto, N. The Bioeconomy in Europe: An Overview. Sustainability 2013 ,5, 2589–2608. [CrossRef] 7. FOREST EUROPE State of Europe’s Forests 2015. Available online: http://www.foresteurope.org/docs/ fullsoef2015.pdf (accessed on 27 June 2016). 8. EUROSTAT Europe in Figures—Eurostat Yearbook. Available online: http://ec.europa.eu/eurostat/ (accessed on 26 February 2016). 9. Zhao, M.; Heinsch, F.A.; Nemani, R.R.; Running, S.W. Improvements of the MODIS terrestrial gross and net primary production global data set. Remote Sens. Environ. 2005,95, 164–176. [CrossRef] 10. Tomppo, E.; Gschwantner, T.; Lawrence, M.; McRoberts, R. National Forest Inventories: Pathways for Common Reporting; Springer: Berlin, Germany, 2010; p. 610. Remote Sens. 2016,8, 554 17 of 18 11. IPCC Good Practice Guidance for Land Use, Land-Use Change and Forestry, CHAPTER 4 FOREST. In 2006 IPCC Guidelines for National Greenhouse Gas Inventories, Prepared by the National Greenhouse Gas Inventories Programme; Eggleston, H.S., Buendia, L., Miwa, K., Ngara, T., Tanabe, K., Eds.; IGES: Kanagawa, Japan, 2006; p. 83. 12. Zianis, D.; Muukkonen, P.; Mäkipää, R.; Mencuccini, M. Biomass and stem volume equations for tree species in Europe. Silva Fenn. Monogr. 2005,4, 63. 13. Baldocchi, D. “Breathing” of the terrestrial biosphere: Lessons learned from a global network of carbon dioxide flux measurement systems. Aust. J. Bot. 2008,56, 1–26. [CrossRef] 14. Jung, M.; Reichstein, M.; Margolis, H.A.; Cescatti, A.; Richardson, A.D.; Arain, M.A.; Arneth, A.; Bernhofer, C.; Bonal, D.; Chen, J.; et al. Global patterns of land-atmosphere fluxes of carbon dioxide, latent heat, and sensible heat derived from eddy covariance, satellite, and meteorological observations. J. Geophys. Res. 2011 , 116, G00J07. [CrossRef] 15. Neumann, M.; Zhao, M.; Kindermann, G.; Hasenauer, H. Comparing MODIS Net Primary Production Estimates with Terrestrial National Forest Inventory Data in Austria. Remote Sens. 2015 ,7, 3878–3906. [CrossRef] 16. Jiang, F.; Chen, J.M.; Zhou, L.; Ju, W.; Zhang, H.; Machida, T.; Ciais, P.; Peters, W.; Wang, H.; Chen, B.; et al. A comprehensive estimate of recent carbon sinks in China using both top-down and bottom-up approaches. Sci. Rep. 2016,6. [CrossRef] [PubMed] 17. Zhao, M.; Running, S.W. Drought-induced reduction in global terrestrial net primary production from 2000 through 2009. Science 2010,329, 940–943. [CrossRef] [PubMed] 18. Hasenauer, H.; Petritsch, R.; Zhao, M.; Boisvenue, C.; Running, S.W. Reconciling satellite with ground data to estimate forest productivity at national scales. For. Ecol. Manag. 2012,276, 196–208. [CrossRef] 19. Thornton, P.E. Description of a Numerical Simulation Model for Predicting the Dynamics of Energy, Water, Carbon and Nitrogen in a Terrestrial Ecosystem; University of Montana: Missoula, MT, USA, 1998; p. 280. 20. Yang, W.; Huang, D.; Tan, B.; Stroeve, J.C.; Shabanov, N.V.; Knyazikhin, Y.; Nemani, R.R.; Myneni, R.B. Analysis of Leaf Area Index and Fraction of PAR absorbed by Vegetation Products from the Terra MODIS sensor: 2000–2005. IEEE Trans. Geosci. Remote Sens. 2006,44, 1829–1842. [CrossRef] 21. Friedl, M.; Sulla-Menashe, D.; Tan, B.; Schneider, A.; Ramankutty, N.; Sibley, A.; Huang, X. MODIS Collection 5 global land cover: Algorithm refinements and characterization of new datasets. Remote Sens. Environ. 2010 , 114, 168–182. [CrossRef] 22. Zhao, M.; Running, S.W.; Nemani, R.R. Sensitivity of Moderate Resolution Imaging Spectroradiometer (MODIS) terrestrial primary production to the accuracy of meteorological reanalyses. J. Geophys. Res. 2006 , 111, G01002. [CrossRef] 23. Kanamitsu, M.; Ebisuzaki, W.; Woollen, J.; Yang, S.-K.; Hnilo, J.J.; Fiorino, M.; Potter, G.L. NCEP–DOE AMIP-II Reanalysis (R-2). Bull. Am. Meteorol. Soc. 2002,83, 1631–1643. [CrossRef] 24. Haylock, M.R.; Hofstra, N.; Klein Tank, A.M.G.; Klok, E.J.; Jones, P.D.; New, M. A European daily high-resolution gridded data set of surface temperature and precipitation for 1950–2006. J. Geophys. Res. 2008,113, D20119. [CrossRef] 25. Moreno, A.; Hasenauer, H. Spatial downscaling of European climate data. Int. J. Climatol. 2015. [CrossRef] 26. Hijmans, R.J.; Cameron, S.E.; Parra, J.L.; Jones, P.G.; Jarvis, A. Very high resolution interpolated climate surfaces for global land areas. Int. J. Climatol. 2005,25, 1965–1978. [CrossRef] 27. Olson, R.J.; Johnson, K.R.; Zheng, D.L.; Scurlock, J.M.O. Global and Regional Ecosystem Modeling: Databases of Model Drivers and Validation Measurements; Oakridge Laboratory: Oak Ridge, TN, USA, 2001; p. 95. 28. Jackson, R.; Mooney, H.; Schulze, E. A global budget for fine root biomass, surface area, and nutrient contents. Proc. Natl. Acad. Sci. USA 1997,94, 7362–7366. [CrossRef] [PubMed] 29. Finér, L.; Ohashi, M.; Noguchi, K.; Hirano, Y. Fine root production and turnover in forest ecosystems in relation to stand and environmental characteristics. For. Ecol. Manag. 2011,262, 2008–2023. [CrossRef] 30. Raich, J.; Nadelhoffer, K. Belowground carbon allocation in forest ecosystems: Global trends. Ecology 1989 , 70, 1346–1354. [CrossRef] 31. Liu, C.; Westman, C.J.; Berg, B.; Kutsch, W.; Wang, G.Z.; Man, R.; Ilvesniemi, H. Variation in litterfall-climate relationships between coniferous and broadleaf forests in Eurasia. Glob. Ecol. Biogeogr. 2004 ,13, 105–114. [CrossRef] Remote Sens. 2016,8, 554 18 of 18 32. Neumann, M.; Moreno, A.; Mues, V.; Härkönen, S.; Mura, M.; Bouriaud, O.; Lang, M.; Achten, W.M.J.; Thivolle-Cazat, A.; Bronisz, K.; et al. Comparison of carbon estimation methods for European forests. For. Ecol. Manag. 2016,361, 397–420. [CrossRef] 33. Hasenauer, H.; Eastaugh, C.S. Assessing Forest Production Using Terrestrial Monitoring Data. Int. J. For. Res. 2012. [CrossRef] 34. Huber, M.O.; Eastaugh, C.S.; Gschwantner, T.; Hasenauer, H.; Kindermann, G.; Ledermann, T.; Lexer, M.J.; Rammer, W.; Schörghuber, S.; Sterba, H. Comparing simulations of three conceptually different forest models with National Forest Inventory data. Environ. Model. Softw. 2013,40, 88–97. [CrossRef] 35. Seidl, R.; Schelhaas, M.-J.; Rammer, W.; Verkerk, P.J. Increasing forest disturbances in Europe and their impact on carbon storage. Nat. Clim. Chang. 2014,4, 806–810. [CrossRef] [PubMed] 36. Reineke, L. Perfecting a stand-density index for even-aged forests. J. Agric. Res. 1933,46, 627–638. 37. Turner, D.P.; Ritts, W.D.; Cohen, W.B.; Gower, S.T.; Running, S.W.; Zhao, M.; Costa, M.H.; Kirschbaum, A.A.; Ham, J.M.; Saleska, S.R.; et al. Evaluation of MODIS NPP and GPP products across multiple biomes. Remote Sens. Environ. 2006,102, 282–292. [CrossRef] 38. MacLean, R.G.; Ducey, M.J.; Hoover, C.M. A Comparison of Carbon Stock Estimates and Projections for the Northeastern United States. For. Sci. 2014,60, 206–213. [CrossRef] 39. Austrias National Inventory Report 2012: Umweltbundesamt Submission under the United Nations Framework Convention on Climate Change and the Kyoto Protocol 2012 National Inventory Report for the German Greenhouse Gas Inventory 1990–2010; Umweltbundesam: Dessau, Germany, 2012; p. 827. 40. Kändler, G.; Bösch, B. Methodenentwicklung für die 3. Bundeswaldinventur: Modul 3 Überprüfung und Neukonzeption einer Biomassefunktion—Abschlussbericht. Tech. Rep. FVA-BW 2013, 1–72. (In German) 41. Neubauer, M.; Demant, B.; Bolte, A. Tree-based estimator for below-ground biomass of Pinus sylvestris L. Forstarchiv 2015,86, 42–47. 42. Umweltbundesamt Berichterstattung unter der Klimarahmenkonvention der Vereinten Nationen und dem KyotoProtokoll 2014 Nationaler Inventarbericht zum Deutschen Treibhausgasinventar 1990–2012. Clim. Chang. 2014,24, 965. (In German) 43. He, L.; Chen, J.M.; Pan, Y.; Birdsey, R.; Kattge, J. Relationships between net primary productivity and forest stand age in U.S. forests. Glob. Biogeochem. Cycles 2012,26, 1–19. 44. Moreno, A.; Neumann, M.; Hasenauer, H. Optimal Resolution for Linking Remotely Sensed and Forest Inventory Data in Europe. Remote Sens. Environ. 2016,183, 109–119. [CrossRef] 45. Härkönen, S.; Lehtonen, A.; Manninen, T.; Tuominen, S.; Peltoniemi, M. Estimating forest leaf area index using satellite images: Comparison of k-NN based Landsat-NFI LAI with MODISRSR based LAI product for Finland. Boreal Environ. Res. 2015,6095, 181–195. © 2016 by the authors; licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC-BY) license (http://creativecommons.org/licenses/by/4.0/).