Study area

The Tibetan Plateau is located in western China and ranges from 73.44° to 103.98°E and 26.87° to 39.81°N. The plateau has a total area of approximately 2,530,000 km² and covers more than one fourth of the land area of China (Figure 1). Known as the “third pole of the earth” (Qiu 2008), the Tibetan Plateau is the highest and largest plateau in the world, with an average elevation of more than 4000 m above sea level. As a result of the high altitude, the climate of the plateau is quite different from that of other regions at nearly the same latitudes. The mean annual air temperature varies between −2.8 and 11.9°C, depending on the topographic conditions. The annual precipitation decreases from more than 800 mm in the southeast to lower than 200 mm in the northwest, with most precipitation occurring from May to September. Influenced by the climate and altitude, land cover on the Tibetan Plateau also shows a clear geographic distribution, following the order of forest, grassland, and desert from the southeast to the northwest (Wang et al 2005).

FIGURE 1

Location and topography of the study area. (Map by Y. Xu)

Data set

The AQUA MODIS LST product (MYD11A2) spanning 2003 to 2010 was obtained from the Warehouse Inventory Search Tool (WIST 2011). MYD11A2 is an 8-day composite daytime and nighttime LST product with a 1-km spatial resolution that is derived from AQUA MODIS data. The LST product is first computed using the generalized split-window algorithm proposed by Wan and Dozier (1996) and then produced by averaging the LST values under clear-sky conditions over an 8-day period (Wan et al 2002) for the purpose of removing cloudy pixels. In this study, the daytime LST was employed with an overpass time of approximately 1:30 PM (local solar time).

First, the MODIS daily LST product was computed using the generalized split-window algorithm (Wan and Dozier 1996):

where T_s is LST; T₃₁ and T₃₂ are the brightness temperatures in MODIS band 31 and 32, respectively; ε₃₁ and ε₃₂ are the surface emissivities in these 2 bands estimated using the classification-based emissivity method (Snyder et al 1998), according to land cover types determined by the MODIS land cover product (MOD12Q1) and snow cover product (MOD10_L2); and C, A₁, A₂, A₃, B₁, B₂, and B₃ are regression coefficients given by multidimensional lookup tables.

Then, MYD11A2 was generated by averaging the daily LST values under clear-sky conditions over an 8-day period (Wan et al 2002) for the purpose of removing cloudy pixels. Each year of the MYD11A2 LST data contained 46 eight-day images.

The MODIS LST product (V005) were evaluated using in situ observations from the Asia–Australia Monsoon Project on Tibetan Plateau (CAMP/Tibet), with an average error of 2.36 K (Wang et al 2012). Compared with the validation results of MODIS LST in other regions (Wan 2008; Wang et al 2008b; Crosman and Horel 2009), the accuracy error in the Tibetan Plateau is somewhat higher, which may be due to the complicated topography and atmosphere conditions, but it is still acceptable for practical applications.

The corresponding MODIS land cover product (MCD12Q1) was also used in this study. The MCD12Q1 product incorporates 5 land cover classification schemes derived from supervised decision-tree classification methods (Friedl et al 2002). We used the International Geosphere–Biosphere Program (IGBP) classification scheme, which contains 17 land cover classes (water, evergreen needleleaf forest, evergreen broadleaf forest, deciduous needleleaf forest, deciduous broadleaf forest, mixed forests, closed shrubland, open shrubland, woody savanna, savanna, grassland, permanent wetland, cropland, urban and built up, cropland/natural vegetation mosaic, snow and ice, and barren or sparsely vegetated). Using MODIS Reprojection Tool software, all MODIS data were mosaicked and reprojected. In addition, elevation data of the plateau were acquired from the Shuttle Radar Topography Mission 90-m digital elevation model (SRTM3 DEM, version 4.1) and were resampled to a 1-km resolution to coincide with the MODIS data.

Daily maximum air temperature (T_a) observation data of 79 weather stations from 2003 to 2010 in the research area (Figure 1) were employed to evaluate the fitting performance of harmonic analysis. These T_a data were averaged over each 8-day period for comparison with the MODIS LST data.

