Snow plays an important role in the global energy and water cycles because of its high albedo and thermal and water storage properties, and can indicate the changes in global climate.12.3.–4 The Tibetan plateau (TP) is the world’s highest region and is often called as the Earth’s third pole.5 It has been characterized as the driving force and the amplifier for global climate change.6 Several studies have proved that there were significant warming trends in the TP,78.–9 and the cryosphere of the TP has been changing rapidly.1011.–12 Snow cover over the TP is a vital water source in western China. The large rivers of China, such as the Yangtze River, Yellow River, etc., have their headwaters there.13 In addition, snow cover is not only sensitive to climate change but also is closely related to hydrological and biological processes in the TP as well as its surrounding areas. Thus, the information of spatial and temporal pattern of snow cover over the TP is very important for scientific studies and management applications of the region.
However, the terrain on the TP is high, steep, and thus difficult to access. Conventional meteorological stations are very rare, and most of them are distributed in lower-altitude river valleys or plains where there is usually less snow. Therefore, the in situ observations are difficult to adequately reflect the spatial and temporal distribution of snow cover over the TP. There are many uncertainties in those previous studies that only applied in situ observations to monitor the long-term snow cover change over the TP.1415.16.–17
Since the launch of television infrared observation satellites (TIROS)-1 in 1960, with capability in monitoring snow cover, dozens of satellites have been used to monitor snow cover and have played an important role in it.18,19 On December 18, 1999, the Terra satellite was launched with a complement of five instruments, including the moderate resolution imaging spectroradiometer (MODIS).20 A suite of MODIS/Terra snow cover products, including binary snow cover map (snow or nonsnow) and fractional snow cover (FSC) map at various spatial and temporal resolutions have been available since February 2000.21 These products have been widely used for regional snow cover monitoring, climatological studies, and hydrological circulation modeling.4,13,2223.–24 However, the extensive and persistent cloud coverage in MODIS snow cover products is the main limitation in many applications. A series of methods have been proposed to mitigate the cloud obscuration in using the binary MODIS snow cover products, such as combining images from two MODIS platforms (Terra and Aqua),2526.27.–28 combining MODIS images with the passive microwave products,29,30 using the experiential knowledge of snowlines or snow-covered period,23 and employing the most recently available cloud-free observations from prior and subsequent days.31,32 Although all of these approaches are particularly useful for reducing the cloud cover pixels, they sacrifice temporal and spatial resolution and introduce some uncertainties to different degrees. Since the FSC data can more clearly reflect the gradual change of the snow cover in each pixel than the binary snow cover data, it can be more accurate when using the FSC data to remove the cloud cover through temporal filtering. Therefore, the MODIS daily FSC data are chosen in this study.
The objective of this study is to provide a reliable and up-to-date analysis of variations in snow cover over the TP by taking advantage of long-term, continuous observations from MODIS daily FSC data, and how they are related to the changing temperature. In this paper, we first introduce the datasets used in this study and evaluate the accuracy of MODIS daily FSC data under clear sky conditions over the TP. Then, we describe a simple technique using original MODIS FSC data to fill the data gaps caused by cloud obscuration, which is used to produce cloud-gap-filled (CGF) daily MODIS FSC datasets from last 11 years (2001 to 2011). Finally, through the use of the MODIS-derived snow-covered area (SCA), snow-covered days (SCD), the characteristics of spatial and temporal variations of snow cover and their association with temperature are analyzed.
In this study, we specifically define the area of the TP as the region with an elevation higher than 2500 m. It is located between 73 to 105 degE and 25 to 40 degN and covers an area of approximately (Fig. 1). The mean elevation of the TP is 4367 m.
Data and Methodology
MODIS daily FSC data
Using an algorithm developed by Salomonson,33 daily FSC is provided in the MOD10L2 and the MOD10A1 snow maps at 500-m resolution in MODIS snow products Collection 5.34 The FSC algorithm is based on a statistical-linear relationship developed between the NDSI from MODIS and the true subpixel fraction of snow cover as determined using Landsat scenes from Alaska, Canada, and Russia. The fraction of snow cover within a MODIS 500-m resolution pixel is provided with a mean absolute error of less than 0.1 over the entire range of FSC from 0.0 to 1.0.21,33
The daily snow product in MOD10A1 is a tile of data gridded in the sinusoidal projection. There are four data fields of snow data; snow cover map (snow or nonsnow), FSC, snow albedo, and quality assessment (QA) in the data product file.34 The MOD10A1 FSC data, for the period 2001 to 2011, are used in this study. Those coded integer values of the MOD10A1 FSC data over the TP include 0 to 100 (fractional snow), 200 (missing data), 201 (no decision), 225 (land), 237 (inland water), 250 (cloud), 254 (detector saturated), and 255 (fill; no data expected for pixel). The images from 2001 to 2011 over the TP are mosaicked and georeferenced into a UTM projection using the MODIS Reprojection Tool.35
In situ temperature and SCD data
In situ daily air temperature data from 102 meteorological stations of the Chinese Meteorological Administration (CMA) (Fig. 1), are used to analyze the relationship between the snow cover and air temperature. In situ SCD data from 2001 to 2006 at 93 stations (which are located in the study area) are applied to validate the MODIS-derived SCD map.
Evaluation of MODIS FSC Data
The quality of FSC data within a MODIS 500-m resolution pixel has been evaluated in some areas of Alaska, Canada, and Russia, with a mean absolute error less than 0.1.33 However, the accuracy of the MODIS FSC data over the TP is unknown. In this study, we validate the quality of the MOD10A1 FSC data in the absence of cloud cover over the TP using the same evaluation method in Ref. 33.
In this work, six Landsat images (selected as reference data) and the same day MOD10A1 FSC data are collected in the evaluation. These images are georeferenced with the same projection for comparison. The specific dates and locations of Landsat images are listed in Table 1 and shown in Fig. 1. For every Landsat images, using the current ‘‘SNOWMAP’’ approach,20 each 30-m pixel is classified as snow or nonsnow (a binary classification). The Landsat bands (band 2 and band 5) corresponding to the MODIS bands 4 and 6 are used to calculate NDSI values.
Information about Landsat scenes used as ground truth in validating MOD10A1 FSC data in the TP.
For each MODIS 500-m grid cell, the true FSC value is determined on the basis of the Landsat 30-m observations by counting the number of Landsat pixels covered by snow versus the total number of Landsat pixels in the cell. In this evaluation, only the grid cells with true FSC value greater than 0.02 are used to compare with MOD10A1 FSC data.
On the TP, the cloud obscuration significantly limits the usefulness of the MODIS FSC data. For instance, the average cloud coverage was about 47.3% in 2008 (Fig. 2), which indicates that clouds cover a wide range of the TP.
Since clouds have the rapidly-changing and daily-shifting features and snow change gradually over time, we can acquire snow information (FSC value) of the cloud-covered pixels through temporal interpolation. In this study, based on the cloud-free days’ observations, a simple interpolation method, cubic spline interpolation algorithm that goes through each data point is used to fill in the data gaps caused by clouds or other reasons. We specially define it as a CGF method.
First, the coded integer values of MODIS FSC images are reclassified into two categories: fractional snow (0 to 100), and cloud (250). The new coded integer values have no change if the original value is 0 to 100 (fractional snow) and 250 (cloud); the 225 (land) and 237 (inland water) are merged into 0 (fractional ); the other values like 200 (missing data), 201 (no decision), 254 (detector saturated), and 255 (fill) are merged into 250 (cloud). Each cloud pixel in these images is then interpolated to FSC value using the cubic spline interpolation algorithm based on the observations (FSC value) of cloud-free days. If the interpolated value less than 0, reassign it with 0; and if the interpolated value greater than 100, reassign it with 100 (because ).
In addition to the filled-in FSC value for each cloud pixel, an associated cloud-persistence-days (CPD) of each cloud pixel is calculated within this CGF method. The CPD represents the number of consecutive days of cloud obscuration from the last view of the surface to the next view of the surface. The CPD map is convenient to assess the accuracy of the CGF method under different CPD conditions.
Accuracy assessment of the CGF method
The performance of the CGF method is evaluated by comparing the original MODIS FSC data against the MODIS CGF FSC results (which are based on cloud assumption). At first, four original MODIS FSC images over the TP with relatively less cloud obscuration (on February 6, April 8, October 12, and November 23 of 2008) are randomly selected as the reference image for validation. Second, the pixels under clear sky conditions of these four original MODIS FSC images are assumed as cloud cover (merged into 250), and fill in these cloud pixels using the CGF method. Finally, we compare the pixels under clear sky conditions of these four original MODIS FSC images against their respective cloud removal results (based on the cloud assumption). The effectiveness of this method under different CPD conditions is evaluated through mean absolute error.
Methodology of Snow Cover Changes
For the period from 2001 to 2011, daily MODIS CGF FSC datasets of TP are derived from MODISFSC data using CGF method. To observe the spatiotemporal variation of snow cover, we obtain the SCA and SCD from daily CGF FSC data and analyze them from four different elevation zones.
In the process of SCA calculation, only the pixels in the MODIS CGF FSC images with the value greater than 50 (i.e., FSC greater than 50%) are seen as snow-covered pixels. The SCD represents the total number of days with snow cover in an annual cycle of the snow. The SCD is calculated using all images within a year from January 1 to December 31 by Eq. (1):
The four elevation zones used are 2500 to 3500 m, 3500 to 4500 m, 4500 to 5500 m, and above 5500 m. The areas of these four elevation zones covered in TP are approximately 541096 (18.2%), 831444 (28.0%), 1462796 (49.3%), and , respectively. The number of meteorological station in these four elevation zones is 54, 38, 10, and 0, respectively. Thus in the later analysis, the temperatures of 10 stations above 4500 m are used instead of the temperatures in the elevation zone higher than 5500 m.
Ordinary least-squares analyses are conducted to estimate linear time trends of the SCD for each pixel over the study period, and their significance levels () are presented by F-test. In addition, we calculate Pearson correlation coefficients between SCA and temperature, and assume that the inter-annual variability in SCA is related to the temporal variability in temperature if the correlation coefficients are statistically significant.
Results and Discussion
Accuracy of MODIS FSC Data
These results of the comparisons between the MOD10A1 FSC data and the true FSC images obtained from Landsat are shown in Table 2. The overall mean absolute error, standard deviation, and correlation coefficient for the MODIS FSC products (in MOD10A1) in TP is 0.098, 0.156, and 0.916, respectively. These results indicate that the MOD10A1 daily FSC data has sufficient accuracy to reflect snow cover information over the TP.
Results of the comparisons between the MOD10A1 FSC data and the true FSC images obtained from Landsat.
|Region||Mean absolute error||RMSE||Correlation coefficient|
Note: the number of points (n) involved in the evaluation is listed in the table.
Effectiveness of CGF Method and SCD Validation
The frequency of the CPD reduces gradually as the increasing of CPD value. The CPD of 92.5% of cloud covered pixels over the TP is less than 15 days [Fig. 3(a)]. The accuracy assessment results show that the mean absolute error of the CGF method increases as the rising of CPD, this method has very high accuracy (with mean absolute error less than 0.1) when CPD within 5 days [Fig. 3(b)]. It can bring relatively greater error (with mean absolute error greater than 0.15) when CPD is larger than 20 days, but this situation is very rare [see from Fig. 3(a)]. The overall mean absolute error of the CGF method using the frequency of CPD [in Fig. 3(a)] as weight is 0.092. These indicate that the CGF method used in this study has a higher accuracy to fill in the missing data due to cloud obscuration or swath gaps.
Figure 4 shows the comparisons of original MODIS FSC map, MODIS CGF FSC map, and CPD map over the TP on 5 March 2008. Using the CGF method, for each cloud-covered pixel on 5 March 2008, the FSC value is effectively retrieved [Fig. 4(b)]. The corresponding CPD map [Fig. 4(c)] probably reflects the confidence of the retrieved FSC data at each pixel. The CPD of most cloud-covered pixels is less than 6 days.
The MODIS-derived SCD in almost all the stations appears lower than in situ observations from 2001 to 2006 (Fig. 5). Compared with the in situ observations, the MODIS-derived SCD are averagely reduced 8.16 days. Although there are some differences between the in situ SCD and MODIS-derived SCD, a strong relationship exists between two data sets. The correlation coefficient between in situ SCD and MODIS-derived SCD is 0.96, which also indicates the higher accuracy of the CGF method. This SCD validation result is somewhat different from that obtained by Wang and Xie.28 In their results,28 the average MODIS-derived SCD at the 20 stations (in northern Xinjiang, China) is 9 days more than the average in situ SCD. But the MODIS-derived SCD (in Ref. 28) are calculated from the multi-day combination of Terra and Aqua MODIS binary snow cover products.
These differences between in situ SCD and MODIS-derived SCD can be explained by the following aspects:
1. The MODIS-derived SCD are greatly affected by mixed pixel, especially the complex terrain in the TP. While the observed SCD only represent the stations.
2. When using MODIS NDSI to obtain snow cover information, the measurement will be severely affected by the presence of dense forest. We analyze the land use of these greatly underestimated stations and find out that these stations are covered or surrounded by dense forest and shrubs.
3. The MODIS-derived SCD are also partly affected by the error in the CGF method, especially the larger error when CPD larger than 20 days (although this situation is very rare).
Annual Cycle of Snow Cover
The values present in Fig. 6 are averages of 11 years from 2001 to 2011. This result leads to the construction of snow depletion curves, which are commonly used as an input for snow-runoff simulating.
The change in the snow cover extent is almost monotonously increasing with the increasing elevation [Fig. 6(a)]. Strong seasonal variations in SCA are found over the whole area of the TP [Fig. 6(b)]. From the end of October to beginning of May, the SCA over the whole area is greater than 15% with relatively large standard deviations reflecting the very high inter-annual variability. While in the July to August period, the SCA over the whole area is less than 5% and with a relatively small inter-annual variability.
The SCA, depending upon the elevation, peaks in January at the lower elevations () and later as in March to April at elevations above 5500 m, and then progressively decreases (melts) and reaches minimum values in the July to August period [Fig. 6(a)]. In June, snow is present above 3500 m in the TP; at the lower elevations (), snow is almost not present after the beginning of May. The accumulation of snow in the TP, also depending upon the elevation, starts in the beginning of September at higher elevations () and as late as in the end of October at elevations below 3500 m. Seen from the whole area of TP [Fig. 6(b)], the obvious snow melting and accumulation periods occur in March to June and September to October, respectively.
At elevations above 5500 m, the SCA is greater than 30% throughout the year, with two maxima in SCA appearing in spring and autumn seasons and a relatively minimum during December to January along with the typical minimum in the July to August period [Fig. 6(a)]. The relative minimum in SCA during winter months at higher elevations can be seen as the unique characteristic of snow cover in TP, which may be principally attributed to the following aspects: (1) very dry weather and less frequent snowfall during the winter months at this elevation in TP under the influence of mid-latitude cold current in winter; and (2) increased sublimation of the snow over this high-elevation area where it is frequently accompanied by high winds. More than half of the snow mass in the TP was lost by sublimation in winter.36
Characteristics of Spatial and Temporal Variations of Snow Cover
The results of SCD from 2001 to 2011 over the TP are shown in Fig. 7. The areas with SCD greater than 60 are usually considered as stable snow-cover areas in China,37 and they are the major containers of snow water resources and the areas where many large rivers originate from.
The areas with higher SCD correspond well with the huge mountains, including Kunlun, Karakoram, Himalaya, Qilian, Tanggula, and Hengduan Mountains. The most persistently SCAs (with SCD greater than 240 days) are concentrated in these huge mountain regions. In contrast, due to strong shielding from these huge mountains, most of the interior of the TP has relatively lower snow-covered persistence (with SCD less than 60 days) although the averaged elevation is over 4000 m.
For each pixel, the linear trend of SCD over the 11 years [Fig. 8(a)] and its statistical significance level [5%, Fig. 8(b)] are estimated. The linear trend analysis of SCD shows that about 34.14% of the study area has a declining trend during 2001 to 2011. However, only 5.56% pixels with a declining trend are statistically significant (). Large ( days per year) and significant () SCD decline are mainly located in the high altitude area, especially in the west of Hengduan mountains and northern Karakoram mountains. In contrast, 24.75% of pixels have an increasing trend (only 3.9% with a significant increase). Most of these pixels are located in the south of Karakoram mountains, northern Qilian mountains, and some area of the interior of the TP (between the Kunlun mountains and Tanggula mountains). About 41.11% of pixels exhibit no or very small changes (i.e., days per years), which are considered to be stable.
Figure 9 depicts time series of the average SCD for different elevation zones and the average SCA (%) for different seasons of the whole TP during 2001 to 2011. It demonstrates that the variability of snow cover in these 11 years is characterized by normal oscillation. Snow cover fluctuates around the mean value. The time series of fluctuations in average SCD and SCA in the whole area show very high inter-annual variability, although no obvious increase or decrease tendency is found in them. For instance, the average snow cover in the TP during the spring of 2004 was only 15%, while in 2005 it was about 26%.
However, with only 11 years of MODIS data, it’s hard to reach any definitive conclusion according to the time trend of the SCD and SCA for this area. Consequently, a longer time series of data are needed to reach any definitive conclusions.
Correlations Between Snow Cover and Temperature
Low temperature is an essential condition to meet the formation of snow cover. Moreover, rising temperature is the main reason for accelerated snow melting. To find a clue to the response of snow cover to temperature change, it is necessary to understand the linkage of variations between snow cover and temperature.
Table 3 shows the Pearson correlation coefficients between the monthly snow cover and the in situ measurements of temperature at different elevation zones. The variability of snow cover and temperature of the months (and elevation zones) with reasonably high correlation coefficients (statistical significance at the 0.01 level) are selected and shown in Fig. 10. The correlation coefficients between snow cover and temperature vary among different months and elevations, showing a large spatial and temporal heterogeneity. The snow cover in most of the months and elevations shows negative correlation with the temperature of the same month, particularly in some elevations of February, April, May, August, and September (Table 3). These high negative correlations between snow cover and temperature indicate that temperature strongly affects the snow cover.
Pearson correlation coefficients between monthly snow cover and temperature at different elevation zones in the TP from 2001 to 2011.
|<3500||3500 to 4500||4500 to 5500||>5500||Whole|
Note: ** and * indicate statistical significance at the 0.01 and 0.05 level, respectively. The periods in which the snow covers are less than 2% are represented by “N.A.” (meaning “not available”).
The impact of the temperature on snow cover can be explained as following two aspects. One is that the rising temperature increases the melting of snow in the snowmelt season, for example, the high negative correlations in April and May. The other is the decreasing temperature increases the snowfall in the early snowfall season, such as the high negative correlations in August and September for the high elevation zones.
Besides, snowfall is also an important influencing factor for the snow cover changes. In the TP, snowfall mainly concentrates in autumn and spring. Because of the very dry weather, there is little snowfall in winter months (especially in February). For instance, the SCA in other snow season months (such as in October, November, and March) shows a weak negative correlation with temperature, which indicates that the snowfall (as a positive factor) also greatly affects the snow cover in these months.
In winter, the influence factors of SCA are more complex, which include the small amount of snowfall and the temperature-related sublimation and melt. In December and January, because the temperature is well below the freezing over the TP, there is little temperature-related snow melting, while the temperature can also reasonably affect the SCA through changing the sublimation of snow. Besides, the small amount of snowfall also positively affects the SCA. Thus, the SCA in December and January shows a weak negative correlation with temperature. While in February, as the temperature rises, the daily maximum temperature is obviously higher than freezing in some low elevation zones, and snow begins to melt in these areas. The temperature-related snow sublimation and melt as well as the rare snowfall lead to the significant negative correlations between SCA and temperature in February.
The MODIS daily snow products (MOD10A1) from the Terra satellite have provided an excellent opportunity to analyze the spatial and temporal distribution of snow cover over the TP. The MODIS FSC data in MOD10A1 have sufficient accuracy to reflect snow cover information over the TP (with mean absolute error is about 0.098).
Based on cubic spline interpolation algorithm, the CGF method can effectively eliminate the cloud obscuration of MODIS FSC images using the observations in cloud-free days for each cloud pixel. The overall mean absolute error of the CGF method is 0.092, indicating that this method has relatively high accuracy on cloud reduction and snow information retrieval.
The most persistently SCAs are concentrated in the southern and western edges of the TP, where there are huge mountains. In the interior of the TP, snow cover is relatively small and less persistent. Moreover, there exist strong intra-seasonal variations in SCA over the TP. But the maximum snow accumulation and melting times over the year vary in different elevation ranges.
The very high inter-annual variability of SCD and SCA during 2001 to 2011 is discovered. About 34.14% (5.56% with a significant decline) and 24.75% (3.9% with a significant increase) of the study area shows declining and increasing trend in SCD, respectively. However, a longer time series of data need to be examined to obtain some definitive conclusions about temporal trends. Some of the inter-annual fluctuation of snow cover can be explained by the high negative correlations which are observed between the snow cover and in situ air temperature. If global warming goes on, the decrease of snow cover caused by the high negative correlation with temperature in February, April, and May, may result in significant changes in the river flows and water resources in the TP, particularly in spring. This will have impacts on ecosystem, irrigation dependent agriculture, and water resources in the densely populated downstream areas.
This study was supported by the National Basic Research Program of China (grant no. 2010CB951403), the National Natural Science Foundation of China (grant no. 41071227), the National Natural Science Foundation of China (Grant No. 91025001), and the FP7 project “CEOP-AEGIS” (Grant No. 212921).