Impacts of land use and land cover change on regional climate: a case study in the agro-pastoral transitional zone of China

The WRF model was configured for fine-scale simulations with three nested domains (figure 1(a)). The coarse outer domain (D01) covers almost all of China and extends to the East and South China Sea to capture the influence of the East Asian Monsoon on the climate of APTZ, and covers a total area of 5310 km × 5490 km with 45 km grid spacing in both horizontal directions. The intermediate domain (D02) covers an area of 3180 km × 3960 km, using 15 km grid spacing. The innermost domain (D03) is centered at 42.3°N and 112.8°E, and covers an area of 2000 km × 1950 km, using 5 km grid spacing. The main physical parameterizations used for all simulations are presented in table A.1. A 10 yr (1 January, 2001–31 December, 2010) averaged global final analysis data obtained from the National Centers for Environmental Prediction were used as initial and lateral boundary conditions. The model was executed using one-way nesting and all simulations were initialized on 1 January at 0000 UTC and run through 31 December at 1800 UTC for a period of 1 yr, with the first 15 d considered spin-up and therefore excluded from the analysis.

Although the coupled WRF-Noah Land Surface Model (LSM) is widely used to simulate regional climate, deficiencies exist with prescribed land surface biophysical parameters. In Noah LSM, the temporal variability of vegetation parameters is determined by changes in vegetation fraction and corresponding annual maximum/minimum values prescribed in look-up tables. LAI and emissivity vary in proportion to vegetation fraction, whereas albedo varies conversely with it. However, there is no evidence showing such linear relationships, and in most cases, albedo is expected to increase as vegetation matures to peak greenness (e.g., Song 1999). In addition, the maximum/minimum LAI, albedo, and emissivity defined in the tables are not suitable for our study area. Thereby, large errors may occur when using the default land surface properties as prescribed by Noah LSM, for longer time-scale simulations (see figures A.1 and A.2).

WRF also allows users to define vegetation fraction, albedo, and LAI from geographical data sets produced based on remote sensing data at a global scale. However, the default data sets utilized also have limitations: (1) the spatial resolution (e.g., 0.144° × 0.144°) is too coarse to describe landscape heterogeneity for high-resolution simulations; (2) interannual changes of vegetation characteristics are omitted. For example, the default vegetation fraction was obtained during the period of 1986–1991, and has not been updated since, despite potentially large changes in the amount of actual vegetation fraction (see figures A.1 and A.2). Therefore, in this study, we used finer spatial resolution with more detailed biogeographic information to improve characterization of landscape properties (see below).