HANTS simulation of the LST time series

To perform the HANTS analysis, 5 parameters should be set: the number of frequencies, a high/low suppression flag, a valid data range, the fit error tolerance, and the degree of overdeterminedness. A detailed description of these 5 HANTS parameters can be found in Roerink et al (2000). In this study, the number of frequencies was set to 2 (frequency 1  =  annual; frequency 2  =  semiannual). The high/low suppression flag was set to low, considering that cloud contamination leads to abnormal low LST values. The valid data range was from −50 to 70°C. The fit error tolerance was set to 6°C, and the degree of overdeterminedness was set to 7.

The HANTS algorithm was applied to the LST data on a per pixel basis in each year between 2003 and 2010 for the entire study area. Two examples of the HANTS-fitted and the original LST time series are given in Figure 2, whose locations are also labeled in Figure 1 as Pt 1 and Pt 2, respectively. In this figure, the light gray bars are the original LST temporal profiles, and the solid line is the HANTS-fitted LST temporal profile. The outliers in the original LST temporal profile that were rejected according to the fit error tolerance in the HANTS iterative procedures are shown as dark gray bars. Figure 2A represents an LST temporal curve that contains only 3 data gaps, and Figure 2B represents an LST temporal curve that contains 9 data gaps. The simulation results show that HANTS can properly and robustly fit the LST temporal curves even if there are conspicuous data gaps.

FIGURE 2

Examples of the HANTS-fitted and the original LST time series in 2003 for 2 pixels. (A) 100.67°E, 36.87°N; (B) 79.28°E, 33.73°N (location of pixels: see Figure 1).

Observed maximum air temperature from 79 weather stations was compared with the original and HANTS-fitted LST to evaluate the fitting performance. Figure 3A shows the relationship between valid values in the original LST and corresponding T_a values. There is an obvious linear correlation between these 2 data, with a correlation coefficient of 0.7885 (P < 0.01). The LST values are generally higher than T_a, but some LST values are lower than T_a. Considering the overpass time of AQUA MODIS, these pixels most likely contain undetected residual clouds, which result in abnormally low surface temperature. Figure 3B shows the relationship between HANTS-fitted LST and T_a. The light gray triangles are the fitted LST values for the invalid values in the original LST caused by cloud cover. These triangles show similar relationship between LST and T_a as the dark gray diamonds in the figure, indicating that no significant inconsistencies were introduced in these fitted values. The correlation coefficient between HANTS-fitted LST and T_a is 0.7751 (P < 0.01), which is slightly lower than that between original LST and T_a and suggests the good performance of HANTS fitting. In addition, the abnormally low values in the original LST do not show very low surface temperatures relative to air temperature, indicating that the influence of residual clouds on the original LST are effectively eliminated.

FIGURE 3

Comparison between observed T_a and remotely sensed LST. (A) T_a versus original LST; (B) T_a versus fitted LST.

The success of HANTS fitting can also be testified by LST spatial distributions. Figure 4 shows the original and HANTS-fitted LST map in the period of 361 days of year (DOY) 2010. In the original LST map, there are several large data gaps in the northeastern, northwestern, and southeastern part of the plateau, which are caused by cloud cover. The HANTS-fitted LST map shows good consistency with the original LST map for the clear-sky pixels, and the missing values in the original LST are effectively filled.

FIGURE 4

LST images over the Tibetan Plateau in the period of 361 DOY, 2010. (A) Original LST; (B) HANTS-fitted LST.

The fitting results were saved as 5 harmonic components: the mean LST and the amplitude and the phase of the first (frequency  =  1) and the second (frequency  =  2) harmonics. The mean LST and amplitudes were represented in units of degrees Celsius, and the phases were represented in units of day of year. An example of the fitted LST curve and its harmonic components is given in Figure 5. The HANTS-fitted LST curve and the mean LST term are shown on the left axis, while the first and second LST harmonic curves are shown on the right axis. The HANTS-fitted curve is achieved according to the summation of the mean LST and the 2 harmonics. The first and second harmonics represent the annual and biannual cycles of the LST time series, respectively, according to their frequencies. Their amplitudes denote the periodic fluctuation range of the LST time series based on the mean value, and their phases illustrate the peak time of the LST periodic cycles, which indicate the temporal character of the LST fluctuation.

