Spatial-Temporal Dynamics of China’s Terrestrial Biodiversity: A Dynamic Habitat Index Diagnostic

Abstract

Biodiversity in China is analyzed based on the components of the Dynamic Habitat Index (DHI). First, observed field survey based spatial patterns of species richness including threatened species are presented to test their linear relationship with remote sensing based DHI (2001–2010 MODIS). Areas with a high cumulative DHI component are associated with relatively high species richness, and threatened species richness increases in regions with frequently varying levels of the cumulative DHI component. The analysis of geographical and statistical distributions yields the following results on interdependence, polarization and change detection: (1) The decadal mean Cumulative Annual Productivity (DHI-$\overline{cum}$ < 4) in Northwest China and (DHI-$\overline{cum}$ > 4) in Southeast China are in a stable (positive) relation to the Minimum Annual Apparent Cover (DHI-$\overline{min}$) and is positively (negatively) related to the Seasonal Variation of Greenness (DHI-$\overline{sea}$); (2) The decadal tendencies show bimodal frequency distributions aligned near DHI-$\overline{min}$~0.05 and DHI-$\overline{sea}$~0.5 which separated by zero slopes; that is, regions with both small DHI-min and DHI-sea are becoming smaller and vice versa; (3) The decadal tendencies identify regions of land-cover change (as revealed in previous research). That is, the relation of strong and significant tendencies of the three DHI components with climatic or anthropogenic induced changes provides useful information for conservation planning. These results suggest that the spatial-temporal dynamics of China’s terrestrial species and threatened species richness needs to be monitored by first and second moments of remote sensing based information of the DHI.

1. Introduction

Global extinction rates have increased due to anthropogenic activities [1,2,3,4,5], which have led to an erosion of ecosystem services [6]. Efforts to reduce this rapid biotic impoverishment have been carried out by many countries in terms of enacting, or planning to enact, legislation for the identification, protection, and recovery of endangered species [5,7,8,9]. The Chinese government has initiated several environmental programs to increase its commitment to environmental protection, including the National Forest Conservation Program (NFCP, 2000; see [10]) and the Grain-to-Green Program (GTGP, 1999; see [11]). However, these programs are carried out statewide so that methods of implementation as well as their effectiveness vary among different areas [12].

Regional evaluation and measurements for habitat conservation, species and endangered species richness still remain poorly integrated into these programs, particularly outside nature reserves [13,14]. Therefore, reliable information about the geographic distribution of species richness and endangered species is necessary. At present, there are three methodologies for measuring species richness and endangered species from sources: they are obtained from field surveys (also called the meta-analysis approach, see [15]); estimated from existing land cover maps by combining information about the known habitat requirements of species [16,17,18], e.g., the USGS National Gap Analysis Program [19]; and derived from remote sensing based indices that consider the relationship between spectral radiance and species richness, e.g., Dynamic Habitat Index DHI [20,21,22]. For applications to extensive and complex land surface areas these methodologies show the following characteristics:

(1) The meta-analysis approach uses the emergent outcome of numerous tests by compiling the results of previously published analyses [15]. It assumes that different areas of ecology share the same sampling system and design properties, although their regional ecologies are virtually unique from study to study [23].

(2) Land cover map analyses extract species occurrences by considering the composition, abundance and distribution requirements of individual species or assemblages with maps of land cover [24]. The pace of land habitat loss and conversion (provided by remotely sensed imagery) allows estimates of the related biodiversity loss [16,17,18]. Compatibility of current state analyses and future projections, and/or of different spatial resolution of land cover systems could be an issue, in particular, when assumptions about missing information have to be made. The project “Geo-Wiki” provides a platform where citizens help in validating global land cover maps [25,26,27], which may lead to improvements in this methodology by using validated land cover maps.

(3) Remote sensing based Dynamic Habitat Index-based analyses propose a relationship between spectral radiance and species richness. The usefulness of DHI has been tested in previous validation research using actual biodiversity data [20] (originally proposed by Mackey et al. [21] and Berry et al. [28]), which show that there is a direct linear relationship between vegetation productivity and the richness of specific species. DHI analysis has been applied to Canada [20,22,29] and North American [30]. Spatial and seasonally induced variation in DHI are closely related to the patterns of avian species richness with highest correlations up to 88% [31]. Furthermore, in combination with long-term climate change, DHI can be used to assess the long-term variability of vegetation biodiversity [32]. However, the index is for the first time being applied in China. Assuming that DHI constitutes a surrogate species richness without testing relationships with actual species richness may be incorrect and could be misleading when expanding the analysis to endangered species richness.

