Detection of Decreasing Vegetation Cover Based on Empirical Orthogonal Function and Temporal Unmixing Analysis

Vegetation plays an important role in the energy exchange of the land surface, biogeochemical cycles, and hydrological cycles. MODIS (MODerate-resolution Imaging Spectroradiometer) EVI (Enhanced Vegetation Index) is considered as a quantitative indicator for examining dynamic vegetation changes.This paper applied a newmethod of integrated empirical orthogonal function (EOF) and temporal unmixing analysis (TUA) to detect the vegetation decreasing cover in Jiangsu Province of China. The empirical orthogonal function (EOF) statistical results provide vegetation decreasing/increasing trend as prior information for temporal unmixing analysis. Temporal unmixing analysis (TUA) results could reveal the dominant spatial distribution of decreasing vegetation.The results showed that decreasing vegetation areas in Jiangsu are distributed in the suburbs andnewly constructed areas. For validation, the vegetation’s decreasing cover is revealed by linear spectral mixture from Landsat data in three selected cities. Vegetation decreasing areas pixels are also calculated from land use maps in 2000 and 2010. The accuracy of integrated empirical orthogonal function and temporal unmixing analysis method is about 83.14%. This method can be applied to detect vegetation change in large rapidly urbanizing areas.


Introduction
Information on vegetation change has practical significance for revealing surface spatial variation and evaluating the regional ecological quality [1][2][3].Vegetation indices are effective quantitative indicators of vegetation health spatial distribution and key parameters to study in landscape ecology, climate change, and soil erosion in various researches of surface processes [4][5][6].MODIS EVI dataset is utilized to examine regional vegetation changes due to its excellent presentation of vegetation information and anti-interference against the soil background and atmosphere [7].
In China, land use change is mainly characterized by urbanization [8,9].Land use and land cover changes are primarily identified based on the repeated acquisition of remote sensing datasets.Proposed approaches for multitemporal analysis include (1) images classification [10], (2) wavelet decomposition [11], (3) a multitemporal dataset which is transformed by principal component (PC) analysis (then the resulting component could reflect various changes [12]), (4) spatial statistical analysis which calculates the quantitative analysis of the changing scope, strength, and trend [13], (5) change vector analysis which can calculate the change type and intensity [14], and (6) temporal unmixing modeling [15].These changing analytical methods have their own characteristics and emphases, but, for multitemporal images, the most important aspect is to remove noise and determine the dominant dimensions [16].In this study, prior information on increasing/decreasing vegetation spatial coverages is calculated by empirical orthogonal function (EOF).
The empirical orthogonal function is usually employed to model the spatial-temporal patterns of the sea surface temperature [17], dynamical atmospheric [18], sea-level rise [19], and shoreline variability [20].The empirical orthogonal function has been applied on the night lights dataset and MODIS EVI for characterization and modeling of the changing extent, intensity, and distribution [21].Temporal unmixing is used to model the spatial distribution of crop types [22,23], forest [24], and sea ice imagery [25].Compared with other approaches [23,26], the empirical orthogonal function can describe vegetation change trend without ancillary information in this research.The integrated of empirical orthogonal function and temporal unmixing method is first provided to model the spatial-temporal patterns of crop types [16].Spectral mixture analysis is used to monitor vegetation change eliminating the background influence, but it is not suitable for a large area [27].Because of the fast urbanization in China, new construction results in widely decreasing vegetation.In this study, the empirical orthogonal function aims to take the decreasing vegetation curves as a prior for temporal unmixing models in Jiangsu Province.The combination method using empirical orthogonal function (EOF) and temporal unmixing analysis (TUA) is introduced to quickly detect decreasing vegetation areas in a large area.
This study aims to evaluate the integrated empirical orthogonal function and temporal unmixing to detect the changing vegetation area and apply the approach in Jiangsu Province, a rapidly urbanized province in southern China [28].The theories of empirical orthogonal function and temporal unmixing and the application results in Jiangsu are presented first.Next, contrasting Landsat data are used to validate the accuracy and consistent spatial distribution of decreasing vegetation with MODIS EVI by empirical orthogonal function and temporal unmixing method.At last, this analysis also identifies strengths and uncertainties of the combined empirical orthogonal function and temporal unmixing method.

