Spatial Representativeness of Gross Primary Productivity from Carbon Flux Sites in the Heihe River Basin, China

Abstract

Studying the spatial representativeness of carbon flux measurement data for typical land cover types can provide important information for benchmarking Earth system models and validating multiple-scale remote sensing products. In our study, daily gross primary productivity (GPP) was firstly derived from eddy covariance observation systems and seasonal variations in field GPP were analyzed at nine flux tower sites for typical land cover types in the Heihe River Basin, China. Then, the real-time footprint distance and climate footprint distance of the field GPP were obtained by using a footprint source area model. Lastly, multiple-scale GPP products were validated at footprint scale, and the impacts (measurement height, surface roughness and turbulent state of the atmosphere) on the footprint distance of field GPP were analyzed. The results of this paper demonstrated that climate footprint distances ranged from about 500 m to 1500 m for different land cover types in the Heihe River Basin. The accuracy was higher when validating MODIS GPP products at footprint scale (R² = 0.56, RMSE = 3.07 g C m⁻² d⁻¹) than at field scale (R² = 0.51, RMSE = 3.34 g C m⁻² d⁻¹), and the same situation occurred in the validation of high-resolution downscaled GPP (R² = 0.85, RMSE = 1.34 g C m⁻² d⁻¹ when validated at footprint scale; R² = 0.82, RMSE = 1.47 g C m⁻² d⁻¹ when validated at field scale). The results of this study provide information about the footprints of field GPP for typical land cover types in arid and semi-arid areas in Northwestern China, and reveal that precision may be higher when validating multiple-scale remote sensing GPP products at the footprint scale than at the field scale.

1. Introduction

Gross primary productivity (GPP) is a fundamental variable in assessing controls on carbon dynamics because it varies with soil water availability, incident solar radiation, temperature, vegetation composition and nutrient availability [1,2]. Therefore, GPP is important in studying terrestrial vegetation carbon cycling and climate change [3]. GPP can be observed using eddy covariance systems at flux towers, which quantify the carbon exchange. However, the observed data can only reflect the carbon exchange to the extent to which the measurement taken in a spatial–temporal domain describes the actual environmental conditions in the space–temporal domain [4]. Remote sensing technology is a good way to simulate regional or global GPP, and to study carbon sources and carbon sinks; however, there is a scale difference between the remote sensing data and footprint of ground-observed GPP. In this context, studying the representativeness of the observed GPP from flux towers, and ways to match ground GPP and satellite-derived GPP is important in validating ecosystem models and remote sensing products. To summarize, the issue of spatial-temporal representativeness is of great significance in model–data benchmarking and remote sensing products’ validation [5].

Several studies have been conducted to analyze flux network-to-region representativeness, i.e., to analyze if the measured data from flux towers reflect the carbon flux in a region. For example, Villarreal et al. [6] analyzed the representativeness of the eddy covariance sites of AmeriFlux based on the ecosystem functional types categories represented by each network. Villarreal and Vargas [7] assessed the representativeness of registered FLUXNET sites across Latin America using GPP and evapotranspiration. Pallandt et al. [8] compared the environmental conditions observed at the tower locations and those within the larger Arctic domain, and mapped the representativeness of these eddy covariance network. Other studies have focused on analysis of the point-to-area representativeness [9], i.e., to what extent do measurements taken at the carbon flux tower reflect the carbon exchange around the measurement instruments. Over the last two decades, a number of footprint models have been proposed to simulate the source area for eddy covariance observations, such as the simple analytical footprint model based on Eulerian coordinates for scalar concentration (SAFE-C) [10], the flux source area model (FSAM) [11,12], the Hsieh model [13], and the Horst–Weil model [14]. Studies of the representativeness of flux tower data for different land cover types have also been conducted. For example, Järvi et al. [15] analyzed the sensitivity of the vertical fluxes at a single measurement point, and estimated the systematic uncertainty of eddy covariance flux measurements in annual cumulative values due to missing data for an urban area based on two towers. Kim et al. [16] analyzed the spatial representativeness of the source area for the vegetation characteristics (density variation and magnitude) within the surrounding area of a flux tower in a mixed forest. Rana et al. [17] studied the representativeness of observed carbon flux data in a Mediterranean urban area with equipment setup restrictions. Chu et al. [4] examined the land cover composition and vegetation characteristics in flux footprints and target areas across AmeriFlux sites, and analyzed the potential biases as a consequence of the footprint-to-target-area mismatch.

