Modeling Net Ecosystem Exchange for Grassland in Central Kazakhstan by Combining Remote Sensing and Field Data

Abstract

Carbon sequestration was estimated in a semi-arid grassland region in Central Kazakhstan using an approach that integrates remote sensing, field measurements and meteorological data. Carbon fluxes for each pixel of 1 × 1 km were calculated as a product of photosynthetically active radiation (PAR) and its fraction absorbed by vegetation (fPAR), the light use efficiency (LUE) and ecosystem respiration (R_e). The PAR is obtained from a mathematical model incorporating Earth-Sun distance, solar inclination, solar elevation angle, geographical position and cloudiness information of localities. The fPAR was measured in field using hemispherical photography and was extrapolated to each pixel by combination with the Normalized Difference Vegetation Index (NDVI) obtained by the Vegetation instrument on board the Satellite Pour l’Observation de la Terra (SPOT) satellite. Gross Primary Production (GPP) of the aboveground and belowground vegetation of 14 sites along a 230 km west-east transect within the study region were determined at the peak of growing season in different land cover types and linearly related to the amount of PAR absorbed by vegetation (APAR). The product of this relationship is LUE = 0.61 and 0.97 g C/MJ APAR for short grassland and steppe, respectively. The R_e is estimated using complex models driven by climatic data. Growing season carbon sequestration was calculated for the modelling year of 2004. Overall, the short grassland was a net carbon sink, whereas the steppe was carbon neutral. The evaluation of the modelled carbon sequestration against independent reference data sets proved high accuracy of the estimations.

1. Introduction

Drylands cover more than 30% of the Earth’s surface and comprise a variety of ecosystems, such as shrublands and grasslands, which are large reservoirs of carbon as well as potential carbon sinks and sources to the atmosphere [1]. In Central Asia, where drylands cover about 85% of the total territory, shrubland and grassland carbon sinks represent the major pool for carbon absorption and are believed to offset a significant proportion of carbon emissions associated with fossil fuel combustion [2]. Thus, the ability to monitor carbon sequestration in drylands is of great interest in relation to understanding the current status of the global carbon cycle and to meeting requirements of the Kyoto Protocol.

The terrestrial carbon cycle is a highly dynamic system that includes several storage pools and different components of the flux. Components of the terrestrial carbon flux are gross primary production (GPP), net primary production (NPP) and net ecosystem production (NEP). The amount of carbon that is fixed from the atmosphere, i.e., converted from CO₂ to carbohydrates during photosynthesis, is called GPP, which is carbon assimilation by photosynthesis ignoring photorespiration. Net ecosystem production, NEP, is the difference between gross primary production (GPP) and total ecosystem respiration (R_e), which determines the amount of carbon lost or gained by the ecosystem without disturbances, such as harvests and fire [3]:

NEP = GPP − R_e

(1)

where R_e is the sum of all autotrophic (R_a) and heterotrophic (R_h) respiration over some time period from instantaneous fluxes to annual total and is defined by the following equation:

$$R_e=R_a+R_h$$

(2)

Net Primary Production (NPP) is an important component of the carbon cycle and key information of ecosystem performance. NPP represents the net new carbon stored as biomass in stems, leaves or roots of plants and defines a balance between GPP and autotrophic respiration R_a [4,5]. The classic model for the net primary production of ecosystem is given by:

GPP (g C/m²/yr) = NPP (g C/m²/yr) + R_a (g C/m²/yr)

(3)

Analysis of biomass inventory data obtained by field survey provides a capability for wide-area monitoring of biomass production and carbon sequestration using the regression models between remote sensing data and measured biomass data [6,7]. Another concept is used for the physical-based models, where the key factors for plant development are extracted from remote sensing data to calculate the production of plant biomass [4,5,8]. These and other related studies used the approach first described by [9] to calculate biomass by linking the incoming radiation to vegetation production through an empirical biophysical conversion factor – the light use efficiency (LUE). According to [9], gross primary productivity can be calculated as a product of LUE and the amount of photosynthetically active radiation absorbed by plants (APAR):

$$GPP=LUE*APAR$$

(4)

This equation includes both above-ground and below-ground compartments of biomass. In the equation above, GPP can be estimated using remotely sensed data. The current methodology for estimating GPP and NPP from satellites combines the amount of photosynthetically active radiation absorbed by the canopy (APAR) with LUE, most commonly as a biome specific constant, modified taking into account the constraints imposed by climatic limitations through using a climatic stress factor, SI, as part of LUE [3,10,11,12,13]. APAR depends upon the geographic and seasonal variability of day length and potential incident radiation, as modified by cloud cover and aerosols, and the geometry of displayed leaf material. For calculation of APAR from the photosynthetically active radiation (PAR), the fraction of photosynthetically active radiation absorbed by vegetation (fPAR) is needed. In remote sensing studies, the Normalized Difference Vegetation Index (NDVI) computed from red and infra-red satellite channels is commonly used to estimate fPAR [14,15]. The final equation for satellite-based GPP estimation will be the following [16,17]:

$$GPP=LUE*fPAR*PAR$$

(5)

The amount of the solar radiation reaching the canopy (PAR) is usually derived from remotely sensed data or computed using common mathematical algorithms [18,19]. An important attribute of remotely sensed data is that vegetation properties are known over the entire area and estimates of GPP/NPP can be obtained for each pixel in the area to be studied. However, satellite estimates of GPP/NPP as well as other ecosystem fluxes may not be very accurate [20].

Operational monitoring of GPP/NPP at the global scale is now underway using imagery from the satellite-borne Moderate Resolution Imaging Spectrometer (MODIS) sensor [3,13]. The MODIS NPP algorithm incorporates separate estimation of gross primary production and autotrophic respiration to account for potential differences in light use efficiency between vegetation types due to differences in respiration costs and climate conditions. These biome- and climate-induced ranges of LUE are taken into account in the MODIS algorithm by incorporating the optimal biome-specific LUE parameter, derived from extensive set of observations, and linear ramp functions of minimum temperature and vapour pressure deficit. The lower and upper limits of these functions are biome-specific constants [13,16,21]. Since the year 2000, by using five different data products from MODIS and two constant datasets, namely a digital elevation model and global soil information, the aboveground biomass in terms of GPP/NPP is operationally being calculated as 8-day time steps and as an annual sum with a spatial resolution of 1-km. However, MODIS GPP/NPP products require extensive validation, which is being in progress using ground truth data from different regions, to be useful for scientific purposes. The ongoing validation of MODIS products reveals significant biome- and region-specific over- or underestimation of GPP/NPP. Thus, the validation study by [21] reported about moderate underestimation in the MODIS products at agricultural land, and strong overestimation at grassland. Despite of the fact that MODIS GPP/NPP products cover the entire land surface and are freely available, the use of them in many cases may not always be reliable. In these cases, estimation of GPP/NPP using remote sensing data from other satellite systems can be preferable.