Study Area and Datasets.
Jiangsu Province is located at 116 ∘ 18  -121 ∘ 57  E, 30 ∘ 45  -35 ∘ 20  N, with a 10.26-million hectare area that accounts for 1.1% of the total terrestrial area in China.The plains area is 7.06 million hectares and the water area is 1.73 million hectares.The elevation of more than 90% of areas in Jiangsu Province is lower than 50 meters.Jiangsu belongs to warm temperate to north subtropical transitional climate (Figure 1).Jiangsu Province's comprehensive economic strength of Jiangsu has been at the forefront in China.After the opendoor policy was issued in 1978 in China, the urban area and the growth rate increased significantly in Jiangsu.In 1990, Jiangsu had an urban population of 14.59 million and a rural population of 53.08 million.As a contrast, there was an urban population of 49.90 million and a rural population of 29.30 million in 2012.In 2012, GDP per capita in Jiangsu reached $11,113.3 compared to the national average of $6251.87[29].The urbanization rate of Jiangsu was 63% in 2012, and more than 80% of the urban growth area occurred outwards from the pregrowth urban fringes at the expense of rural lands [28].Due to urbanization, the arable land area per farmer decreased to less than 335 m 2 [8].Urban sprawl has environmental impacts, such as enhancing urban heat island and increasing carbon emissions, affecting the quality of life in urban areas [30].Therefore, during the growing process, the timely and effective supervision of vegetation is of importance.