With the implementation of the Heihe Watershed Allied Telemetry Experiment Research (HiWATER) in the Heihe River Basin in Northwestern China, many algorithms have been developed and several remote sensing products have been generated in this area [18,19,20,21]. However, it is difficult to verify the accuracy of the proposed models and to validate multiple-scale remote sensing products due to the inadequacy of the footprint information for the carbon flux towers. Challenges still exist in interpreting and incorporating footprint information into these Earth system models due to the lack of footprint information for carbon flux sites [22,23]. Besides, the biases and uncertainty are still not clear when validating these Earth system models and multiple-scale remote sensing products due to a lack of footprint information and because most flux tower sites are located in more-or-less heterogeneous landscapes. In this context, studying the representativeness of the carbon flux data for typical land cover types in the Heihe River Basin, is important for benchmarking these models and multiple-scale remote sensing products in arid and semi-arid areas, which cover more than one-third of China [24].

The aims of this paper were to study the spatial representativeness of field GPP for several typical land cover types in the Heihe River Basin in China, and to analyze the influence factors related to the footprint of field GPP. The results of this paper provide information about the footprint of field GPPs in arid and semi-arid areas in Northwestern China, and reveal a number of impact factors related to the footprint distances of GPP from field eddy covariance systems at flux tower sites in the Heihe River Basin, China.

2. Materials and Methods

2.1. Study Area

The Heihe River Basin is located in Northwestern China (37°41′~42°42′N, 96°42′~102°00′E) (Figure 1). The elevation in this area ranges from 500 m in the Gobi in the northern areas to 5500 m in the mountains in the southern areas. The annual average temperature in the Heihe River Basin is about 6.0 °C to 8.0 °C, annual average precipitation is about 100 mm~400 mm [25], and annual average pan evaporation is about 1200 mm~1800 mm [26]. Croplands are mainly located in the upstream and middle stream in the Heihe River Basin with an area of about 8.5 × 10⁵ ha [27]. The Heihe River Basin is suitable for wheat, corn and Polish canola, with an annual output of 9.9 × 10⁵ tons [28]. The total area of forest is about 1.8 × 10⁵ ha [28], which is mainly located in the southwest of this region. The Heihe River Basin is a typical inland river basin in China, where water-stressed ecosystems are widely distributed, the eco-hydrological processes are complex, and the environment is more sensitive to climate change. Therefore, studying the representativeness of the observed carbon flux in this region is meaningful in revealing the applicability of Earth system models and validating remote sensing products in arid and semi-arid areas in Northwestern China.

2.2. Data and Data Processing

2.2.1. Carbon Flux Observed Data

The carbon flux data used in this paper were from the Multi-Scale Observation Experiment on Evapotranspiration for heterogeneous land surfaces, 2012 of the Heihe Watershed Allied Telemetry Experiment Research (HiWATER-MUSOEXE) [29,30,31]. Carbon flux data with a time interval of 30 min were obtained from nine flux tower sites from May to September in 2014, which is considered to be the growing season in the study area, and were used to study the spatial representativeness of GPP in this paper. The locations of these carbon flux towers are shown in Figure 1. The underlying surface types of the eddy covariance systems include Populus, crops (maize), grass, shrubs, wetland, Gobi and desert. The steps required to process the datasets included time delay corrections, density fluctuations, spectral loss, secondary coordinate rotation, sonic virtual temperature conversion, gaps filled, and the selection of high-quality carbon flux data [3,30,31,32,33,34,35]. The daily carbon flux was the summation of the 30-min data for one day. We calculated the difference in the estimated ecosystem respiration ($R_{eco}$) and observed net ecosystem carbon dioxide exchange (NEE) [36,37,38] as follows:

$$GPP=R_{eco}-NEE$$

(1)

2.2.2. Meteorological Data

The air pressure, wind speed and direction, air temperature and humidity, precipitation, and radiation components (upward short-wave radiation, downward short-wave radiation, upward long-wave radiation, downward long-wave radiation) from May to September in 2014 were collected from meteorological stations at the flux towers to study the source area of field GPP. Half-hourly data were derived from the observed 10-min data. Missing data were filled-in by using the linear interpolation.