The relationship between DHI components and species richness has been demonstrated in a number of studies to estimate the richness of breeding bird species in the United States and Canada [22,30,31,33]. In this study we focus on estimating the species distribution, its change and variability using DHI information. In addition, we expand the original DHI components (related with avian species) to: (1) a relationship that combines mammals, birds, and amphibians; (2) a grid-scale (10 × 10 km²) based analysis for further enhancing the spatial resolution; and (3) a dynamic change and variability analysis associated with the static species richness information obtained from annual climate means.

Therefore, we first revisit the DHI components (in Section 2) and then test the relationship between the spatial patterns of indices derived from remote sensing information as possible surrogates of biodiversity and of the actual overall species richness and threatened species richness in China (in Section 3). This includes the distributions of species richness and threatened species richness and the distributions of change tendencies, from 2001 to 2010. A discussion and summary conclude the analysis (Section 4).

2. Data and Methodology

2.1. Data

The Fraction of Absorbed Photosynthetically Active Radiation (fPAR) time series from MODIS (MOD15A2, 2001–2010) is used in this study. The fPAR data is provided every 8 days at 1 × 1 km² spatial resolution as a gridded level-4 product in the sinusoidal projection. Data processing involves the calibration of sensor changeover and corrections to remove the effects of sensor degradation, satellite orbit drift, solar zenith angles, and other factors [34]. Species and threatened species census (2013, types/10 × 10 km²) of mammals and amphibians from the International Union for Conservation of Nature, and of birds from BirdLife International and NatureServe [35,36] are added to account for the relationship between field survey datasets and the remote sensing derived surrogate of species richness and threatened species richness (see Figure 1, data available from the website [37]).

2.2. Dynamic Habitat Index

The dynamic habitat index (DHI-(cum,min,sea)) is a composite vector deduced from fPAR time series, representing the underlying vegetation dynamics [20], where fPAR is expressed as a non-dimensional fraction of incoming radiation received by the land surface. Monthly maximum fPAR-values on 1 × 1 km² resolution are the basic input data set to compute the three relevant annual indices for the subsequent habitat analysis (Section 3).

(1) Cumulative Annual Productivity (DHI-cum), which is estimated by summing all monthly fPAR values for each year, indirectly explains species distribution and abundance, because it characterizes the supply of resources like food supply [22,38,39]. This is an annual total of all monthly productivity contributions:

$$DHI-cum={{\displaystyle \sum }}_{month}MA{X}_{layer,fPAR}$$

(1)

(2) Minimum Annual Apparent Cover (DHI-min) represents the lowest (minimum value) level of vegetative cover in a year, it indicates the capacity of the landscape to sustain an adequate level of food and habitat resources during this period [22,40]. This annual value represents the monthly minimum productivity of the year:

$$DHI-min=MIN\left\{{\left(MA{X}_{layer,fPAR}\right)}_{month},\dots \right\}$$

(2)

(3) Seasonal Variation of Greenness (DHI-sea) refers to a surrogate of habitat quality associated with its relationship to natural resources (e.g., food, water and nutrients), it is expected to exert selective pressure on life history traits [41,42]. This annual index represents a coefficient of variability to characterize the month-to-month standard deviation in relation to its annual mean:

$$DHI-sea=\frac{STD\left\{{\left(MA{X}_{layer,fPAR}\right)}_{month},\dots \right\}}{MEAN\left\{{\left(MA{X}_{layer,fPAR}\right)}_{month},\dots \right\}}$$

(3)

where layer is an index representing individual fPAR dataset of 8 days temporal resolution. For each month, MAX_layer,fPAR represents a monthly maximum value calculated by the Maximum Value Composite methodology [43]. MIN, MEAN and STD are the minimum, mean and standard deviation of monthly maximum values in a particular year. A ten-year mean climatology (DHI-($\overline{cum},\overline{min},\overline{sea}$)) will be calculated to represent a long-term average state (or climate mean).