We obtained land use/cover data from the Data Sharing Infrastructure of Earth System Science (http://www.geodata.cn/). The data are produced by the Chinese Academy of Sciences (CAS) for 2000 (hereafter CAS2000) and 2010 (hereafter CAS2010), with a spatial resolution of 1 km × 1 km. The data have been converted to the International Geosphere–Biosphere Program land use classification scheme and the overall accuracy of the classified data is 83.14% (Wu et al 2013). The land use/cover map corresponding to 2000 used in the representation of our innermost domain is shown in figure 1(b).

We computed vegetation fraction based on 16 d, 1 km MODIS vegetation indices data sets (MOD13A2; https://lpdaac.usgs.gov/) following the method proposed by Gutman and Ignatov (1998):

$$\begin{eqnarray}&&{\rm{FVC}}=(N-{N}_{{\rm{S}}})/({N}_{{\rm{V}}}-{N}_{{\rm{S}}}),\end{eqnarray} \tag{ 1 }$$

where ${\rm{FVC}}$ denotes the vegetation fraction, $N$ the NDVI at each pixel, ${N}_{{\rm{S}}}$ the bare soil NDVI, and ${N}_{{\rm{V}}}$ the dense vegetation NDVI. Here, we defined ${N}_{{\rm{S}}}$ and ${N}_{{\rm{V}}}$ as the lower and upper 5% NDVI of the study area (Sellers et al 1996), respectively. LAI, albedo, and emissivity were obtained from the global land surface satellite (GLASS; http://glass-product.bnu.edu.cn/), with a temporal resolution of 8 d and spatial resolution of 1 km × 1 km. It should be noted that the GLASS albedo products provide both white-sky and black-sky albedo, which are numerically related to each other as follows:

$$\begin{eqnarray}&&{\rm{BHR}}=({\rm{1}}-s){\rm{DRF}}+s\mathrm{BHR}\_\mathrm{diff},\end{eqnarray} \tag{ 2 }$$

where $s$ is sky light ratio factor determined by proportions of direct beam and diffuse illumination, ${\rm{BHR}}$ the blue-sky albedo, ${\rm{DRF}}$ the black-sky albedo, and $\mathrm{BHR}\_\mathrm{diff}$ the white-sky albedo (Liang et al 2005). In our study, the weight coefficient, $s,$ was calculated as a function of solar angle as per Long and Gaustad (2004).

The semi-monthly and 8 d data were then linearly interpolated to a daily interval and aggregated to the innermost domain resolution of 5 km × 5 km. The four land surface biophysical parameters along with the CAS land use/cover data were incorporated within the innermost domain. For the two outer domains in all simulations, we used vegetation fraction, LAI, albedo, and land use/cover data from the geographical data set provided by the WRF model, whereas emissivity was prescribed by Noah LSM.

We designed two simulation experiments to explore the effects of human-induced LUCC on regional climate. Experiment2001, corresponding to the 2001 landscape pattern, used CAS2000 (in lieu of 2001 land use/cover data, which were not available) and 2001 land surface biophysical parameters (vegetation fraction, LAI, albedo, and emissivity). Experiment2010, corresponding to the 2010 landscape pattern, used CAS2010 and 2010 land surface biophysical parameters (vegetation fraction, LAI, albedo, and emissivity). In addition, two control experiments were carried out to examine whether the use of non-default land surface biophysical parameters that we derived improved model performance. Control2001 used CAS2000 while Control2010 used CAS2010, and both utilized the default land surface biophysical parameters prescribed by WRF-Noah LSM based on tabular values.

The observation data used for evaluating model performance were obtained from 'SURF_CLI_CHN_MUL_MON', a data set provided by China meteorological data sharing service system (http://data.cma.gov.cn/). The data set contain monthly climatological information (pressure, temperature, wind, precipitation, etc) for the meteorological stations across China. For model evaluation, we used observations from 85 stations, which were located in APTZ. To examine the spatial skill of the WRF model and directly compare model simulation results performed at 5 km grid spacing, mean air temperature during summer (i.e., June, July, and August) and winter (i.e., December, January, and February) from 2001 to 2010 were calculated for each station. The 10 yr mean observations were then spatially interpolated to a 5 km × 5 km resolution using the Kriging method (Kalnay 2003).

In addition, the observation minus reanalysis (OMR; Kalnay and Cai 2003) was used to evaluate the temperature changes due to LUCC as simulated by the WRF model. The data used for OMR analyses include the gridded reanalysis data described by Kanamitsu et al (2002) and station-based observations obtained from the data set described above. Prior to use of the OMR analyses, the station-based observations were homogenized using the method proposed by Li and Yan (2009). Temperature trends during summer and winter seasons based on time series of 13 years from 1999 to 2011 were calculated for each station within APTZ by using observations and reanalysis, respectively. The differences in trends of temperatures between those derived from observations and those derived from reanalysis can then evaluate temperature changes induced by LUCC.

Statistical tests were carried out for summer and winter temperature differences simulated by the two experiments (Experiment2010—Experiment2001). As the Experiment2001 and Experiment2010 simulations were carried out for a period of 1 yr respectively, all the grid cells representing temperature differences in summer (and winter) for the Experiment2001 and Experiment2010 comparison were considered as one sample set. Only grid cell values within two standard deviations of the mean of the sample set were considered as confident differences and thus retained for the subsequent analyses (Kalnay 2003). We further enhanced the criteria defined here by additionally requiring temperature change exceeding 0.2 °C for the Experiment2001 and Experiment2010 comparison.

Model evaluation was carried out during winter and summer seasons, and the results indicate improved model performance with our newly developed landscape representation relative to the WRF default representation (figure 2). During winter, the Control2001 simulation produces a warmer regime in Northern APTZ and a cooler regime in central areas compared to the observations, and the spatial correlation coefficient of air temperature between the observations and Control2001 simulation is 0.824. Although deficiencies still exist, the Experiment2001 simulation largely corrects the deviations of the Control2001 simulation, with a spatial correlation coefficient rising to 0.941. During summer, the Control2001 simulation produces warmer air temperature in central North and Southern APTZ, and a cooler regime in the Northeast, with a spatial correlation coefficient equal to 0.627. In contrast, the Experiment2001 simulation better captures the general pattern of the observed temperature, with a spatial correlation coefficient increasing to 0.732.

Zoom In Zoom Out Reset image size

Because the limited number of stations may not adequately represent the spatial variability of temperature, particularly during summer, additional model evaluation was carried out by calculating bias and root mean square error (RMSE) between the observations and simulations (table 1). Model bias was estimated by comparing the averaged station-based observations within APTZ to the corresponding average of simulated grid cells, which were nearest to the station locations, whereas RMSE was computed to measure the discreteness between the observations and simulations. Results show that the model bias and RMSE revealed by the Experiment2001 are smaller than those revealed by the Control2001 in both seasons. Model evaluation for the Experiment2010 and Control2010 simulations further indicates improved model performance with our landscape data (see figure A.3 and table A.2).

LUCC includes changes in both land use/cover types and land surface biophysical properties. The major land use/cover changes in APTZ during 2000–2010 were conversions between cropland and grassland, and conversions between grassland and barren land (Wu et al 2013; figure 1(c) and table 2). Compared to changes in land use/cover types, changes in land surface biophysical parameters were much more remarkable (figure 3).

Zoom In Zoom Out Reset image size

During summer, vegetation fraction and LAI increased over most regions, although this amplification was most evident for Southern portions, with vegetation fraction rising by more than 20% and LAI by more than 0.7. In contrast, both summer vegetation fraction and LAI decreased substantially in Northeastern Inner Mongolia, with vegetation fraction decreasing by more than 10% and LAI by 0.5 on average. During winter, because of vegetation defoliation, variations in vegetation fraction and LAI were generally modest. On the contrary, changes in albedo and emissivity were greater in winter than in summer. During winter, albedo decreased considerably for most parts of APTZ, with the maximum values exceeding 0.1, while patchy increases, in excess of 0.08, took place in Inner Mongolia. The decreases in albedo over the study area are mainly due to the conversion of cropland to grassland in the past decade. Decreased summer season albedo mainly occurred in the South with a range of 0.02–0.04. Unlike albedo, emissivity increased noticeably for most areas with maximum increases exceeding 0.025, while patchy decreases of more than 0.015 occurred around central and Northern areas.

We use differences in summer and winter temperature between Experiment2001 and Experiment2010 (Experiment2010–Experiment2001) to estimate the potential impacts of LUCC on regional climate of APTZ. Results show that differences in land use and land cover patterns led to widespread decreases in summer temperature (figure 4(a)) and extensive increases in winter temperature (figure 4(b)). During summer, the most noticeable decreases in 2 m air temperature occurred in Southern and Northernmost regions, with a magnitude of change on the order of 1 °C, while the most noticeable increases were located in Northeastern Inner Mongolia, with local warming of lesser magnitude (generally less than 0.4 °C). Unlike the variations in summer temperature, winter temperature increased by more than 0.8 °C in parts of the South and the North, but decreased by about 0.6 °C in central and Northern APTZ.

Zoom In Zoom Out Reset image size

Differences in summer and winter temperature between the two controls (Control2010–Control2001) were comparatively smaller. Although differences in summer temperature (figure 4(c)) were more pronounced relative to winter temperature (figure 4(d)), the spatial pattern of temperature variations, showing widespread warming in Northern and central areas and patchy decreases across APTZ, differed from that revealed by the two experiments. The APTZ-averaged changes in summer, winter, and annual mean temperatures from the two controls were close to zero with insignificant variability, whereas those from the two experiments indicated evident summer and annual cooling effects (figure 5). Overall, the net effect of simulated regional climate change due to LUCC was to decrease the annual range of temperature by lowering maximum (i.e., summer season) but increasing minimum (i.e., winter season) 2 m air temperature.

Zoom In Zoom Out Reset image size

The simulated results coincide with the OMR analyses, which also indicate non-negligible summer cooling (figure 4(e)) and widespread winter warming (figure 4(f)). During winter, 41 out of 51 stations show the same signal of temperature changes with Experiment2010–Experiment2001 over APTZ (only those exceeding 0.2 °C were retained). For the 41 stations, the correlation of temperature change between the OMR and Experiment2010–Experiment2001 is statistically significant at 99% confidence level, and the Pearson correlation coefficient, r, between the two is 0.743. Similarly, 35 out of 43 stations show the same signal of temperature changes during summer, and the r value between the two is 0.658 at 99% confidence level among the 35 stations.

As illustrated by figure 6, simulated differences in spatial variability of net shortwave radiation and albedo displayed considerable similarity. In regions where albedo increased, the net shortwave radiation flux decreased accordingly, and vice versa. Net longwave radiation increased in the summer except over Northeastern Inner Mongolia, and decreased in winter, particularly in Southern areas. The total net radiation flux of the land surface is defined as the sum of net shortwave and longwave radiation flux. During summer, increased net radiation flux in excess of 15 W m⁻² was prevalent over Southern regions, while patchy decreases of relatively smaller magnitude occurred over Northeastern Inner Mongolia. During winter, the net radiation flux decreased by about 10 W m⁻² over central APTZ and increased by the same magnitude for other regions.

Zoom In Zoom Out Reset image size

Variations in surface radiation flux resulted in differences in surface energy budget (figure 7), which drive the changes in the 2 m air temperature presented earlier. During summer, in Southwestern APTZ, sensible heat flux decreased by about 15 W m⁻² on average with local maximum decreases reaching 25 W m⁻², whereas latent heat flux increased by more than 25 W m⁻² with local maximum increases approaching 35 W m⁻². Opposite changes were found in Northeastern Inner Mongolia: increases in sensible heat flux and decreases in latent heat flux exceeding, in absolute value terms, 10 W m⁻². During winter, changes in sensible heat flux were reduced relative to summer, whereas latent heat flux impacts were too small to be detected. The spatial variations of sensible heat flux match those of net radiation quite well (figures 6(f) and 7(b)).

Zoom In Zoom Out Reset image size

Our results show that LUCC between 2001 and 2010 considerably affected the 2 m air temperature of APTZ, via alterations of land surface biophysical parameters (the parameters we examined include: vegetation fraction, LAI, albedo, and emissivity). Results demonstrate changes in landscape patterns lowered the summer temperature but raised the winter temperature, through modification of the net radiation and surface energy budget. During summer, although the reduced albedo led to rising net radiation, increases in vegetation fraction and LAI substantially affected the partitioning of available radiation into latent heat flux, thus decreasing the Bowen ratio (i.e., the ratio of sensible to latent heating) with additional consequences for near-surface cooling. However, the increases in wintertime net radiation, while smaller in magnitude compared to summer, contributed mainly to changes in sensible heat flux because of vegetation defoliation. This effect led to an increasing Bowen ratio, and therefore, near-surface warming.

Our results provide valuable information for landscape planning to mitigate and adapt to climate change, particularly in arid and semiarid regions such as APTZ, whose environmental and socioeconomic conditions are sensitive to climate change. Generally, vegetation recovery due to the Grain to Green Program in Southern APTZ led to overall cooling effects, whereas grassland degradation in Inner Mongolia, largely due to overgrazing in the past decade (Wu et al 2015), resulted in an overall near-surface warming. Compared to grassland and cropland, however, changes in vegetation fraction and albedo are minor over the woody areas, thus leading to less pronounced temperature changes in both seasons. This suggests that replacement of grassland by woody plants in arid and semiarid regions where the environmental conditions restrict the growth of woody plants may not be a rational choice for land system planning in terms of climate mitigation and adaptation.

Our study has demonstrated that using the default biophysical parameters prescribed by WRF-Noah LSM, as shown by the two controls, failed to accurately capture the actual effects of LUCC on regional climate. Obviously, the output of any modeling system will be sensitive to changes in biophysical parameters, such as vegetation fraction, LAI, albedo, and emissivity. Because of strong effects of seasonality, vegetation phenology changes in both time and space, and these changes affect vegetation fraction, LAI, albedo, and emissivity. The default parameter set did not consider these changes, and was based on a coarse-grained land use classification in which each land use type may include multiple distinct land cover types with substantially different biophysical properties.

As a result, the model underestimated changes in the amount and spatial heterogeneity of land use and land cover and associated land surface properties, further leading to incorrect partitioning of available energy between sensible and latent heat fluxes, and consequently prediction of changes in temperature of reduced skill in our study region. The discrepancies between the control and experimental simulations indicate that the composition and configuration of local landscapes must be explicitly considered in assessing the impacts of LUCC on climate at regional and larger spatial scales (Pielke and Avissar 1990, Raupach 1991, Pielke and Niyogi 2010) for seasonal or longer time scale simulations.

Our study provides the first comprehensive assessment of impacts of LUCC on the spatiotemporal pattern of temperature in the extensive agricultural-pastoral transitional zone of North China. However, several aspects of the study can, and will be improved in future work. First, although our study demonstrates that spatial patterns of local landscapes may substantially affect the simulation results of WRF, it is not clear how much spatial detail of the landscape should be incorporated for studies of complex dynamic land mosaics such as APTZ. In other words, how small should the grid cell be for execution of the WRF model—or other similar modeling systems—in such studies and what value can be attained via the utility of, for example, a mosaic approach that accounts for sub-grid variability (e.g., Li et al 2013)? Second, comprehensive sensitivity analysis of WRF is needed to quantify changes in model output with respect to changes in the spatial and thematic resolutions of the input land use and land cover maps. Due to limitations of data and computational resources, our study only considered landscape patterns in two separate years (2001 and 2010), and our findings need to be further confirmed by future research that considers more years and longer time horizons (e.g., decadal timescales or longer). With improving computing capacities and data sets, efforts are underway to address these issues.