FIGURE 5

Examples of the HANTS components and fitted curve in 2003 for the pixel located at 79.28°E, 33.73°N (the mean term and HANTS-fitted LST are shown on the left axis; the first and second LST harmonics are shown on the right axis).

Figure 5 proves that the first harmonic shows a much larger amplitude than the second harmonic. Moreover, the mean amplitude of the first harmonic of all pixels is 12.14°C, which is also much higher than that of the second harmonic (2.86°C). This result suggests that the first harmonic (annual term with frequency  =  1) contributes far more to the overall LST time series than the second harmonic (semiannual term with frequency  =  2), indicating a strong unimodal periodic pattern of the LST temporal series on the Tibetan Plateau. The second harmonic contributes less to the overall fit but helps improve the HANTS fit performance.

Spatiotemporal variations of LST

The HANTS simulation result of the LST time series is composed of the mean value of the LST temporal profile, amplitudes, and phases of the 2 harmonics. According to the preceding analysis, the mean and annual harmonic contribute more to the overall LST time series than the semiannual harmonic. Therefore, the 8-year (2003–2010) averaged mean LST and annual harmonic were used to represent the spatiotemporal variations of LST in the Tibetan Plateau. The mean term represents the annual mean LST value. The annual harmonic depicts the annual cycle of LST, its amplitude illustrates the seasonal fluctuation magnitude, and its phase illustrates the peak time of the LST annual cycle.

Figure 6 provides the 8-year (2003–2010) averaged mean LST on the Tibetan Plateau. The study area exhibits a highly contrasted spatial pattern of annual mean LST values, suggesting significant regional variabilities. The annual mean LSTs in most areas were generally within the range of 0 to 30°C, but some pixels were even larger than 30°C or lower than 0°C. Consulting the DEM map of the plateau (Figure 1), the spatial distribution of annual mean LST is influenced by the altitude. Areas with lower altitude generally have higher annual mean LST, indicating that higher altitude leads to lower temperature. The Qaidam Basin located in the northeast shows higher annual mean temperatures above 28°C due to its lower altitude and desert landscape. However, the annual mean LSTs of the South Tibetan Valley in the southwest, with relatively high altitudes ranging from 4500 to 5000 m, are also among the highest in the study area. This can be attributed to the area's lower latitude and valley terrain and to the influence of the warm Southwest Monsoon. The Kunlun Mountains, located in the northwest, and the Nyainqentanglha Mountains, located in the south, show lower annual mean temperatures below 3°C. Large lakes in the plateau, such as Qinghai Lake and Nam Lake, also show relatively lower annual temperatures.

FIGURE 6

Spatial distributions of the 8-year averaged annual mean LST of the Tibetan Plateau. (Map by Y. Xu)

Figure 7A and B represents the 8-year (2003–2010) averaged amplitude and phase, respectively, of the annual cycle of the LST series on the Tibetan Plateau. Figure 7A clearly reveals that the annual amplitude shows a spatial tendency that increases from the southeast to the northwest, which shows a similar spatial pattern with precipitation (Figure 7C). According to the surface heat balance, LST is controlled by solar radiation and surface thermal properties. High precipitation caused by the monsoon indicates heavy clouds, mainly in summer. The occurrence of clouds reduces solar radiation and in turn reduces the surface temperature, especially in summer. Lower summer temperatures decrease lower seasonal fluctuations. Moreover, in this arid and semiarid area, higher precipitation usually results in higher vegetation cover and higher soil moisture, which are more resistant to seasonal oscillations of LST.

FIGURE 7