2.3. Regression Tendency Analysis

To estimate the long-term vegetation dynamics (2001–2010), a simple linear regression analysis by least squares [44,45] is applied to yearly DHI-(cum,min,sea).

$$Slop{e}_i=\frac{N\times {\sum }_{n=1}^{N}n\times DH{I}_i-\left({\sum }_{n=1}^{N}n\right)\left({\sum }_{n=1}^{N}DH{I}_i\right)}{N\times {\sum }_{n=1}^N{n}^2-{\left({\sum }_{n=1}^{N}n\right)}^2}$$

(4)

where Slope_i is the slope of the regression line of DHI_i ~ DHI-(cum,min,sea), $i\in \left\{1,2,3\right\}$. $n\in \left\{1,2,\dots ,10\right\}$ represents the number of year. The trend of DHI_i increases (decreases) when Slope_i is positive (negative).

2.4. Decadal Mean and Inter-Annual Variability

Transforming the input time series data into the DHI vector is a critical step for testing the usefulness whether the DHI constitutes a useful surrogate about species richness. A direct linear relationship between the DHI and particular species richness has been documented [31,33]. However, the inter-annual variability of long-term DHI components, namely the standard deviation of DHI-cum has, to our knowledge, not been implemented into this dynamic analysis nor have any suggestions been made on how to introduce it. Therefore, in this study, the following three-step procedure was developed to derive not only a remote sensing based surrogate of species but also to implement a novel procedure using the inter-annual coefficient of variability to obtain a surrogate of threatened species richness for global/regional application (application in China, see Section 3):

Step-1: Data-preprocessing: ArcGIS© 10.1 software was used to relate/combine and resample MODIS based DHI components (1 × 1 km²) in order to obtain a comparable spatial resolution with 10 × 10 km² of actual field survey based species and threatened species richness.

Step-2: Remote sensing surrogate for species richness: The decadal mean of DHI-cum is used to test the assumption that the DHI-cum constitutes a surrogate of species richness at the regional/continental scale for further classifying species richness in large areas and their changes through time based on a simple linear regression model: DHI-$\overline{cum}$ = a₁ × actual species richness + a₀. Linear regression analysis is employed to relate remote sensing-based information (DHI-$\overline{cum}$) to field survey based species richness.

Step-3: Remote sensing surrogate for endangered species richness: The standard deviation SD of DHI-cum is compared with the field survey data of richness of threatened species to analyze its potential as a remote sensing based surrogate of threatened species richness by employing a simple linear regression: SD(DHI-cum) = b₁ × actual endangered species richness + b₀. The confirmation of the linear relationship is achieved through the R² and p-value of the regression, which reveals that SD(DHI-cum) is a possible remote sensing surrogate for endangered species richness.

3. Application of the Dynamic Habitat Analysis: Results for China

3.1. Ecological Setting

China (see Figure 2) is the third largest country (~9.6 million km²), extending over 37 (62) degrees of latitude (longitude), and the elevation ranges from −154 to 8844 m. It comprises climatically, topographically, and biologically complex areas [46]: five different climates cover China from the cold-temperate, temperate, warm-temperate, subtropical to the tropical regime. The highest mountain range is the “Third Pole” or Tibetan Plateau in western China, from where several major rivers of Asia originate, such as Mekong, Brahmaputra, Yangtze and Yellow rivers [47,48]. Species richness in China contributes about one-eighth of all species on Earth [49]. Now, applying the tools introduced in the previous section provides further information [50] on productivity and species richness and the changes and trends of biodiversity in China (see also [51,52]).

3.2. Remote Sensing Surrogates for Species and Threatened Species Richness

First, geographical distributions of field survey-based species richness are shown (in Figure 1), which, including mammals, birds and amphibians, are the focus in this study. To target biodiversity monitoring strategies for China, we apply remote sensing based DHI statistics (first and second moments, and trends) to China to serve as surrogates for species richness and threatened species richness. That is, to provide biodiversity classification in large areas, temporal measures of tendencies and variability are tested with actual field survey information (Figure 3 and Figure 4).