Many satellite-based products only yield GPP [22,23,24], which, however, is only one of the carbon cycle processes that affect carbon sequestration and by its definition does not include respiration. Respiration is often incorporated as a climate-dependent factor in carbon sequestration models, where it is commonly calculated using complex mathematical models driven by climatic data [3,17]. This approach, being used creatively, may supply a feasible alternative to the estimation of ecosystem respiration based on cost-expensive and time-consuming ecosystem carbon flux measurements. In a number of recent international projects such as AMERIFLUX, EUROFLUX, FLUXNET, HAPEX, ecosystem carbon fluxes are investigated in different ecosystems [25,26]. The main deficiency of the ecosystem carbon flux measurements is that their results are local. The challenge for modellers, therefore, is to find bounds to the application of results from small-scale studies in the elucidation of large-scale processes. Spatial extrapolation of small-scale results of flux measurements and monitoring carbon fluxes at larger spatial scales is feasible through combining ecosystem flux measurements and remote sensed data. The approach combining remotely sensed and ecosystem-flux data has been discussed several times [16,17,20,26,27]. However, most of the study sites concentrate in the western countries except a few remote areas like the Arctic tundra, Siberian taiga, the Amazonian forest and Sahelian savannah. For the drylands ecosystems of Central Asia, knowledge on carbon fluxes is still very poor. The challenge for the modellers working in this region, therefore, is to find appropriate methods for estimating carbon sequestration without undertaking the expensive and often not realizable ecosystem flux measurements. In grasslands of Central Asia and Kazakhstan, ecologists have measured annual NPP through direct sampling of biomass growth increment for decades. However, most of these studies only provide results of local small-scale measurements from different decades that are not feasible for producing carbon stocks maps [28,29,30,31,32,33].

This study describes and tests a satellite data-based model for mapping total carbon sequestration over a semi-arid grassland-dominated region in Central Kazakhstan. This model employs the algorithm of Monteith and combines remote sensing data from SPOT-Vegetation satellite and ground inventory data on biomass. The approach presented in the study expands the algorithm of Monteith incorporating ecosystem respiration. Light use efficiency was estimated from the ground-based data using a simple empirical approach. Results were considered in relation to design of algorithms for satellite-based monitoring of the Kazakhstan-wide carbon cycle.

2. Study Area

The study area is located in the middle part of Kazakhstan between 48°20′ and 49°30′ northern latitude and 72° and 74°10′ eastern longitude and encompasses the southern margin of the Kazakh Hills. It comprises the northern area of the Shetsky raion (district) in Karaganda oblast’ (province). The climate of the region is dry, cold and high continental. Average annual precipitation is about 250-300 mm. The most part of precipitation falls during warm period from March to October. Inter-annual rainfall variation has a coefficient of variation of 20-35%. The region is often affected by drought hazard. The temperature amplitude is relative high, average January temperature is below –12°C and average July temperature is about 26-28 °C. The growing season starts in April and continues till October.

In the vegetation cover, xerophylous bunch grasses are dominant in the northern part of the study area and dwarf shrubs in the south. The steppe grassland is dominated by the genera Festuca and Stipa. Few euryxerophilous forbs occur; the co-dominants are dwarf shrubs of the genus Artemisia and sometimes of other genera, particularly Anabasis and Salsola. Species diversity is about 12-15 species per square metre. The height of the canopy decreases from 30-40 in the north to 15-20 in the south, while vegetation cover decreases from 50-70% to 20-30%, and even less. The vegetation growth in the study area is strongly dependent on precipitation dynamics. Grasses and shrubs in the vegetation cover grow during the whole vegetative period, but the vegetation growth is most rapid during May and early June (the period of greatest precipitation) in the southern part, and during June in the northern part of the study area. During droughty months in summer (July and early August) their rate of development is hindered. This period of semi-dormancy occurs throughout the study region [34,35].

Despite of the fact that the study area shows a wide spectrum of vegetation communities due to variations in topography and soils, there are only two main vegetation classes, short grassland and steppe grassland (Figure 1). The land cover data in the study area were taken from the digital land-cover map derived from Moderate Resolution Imaging Spectroradiometer (MODIS) by [36] that has been accessible in the United States Geological Survey (USGS) archive centre. According to this map, short grassland covers 74.22% of the whole territory, while steppe grassland covers 25.77%.

Figure 1. Location of the study region on the map of Kazakhstan and distributions of land cover classes in the study region based on the MODIS land-cover map [36].

Figure 1. Location of the study region on the map of Kazakhstan and distributions of land cover classes in the study region based on the MODIS land-cover map [36].

3. Data

3.1. Field data

3.1.1. Carbon data

Carbon stocks residing in the vegetation and soils were calculated from field records collected along a 230-km transect in June 2004. The fieldwork was carried out by colleges from the Institutes of Botany and Soils science (Almaty, Kazakhstan). 14 test sites were selected for analysis with four replicates (sub-plot) per site. Each replicate occupied 0.0625 ha. The test sites were classified by the land cover categories in accordance with the land cover map in Figure 1. Of the test sites, 6 were located within the steppe grassland and 8 within the short grassland. The data consisted of grass biomass (g DM/m²), root biomass (g DM/m²), litter mass (g DM/m²), vegetation cover estimations (%), and mass of soil carbon at the depths from 10 to 50 cm (g C/m²). The peak-season living aboveground biomass (AGB₁), root biomass (belowground biomass, BGB) and litter was measured by destructive sampling of 1.0 m², with samples weighted, sub-sampled, and dried at 65 °C to constant weight to correct for moisture content. Four samples per each replicate (i.e., 16 per plot) were taken. Roots were collected by excavating a square of 1 × 1 m² to a depth of 50 cm with a narrow, flat-bladed shovel and handsaw. Coarse roots were hand sorted and washed. The remaining sample was dispersed in tap water, passed through a 2-mm sieve, and fine roots were collected without differentiation between live and dead roots. Roots were washed of gross mineral contamination, dried at 65 °C to constant weight, and weighted. The proportions of litter and roots to aboveground biomass were calculated. The original values of dry matter (g DM) were converted to carbon (g C). The proportion of C in all biomass pools was assumed to be 0.47.

We reinterpreted the field data with regard to the land cover classes distributed in the study area. From the field records of the individual test sites average values for the two land cover categories, short grassland and steppe grassland, were determined. The results of these calculations are given in Table 1.