Datasets. The 16-day MODIS EVI (MOD 13Q1
) composites with a 250 m spatial resolution were downloaded from the USGS website (http://glovis.usgs.gov/).EVI temporal profiles span from February 2000 to December 2012 (296 images).EVI is less susceptible to cloud and haze contamination than NDVI [31].The EVI time series are mosaicked, reprojected, and resampled to 1000 m for displaying dominant vegetation change trend.
One validation dataset is Landsat TM images in 2000, 2002, 2006, and 2009 by linear spectral unmixing method.The three net spectral endmembers (substrate, vegetation, and dark) are provided and validated from this research results [32].We use spectral unmixing model of ENVI software.
Another dataset to validate the temporal unmixing analysis method accuracy is land use maps in 2000 and 2010.The land use map production is from Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences.We also use ArcGIS software to generate 60 random points and compare the land surface type (Agriculture, Grass, Forest, Water, and impervious surface) of google earth and the land use map in 2000 and 2010.The accuracy of the land use map is about 91%.

Empirical Orthogonal Function.
The empirical orthogonal function method decomposes the original data into the product of temporal function and spatial function [18,33].The curves of empirical orthogonal function represent temporal patterns, which are the eigenvectors of the covariance matrix from the principal transform of the original data.PCs represent the spatial weight of the corresponding curves from empirical orthogonal function.In this research, the curves of empirical orthogonal function are the vegetation increasing/decreasing change curves.PCs display the spatial distribution of corresponding curves of empirical orthogonal function.There is principal components function in the ENVI software.The eigenvectors from the statistic file are EOFs.
In the empirical orthogonal function method, the original data () is divided into the product of a temporal function (EOFs-) and spatial function (PCs-) [34,35]: Suppose that   has large projection on the first  vectors in the  spatial field: where  is the residual error when  is expressed by  vectors.where   () represents error. is the spatial dimension number, and  is the different time.() is the variance for the total error.It is the sum of variance of each pixel error: where () is the exception of error equation.The constraint condition is where  is Lagrangian constant.
Moreover, the empirical orthogonal function method aims to reduce the dimensionality with a minimum loss of information while maintaining the majority of the variation affected by independent processes and capturing the essential features [36,37].
In this analysis, the curves from empirical orthogonal function with decreasing trends are temporal patterns and represent the vegetation cover reduction.In conventional empirical orthogonal function, the PCs and empirical orthogonal function are separately interpreted in terms of spatiotemporal processes, and the empirical orthogonal function only represents statistically unrelated modes of variance.In the context of this study, much more attention  is paid to the decreasing curves of empirical orthogonal function related to vegetation reduction, which provides prior information for temporal unmixing model.The curves of empirical orthogonal function and PCs can be obtained by ENVI software.

Temporal Unmixing Analysis.
Temporal unmixing analysis is an extension of the linear spectral unmixing.The concept of the temporal unmixing model is that each pixel is the linear combination of temporal endmembers and corresponding fractions [22].Fractions of endmember should be equal or greater than zero and the sum of fractions in one pixel should equal to one.
Accurate endmembers and temporal dimensions are the keys to the temporal unmixing model.Endmembers are in the extreme position of the feature space and represent different fundamental processes.The selection of endmembers is crucial for the temporal unmixing model; here, endmembers are selected by the geometric vertex method [38,39].The curves of empirical orthogonal function provide vegetation increasing/decreasing trend as prior information.From the EOF curves with decreasing trend, the pixels with decreasing trend can be found in the corresponding PC.Temporal vegetation decreasing endmembers can be extracted from the corresponding PC scatter plot.
The temporal unmixing used here has two differences with the traditional temporal unmixing [40].First, the approach only selects one decreasing vegetation endmember to model the temporal unmixing analysis.Second, here the curves of empirical orthogonal function provide vegetation decreasing prior information for temporal unmixing [16].In this research, the integration of empirical orthogonal function and temporal unmixing is useful for identifying the processes of decreasing and increasing vegetation cover.The temporal unmixing analysis is completed in ENVI software.
The technical flow chart for the research to detect decreasing vegetation trend and validate the accuracy is as follows (Figure 2).

Results
The first curve of empirical orthogonal function has primary eigenvalues, which contributes to approximately 89.25% of the variance (Figure 3(a)).Other curves' variances continuously decrease.The variance of the second to the seventh eigenvector of empirical orthogonal function account for 2.04%, 1.09%, 0.97%, 0.44%, 0.41%, and 0.28%, respectively.It is important to acquire the vegetation decreasing/increasing trends as prior information from the eigenvectors which account for large variance.The amplitudes of the first ten curves of empirical orthogonal function could be quantified in the time domain (Figure 4).The first ten curves of empirical orthogonal function are temporal eigenvectors of the EVI variance structure.The first curve has relatively low amplitude because it is the mean value of EVI with no variance.The second and third curves have annual and biannual peaks.Distinctly, the fourth curve shows an increasing trend.The fifth curve also displays periodic cycles with biannual peaks.The sixth and seventh curves have decreasing trends before 2006 and gradually increase afterward.
The temporal curves of empirical orthogonal function provide prior information for the temporal unmixing model.The fourth, sixth, and seventh curves are related to vegetation change trend, but the fourth curve accounts for more variance than the sixth and seventh curves.The fourth curve could reveal increasing vegetation cover changes, so the opposite pixels of the fourth curve are related to vegetation decreasing.
In this analysis, the third PC and the fourth PC are defined as  and  apexes to obtain endmembers that can describe the details for vegetation trend change (Figure 5(a)).The EOF provides prior information for temporal unmixing.The fourth EOF has the increasing trend, so there are pixels with increasing trend in the fourth PC; further, the vegetation decreasing pixels are in the opposite direction of vegetation increasing endmembers in the fourth PC.In the scatter plot, the pixels in the top position are related to vegetation increasing as the fourth EOF curve provides prior information.Oppositely, the vegetation decreasing endmembers can be found in the bottom vertex pixels.
Scatter plots of PC3 and PC4 were used to select endmembers (Figure 5 Spatial distribution of decreasing vegetation endmember by temporal unmixing model is shown in Figure 5(c).It can be observed that the decreasing vegetation endmember is mainly located in the suburbs.In Suqian City, the decreasing vegetation cover is displayed in the suburbs.In Nanjing City, the decreasing vegetation is shown in the suburbs and in the south new area.In Taizhou City, the decreasing vegetation is in the suburbs and in the southern part.In the middle part of Jiangsu Province, the decreasing vegetation endmember is along the Yangtze River.In the southern part of Jiangsu, Suzhou City, the decreasing vegetation endmember not only is in the urban edge but also has scattered distribution, because there is a high speed economic development in Suzhou.In the development of Jiangsu Province during 13 years, the vegetation decreasing speed in the south is much faster than that in the north.The typical urbanization processes around the old city center, leading to vegetation decreasing in suburbs of Suqian City (Figure 6(a)).Suqian is in the north of Jiangsu Province.Comparing Figure 6(a) with Figure 5(c), the linear spectral unmixing and temporal unmixing methods both display same vegetation decreasing area in the suburbs of Suqian City.
The decreasing vegetation cover area in the middle of Nanjing is located in the suburbs due to urban expansion (Figure 6(b)).At the same time, the vegetation reduction in the southern part is due to new construction.Contrasting Figure 6(b) with Figure 5(c), the empirical orthogonal function and temporal unmixing methods detect the same decreasing vegetation area with linear spectral unmixing method.
Vegetation reduction in Suzhou presents a star-scattered pattern.Suzhou City is in the southern region of Jiangsu Province.Urbanization and economic development are the main reason for vegetation decreasing.The spatial changes of the vegetation fractions (Figure 6(c)) are consistent with the decreasing vegetation distribution of the empirical orthogonal function and temporal unmixing method (Figure 5(c)).Both display decreasing vegetation of star-scattered patterns and similar spatial distribution in the suburbs of Suzhou City.

Validation Based on Land Use Map.
Decreasing vegetation areas are calculated from land use map from 2000 to 2010 (Figure 7).The blue part is vegetation decreasing area from empirical orthogonal function and temporal unmixing method.The red part is vegetation decreasing area from the land use map in 2000 and 2010.Vegetation decreasing area is larger in the south due to faster economic development than that in the north of Jiangsu.
According to decreasing vegetation pixels coincidence from the land use map and temporal unmixing analysis, the accuracy of empirical orthogonal function and temporal unmixing analysis is 83.14% (Table 1).Vegetation decreasing area from empirical orthogonal function and temporal unmixing analysis is 6956 km 2 .Vegetation decreasing area from the land use map during 2000 and 2010 is 7111 km 2 .As a result, the spatial coincidence is 5912 km 2 .

Strengths and Uncertainties.
The combination method of empirical orthogonal function and temporal unmixing was first mentioned [16] for identifying and representing phenology spatiotemporal patterns.The approach used the number of phenology dimensions based on empirical orthogonal function and modeled the vegetation phenology distribution by temporal unmixing analysis.Here, we pay attention   to vegetation decreasing and increasing eigenvectors from empirical orthogonal function and select the vegetation decreasing endmembers based on temporal feature space.The two approaches both take statistical empirical orthogonal function results as prior information but emphasize different vegetation changes.The approaches described in [16] aim to describe the temporal phenology endmembers and spatial distribution, whereas the approach here aims to display the spatial distribution of decreasing vegetation.The advantage of this combined method is using the EOF as prior information for temporal unmixing analysis and using the vegetation decreasing endmembers to unmix the spatial distribution of vegetation decreasing area.In another research [41], the approach enables the detection of different types of changes occurring in time series, including the dates of changes occurring within seasonal and trend components.This research here only pays attention to the decreasing and increasing vegetation trend and does not emphasize the accurate phenology dates of vegetation changes.
In this analysis, the unmixing processes corresponding to the spatial distribution are selected manually.Manual selection could allow for the consideration of stable endmembers.Compared with other methods, this method highlights the important benefit to quickly detect decreasing vegetation over large areas without classification and auxiliary [21].

Further Application.
Further research is necessary to apply the empirical orthogonal function and temporal unmixing method to different study areas to detect boundary sensitivity of the endmembers when the decreasing vegetation endmembers are selected.In a research [16], pixels with a strong trend of vegetation increase and decrease are identified due to the annual cycle of rising and falling of water.In this study, we have focused on the spatial distribution of decreasing vegetation.Future work may improve the accuracy of decreasing vegetation endmembers.Here, Landsat data at a 30 m spatial resolution can serve to illustrate the spatial mapping accuracy.
Decision makers could use MODIS EVI by empirical orthogonal function and temporal unmixing to quickly detect the spatial extent of decreasing vegetation and it could help in land use planning.Vegetation plays an important part in the land surface characterization, climate change modeling, and biogeochemical cycles.During the processes of urbanization in China, urban expansion is the significant driver for changes toward decreasing vegetation [42].However, monitoring vegetation over large areas at regular intervals is expensive.The combination approach of empirical orthogonal function and temporal unmixing analysis to detect decreasing vegetation could be seen as a preliminary tool.Furthermore, the vegetation fraction by linear spectral unmixing could be utilized to focus on the plots.The new combination method of empirical orthogonal function and temporal unmixing analysis contributes to vegetation cover mapping and monitoring.

Conclusions
Empirical orthogonal function analysis uses principal temporal and spatial patterns to present the original dataset.In this study, much more attention is paid to the increasing and decreasing vegetation eigenvectors which provide prior information for the temporal unmixing analysis.
Here temporal unmixing analysis identifies the spatial distribution of decreasing vegetation endmembers.This approach extracts decreasing vegetation endmembers from temporal principal components to model the spatial distribution.In Jiangsu Province, the decreasing vegetation mainly is distributed in the suburbs due to urbanization.
The Landsat dataset by linear spectral mixture is used for analysis in consistency of decreasing vegetation distribution with the integrated empirical orthogonal function and temporal unmixing analysis.The three components linear spectral unmixing provide estimates of vegetation fraction and the vegetation decreasing patterns.The decreasing vegetation in Suzhou displayed star-scattered pattern around the old city.The decreasing vegetation in Suqian is located in the suburbs.The decreasing vegetation in Nanjing is in the suburbs and new constructed area in the south.The empirical orthogonal function and temporal unmixing method display the same spatial extent of decreasing vegetation with linear spectral unmixing based on the Landsat dataset.Compared with vegetation changes from land use map in 2000 and 2010, the accuracy of the integrated empirical orthogonal function and temporal unmixing method is about 83.14%.

Figure 1 :
Figure 1: Location of Jiangsu Province and land use classification in 2010 (the land use production is from Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences).

Figure 4 :
Figure 4: Curves amplitude of the first ten of empirical orthogonal function could be quantified in the time domain.Importantly, the fourth curve has increasing trend which is opposite to decreasing vegetation trend.
(a)) representing decreasing vegetation cover between 2000 and 2012 (Figure 5(b)).The decreasing vegetation endmember reveals the vegetation reduction processes from 2000 to 2012.It is taken as the average vegetation decreasing representation in Jiangsu for temporal unmixing model, because the third and fourth PCs are used to select the decreasing vegetation endmembers.

1 .
Validation Based on Landsat Dataset by Linear Spectra Unmixing.The vegetation fractions in 2000, 2006/2002, and 2009 in Suqian, Nanjing, and Suzhou correspond to the blue, green, and red channels in Figure 6.The dark areas mean no vegetation from 2000 to 2009.The blue areas mean vegetation cover in 2000 but with no vegetation in 2006/2002 and 2009, which clearly exhibits vegetation change processes.

Figure 5 :
Figure 5: Feature space representation of the third and fourth PCs for temporal unmixing analysis.The decreasing vegetation endmembers are found in the south part of PC space (a).Decreasing vegetation endmember (b).Spatial distribution of decreasing vegetation area estimated by temporal unmixing analysis in Jiangsu (c).

Figure 6 :
Figure 6: Vegetation fraction map in Suqian City (a), Nanjing City (b), and Suzhou City (c) in Jiangsu Province.Vegetation endmember fractions in 2000, 2006/2002, and 2009 are shown as blue, green, and red channels.The dark areas mean no vegetation from 2000 to 2009.The blue areas mean vegetation cover in 2000, but with no vegetation in 2006/2002 and 2009, which clearly exhibits vegetation change processes.The whole area is the city and the red rectangles refer to central urban.

Figure 7 :
Figure 7: The blue part is vegetation decreasing area from empirical orthogonal function and temporal unmixing analysis.The red part is vegetation decreasing area from land use map in 2000 and 2010.Vegetation areas are mainly distributed in the suburbs.

Table 1 :
EOF and TUA method accuracy analysis.