Frequency distribution and box plot based statistics are presented to relate richness of species and threatened species from field survey and decadal averaged Cumulative Productivity DHI-$\overline{cum}$ (and inter-annual standard deviation of DHI-$cum$). Figure 3 and Figure 4 show a linear relationship using the mean value of each actual field survey interval for mammals, birds, amphibians individually, and sum of three, with high correlation R² ~ 0.86 (R² ~ 0.85, for more details see

Table 1. The correlation coefficients between field survey species richness and remote sensing based DHI-$\overline{cum}$ using the simple linear regression.

Species Richness	Slopes	R2	p Value	DF	Pixels
Mammals	0.4516	0.8538	5.34e-07	13	64817
Birds	0.7980	0.9289	1.125e-06	9	64803
Amphibians	0.5586	0.8124	0.0002	8	56219
ALL	0.5393	0.8611	1.062e-06	12	71380

and

Table 2. The correlation coefficients between field survey based threatened species richness and the standard deviation of DHI-cum using the simple linear regression.

Threatened Species	Slopes	R2	p Value	DF	Pixels
Mammals	0.0142	0.8766	7.989e-10	18	62598
Birds	0.0220	0.8372	9.895e-09	18	67034
Amphibians	0.0115	0.8491	9.364e-05	8	24433
ALL	0.0227	0.8448	4.367e-08	16	67891

). This linear relationship supports the assumption that DHI components constitute surrogates of species richness and threatened species richness, at least for mammals, birds, and amphibians in China. The geographical distribution of remote sensing and field survey data are displayed in Figure 5 and the following results are noted:

(1) Productivity and species richness are positively correlated with areas of high productivity supporting a relatively high level of species richness [54,55,56,57,58,59]. Remote sensing based productivity, more specifically DHI-cum shows a linear relationship with estimated species richness. Lack of vegetation productivity, or DHI-cum, limits the species spatial distribution due to restricted food resources and habitat.

(2) Similarly, remote sensing based variability of land surface productivity, SD(DHI-cum), shows a linear relationship with estimates of endangered species richness. In addition, variability of land surface productivity occurs in regions with relatively high species richness. As regions with high species richness are more resilient, they recover from perturbations [60] and thus provide a relatively stable environment for their biomes and vice versa.

(3) Species richness and human population density in China tend to concentrate in regions with high species richness, which decreases from the Southeast towards the Northwest. The potential species richness and human population density are positively correlated (see Figure 2 and Figure 6a), because both respond positively to increasing levels of primary productivity [61].

3.3. Dynamic Habitat Analysis in China

The linear relationship obtained in the previous section provides the underpinning of the subsequent dynamic habitat analysis. The geographical distributions of DHI-($\overline{cum},\overline{min},\overline{sea}$) components are displayed in Figure 6a–c:

(1) DHI-($\overline{cum},\overline{min}$) decrease from the warm and wet Southeast to the cold and arid Northwest, except for the mixed forest cover of the Greater and Lesser Khingan Mountains [62] and for the rainforest covered low-altitude regions on the Tibetan Plateau in the Shannan and Linzhi Prefectures [63].

(2) Megacities in eastern China (latitude 30 to 40°N) could not be distinguished from planted regions in the Northwest Plain when analyzing DHI-$\overline{cum}$ only (see light green regions). But DHI-$\overline{min}$ shows an obvious change (that is, a decrease, as confirmed by the scattered brown regions) separating megacities from planted regions and leaving only evergreen biomes at high values (dark green).

(3) DHI-$\overline{sea}$ associated with DHI-$\overline{cum}$, show low value is low in regions characterized by both high (southern evergreen forest regions) and low (northwestern sparsely vegetated regions) cumulative productivity. Furthermore, DHI-$\overline{sea}$ is positively related to the number of cropping times per year (Northeast > North > South China). For example, a relatively short growing season (May to September) and low winter temperature (−24 to −9 °C in January) in northeastern China allows only one crop per year [64,65,66,67]. But two crops per year occur in the North China Plain (−5 to 0 °C in January) governed by a long growing season from March to September [66,68,69]. Furthermore, South China has a longer growing season (February to December) and higher winter temperatures (0 to 16 °C in January) which allow two or three crops per year [70].

Frequency distributions (Figure 6d–i) allow a comparison of all DHI components:

(1) DHI-$\overline{cum}$~4 separates the relationship with DHI-($\overline{min}$, $\overline{sea}$). DHI-$\overline{min}$ (or $\overline{sea}$) is independent of (or positively related to) DHI-$\overline{cum}$ where DHI-$\overline{cum}$< ~4 but positively (or negatively) related where DHI-$\overline{cum}$ > ~4. The change from low to high DHI-$\overline{cum}$ follows the direction from the Northwest (red) to the Southeast (dark green, see Figure 6a). Exceptions are the first two bins in Figure 6h with 3.04% (227 + 791) pixels appear to be not aligned with the trend. Due to the small number and the discrete distribution, no exceptions appear in the zero lines of the density plot (Figure 6e).

(2) Likewise, DHI-$\overline{min}$~0.1 separates the relationship with DHI-$\overline{sea}$. DHI-$\overline{sea}$ is independent where DHI-$\overline{min}$ < ~0.1 but negatively related where DHI-$\overline{min}$ > ~0.1. The change from low to high of DHI-$\overline{min}$ is also from the Northwest (red) to the Southeast (dark green, see Figure 6b).

(3) The DHI-$\overline{cum}$ increase is related to vegetation types changing from grass to deciduous forest and to evergreen forest. Low DHI-$\overline{min}$ < ~0.1 corresponds to a wide range of DHI-$\overline{sea}$ in the geographically less-vegetated Tibetan Plateau and widespread cropped Northeast Plain (Figure 6b–c).

3.4. Regression and Change

Linear regression analyses applied to the long-term (2001–2010) DHI-(cum,min,sea) components are displayed geographically in Figure 7a–f including slope and significance (see also

Table 3. Significance of slope of Dynamic Habitat Index components: Cumulative Annual Productivity DHI-cum, Minimum Annual Apparent Cover DHI-min and Seasonal Variation of Greenness DHI-sea.

Significant	Slope of DHI-cum	Slope of DHI-min	Slope of DHI-sea
p ≤ 0.05	15.76%	8.70%	9.85%
p > 0.05	84.24%	91.30%	90.15%

). Four sub-regions, North China Plain, Nanling Mountains, Yunnan–Guizhou Plateau, and eastern Inner Mongolia, are highlighted due to the strong tendencies of the three DHI components, as well as the regional differences of climate, population density, species and endangered species richness. The following results are noted:

(1) The maximum slopes of decreasing DHI-sea and increasing DHI-(cum,min) occur in the cropland dominated (cropland accounts for 78.25%) North China Plain (Frame A), which is related to increasing temperatures in winter [71] and the so-called agricultural double-delay technology [69] delaying both the harvesting of maize and the sowing of wheat by about seven days.

(2) These phenomena are also observed in the Yunnan–Guizhou Plateau (Frame B) due to the implementation of China’s National Forest Conservation Program since 2000. Reforestation in Yunnan-Guizhou Plateau has shown forest coverage increasing from 40.77% (23.83%) in 1999 to 47.50% (31.61%) in 2008 [72].

(3) In contrast, the spatial–temporal distribution of the slope of DHI-sea can be associated with extreme climate events, varying from the well-vegetated south to less-vegetated north: the maximum increasing (decreasing) slopes of DHI-sea (cum,min) occur in Nanling Mountains (Frame C), where an ice storm occurred in 2008 [73,74]. The maximum slope decrease (decrease) of DHI-sea (cum,min) occurs in Inner-Mongolia (Frame D) perhaps due to a decline of accumulated precipitation over northern China during the summer season of 2009 [75].

Comparison between the three DHI components and their relationship tendencies lead, for most regions of significant change, to the following results (Figure 7g–i): (1) Significant slopes of DHI-cum, DHI-min, and DHI-sea occur in 15.76%, 8.70% and 9.85% of the total area, respectively. (2) DHI-cum tendencies increase if their values are high. (3) The slopes of DHI-(min,sea) decrease (increase) if their values are low (high), which is related to the bi-modality of the distribution separating DHI-$\overline{min}$ (DHI-$\overline{sea}$) near 0.05 (0.5).

4. Conclusions