2.2.3. Land Cover Data

A 30 m land cover map was used to analyze the composition of land cover types on the footprint of field GPP in our study. This map was classified by using Chinese HJ-1 data [39]. Studies have shown that this land cover map has high accuracy (the overall classification accuracy of the map is 92.19%) [40,41].

2.2.4. Multi-Scale GPP Products

Moderate Resolution Imaging Spectroradiometer (MODIS) GPP products were used in this paper to compare the accuracy when validating field GPP at the footprint scale. MODIS GPP products (MOD17 A2H, Version 6) with a spatial resolution of 500 m and temporal resolution of 8 days were produced based on a light use efficiency model [42].

A 30 m GPP product with a temporal resolution of 8-day in the study area was also used in this paper to compare the accuracy when validating at field scale and at footprint scale. This GPP product was generated by using a downscaled model based on MODIS and Landsat data [43]. Specifically, fraction photosynthetically active radiation (FPAR) and leaf area index (LAI) with a spatial resolution of 30 m and a temporal resolution of 8 days were firstly generated from a downscaled method, and then 30 m GPP were obtained by using a light use efficiency model [43].

2.3. Methods

A flowchart of the process for analyzing the spatial representativeness of field GPP is shown in Figure 2. Firstly, the field GPP was derived from field carbon flux data and the seasonal trends of the GPP were analyzed. Then, a footprint source area model was used to reveal the real-time footprint and climate footprint at the flux tower sites for different land cover types. Lastly, MODIS GPP and downscaled high-resolution GPP products were validated at the footprint scale, and the impact of measurement height, surface roughness and atmospheric stability on the footprint distance of the field GPP were assessed, the composition of land cover types on the footprint of the field GPP was analyzed, and the uncertainly and limitations of this work were discussed.

2.3.1. Real-Time Footprint of Field GPP

The real-time flux footprint depends on the wind direction, measurement height, underlying ground surface characteristics and turbulent state of the atmosphere [44]. In this paper, the FSAM [11,12,45,46] was used to study the representativeness of field GPP for different land cover types in the Heihe River Basin. The FSAM is an analytical method model used to simulate the source area based on the K-theory and advection diffusion model. The carbon flux observed by this instrument can be described as [11,12]:

$$\eta ={\displaystyle \iint }Q\left(x,y\right)f\left(x,y,z_m\right)dxdy$$

(2)

where$\eta$ is the measured flux, $Q\left(x,y\right)$ is the surface source distribution function, $f\left(x,y,z_m\right)$ is the footprint function, $\left(x,y,z_m\right)$ is the location of the instrument.

Based on the K theory, and a two-dimensional advection diffusion equation, the FSAM can obtain the source area under different contribution rates [11,12]:

$$P=\frac{1}{F_u}{\displaystyle \int }{\Omega }_p{\displaystyle \int }[-{{\displaystyle \int }}_{z_0}^{z_m}\overline{u}\left(z\right)\frac{\partial }{\partial x}\overline{C^y}\left(x,z\right)dz{D}_y\left(x,y\right)]dxdy$$

(3)

where P is the contribution rate, $F_u$ is the surface source–sink, ${\Omega }_p$ is the flux source area, $z_m$ is the measurement height, $z_0$ is the surface roughness length, $\overline{u}$ is the mean wind speed profile, $\overline{C^y}$ is the crosswind integrated concentration, and $D_y\left(x,y\right)$ is the crosswind integrated concentration distribution function. The key parameters used to run the FSAM, such as surface roughness length ($z_0$), were obtained from experimental parameters [47,48]. The meteorological datasets used to run the FSAM were obtained from the automatic meteorological stations at the flux towers. The outputs of the FSAM include the length and location of the source area at different contribution rates.

2.3.2. Climate Footprint of Field GPP

The flux climate footprint, which is calculated by aggregating the real-time footprint, can be used to identify the spatial extents and temporal dynamics of the areas contributing to the observed fluxes at a tower site [44]. In this paper, a weighted footprint model was adopted to describe the climate footprint of field GPP [3]:

$$f_{climatology}\left(x,y,z_m\right)={\displaystyle \sum }_{i=1}^{N}f\left(x,y,z_m\right)\frac{Flux\left(i\right)}{{\displaystyle \sum }Flux\left(i\right)}$$