Table 1. Description of field data used for the modelling (averaged by land cover type).

Land cover	Dominate vegetation	Carbon, g/m²
		BGB	AGBl	Litter	Total
Short grassland (8 sample plots)	Convolvulus arvensis, Artemisia sublessingiana	38.1	54.2	80.5	172.8
Steppe grassland (6 sample plots)	Festuca valesiaca, Stipa capillata, Coeleria cristata	58.0	56.4	61.3	195.7
Average of 14 sample plots		48.0	61.6	70.9	179.5

Table 1. Description of field data used for the modelling (averaged by land cover type).

Land cover	Dominate vegetation	Carbon, g/m²
		BGB	AGBl	Litter	Total
Short grassland (8 sample plots)	Convolvulus arvensis, Artemisia sublessingiana	38.1	54.2	80.5	172.8
Steppe grassland (6 sample plots)	Festuca valesiaca, Stipa capillata, Coeleria cristata	58.0	56.4	61.3	195.7
Average of 14 sample plots		48.0	61.6	70.9	179.5

3.1.2. Vegetation structure variables

Leaf area index (LAI) and the fraction of absorbed photosynthetically active radiation (fPAR) are two key variables of vegetation structure used as input variables to models that predict NPP over extensive terrestrial areas [5,27,37]. Either direct or indirect methods have been used to estimate in situ LAI and fPAR. In the present study, LAI and fPAR were measured at 25 sample plots across the study area using gap fraction analysis of hemispherical photos. The gap fraction is an integrated canopy structural quantity. The vertical gap fraction corresponds directly to the complement to unity of the cover fraction and allows estimating both the LAI and fPAR of canopies with an assumption about the spatial organisation of the vegetation elements [38,39]. Recent advances in instrumentation and the improvement of gap fraction inversion models have significantly increased the attractiveness of this technique [40].