Spatial distribution of the 8-year averaged (A) amplitude and (B) phase of the first harmonic in the Tibetan Plateau. (C) Precipitation map and (D) land cover map of the Tibetan Plateau. (Maps by Y. Xu)

Despite this general southeast–northwest difference, the most striking seasonal oscillation occurs at the Qaidam Basin, which is ascribed to its sparse vegetation cover, low soil moisture, and low altitude. The annual amplitude in this basin is usually higher than 15°C. In the Yarlung Zangbo Grand Canyon located in the southeast, the annual amplitude is mostly lower than 5°C, indicating slight seasonality. As the main water vapor channel of the Tibetan Plateau (Lin and Wu 1990), this area is highly influenced by the warm–humid Southwest Monsoon and shows a warm–temperate monsoon climate. Figure 7B demonstrates the spatial distributions of peak time for the LST annual cycle. Most of the study area has annual phase values ranging from 148 to 204 DOY, suggesting they achieve their peak value of the LST annual cycle at June and July.

In contrast to the amplitude, the spatial distribution of the annual phase does not present an obvious regular pattern, so it may be controlled by complex interactions of environmental variables, including solar radiance, underlying surface type, topography, vegetation abundance, climate, and soil moisture conditions. On the whole, the annual phase exhibits an inconspicuous southeast–northwest spatial tendency. The largest annual phases were found over large water bodies, suggesting that these areas have the latest peak time of LST due to their high heat capacities. High mountains, such as the Karakoram Mountains and Nyainqentanglha Mountains, also exhibit delayed peak times of LST for their permanent snow cover. Compared with the DEM map, there is an inconspicuous analogy between the spatial distributions of altitude and the annual harmonic, suggesting that altitude is not the decisive factor for the annual variations of the LST peak time.

LST variations for different land cover types

To analyze the LST dynamic characters of the different land cover types in the study area, the original MODIS land cover map was reclassified by aggregating the 17 IGBP categories into 5 major land cover classes: forest, grassland, water, snow/ice, and bare land (Figure 7D). Based on the land cover map, the average values of the annual mean LST and the amplitude and phase of the first harmonic for each land cover type were calculated (Table 1). Table 1 shows that bare land exhibits the highest annual mean LST (18.00°C), followed by grassland (17.24°C), and forest (16.62°C).

TABLE 1

General statistical results of the HANTS harmonic characteristics for each land cover type.

Vegetation cover has a dampening effect on LST due to the cooling effect of shading and evapotranspiration. For example, bare land displays a high temperature that can be related to its high insolation (caused by low cloud cover), low heat capacity, and lack of evaporative cooling (caused by low vegetation cover). In contrast, water displays a relatively low temperature because of its rather high thermal capacity and evapotranspiration. Thus, water shows a much lower mean temperature of 6.75°C, and snow/ice shows the lowest mean temperature of 0.53°C.

For the amplitude of the first harmonic, bare land also shows the most distinct annual fluctuation of LST, with an amplitude of 15.63°C, due to its low heat capacity and low evapotranspiration. Grassland, water, and snow/ice show lower seasonalities, with amplitudes of 10.40, 10.02, and 9.81°C, respectively. However, forest exhibits the lowest seasonal oscillation, with an amplitude of 9.48°C. The lowest amplitude of forest is probably due to the evaporation and shelter effect of the vegetation canopy, which efficiently reduce temperature in growing season, especially in summer. For the phase of the first harmonic, grassland and forest show relatively lower respective values of 170.88 and 174.40 DOY, indicating that their maximum LSTs occur on June 19 and June 23, respectively. Forest and grassland show earlier peak times than other land types, which can be explained by the high cloud cover and high vegetation abundance in the growing season (from July to September) that reduce the summer temperature. Bare land reaches its peak LST on 181.12 DOY, approximately 10 days later than grassland. Snow/ice reaches its highest LST on 191.74 DOY, approximately 10 days later than bare land. Water shows the largest amplitude of 200.95 DOY, suggesting that its peak LST occurs on July 19. The lag effect of water surface temperature results from its high heat capacity.