(4)

where $f_{climatology}\left(x,y,z_m\right)$ is the climate footprint, i is the time intervals, N is the number of observed data in the time intervals, (x,y) is the locations, x and y is the direction of the up-wind and cross wind, $z_m$ is the measurement height, $f\left(x,y,z_m\right)$ is the real-time footprint, and $Flux\left(i\right)$ is the measured carbon flux.

2.3.3. Validation of Multiple GPP Products at Footprint Scale

Firstly, we directly compared the field GPP with MODIS GPP products and the downscaled GPP. Considering the co-registration errors between downscaled high resolution t images and field flux towers, the average GPP from a 3 × 3 pixels window around the flux towers was used when validating the downscaled GPP. Then, the weighted MODIS GPP or downscaled GPP in the footprint was compared with the measured data. Specifically, the weight of each pixel in the footprint was obtained according to the contribution rate to the total source area [3]:

$$R=\frac{f{\left(x,y,z_m\right)}_i}{{{\displaystyle \sum }}_{i=1}^{N}f{\left(x,y,z_m\right)}_i}\times 100\%$$

(5)

where R is the weight, N is the total number of the pixels of GPP products in the footprint, and $f{\left(x,y,z_m\right)}_i$ represents the footprint of i-th pixel.

The weighted average of the MODIS GPP or downscaled GPP was calculated and compared with measured GPP, and the determination coefficient (R²) and root mean square error (RMSE) were used to quantify the accuracy of GPP products.

3. Results

3.1. Seasonal Variation in Field GPP

We obtained time series for GPP, NEE and $R_{eco}$ from May to September in 2014 (Figure 3). We found that GPP increased visibly from May to July, and then decreased from July to September in the study area. Generally, GPP was 20.60 g C m⁻² d⁻¹ in the peak in Daman (cropland, mainly maize), which was the highest among these land cover types. This was followed by Shidi (wetland), Huyanglin (deciduous broadleaved forest, DBF) and Hunhelin (mixed forest, MF), where the GPP was 10.28 g C m⁻² d⁻¹, 7.63 g C m⁻² d⁻¹, and 7.08 g C m⁻² d⁻¹, respectively. At Bajitan (Gobi) and Huazhaizi (desert), the GPP was 2.02 g C m⁻² d⁻¹ and 3.38 g C m⁻² d⁻¹, respectively, demonstrating that GPP was lowest in bare land among these land cover types. The average GPP in July was also highest in Daman (cropland), where the average GPP was 14.27 g C m⁻² d⁻¹. This was followed by Shidi (wetland), Huyanglin (DBF) and Hunhelin (MF), where the average GPP in July was 7.65 g C m⁻² d⁻¹, 5.38 g C m⁻² d⁻¹ and 5.47 g C m⁻² d⁻¹, respectively. The average GPP was 0.81 g C m⁻² d⁻¹ in Bajitan (Gobi), which was the lowest among these land cover types.

3.2. Footprint of Field GPP

3.2.1. Real-Time Footprint of Field GPP

A map of the real-time footprint of field GPP on 24 July, 2014, in which the atmosphere is relatively stable on a clear day, is shown in Figure 4. In general, the real-time footprint was consistent with the dominant wind directions, i.e., west and northwest wind. Footprint distances in the east to west direction were more than 500 m in Huyanglin (DBF) and Hunhelin (MF), demonstrating that footprint distances were the furthest in forest. In Sidaoqiao (shrub) and Shidi (wetland), the footprint distances were about 150~250 m. Footprint distances were less than 50 m in Huazhaizi (desert) and Bajitan (Gobi), which were the nearest among these land cover types. The wind direction and height of the instruments were the main factors influencing the real-time footprint, which will be discussed in Section 4.2.

3.2.2. Seasonal Variation in Climate Footprint of Field GPP

The monthly climate footprints of field GPP in the growing season in the study area are shown in Figure 5. The footprint distances were about 500 m to 1000 m at Arou (grass) in the northwest–southeast direction. Maximum footprint distances were about 600 m at Dashalong (grassland) in the northwest–southeast direction. In the middle stream area, the difference in the footprint was caused by the difference in the prevailing wind direction. The prevailing winds were west and northwest at Daman (cropland), northwest and southeast at Shidi (wetland), southwest at Huazhaizi (desert), north and southwest at Bajitan (Gobi). At Huyanglin (DBF) and Hunhelin (MF), the footprint distances were about 1200 m to 1500 m. At Sidaoqiao (shrub), the footprint distances were about 800 m to 1200 m.