For hemispherical photography, a WinScanopy Image Acquisition instrument from “Regent Instruments” (http//www.regentinstruments.com) was used. Field survey was carried out at the peak of vegetation activity in June 2008. The plot size was chosen to correspond to an area observed by nine 30-m Landsat ETM+ pixels. Plot locations and coordinates were determined using global positioning system (GPS) and hemispherical photos were taken at each of the sub-plots. Each plot has a size of 90 × 90 m with multiple sub-plots over which measurements were averaged to produce a unique value per plot. The measurements were made in 30-m transect spacing within each plot. This spacing was selected on the basis of test trials conducted in reference plots and analysis of fine-scale vegetation patterns at Landsat ETM+ imagery. In order to take pictures underneath the canopy, the hemispherical camera was placed on ground surface on such manner that the objective surface was at the similar level with the ground surface. On such way, the portion of the light absorbing canopy, missed by the measurements, was minimal.

The “Can Eye” software (INRA, France, http//www.avignon.inra.fr/can_eye/) was used for the processing of hemispherical photographs and derivation of LAI and fPAR values for each measurement site and than for each plot. “Can Eye” uses a look-up-table approach which is composed of gap fraction measurements in different view zenith angles and the corresponding LAI and fPAR parameters.

The hemispherical measurements (non-contact method) were supported by calculation of LAI from harvesting biomass by destructive sampling (contact method) at the test sites described above in Section 3.1.1. A value of leaf area index for samples of contact measurements was calculated as:

LAI = AGB (gC/m²)*SLA(m²/gC)

(6)

where AGB is aboveground biomass, and SLA is the specific leaf area. The SLA is an ecophysiological biome-dependent constant with values taken from [37]. The LAI calculated from contact measurements served as reference for LAI estimations by hemispherical photography. The study by [41] compared the LAI values from contact and non-contact measurements in the data set used in the present study. They found a high consistency of the LAI values from different sources.

3.2. Satellite data

The Normalized Difference Vegetation Index (NDVI) has been the most common vegetation index used throughout the history of satellite data applications. NDVI represents the absorption of photosynthetic active radiation and hence is a measurement of the photosynthetic capacity of the canopy. NDVI is calculated from reflectance in the red (R) and near-infrared (NIR) bands using the following equation:

$$NDVI=\frac{(NIR-R)}{(NIR+R)}$$

(7)

NDVI has been established to be highly correlated to green-leaf density, absorbed fraction of photosynthetically active radiation and above-ground biomass and can be viewed as a surrogate for photosynthetic capacity [42,43].

10-day 1-km NDVI data set for the growing season (April-October) 2004 was obtained from SPOT-Vegetation satellite data. The data set had been generated using a maximum value composite procedure, which selects the maximum NDVI value within a 10-day period for every pixel. On this way, non-vegetated noise originally presented in the row data was significantly reduced [44]. Nonetheless, some noised pixels remained in the NDVI data set and were removed by employing a spatial filter, which calculated a new value for each noised pixel from values of neighbouring non-noised pixels.

3.3. Climate data

10-day maximum, mean and minimum temperature, air humidity and the cloud cover of sky data for 2004 at nine climate stations located in the study area were obtained from the National Hydro-Meteorological Centre of Kazakhstan. Gridded maps of 10-day values for each variable were constructed using the interpolation method kriging with external drift (KED) where a digital elevation model, scaled in meters, was used as an external drift. In order to assess the accuracy of the data preparation, we randomly removed three weather stations from the interpolation for one of the 10-day from the recorded values of mean temperature. Average error resulted in a value less than 6%. It means that the KED approach worked very effectively.

4. Methods

Figure 2 illustrates the processing stream. The parameters inside the ellipsoids, represents the raw data (including constants derived from the literature) used in the model, whereas the parameters inside the boxes are derived during the execution of the model. Note that Equations 1-4 represent a highly generalized form of this flowchart, and are depicted in the white boxes in Figure 2. Derivation of individual parameters at per-pixel level will be described below. Parameters used as constants are joined together in Table 2.

Figure 2. Flowchart of the NEE model.

Figure 2. Flowchart of the NEE model.

Table 2. The constant parameters, their descriptions, and values used for the light use efficiency model in this study.

Parameter	Units	Description	Land cover class		Source
			Short grassland	Steppe grassland
LUEopt	gC/MJ	Optimum radiation conversion efficiency	0.61	0.97	This study
$Q_{10}$	none	Temperature sensitivity factor	2.0	2.0	White et al. [37]
$T_b$	°C	Base temperature	20.0	20.0	White et al. [37]
SLA	m²/kg C	Specific leaf area	40.0	40.0	White et al. [37]
$r_{m,root}$	gC/gC/day	Maintenance respiration coefficient for fine roots at 20°C	0.005	0.006	White et al. [37]
$r_{m,leaf}$	gC/gC/day	Maintenance respiration coefficient for above-ground plant parts at 20°C	0.009	0.01	White et al. [37]
$r_{g,root}$	gC/gC	Growth respiration coefficient for fine roots	0.3	0.3	White et al. [37]
$r_{g,leaf}$	gC/gC	Growth respiration coefficient for above-ground plant parts	0.3	0.3	White et al. [37]
$r_{a,root}$	none	Carbon allocation fraction for fine roots	1.0	1.0	White et al. [37]
$r_{a,leaf}$	none	Carbon allocation fraction for above-ground plant parts	1.0	1.0	White et al. [37]
BGB/AGB1	none	Ratio of below ground biomass to above ground biomass	0.63	1.08	This study
VPDstop	Pa	Vapour pressure deficit at photosynthesis stop	5000	4100	White et al. [37]
VPDstart	Pa	Vapour pressure deficit at photosynthesis start	1000	970	White et al. [37]
$T_{opt}$	°C	Temperature at which LUE = LUEopt	25.0	25.0	White et al. [37]
$T_{range}$	°C	Maximum temperature range	40.0	40.0	White et al. [37]
F	gC/m²/day	Soil respiration rate at 0°C	1.250	1.480	Raich et al., [55]; Tesarova & Gloser, [62]

Table 2. The constant parameters, their descriptions, and values used for the light use efficiency model in this study.

Parameter	Units	Description	Land cover class		Source
			Short grassland	Steppe grassland
LUEopt	gC/MJ	Optimum radiation conversion efficiency	0.61	0.97	This study
$Q_{10}$	none	Temperature sensitivity factor	2.0	2.0	White et al. [37]
$T_b$	°C	Base temperature	20.0	20.0	White et al. [37]
SLA	m²/kg C	Specific leaf area	40.0	40.0	White et al. [37]
$r_{m,root}$	gC/gC/day	Maintenance respiration coefficient for fine roots at 20°C	0.005	0.006	White et al. [37]
$r_{m,leaf}$	gC/gC/day	Maintenance respiration coefficient for above-ground plant parts at 20°C	0.009	0.01	White et al. [37]
$r_{g,root}$	gC/gC	Growth respiration coefficient for fine roots	0.3	0.3	White et al. [37]
$r_{g,leaf}$	gC/gC	Growth respiration coefficient for above-ground plant parts	0.3	0.3	White et al. [37]
$r_{a,root}$	none	Carbon allocation fraction for fine roots	1.0	1.0	White et al. [37]
$r_{a,leaf}$	none	Carbon allocation fraction for above-ground plant parts	1.0	1.0	White et al. [37]
BGB/AGB1	none	Ratio of below ground biomass to above ground biomass	0.63	1.08	This study
VPDstop	Pa	Vapour pressure deficit at photosynthesis stop	5000	4100	White et al. [37]
VPDstart	Pa	Vapour pressure deficit at photosynthesis start	1000	970	White et al. [37]
$T_{opt}$	°C	Temperature at which LUE = LUEopt	25.0	25.0	White et al. [37]
$T_{range}$	°C	Maximum temperature range	40.0	40.0	White et al. [37]
F	gC/m²/day	Soil respiration rate at 0°C	1.250	1.480	Raich et al., [55]; Tesarova & Gloser, [62]

4.1. PAR and APAR

The photosynthetically active radiation (PAR) is defined as the domain of incoming solar radiation exploited by green vegetation for photosynthesis (400 – 700 nm). PAR can be derived from a budget modelling approach, with the potential solar radiation calculated from the geographical location [45,22], or modelled from remote sensing data [18,19]. In this study, PAR was obtained on pixel-by-pixel basis from the budget modelling approach, which computes the solar irradiance at the top of the atmosphere and transforms it to the amount of solar radiance coming to the Earth’s surface. For obtaining the solar radiance at the top of the atmosphere this algorithm uses the following variables, Earth-Sun distance, solar inclination, the angle between the Earth’s orbital and equatorial planes, solar elevation angle, geographical position, and day of year [45]. The variables of surface elevation, day length and cloud cover information at localities were used to compute optical depth of the atmosphere and to estimate solar irradiance reaching the ground. Spatial distribution of the radiation reaching the ground strongly depends on the geometry of terrain (relief slope, exposition, aspect). These variables were obtained from a digital terrain model and inputted into the equation given by [46]. The incoming solar energy is reduced to the photosynthetically active radiation by a factor of 0.48 [18]. Amount of the absorbed photosynthetically active radiation (APAR) was calculated through multiplication of PAR with fPAR.

4.2. fPAR and LAI

Recent studies have shown that both fPAR and LAI are strongly related to multispectral vegetation indices and can be determined from satellite derived NDVI data [14,15,43,47]. Relationships between these both vegetation structure variables and NDVI have been found to be remarkably consistent over a range of non-woody vegetation types. In the present study, we estimated fPAR on pixel-by-pixel basis using a linear regression model between the ground-based fPAR values and corresponding NDVI values from SPOT VGT recorded on the sample plots at the time of the in situ measurements (3rd decade of June 2008):

$$fPAR=1.25*NDVI-0.179$$

(8)

For calculation of LAI we employed an empirical model described by [41,48]. This model was established using the data set of ground-based LAI from the study region with combination of Landsat ETM+ and SPOT-VGT data. The formulation of the model is given by the following equation:

$$LAI=0.4377*\ln (NDVI)+1.0205$$

(9)

Previous results have shown that both models were statistically significant at the level of p-value < 0.01 with values of R² = 0.59 and R² = 0.67, respectively [41,48]. The models were considered to have sufficient prediction accuracy and were used for modelling 10-day time series of fPAR and LAI from SPOT VGT NDVI data. Thereafter, the derived fPAR data set was used for computing APAR as the product of PAR and fPAR, whereas the derived LAI data set was used for estimating leaf growth and maintenance respiration (see Equations 16 and 17).

4.3. LUE

The most common approach to determining LUE is through its relation with various environmental and vegetation-type specific factors that describe the limits to the carbon uptake [20,50,51,52]. For remote sensing based estimations of GPP, the efficiency of light conversion has been commonly assumed to have a constant value [22,24]. [52] used estimates of GPP and APAR from micrometeorological flux measurements at test sites for calculation of daily light use efficiency:

$$LUE(gC/MJ)=\frac{GPP}{APAR}$$

(10)

Unfortunately, we did not have micrometeorological flux measurements of net ecosystem exchange from the study area. However, we had the peak growing season estimations of carbon from 14 test sites. For that reason, we modified Equation (10) for the use with the available ground data. The modified equation is presented by:

$$LUE(gC/MJ)=\frac{AG{B}_l+L+BGB/3}{{\displaystyle \sum _{i=1}^{9}APA{R}_i}}$$

(11)

where $AG{B}_l$ is the peak-season living aboveground biomass (g C/m²), L is the peak-season litter (g C/m²), BGB is belowground biomass divided by 3 (turnover time of roots for grassland biome) [29], and APAR₁ is amount of photosynthetically active radiation absorbed by vegetation during a 10-day period. The equation sums APAR during nine decades over a period from the beginning of the growing season (1^st decade of April) to the time-point of the field measurements (3^rd decade of June). A LUE value for each test site was estimated using Equation (11).

4.4. Modelling factors limiting light use efficiency

In general, a variation in LUE over time for an ecosystem depends mainly on temperature and moisture deficit. To simulate the reduction of the light use efficiency coefficient due to climate conditions departed from the optimum, the SI was modelled as a product of two limiting indices based on air temperature (T) and vapour pressure deficit (VPD):

$$SI=f(T)*f(VPD)$$

(12)

where the environmental functions T and VPD are scalars between 0 and 1 which are formed in the same way as in BIOME-BGC [37,53]. The temperature function is based on an equation for a Gaussian temperature response that scales the temperature variable to an index between 0 and 1. At the temperature minimum the plant growth amounts to 0; the photosynthesis rate achieves its maximum at the temperature optimum; a further increase of temperature reduces the plant growth rate; and at the temperature maximum photosynthesis stops. The VPD index scales the availability of moisture in air between the values at which the photosynthesis process stops and the optimum value corresponding to each vegetation type. These functions are expressed as:

$$f(T)=\{\begin{array}{ll}0& T<1\\ \cos (\frac{\pi }{2}(T-T_{opt})/(T_{range}-T_{opt}))& T<T_{opt}\\ \ln (T)/\ln (T_{opt})& T<T_{opt}\end{array}$$

(13)

$$f(VPD)=\{\begin{array}{ll}0& VPD>VP{D}_{stop}\\ (VP{D}_{stop}-VPD)/(VP{D}_{stop}-VP{D}_{start})& VP{D}_{start}<VPD<VP{D}_{stop}\\ 1& VPD<VP{D}_{stop}\end{array}$$

(14)

The values for minimum, maximum and optimum of T and VPD parameters were taken from BIOME-BGC [37].

4.5. Modelling ecosystem respiration

The annual sum of autotrophic and heterotrophic respiration was estimated using complex models driven by climatic data. These models combine a number of biome-specific constants, which were borrowed from recent literature, and variables, which were empirically estimated in this study (see Table 2). Conventionally, autotrophic respiration is separated into maintenance respiration R_m and growth respiration R_g [3,54]:

$$R_a=R_m+R_g$$

(15)

Maintenance respiration is strongly dependent on temperature and can be modelled using the following equation [3,54]:

$$R_{m,i}=M_i{r}_{m,i}{Q}_{10}^{(T-T_b)/10}$$

(16)

where R_m,i is maintenance respiration of plant component $i$, $M_i$ is the biomass of plant component $i$, r_m,i is maintenance respiration coefficient for component $i$, $Q_{10}$ is the temperature sensitivity factor, and T_b is the base temperature. R_m,i was calculated separately for plant leaves (above-ground biomass) and roots (below-ground biomass). The M_i for plant leaves was estimated as LAI/SLA, where LAI, the leaf area index, is obtained from remote sensing data (see description above), and the specific leaf area is a biome-specific constant (Table 2). The roots mass (M_i) was calculated as a product of leaf mass multiplied with shoot-root ratio (BGB/AGB₁ in Table 2).

Growth respiration is generally considered to be independent of temperature and is proportional to GPP [3,54]:

$$R_{g,i}=r_{g,i}{r}_{a,i}GPP$$

(17)

where r_g,i is a growth respiration coefficient for plant component i, and r_a,i is the carbon allocation fraction for plant component i. The values of r_g,i and r_a,i are biome-specific constants (Table 2). R_g,i was computed separately for above-ground and below-ground plant components. The separation of annual GPP into below-ground and above-ground components was achieved using the variable BGB/AGB_l empirically derived from the ground-based data.

Heterotrophic respiration is the sum of respiration by decomposers (animals, bacteria and fungi feeding on dead tissue or soil organic matter) and is most closely associated with the decomposition of residue and soil organic matter. The relative contributions of autotrophic and heterotrophic respiration to the total are important, whereas heterotrophic respiration includes soil respiration and respiration of dead biomass on the ground surface (litter). In the present study R_h was modelled using equation given by [55]:

$${\mathrm{R}}_{\mathrm{h}}=\mathrm{F}*\exp (0.05452*{\mathrm{T}}_{\mathrm{a}})*\mathrm{P}/(4.259+\mathrm{P})$$

(18)

In this equation R_s refers to the mean monthly soil efflux in g C/m²/day, F is the base soil respiration rate, T_a refers to the mean monthly air temperature (°C), and P is the mean monthly precipitation (cm). The base soil respiration rates for the soil types dominating steppe grassland and short grassland were borrowed from literature.

5. Results

5.1. Light use efficiency

The derived values of the LUE coefficient for the individual test sites ranged from 0.34 to 1.59 g C/MJ, with a mean value of 0.86 g C/MJ and a standard deviation of 0.42 g C/MJ. LUE value was lowest at a degraded test site, highest at a steppe test site covered by perennial gramineous grasses used as pasture, and intermediate at test sites in short grassland and slightly degraded steppe grassland. The values of LUE for the two land cover categories, short grassland and steppe grassland, were 0.61 gC/MJ and 0.97 gC/MJ, respectively. This pattern is generally consistent with differences among the biomes in net primary production and in light use efficiency for net primary production [37,56,57].

5.2. Carbon sequestration

The model is used to calculate GPP, NPP and NEE at a spatial resolution of 1 km² and a temporal resolution of ten days. The basis for this modelling is a time series of the described SPOT-Vegetation data set. From 10-day values we calculated growing season products of GPP, NPP, R_e and NEE for the year 2004. Figure 3 shows final maps of these products. We also disaggregated our estimates on the basis of land cover classes (Table 3). Class totals of GPP, NPP, R_e and NEE are given, as well as their means per category and standard deviations.

Table 3. Modelled GPP, NPP and NEE* by cover type.

Cover type	Area (km²)	Modelled variables	Mean (g C/m²/year)	Standard deviation	Total (g C×106)
Short grassland	20,028	GPP NPP Re NEE	211.03 131.30 187.05 13.17	51.42 44.29 35.12 32.03	4,226,509 2,629,877 3,446,237 263,768
Steppe grassland	6,955	GPP NPP Re NEE	243.70 145.61 243.79 -0.09	59.50 54.22 44.47 39.56	1,694,934 1,012,648 1,695,351 -625
Total	26,983	GPP NPP Re NEE	227.35 138.9 215.4 6.54	55.46 49.25 44.9 35.80	5,921,442 3,642,525 5,658,300 263,142

*Positive NEE indicates that site was a sink for carbon, and negative NEP indicates the sites were a source for carbon to the atmosphere.

Table 3. Modelled GPP, NPP and NEE* by cover type.

Cover type	Area (km²)	Modelled variables	Mean (g C/m²/year)	Standard deviation	Total (g C×106)
Short grassland	20,028	GPP NPP Re NEE	211.03 131.30 187.05 13.17	51.42 44.29 35.12 32.03	4,226,509 2,629,877 3,446,237 263,768
Steppe grassland	6,955	GPP NPP Re NEE	243.70 145.61 243.79 -0.09	59.50 54.22 44.47 39.56	1,694,934 1,012,648 1,695,351 -625
Total	26,983	GPP NPP Re NEE	227.35 138.9 215.4 6.54	55.46 49.25 44.9 35.80	5,921,442 3,642,525 5,658,300 263,142

*Positive NEE indicates that site was a sink for carbon, and negative NEP indicates the sites were a source for carbon to the atmosphere.

The maps of the growing season GPP and NPP reveal clear spatial patterns in their distribution over the study area (Figure 3a,b). The growing season GPP and NPP was significantly higher in the areas covered by the steppe grassland. This difference was primarily a reflection of differences in the productive potential of the vegetation in the two land cover types. As the steppe grassland occupies territories with significantly higher precipitation amount, the conditions for vegetation growth are much better. It is well known that, because precipitation is the major limiting factor for vegetation growth in drylands, it controls strongly the spatial patterns in vegetation [35,58,59]. Within the short grassland, there is a distinct gradient in the distribution of GPP and NPP from South to North. The southern areas show low biomass production with values of 12-150 g C/m² for GPP and 100-120 g C/m² for NPP, while the northern areas of the short grassland demonstrate much higher values of GPP and NPP, 200-250 g C/m² and 130-160 g C/m², respectively. This pattern was driven primarily by the South-North precipitation gradient clearly presented in the flatland areas, which are occupied mainly by short grassland. The mean GPP was 243.7 g C/m²/year for the steppe grassland and 211.03 g C/m²/year for the short grassland (Table 2). The mean NPP was 145.6 g C/m²/year for the steppe grassland compared to 131.31 for the short grassland.

Figure 3. Spatial distribution of the modelling results, (a) Gross primary production, GPP, (b) Net primary production, NPP, (c) Total ecosystem respiration, R_e, and (d) Net ecosystem exchange, NEE. The units are g C/m²/year.

Figure 3. Spatial distribution of the modelling results, (a) Gross primary production, GPP, (b) Net primary production, NPP, (c) Total ecosystem respiration, R_e, and (d) Net ecosystem exchange, NEE. The units are g C/m²/year.

The modelled NEE ranges throughout the study area from highly negative to positive values (about –100 to more than 100 g C/m²/year) (Figure 3d). Spatial pattern in the NEE is different from those of GPP and NPP. The spatial distribution of the net ecosystem production, NEP, was strongly imposed by the patterns in the ecosystem respiration, R_e (Figure 3c). There is a strong gradient in the distribution of NEE from SE to NW. Most of the areas with positive values of NEE are located in the western and southern parts of the study area, while the north-eastern part is dominated by low or negative values of NEE. This gradient is more consistent with the temperature gradient than with the precipitation gradient. This result agrees with the results reported by [20] who found no simple relationship between annual NEE with annual precipitation for grass prairie in North America. The association of NEE distribution with the land cover categories was not as strongly as that of GPP and NPP. As a whole, the short grassland was a carbon sink, while the steppe grassland was nearly in carbon balance (Table 3). The NEE values totalled by the cover type showed that the steppe grassland as a whole was characterized by a slight escape of carbon into the atmosphere (-625.9 tons), while the short grassland assimilated a total of 263768 tons carbon (about 13.7 tons C/km²/year). A positive balance of carbon characterizes the whole study region in the year 2004.

5.3. Evaluation of the modeling results

We used three techniques to assess model performance. First, a comparison of modelling results with biomass and GPP values reported in recent literature for the similar biomes was made. Secondly, the modelling results were directly compared with the ground-based field data on biomass at 14 test sites. And finally, the results of the modelling were compared with the MODIS GPP/NPP product for the year 2004. In this section we describe the results of model validation.

5.3.1. Qualitative comparison with previous studies

Intensive efforts have been made to collect all credible published estimates of biomass, GPP and NEP from grassland biomes, particularly in the study region. The method used to predict NPP/GPP in previous studies range from biomass harvesting, simple above-ground NPP-NDVI correlations to complex ecophysiological simulation models, and cover a number of scales from point to regional and global. For qualitative comparison with our results, we compiled a number of published studies carried out in grasslands of Kazakhstan [28,29,30,31,32,33,34]. We also took some studies providing averages of vegetation productivity with regard to vegetation biome and climatic zone, which are calculated globally [60] or regionally [61].

Table 3 gives a summary of the comparison, which we based on the literature survey. The table also contains some studies carried out in grasslands in other regions (USA, Sahel, former Czechoslovakia). Qualitative agreement of productivity estimations in previous studies, particularly in the studies from Central Asia, with our LUE model is good. The compilations of the studies from Central Asia have defined the known range of NPP in the dry steppe biome, from 126 g C/m²/year to 326 g C/m²/year [28,29,30] and the semidesert short grassland from 114 gC/m²/year to 220 gC/m²/year [31,32]. These values are highly consistent with values obtained in the present study, 145 gC/m²/year and 131 gC/m²/year for dry steppe and semidesert grassland, respectively. The NPP estimates for semidesert in Central Kazakhstan by [33] has in some extent higher values, from 90 to 310 gC/m²/year. This may be explained by more favourable precipitation conditions, which were observed in the years of their field investigations.

The results of the present study go well with the scope of NPP values averaged for grassland biomes by [60] and [61]. Our results are also consistent with NPP and GPP estimations in other regions such as South of the USA [20] and Sahel in Niger [22]. Of particular interest is a quantitative comparison of NEP estimates. However, as NEP estimations need collection of much more field data than NPP/GPP estimations, they have been carried out relatively not often. Recent literature is characterized by a scarcity of NEP estimates from grasslands in Central Asia. A comparison with other regions revealed that the magnitude of NEP obtained in the present study was, in general, consistent with the NEP values in grasslands of Wyoming, USA, and former Czechoslovakia [20,62,63]. These results boost confidence in the NPP, GPP and NEP values obtained in our LUE model with values from the previous studies in the study region and in other regions with similar vegetation classes. However, the comparison of the present results with the results of previous studies should be taken with a caution, because there is a lot of variability within the results of previous studies. Even though, that the mean values of NPP, GPP and NEE from this study are well within the range of previous grassland studies, the large variability within the previous studies puts the comparison into perspective.

Table 4. Previous estimates of primary production (GPP and NPP) and net ecosystem production (NEP) in arid and semi-arid zones.

Location and vegetation type	Estimating technique	GPP (gC/m²)	NPP(gC/m²)	NEP*(gC/m²)	Reference
Central Asia Dry steppe Dry steppe Dry steppe Semidesert Semidesert Semidesert Semidesert	Field observations		326 126 148 90-310 117-189 220 114		Pershina & Yakovlewa [28] Makarova [29] Tyurmenco [30] Robinson et al. [33] Tyurmenco [30] Fartushina [31] Gristchenco [32]
Former Czechoslovakia Grassland Grassland	Field observations	664 492	481 318	18 21	Tesarova & Gloser [62] Rychnovska et al. [63]
Global Grassland Grassland	Field observations		70-410 91-385		Zheng et al. [60] Rodin et al. [61]
Wyoming, USA Mixed-grass prairie Sagebrush steppe	Measurements of net ecosystem exchange	321 239		14 42	Hunt et al. [20]
Sahel, Niger Grassland	Satellite-based LUE model	352.1	169.3		Seaquist et al. [19]
Central Kazakhstan Dry steppe Semidesert	Satellite-based LUE model	243.70 211.03	145.61 131.30	-0.09 13.17	The present study

*Positive NEP indicates that site was a sink for carbon, and negative NEP indicates the sites was a source for carbon to the atmosphere.

Table 4. Previous estimates of primary production (GPP and NPP) and net ecosystem production (NEP) in arid and semi-arid zones.

Location and vegetation type	Estimating technique	GPP (gC/m²)	NPP(gC/m²)	NEP*(gC/m²)	Reference
Central Asia Dry steppe Dry steppe Dry steppe Semidesert Semidesert Semidesert Semidesert	Field observations		326 126 148 90-310 117-189 220 114		Pershina & Yakovlewa [28] Makarova [29] Tyurmenco [30] Robinson et al. [33] Tyurmenco [30] Fartushina [31] Gristchenco [32]
Former Czechoslovakia Grassland Grassland	Field observations	664 492	481 318	18 21	Tesarova & Gloser [62] Rychnovska et al. [63]
Global Grassland Grassland	Field observations		70-410 91-385		Zheng et al. [60] Rodin et al. [61]
Wyoming, USA Mixed-grass prairie Sagebrush steppe	Measurements of net ecosystem exchange	321 239		14 42	Hunt et al. [20]
Sahel, Niger Grassland	Satellite-based LUE model	352.1	169.3		Seaquist et al. [19]
Central Kazakhstan Dry steppe Semidesert	Satellite-based LUE model	243.70 211.03	145.61 131.30	-0.09 13.17	The present study

*Positive NEP indicates that site was a sink for carbon, and negative NEP indicates the sites was a source for carbon to the atmosphere.

5.3.2. Direct comparison with ground-based data

Next, we compared the modelled GPP with ground-based GPP at the 14 test sites. The main problem, which limits a comparing the field data with the modelling results, is a significant mismatch in the spatial resolution of the SPOT-VGT data (1-km) and the test sites. However, we suggested that the problem of the scale difference could not hinder the assessment of the model through the direct comparison. The reasons for this suggestion are the following:

The scale problem had been partly solved by a collection of several replicates at each test site. An average of these replicates for each test site was than computed.

The ground surface at the most test sites and around them was relatively homogeny. We examined the homogeneity of the surface using a fine-resolution satellite image (Landsat ETM+) and considered that the variability of spectral reflection of the ground surface around the test sites was relatively low.

In Figure 4, the scatter plot between the peak aboveground biomass and the growing season GPP is presented. The accumulated GPP accounts for more than 82% of all variations in the biomass (p < 0.001). Taken into account the mismatch between spatial resolutions of the source data, this result is very positive for a regional modelling approach.

Figure 4. Relationship between the growing seasons GPP modelled in this study and ground-based GPP. Points represent the observed GPP collected at the 14 test sites across the study area.

Figure 4. Relationship between the growing seasons GPP modelled in this study and ground-based GPP. Points represent the observed GPP collected at the 14 test sites across the study area.

5.3.3. Comparison with MODIS GPP product

The final evaluation of the modelling results was done using a pixel-by-pixel comparison between the growing season GPP modelled in this study and the MODIS GPP product for the growing season 2004. A scatter plot between these both data sets is shown in Figure 5. The figure proofs strong relationship between the both GPP products, even though individual pixels show relatively high scattering in the plot, particularly at high values of GPP. The reason for the scattering may be a number of distinct artefacts within the MODIS GPP product, which cause a mismatch between the MODIS GPP and the GPP from this study. The global approach used for construction of the MODIS GPP cannot model the conditions of specific regions in details. The MODIS approach uses global climate data with a spatial resolution of 1° × 1° to parameterise the growing conditions. All other input variables are considered as biome-specific constants, which ignore local variations in environmental factors and vegetation conditions [37].

A visual analysis of the scatter plot in Figure 5 proved that the relationship between the GPP modelled in this study and the MODIS GPP deviated from the 1:1 line. As a consequence, the model presented in this study underestimates lower values in the MODIS data set and overestimates higher values of MODIS GPP. Our model has a wider value range than the MODIS model. This is also reflected in descriptive statistics calculated for the data sets. The MODIS GPP ranged from 87.07 to 474.3 g C/m² with a mean value of 201.4 and a standard deviation of 38.1 g C/m², while the GPP from this study varied between 43.86 and 494.1 g C/m² with a mean value of 220.4 g C/m² and a standard deviation of 55.5 g C/m². With regard to the high consistency between the GPP predicted by our model and the ground-based GPP (Section 5.3.2 and Figure 4), the difference between the mean values may suggest that the MODIS algorithm appears to slightly underestimate the growing season GPP in the grassland biome of Central Asia. This suggestion agrees with the results reported by [64] who also found an underestimation of GPP in the grassland biome in the semi-arid Sahel by the MODIS algorithm.

Figure 5. Pixel-by-pixel relationship between MODIS GPP and the GPP modelled in this study. Points represent individual pixels (1 × 1 km, n = 31500). The thin line represents 1:1 line, while the thick one is the regression line between the data sets.

Figure 5. Pixel-by-pixel relationship between MODIS GPP and the GPP modelled in this study. Points represent individual pixels (1 × 1 km, n = 31500). The thin line represents 1:1 line, while the thick one is the regression line between the data sets.

6. Discussion and Conclusions

Based on SPOT-VGT, in situ measurements of carbon residing in the aboveground and underground vegetation and climate data, the values of gross primary production (GPP), net primary production (NPP), total respiration of ecosystem (R_e) and net ecosystem exchange (NEE) were calculated for a large grassland region in Central Kazakhstan. The presented study used the well-known Monteinth’s approach, but the most important advantage of the model presented is the exclusive use of the variables obtained from field survey data for model’s calibration and validation. The model is fixed in its regional scale with the spatial resolution of 1 × 1 km and temporal resolution of 10 days.

Our analysis showed that for the investigated year 2004 the short grassland was a net carbon sink, while the steppe grassland was a net carbon float. The average accumulation of carbon in the short grassland was 13.17 g C/m²; this resulted in 263,768 tons over the whole area of the short grassland. The steppe grassland was losing carbon with a low rate of –0.09 g C/m². One of the reasons for the carbon loses by the steppe grassland was a higher rate of respiration and lower NPP modelled in this zone. Other reasons were not investigated in this study and are objective of an ongoing research. However, basing on the conclusions by [65], we can suggest that climate, particularly precipitation amount, should be the major driving force for dynamics of carbon flux. [65] concluded that globally about 44% of interannual variability in NEP resulted from precipitation, about 16% from temperature and about 12% from soil carbon.

In general, the evaluation of the modelling results reveals a high reliability. On the one hand, the comparison of the results with directly measured field data suggested a strong relationship between the values of GPP derived by the model and the peak aboveground carbon amount at 14 test sites. No attempts have been made to quantify the uncertainty of the direct comparison, but as in all cases of evaluating coarse resolution satellite based primary production models by harvest measurements some uncertainty must be expected due to differences between the spatial unit size of the ground measurements and the coarse resolution of the satellite sensor [66,67,68,69,70,71]. In the present study, the total size of the samples taken in each of the 14 test plots was only 256 m² (four samples of 1 x 1 m per replicate and four replicates per plot). It encompasses 0.0256% of a MODIS pixel. Without doubt, the surface of an area covered by one MODIS pixel has more variations than the surface covered by one test site. Fortunately, the heterogeneity of surface for grassland is much less than it is for forested ecosystems. The study by [72] estimated magnitude of surface heterogeneity and its influence on the results of LAI modelling with MODIS data. They found that the surface heterogeneity in grassland of the Sahel within one MODIS pixel is low, and its influence can be neglected by direct comparison of MODIS data with field measurements. Taking some uncertainty related to scale differences between observed biomass and satellite data into account the model presented in this study performs well for the grassland of Central Kazakhstan.

On the other hand, the model bears comparison with the global MODIS GPP dataset. This comparison reveals a general similarity in the spatial pattern of the MODIS GPP and GPP from this study, even though the first dataset includes a number of block artefacts and some data gaps at the level of individual pixels. The average value of GPP for the MODIS dataset was slightly underestimated (about 8%, considering mean values of the two data sets) in comparison to the dataset in this study. The comparison with data from recent literature on GPP/NPP in similar vegetation types also proofed the reliability of the modelling results. Estimates of the NEE are more difficult to evaluate than the GPP or NPP because of the greater uncertainty about the modelled NEE. Quantifying the NEE requires estimates of carbon budget components that each has associated errors. NEE can range from 10-50% of NPP depending on the soil carbon pools available for decomposition [73]. The results of the present study revealed the value of NEE for short grassland accounts to about 10% of NPP, while steppe grassland was approximately in balance with a slightly negative value. For the modelled values, one of the greatest uncertainties is the amount of soil respiration. In this study, soil respiration was calculated using empirical equation for estimation soil respiration at the global scale [55]. Obviously, more accurate results would be achieved in the presented model if we had had in situ measurements of soil respiration from the study area, for example, using the chamber approach. Such measurements are planned for the prospected research work in the study area. Certainly, they will greatly improve the model’s accuracy.

Remotely sensed data provide necessary inputs, such as land cover categories and leaf area index, but other required spatial data such as root-leaves ratio, litter-leaves ratio, base rate of maintenance or growth respiration etc. are poorly known and considered as biome-dependent constants, which limit the accuracy of model predictions. In the model presented here, most of these constants were calculated from the in situ measurements of carbon residing in the vegetation and soil at a number of test sites. This improves the general accuracy of the model. However, the values of these constants were averaged over the two land cover categories without taking into account spatial variance of them between localities. Obviously, values of all input variables might demonstrate significant spatial variability within each land cover category. This variability was not incorporated into the model presented in the study. Additional studies are required in order to investigate possibilities of mapping the input variables at the per-pixel scale. The use of such information would significantly improve modelling results at localities.

Finally, considering the applicability of the model, we should also take the different time scales of the application into account. In this study, in situ measurements of biomass were carried out at the end of June. This period is considered to represent the peak of the growing season in the study area. Normally, the middle part and the end of June in this region are characterized by the peak of precipitation [30,74]. However, there is a lot of variability in distribution of precipitation throughout the growing season, and this intra-annual variability is overlapped with inter-annual variability of precipitation. An additional source for uncertainty in localization of the growing season peak is a different reaction time of the vegetation to precipitation. Previous studies proved that vegetation types in Central Asia and the study region are characterized by different time-lags in their response to precipitation. Thus, the short grassland in the region of the study shows the peak of plant growth about 10-20 days after the precipitation peak, while the time-lag for the steppe grassland accounts to about 20-40 days [74]. For this reason, it is difficult to forecast exactly the peak of the growing season for a certain year before making in situ measurements of the peak biomass. This may limit to some extent the applicability of the LUE model based the one-time measurements of biomass. To reduce the uncertainty in the model’s applicability, data on biomass from different periods within the growing season are needed.

The findings of the study serve to a better understanding of carbon cycle in dry lands of the interior Eurasia and should play an important role in the establishing of an appropriate model for calculation of carbon assimilation in grasslands of Kazakhstan for annual reports for the Kyoto Protocol.