In summarizing, although DHI is not a new index and its usefulness has been tested using actual biodiversity data in Australia, Canada and North America, we have, for the first time, applied this index in China. Assuming that DHI constitutes a viable surrogate of species richness without testing any relationship with actual species richness may be incorrect and could be misleading when extending the analysis to endangered species richness. Therefore, we first used field survey-based species richness data to test MODIS-based DHI for its applicability in China employing the DHI linear relationship. Similarly, we use field survey-based endangered species richness data to test MODIS-based standard deviation of DHI, while also using their respective linear relationships. Note that estimating the uncertainty of listed species as “potential” compared to actual species is an important issue for practical application on the regional scale (Guest Editor’s personal communication) but beyond the scope of this analysis. Here, the remote sensing derived statistic of annual productivity (i.e. first and second moments) is the only attributor used for estimating species and endangered species richness. Variables characterizing external climate and internal anthropogenic impact also contribute to the biodiversity dynamics as well [50,58,76,77]; for details in an ecohydrological, framework see Cai et al. [51,52] and Fraedrich et al. [78].

Geographic distributions of the species and threatened species census (2013, types/10 × 10 km²) of mammals, amphibians and birds are globally available, which thus makes our approach applicable for other continents or on a global scale. Furthermore, a large number of processes or events can cause biodiversity change which, for sake of simplicity, can be categorized roughly as either abrupt or gradual. Abrupt changes are generally referred to as disturbances, including wildfires, high winds, landslides, floods, avalanches and alike. Taking long term means (decadal averaged DHI-$\overline{cum}$) and a dynamic metric (standard deviation of inter-annual DHI-cum) to compare them with static metrics (field surveys based species and endangered species richness) may lead to temporal mismatching, particularly when spatial resolution is high and abrupt changes occur frequently. The linear regression comparisons do not exclude the identification of abrupt biodiversity changes or the distinction between subsequent time series, which can be very important for the interpretation of biodiversity processes. Gradual biodiversity change comprises an ecological succession, and these successional processes are subject to many variables that may or may not lead to some static or final state, and they set limits on the predictability of future living conditions of individual species.

This analysis is focused on China as the study region, which, consisting of 71,380 pixels with 10 × 10 km² spatial resolution, allows the following estimates of the underlying dynamics of species distribution and abundance, the adequate level of resources, and the pressure on biomes from the external environment.

(1) The decadal mean Cumulative Annual Productivity, DHI-$\overline{cum}$ < 4, in Northwest China (DHI-$\overline{cum}$ > 4 in Southeast China) shows a stable (positive) relation with the Minimum Annual Apparent Cover DHI-$\overline{min}$~0.1 and is positively (negatively) related to Seasonal Variation of Greenness DHI-$\overline{sea}$.

(2) The significant decadal tendencies reveal bimodal frequency distributions aligned near DHI-min~0.05 and DHI-sea~0.5 separated by zero slopes; that is, regions with small DHI-min (and DHI-sea) values are becoming smaller (and smaller), and regions with high DHI-min (and DHI-sea) are becoming higher (and higher).

(3) Land surface vegetation productivity and its variability have an effective impact on the richness of species and threatened species. Areas with high productivity support a relatively high level of species richness, and threatened species richness increases in frequently varied productivity regions. The restricted food resources and habitat, which lead to a reduced vegetation productivity, confine the species’ spatial distributions and also the human population density.

(4) For conservation planning, regions with high variability require more attention due to the high positive relationships with threatened species richness. Therefore, monitoring spatial-temporal dynamics of China’s terrestrial species and threatened species richness (by first/second moments and tendencies of remote sensing based DHI) may also provide the data necessary to support decision making processes.

These results promote our current understanding of the spatial-temporal dynamics of China’s terrestrial biodiversity, especially of effects induced by climate and anthropogenic changes on habitat dynamics based biome accumulation, minimum coverage and vegetation seasonality. Areas with anthropogenic activities, including cropping and reforestation contribute mainly to the DHI components changing in the North China Plain and Yunnan-Guizhou Plateau. Impacts by climate extremes, such as severe droughts and or ice storms, are detected in the Nanling Mountain and Inner Mongolia. In conclusion, this research contributes to a remote sensing based method to support the planning of natural resources conservation, and in turn, to protect biodiversity.