We also found that iat Daman (cropland), Huazhaizi (desert) and Bajitan (Gobi), the monthly climate footprint varied significantly. The reasons may be that the surface roughness length changed with the state of growth of crops, and the atmospheric stability may change as these sites are located in the oasis regions. At Huazhaizi (desert) and Bajitan (Gobi), the land cover types were Gobi and desert. The wind speed was high and the wind direction varied, leading to varied footprints in these sites.

3.3. Comparison of Validation Multi-Scale GPP Products at Field Scale and at Footprint Scale

The MODIS GPP and downscaled GPP was validated at field scale and at footprint scale, as shown in Figure 6. Although studies have demonstrated that MODIS GPP was underestimated in the Heihe River Basin [49], we found that the accuracy of MODIS GPP was higher when validated at footprint scale (R² = 0.56, RMSE = 3.07 g C m⁻² d⁻¹) than at field scale (R² = 0.51, RMSE = 3.34 g C m⁻² d⁻¹). The same situation occurred in the validation of downscaled GPP. The precision of downscaled GPP was higher when validating in the footprint scale (R² = 0.85, RMSE = 1.34 g C m⁻² d⁻¹) than validating GPP using field data directly (R² = 0.82, RMSE = 1.47 g C m⁻² d⁻¹). This indicated that the representativeness of the observed carbon flux data provided some important information for the validation of remote sensing GPP products at different scales. This may be useful when taking into account the footprint of field data in the validation of multiple-scale remote sensing products [50,51].

3.4. Land Cover Types in the Footprint

Studies have shown that the composition of the underlying vegetation and its spatial distribution around the eddy covariance observation systems certainly has an influence on the footprint [4,52,53]. The dominant land cover type within the footprints has the highest footprint-weighted percentage. To illustrate the contribution of different land cover types to the footprint, we analyzed the proportions of each land cover types in the footprint based on the 30 m land cover data, and the percentage of land cover types in the footprint is shown in Figure 7. The land cover type in the footprint was homogeneous at Huazhaizi, Bajitan, Dashalong and Arou. Specifically, bare land accounted for 100% of the footprint at Huazhaizi and Bajitan, grass land accounted for 100% in the footprint at Dashalong and Arou from May to September in 2014. However, multiple land cover types were located in the footprints of the other sites (Daman, Shidi, Huyanglin, Hunhelin, Sidaoqiao). The area of the dominant land cover type accounted for more than 90% of the footprint at Daman (cropland), and more than 80% at Shidi (wetland). At Hunhelin, Huyanglin and Sidaoqiao, the dominant land cover types only accounted for less than 50% of the footprint. This demonstrated that the underlying surface was related to the distance of the footprint.

3.5. Main Impact Factors of the Field GPP Footprint

3.5.1. Influence of Measurement Height on Footprint of Field GPP

Measurement height was one of the main factors affecting the footprint distance of field GPP [44,54]. In general, the GPP footprint distance increased with the increase in the height of the eddy covariance systems. In the study area, the footprint distances were the furthest (more than 800 m) at the forest sites (Huyanglin and Hunhelin, where the measurement heights were also the highest (22.00 m). At Huazhaizi and Arou, the footprint distances (less than 200 m) were the nearest, and the measurement height were also the lowest (2.85 m at Huazhaizi, 3.50 m at Arou). In order to illustrate this, we simulated the footprint distance as the measurement height increased by 25% and 50% or decreased by 25% and 50%. As shown in Figure 8, the footprint distance decreased by 23.87%~30.2% when the measurement height decreased by 25%, and decreased by 43.90%~65.20% when the measurement height decreased 50%. The footprint distance increased by 27.72%~36.21% with a 25% increase in the height of the eddy covariance systems, and increased 56.51%~74.90% with a 50% increase in the measurement height.

3.5.2. Influence of Surface Roughness on Footprint of Field GPP

In general, field GPP footprint distance decreased with the increase in surface roughness [38]. When the surface roughness increased from 0.01 m to 0.10 m, footprint distance of field GPP decreased about 33.3%~49.2% (footprint distances decreased more than 40% at Huazhaizi, Huyanglin, Bajitan, Sidaoqiao and Daman, and decreased about 33.0%~38% at Hunhelin, Arou, Shidi and Dashalong) (Figure 9).

3.5.3. Influence of Atmospheric Stability on Footprint of Field GPP

Field GPP footprint distance increased with the increase in atmospheric stability. In this paper, we analyzed the footprint distance when the atmospheric stability was at an unstable, stable and intermediate state (Figure 10). We found that the footprint distance became very small when the atmospheric state was unstable (footprint distances ranged from 15 m to 117 m). When the atmospheric stability was at the intermediate state, the footprint distance was further (footprint distances ranged from 96 m to 1122 m). When the atmospheric stability was at a stable state, the footprint distances were the furthest (footprint distances ranged from 154 m to 2468 m).

4. Discussion

4.1. Scale Mismatching between Satellte GPP and Ground Observed GPP

Studies have shown that GPP estimation from remote sensing are usually robust at most regional or global scales, but that the methods used from benchmarking and validation will impact the accuracy of the GPP estimation models [55]. The spatial resolution of GPP products from satellite observations range from a few hundred meters to tens of meters. Therefore, differences always exist due to the scale mismatch between the grid cells of remotely-sensed data and the scale of the ecological variable on the ground [56], which makes it difficult to directly evaluate remotely-sensed GPP using field carbon flux observations on the ground. Taking the footprint of ground-observed carbon flux into consideration is a good way to eliminate the scale difference between remote-sensing GPP and field-observed GPP. In this paper, the weighted MODIS GPP or downscaled GPP in the footprint were compared with the measured GPP, and higher precision was obtained when validating the multiple GPP products at footprint scale.

Several studies [3,57,58] demonstrated that upscaling the field GPP through machine learning methods was also a good way to make the field GPP comparable with remote sensing products. However, the uncertainties and error transfers in the upscaling process are still not very clear [59]. With the development of airborne and unmanned aerial vehicles (UAV) data, the gap between satellite data and field-measured carbon flux data may be filled, which may be helpful in satellite-derived GPP estimation and validation.

4.2. Limitations and Future Work

Uncertainties in the footprint calculation of observed GPP can result from the footprint models or errors in the input parameters [60]. Specifically, the key parameters to run FSAM, such as surface roughness length ($z_0$), are related to the accuracy of footprint model. In this paper, the key parameters were obtained from the experimental parameter, which could bring some uncertainties in the estimation of the source area of the observed carbon flux data. Besides, footprints from nine flux tower sites covered most of the land cover types in the Heihe River Basin, but they are still insufficient to reflect all the land cover types in this region. In the future, we may conduct further research to assess the representativeness of the locations of the flux tower sites in this area. Furthermore, uncertainties related to the FSAM could also bring some errors. In the future, we may compare the difference between several footprint models, such as the flux footprint prediction (FFP) model [46] and the footprint symmetry index [4].

5. Conclusions

Studying the spatial representativeness of the field GPP for typical land cover types may provide important information for the measurement of carbon flux, as well as the validation of GPP at multiple scales. In this paper, field GPP was firstly derived from field eddy covariance systems. Seasonal trends and daily trends of field GPP were analyzed at nine flux tower sites over seven land cover types in the Heihe River Basin, China. Then, the real-time footprint distance and climate footprint distance of field GPP were described by using a footprint source area model. Lastly, MODIS GPP and downscaled GPP products were validated at the footprint scale, and the impacts (measurement height, surface roughness and atmospheric stability) on the footprint distance of field GPP were analyzed. The results of this paper provided information about the spatial representativeness of field GPP for the typical land cover types in the Heihe River Basin, and also demonstrated the influence factors on the footprint of field GPP in arid and semi-arid areas in Northwestern China.

Generally, significant diurnal and seasonal variation trends could be found at all land cover types, except bare land. Climate footprint distances ranged from about 500 m to 1500 m for different land cover types in the Heihe River Basin. When validating at footprint scale, the accuracy of the MODIS GPP and downscaled GPP was higher compared with validation at field scale. The field GPP footprint distance increased with the increase in the height of the eddy covariance systems and atmospheric stability, while the GPP footprint distance of field GPP decreased with the increase in surface roughness.