Estimating Temperature Fields from MODIS Land Surface Temperature and Air Temperature Observations in a Sub-Arctic Alpine Environment

Abstract

Spatially continuous satellite infrared temperature measurements are essential for understanding the consequences and drivers of change, at local and regional scales, especially in northern and alpine environments dominated by a complex cryosphere where in situ observations are scarce. We describe two methods for producing daily temperature fields using MODIS “clear-sky” day-time Land Surface Temperatures (LST). The Interpolated Curve Mean Daily Surface Temperature (ICM) method, interpolates single daytime Terra LST values to daily means using the coincident diurnal air temperature curves. The second method calculates daily mean LST from daily maximum and minimum LST (MMM) values from MODIS Aqua and Terra. These ICM and MMM models were compared to daily mean air temperatures recorded between April and October at seven locations in southwest Yukon, Canada, covering characteristic alpine land cover types (tundra, barren, glacier) at elevations between 1,408 m and 2,319 m. Both methods for producing mean daily surface temperatures have advantages and disadvantages. ICM signals are strongly correlated with air temperature (R² = 0.72 to 0.86), but have relatively large variability (RMSE = 4.09 to 4.90 K), while MMM values had a stronger correlation to air temperature (R² = 0.90) and smaller variability (RMSE = 2.67 K). Finally, when comparing 8-day LST averages, aggregated from the MMM method, to air temperature, we found a high correlation (R² = 0.84) with less variability (RMSE = 1.54 K). Where the trend was less steep and the y-intercept increased by 1.6 °C compared to the daily correlations. This effect is likely a consequence of LST temperature averages being differentially affected by cloud cover over warm and cold surfaces. We conclude that satellite infrared skin temperature (e.g., MODIS LST), which is often aggregated into multi-day composites to mitigate data reductions caused by cloud cover, changes in its relationship to air temperature depending on the period of aggregation.

1. Introduction

High latitudes in the Northern Hemisphere have experienced significant recent warming, with Yukon and parts of Alaska experiencing the greatest warming of sub-arctic environments over the last 50 years [1,2]. This warming trend is expected to continue throughout the Arctic [3]. Temperature changes will affect many aspects of northern alpine ecosystems and the associated cryosphere, including snow extent, tundra land cover composition and distribution, permafrost, net ecosystem productivity, and population dynamics of animals and plants [4].

Mean surface temperature provides a fundamental measure for understanding change occurring in Arctic, sub-Arctic and alpine land surface processes. However, the absence of fine scale, continuous temperature monitoring over large geographical areas makes identifying climate induced changes difficult. Air temperature is commonly measured hourly at ground-based monitoring stations which are usually sparsely located in valley bottoms, and thus the lower end of the elevation gradient in alpine regions. For example, in Yukon, Canada, which has an area of 483,450 km², there are only eleven meteorological stations maintained by Environment Canada ( http://www.climate.weatheroffice.gc.ca). Furthermore, in the mountainous southwest Yukon the highest elevation station is 807 m above sea level. Spatial interpolation of air temperature data can lead to considerable uncertainties in the resulting temperature fields [5], especially at higher elevation.

An increasingly common method for tracking surface temperature trends in the Arctic involves the use of infrared satellite measurements of surface temperature [6]. Polar and near polar orbiting satellites exhibit progressively overlapping swaths, thus at high latitudes many ground surface observations are available each day. However, the swath overlap of near polar orbiting satellites over sub-polar regions is often insufficient to construct diurnal temperature curves because of cloud cover, thus necessitating a robust gap filling method to convert single daily observations to daily averages.

Land surface temperature (or skin temperature) is likely a better descriptor than air temperature for processes that are strongly linked to the ground surface such as low stature arctic vegetation growth, permafrost dynamics, and gas fluxes. Terrestrial ecosystems are often described in terms of their characteristic annual and seasonal temperature and precipitation patterns [7,8]. The distribution and abundance of arctic vegetation is dictated to a large extent by summer temperature, which can be characterised by satellite derived land surface temperature metrics [9,10].

MODIS Terra and Aqua satellites are in sun-synchronous near polar orbits, with Aqua in an ascending orbit and Terra in a descending orbit with equatorial crossings at 10:30 AM for Terra and 1:30 PM for Aqua (local solar time). Both satellites require approximately 90 min to complete an Earth orbit. MODIS Terra data became available in February 2000 and Aqua data became available in July 2002. At the latitude of our study site (60°N), Aqua coverage coincided with dawn and predawn. MODIS is in a low altitude orbit (705 km), and has 36 spectral bands thus improving the cloud detection ability of MODIS over previous satellites.

The daily record of LST observations is limited to those from the MODIS Terra platform between 2000 and 2002, which halves the amount of daily observations provided by the Terra and Aqua ensemble. Validation of MODIS land surface temperature products has emphasised the night-time product [11,12], which is easier to validate. The effects of relative humidity, wind speed, soil moisture, air temperature and sensor view zenith angle on night time land surface temperatures have been investigated, of which only sensor view zenith angle showed a weak influence on LST error propagation [11].

Analysis of thermal images and up to 18–20 daily MODIS LST returns collected over wet polygonal tundra in Siberia [13] and high Arctic tundra in Svalbard [14] indicated that several improvements regarding the performance of MODIS LST over tundra should be considered: (i) improved cloud cover masking and gap filling techniques; (ii) accounting for water bodies; (iii) accounting for snow cover and soil properties. The acquisition of daily satellite observations is complicated by extensive cloud cover commonly experienced in the arctic, and summer cloud cover has been shown to be increasing in the circumpolar arctic [15]. Cloud contamination, due to a failure of the cloud detection algorithm, is a known cause of LST error in the split window temperature extraction method [16]. However cloud contamination continues to be an issue with MODIS LST data where approximately 15% of data contain unidentified cloud contamination [17,18].

Due to the spatial and temporal autocorrelation of thermal infrared satellite data, spatial distances of 100–300 km and temporal intervals of up to two days were sufficient to represent missing data [19]. Furthermore, the diurnal cycle of LST could be modeled with solar geometry and two daily LST data points [20]. Maximum air temperature can be modelled with the diurnal cycle, cloud fraction and minimum air temperature derived from night-time LST [21]. The use of a generic modelled diurnal cycle in the interpolation of surface temperature to a diurnal cycle has also been proposed [22], where maximum daily LST, under clear sky conditions will correspond to peak solar insolation, with a small adjustment for phase lag [22] and minimum LST will correspond to time of sunrise [23]. It is important when aggregating LST observations in high latitudes to pay close attention to daytime length, satellite tile area and over pass time. These factors will contribute to the magnitude of within scene and between scene measurements, especially for observations on the continental scale.

A strong linear correlation between AVHRR LST and air temperature has been reported for temperatures < 0 °C [6]. Many recent studies have found strong linear correlations between MODIS LST (maximum, minimum and average) and air temperatures for many land cover types in Africa [21], in Portugal [24], on the Tibetan Plateau [25] and over the conterminous United States [26]. Steps have been taken to compare MODIS LST with much coarser spatial resolution passive microwave surface temperatures [27] and climate reanalysis products [28]. The typical range of errors when correlating LST to air temperature is approximately 2–3 °C [29] irrespective of the methodology, spatial or temporal resolutions.

This study was conducted to assess the viability of using an interpolated air temperature curve and a single day time LST value for the purpose of extending the daily average LST in a mountainous sub-Arctic region when (1) persistent cloud cover often reduces LST acquisition to a single daily value and (2) for early years of Terra data acquisition, when MODIS Aqua was not operational. We produced a Interpolated Curve Mean Daily Surface Temperature (ICM) product by combining single daily tiled day-time MODIS LST observations (regardless of acquisition time during the day) with the daily average air temperature and daily air temperature curve, using data from 2008. We then compared daily average temperatures from these seven sites to daily average LST created from minimum and maximum LST values (MMM) produced from both MODIS Aqua and Terra swath data. We compared both models to daily average air temperature from seven independent meteorology stations located on glacier, barren and tundra land covers. Lastly, we aggregated the MMM model to 8-days and compared the result to similarly aggregated air temperature to assess the consequences of aggregating LST data. The methods outline above identify limitations in previous approaches by (i) incorporating air temperature observations with LST, to bolster limited LST observations, instead of trying to convert LST to air temperature; and (ii) requiring the use of information about the diurnal curve and data availability for improving LST aggregation.

2. Study Area

The southwest Yukon study site (Figure 1) is bounded by the Yukon–Alaska border (141°W) on the west and southwest, the Yukon–British Columbia border on the south (60°N), 134°W on the east and 62.5°N on the north. Eastern portions of the study area are situated in the rain shadow of the St. Elias Mountain range caused by the westerly atmospheric circulation of air masses originating in the Arctic Ocean, Gulf of Alaska and the Bering Sea all of which are controlled by the position of the Aleutian Low [30]. However, this area receives enough annual precipitation (approximately 30–80 cm) that the semi-arid and arid landscape impediment to the function of the MODIS LST algorithm [31] is likely negligible. This region of the Yukon contains four Environment Canada meteorology monitoring stations situated near airports or airfields and located at elevations between 599 m and 807 m (

Table 1. Description and location of meteorology stations in southwest Yukon.

Station	Elevation (m)	Lat./Long.	Description	Measurement Frequency	Distance from Burwash A (km)
Whitehorse (A)	706	60.71°N; 135.07°W	Environment Canada Monitoring Station	On the Hour	NA
Haines Junction	599	60.77°N; 137.58°W	Environment Canada Monitoring Station	On the Hour	NA
Carmacks	542	62.12°N; 136.19°W	Environment Canada Monitoring Station	On the Hour	NA
Burwash (A)	807	61.37°N; 139.05°W	Environment Canada Monitoring Station	On the Hour	NA
John Creek	1,408	61.20°N; 138.25°W	Grass meadow surrounded by tall shrubs in large patches, little exposed bedrock	Hourly average of 5 min readings	47.1
Ruby Range–North	1,926	61.25°N; 138.19°W	Sparse vegetation, much exposed bedrock, no shrubs	Hourly average of 5 min readings	48.5
Pika Camp	1,635	61.21°N; 138.28°W	Grass and sparse low shrub, little exposed rock	On the Hour	45.0
Transect Canada Creek	2,184	60.88°N; 138.97°W	Sparse vegetation, much exposed bedrock, no shrubs	Hourly average of 5 min readings	55.4
Transect Duke River	2,214	60.94°N; 138.90°W	Sparse vegetation, much exposed bedrock, no shrubs	Hourly average of 5 min readings	49.0
South Glacier	2,280	60.82°N; 139.13°W	Glacier	Hourly average of 5 min readings	61.8
North Glacier	2,319	60.91°N; 139.16°W	Glacier	Hourly average of 5 min readings	51.2

). These stations provided diurnal temperature curves, which when used in conjunction with MODIS LST, were used to predict mean daily LST. Seven validation meteorological monitoring stations are located on barren, tundra and glaciated land cover at elevations between 1,408 m and 2,319 m (Table 1).

3. Methods

3.1. MODIS Version 5 Clear Sky Daytime Land Surface Temperatures (LST)

The MODIS LST data used in this study are the 1 km gridded clear-sky reprocessing version 5 MOD11A1 h11 v02 tile (ICM method) and MODIS Terra and Aqua (MOD11L2 and MYD11L2) LST swath data (MMM method). The 1 km gridded MODIS land surface temperature is produced using the split window technique [31] that uses MODIS bands 31 and 32 (10.78–11.28 μm and 11.77–12.27 μm respectively). This technique uses a global land surface emissivity derived look-up table [32] to correct emissivity variations. MODIS data were downloaded from the Land Processes Distributed Active Archive Center—LPDAAC ( https://lpdaac.usgs.gov/). The MOD11A1 data are produced from Land Surface Temperature (MOD11L2) swath data, which are dependent on several parameters related to the satellite platform and sensor and is generated using the Product Generation Executive (PGE16) code. The MOD11L2 data are comprised MODIS L1B calibrated and geolocated radiances, geolocation, cloud mask, atmospheric profiles, land and snow cover. The MOD11A1 product used in this study is produced by mapping daily single clear sky observation MOD11L2 swath data onto 1 km tiled grids in sinusoidal projection with an accuracy of 0.02 K.

3.2. Pre-Processing of MODIS Terra Tile and Terra and Aqua Swath Data

The version 5 MOD11A1 reprocessing data contain single observations in each grid cell, rather than averaged swath observations (averaged observations were contained in earlier reprocessing versions) where multiple inputs to grid cells from overlapping swaths occur. MODIS LST (MOD11A1) tile data can contain non-sequential and temporally disparate surface temperature data. Spatial subsets matching the study area were extracted from the 2008 MOD11A1 tile data. MOD11A1 LST and the coincident local solar View Time (VT) were reprojected from sinusoidal projection to 1 km gridded geotiff Albers Equal Area Projection, (WGS84) using nearest neighbour resampling [33]. View time is extracted because surface temperature (and thus LST) is critically dependent on solar azimuth, and thus the time of day. These LST product layers were converted to local solar time and temperature using the conversion factors provided at http://lpdaac.usgs.gov/modis. The conversion of spatial subsets of the MODIS native local solar VT to local standard time requires longitude information for each grid cell. To facilitate conversion of view time from local solar time to local time the view time layers were subdivided into 1° longitude bands, across which uniform time conversion was assigned. This procedure introduces error in local time, which is at a maximum of 2 min at whole number lines of longitude and decreases to zero at the half way line between whole number lines. Since the local solar time to local time conversion is dependent upon the day of year, a new conversion is required for each day of the year. Thus a different time conversion is applied to each daily one degree of longitude subset.

Lastly, LST was filtered for cloud contamination. We adopted an approach similar to that used by Raynolds et al. [9] on AVHRR LST data for cloud masking of MODIS LST data. However we made several changes to complement the more robust cloud masking and better spatial resolution of MODIS LST data. If the MODIS LST grid cell in question was >4 K lower than the surrounding grid cells, it was eliminated from the analysis. Daily differencing was also applied to the ICM product. A difference of >5 K between the day in question and surrounding days determined elimination. This masking was implemented beyond the normal cloud mask to eliminate any anomalous values, such as those resulting from optical leaks [24].

LST view time in the MOD11A1 tiles forms irregular temporal patterns (an example is presented in Figure 2) resulting from multiple daily imaging of the study area. Cloud cover is a frequent phenomenon in this part of the Yukon and results in areas with no data in the MODIS products (i.e., the NoData values in Figure 2 largely result from cloud cover). For our study area between 1 April and 31 October 2008 the number of daily grid values containing no data compared to the total number of grids averaged 61% ± 26% (1 standard deviation).

The 2008 MYD11L2 and MOD11L2 swath data were reprojected from sinusoidal projection to 1 km gridded geotiff Albers Equal Area Projection (WGS 84) using nearest neighbour resampling from separate geolocation and LST data files. This processing step was batch run using MODIS Swath Reprojection Tool available from the LPDAAC.

3.3. Meteorology Station Air Temperature Data

The 2008 air temperature data used to provide diurnal curves for the ICM model were recorded on the hour at four long-term Environment Canada meteorology monitoring stations (Figure 1) and were obtained from Environment Canada’s National Climate Archive ( http://www.climate.weatheroffice.gc.ca/Welcome_e.html). Seven meteorology monitoring stations located in the southwest Yukon (Figure 1) provided the validation data set and collect air temperature as hourly averages of measurements recorded on five minute intervals, except the Pika Camp station which collected air temperature on the hour. Three stations were located in altitudinal tundra, two stations on non-vegetated rocky outcrops, and two stations on glaciers (Figure 1). Air temperature was measured with Campbell Scientific HMP45 sensors (Campbell Scientific Inc., Logan, UT, USA). Air temperature was typically recorded at 2 m above ground to an accuracy of ±0.1 °C. The validation data covered the 1 April to 31 October 2008 study period, except the Ruby Range North Station and John Creek Station data, which covered the period from 15 July to 31 October 2008. Table 1 provides information regarding the location, elevation and other characteristics of the four reference and seven validation meteorology station sites.

3.4. Interpolated Curve Mean Daily Surface Temperature (ICM)

Air temperatures recorded at Environment Canada stations in the southwest Yukon (2000 to 2008) were constrained in previous work to within 12 min of MODIS Terra LST observations from 8 AM to 2 PM local time; these produced strong linear correlations between LST (MOD11A1) and measured air temperature relationships (0.94 ≤ R² ≤ 0.96) [18]. The time of daytime LST capture did not affect the LST to air temperature relationship, nor did the season in which it was collected. The strong linear correlations between daytime LST and air temperatures provide the theoretical basis for production of the Interpolated Curve Mean Daily Surface Temperature (ICM) model through the spatial interpolation of diurnal air temperature curve shape.

The ICM was constructed using daily air temperature curves from each of the four Environment Canada meteorology stations and spatially intersected MODIS LST and view time between 1 April and 31 October 2008. This procedure was comprised of three steps (Figure 3), outlined in detail below. The first step created daily daytime air temperature curves that facilitated the conversion of daily LST to daily average LST values. The second step created a spatial interpolation of the four sets of curves and converted time of LST capture to a ratio of what the air temperature should be at that time to the daily air temperature average. The third step created the ICM product.

Step 1: Air temperature curves for each of the four Environment Canada meteorology station were produced from air temperatures, recorded on the hour, between 8 h and 14 h (Pacific Standard Time) for the days between 1 April and 31 October. This interval includes multiple MODIS Terra swaths which results in a wide range of LST capture times found within the MOD11A1 tile. Only those air temperature curves displaying a polynomial fit of R² = 0.95 or greater were retained (see Figure 4 for an example). These curves were fit with a polynomial curve between 2nd and 5th order, in a successive manner starting at 2nd order and working towards 5th order. The range of polynomial curves was required because the 8 h to 14 h time interval produced curves with a range of inflection points and local troughs and peaks. The process of curve fitting stopped when the curve showed an R² value greater than or equal to 0.95. Due to abnormal atmospheric conditions (e.g., inversions), approximately 10% of the period indicated previously did not have curves greater than or equal to R² = 0.95, and so were omitted from this study.

Step 2: Four Thiessen polygons (polygons that are defined by one internal point each, where the sides of each polygon are closest to each internal point in relation to the other points) were created around the four Environment Canada stations (Figure 1). In each Thiessen polygon, the corresponding daily air temperature curve was used to determine what the air temperature should be at the time of LST capture for each MODIS grid cell contained within the Thiessen polygon. This value was then used to create a normalized dimensionless value by dividing into average daily air temperature (see Figure 4 for an example).

Step 3: The daily unitless values created in Step 2 were multiplied by the corresponding MODIS LST grid value to produce the ICM, which retained the same spatial resolution as the input MODIS LST data. MODIS LST over the study area was collected during the overpass time, not specifically on the hour as the Environment Canada station air temperature was. Therefore, to fully exploit the MODIS LST data, the curves in step one are required.

3.5. Mean Daily Surface Temperature Calculated from Maximum and Minimum LST (MMM)

Daily mean LST was calculated as the average of maximum and minimum LST and retained the same spatial resolution as the input MODIS LST data. Maximum LST was determined through the following steps: (1) MODIS Aqua and Terra swaths that corresponded to within 1 h before and after solar noon at 61°N and 139.5°W were identified; (2) Images were processed as outlined in Section 3.2; (3) The maximum value from multiple daily images (2 to 4 depending on how the 90 min repeat orbital corresponds to the 2 h interval around solar noon) was selected for each 1 km grid cell. Minimum LST was determined through the following steps: (1) MODIS Aqua swaths that correspond to 2 h before dawn to dawn were selected. Due to orbital configurations only MODIS Aqua imaged the study area for the two predawn hours; (2) Image processing was undertaken as outlined in Section 3.2; (3) The minimum value for daily images was selected for each 1 km grid cell. Depending on how the satellite repeat orbit of ∼90 min corresponds to the two hour dawn interval, 1 to 2 Aqua images were available with which to determine the minimum value. Minimum values were calculated for the period of 1 May to 31 August, as this is the period when two hours before dawn to dawn intersect with satellite coverage. During the remainder of the year at this latitude, dawn occurred too late in the day to coincide with MODIS Aqua coverage. Cloud contamination masking was performed in a similar manner to that for the ICM model as detailed in Section 3.2. A daily average was produced if both maximum and minimum values occurred on the same day. Lastly, the daily MMM product was aggregated to 8-day averages and compared to the coincident 8 day averages of air temperature. An 8-day average was produced if at least one minimum and one maximum value occurred during the 8-day period, but not necessarily on the same day.

All GIS analysis and data manipulation was conducted with ESRI's ArcGIS 9.2. Statistical analysis was conducted using the R statistical software package [34].

4. Results

ICM validation showed statistically significant strong linear correlations (0.72 ≥ R² ≥ 0.85) with mean daily air temperature (Figure 5;

Table 2. Coefficient of determination (R²), root mean square error (RMSE) and skew between daily air temperature averages and modeled mean daily surface temperature (ICM and MMM).

Station-Model	R2	RMSE (K)
South Glacier–ICM	0.86	4.09
North Glacier–ICM	0.82	4.13
Transect Duke River–ICM	0.72	4.62
Transect Canada Creek–ICM	0.76	4.90
John Creek–ICM	0.72	4.36
Ruby Range North–ICM	0.77	4.15
Pika Camp–ICM	0.80	4.53
All Stations–MMM	0.90	2.67
All Stations–MMM (8 day Composite)	0.84	1.54

). RMS errors ranged between 4.09 and 4.90 K (Table 2), which were approximately 1 to 2 K larger than typical when compared to air temperatures [29]. Variations in R² values did not correlate with distance of the interpolated diurnal curve form (Table 1) indicating factors other than spatial interpolation were playing a dominant role in correlation variability. The strength of the R² values correlated to land surface structure, with the largest R² values corresponding with glacier and sparsely vegetated surfaces, and the lowest R² values corresponding with sites with large relief in topography, tall shrubs and complex vegetation canopies (Tables 1 and 2). The steepest slopes and largest y-intercepts most consistently occurred on glaciers, and to a lesser degree on exposed rock (barren) that displayed high topographic relief. All of the validation data sets were contained in the Burwash Landing Thiessen polygon (Figure 1) for the ICM method and all are above tree-line.

The MMM show statistically significant strong linear correlations (R² = 0.90, RSME = 2.67 K) with mean daily air temperature (Figure 6; Table 2). However, only 44 out of a possible 1,575 LST observations intersecting the seven air temperature monitoring stations, had maximum and minimum values originating on the same day (

Table 3. Number of daily minimum and maximum values observed and calculated (mean) for seven meteorology stations for period between 1 May and 31 August 2008.

Station	Maximum Count	Minimum Count	Mean Count
South Glacier	4	42	2
North Glacier	7	50	3
Transect Duke River	12	50	7
Transect Canada Creek	7	60	5
John Creek	25	27	12
Ruby Range North	25	23	6
Pika Camp	20	31	9

). The majority of the observations were from the predawn over the icefields, whereas the number of observations were roughly equal further to the east over the lower elevation tundra. Aggregating to 8-day averages had little effect on the number of LST values correlated to air temperature (

Table 4. Number of 8-day composite values observed and calculated (mean) for seven meteorology stations for period between 1 May and 31 August 2008.

Station	Maximum Count	Minimum Count	Mean Count
South Glacier	2	14	2
North Glacier	4	16	4
Transect Duke River	6	16	6
Transect Canada Creek	5	16	5
John Creek	12	13	11
Ruby Range North	11	13	9
Pika Camp	15	12	11

). However, the aggregation did have a noticeable effect on the linear regression (R² = 0.84, RSME = 1.54 K), where the slope become flatter than the 1:1 line and the y-intercept became positive and increased by 1.6 °C compared to the daily MMM result. The RMS error reported here for daily and aggregated MMM method is similar to that reported elsewhere [29].

5. Discussion

The simple process of temperature curve shape interpolation applied across large geographic areas defined by Thiessen polygon allocation produced Mean Daily Temperatures. The interpolation of diurnal air temperature curve at distances of 45–62 km from Burwash (A) station, is supported by spatial auto-correlation of LST values across distances of 100–300 km [19]. The reduced requirement of LST data (e.g., one point instead of two or more daily points) in the ICM is important for regions where extensive cloud cover is a perennial concern and where the number of daily overpasses is reduced compared to high latitude locations. Although the average LST calculated from minimum and maximum LST values (MMM method) provided better agreement with daily air temperature averages, the very small number of daily averages indicated that direct calculation of average LST was not always feasible for the study area. The ICM model should find application in many instances of land surface monitoring where available daily LST observations are insufficient represent measurements at dawn and solar noon. Regardless of which of the two above methods are considered, both methods require independent observations of air temperature and surface thermal infrared temperature, which cannot replace each other, even though they are strongly correlated.

The LST diurnal temperature curve is mostly, but not wholly, comprised of absorbed radiation and atmospheric interactions [22,35]. The high measurement accuracy of the air temperature and MODIS LST is very small compared to RMSE values presented here, which indicates that the variability in the results is a function of a time lag difference in processes that affect LST and air temperature differently and not measurement error. The difference between ICM and daily average air temperature was likely the combination of several factors, including:

(1)

Cloud contamination or surface cloud shadow likely contributed to variation in the ICM values because this method uses a single daily input compared to two inputs required for the MMM method, which has been shown to increase correlation to air temperature [24]. Variation in LST was likely influenced by cloud contamination, which likely disproportionately affected the ICM product because of its dependence on observations of maximum LST [36]. Furthermore, cloud shadow can cause differences in LST [11] across spatial and temporal scales and might also contribute to the influence of cloud contamination in LST.

(2)

Decoupling of air temperature and LST at higher temperatures, caused by low albedo vegetation cover, was likely playing a role in the high RMS error reported here, especially for the ICM product because its sole input originates from the warmest part of the day. This interpretation was supported by the two non-vegetated sites (glaciers) displaying the lowest RMSE values for the ICM. Meltwater ponds at the glacier sites are small compared to the LST grid cell. However the air temperature which was recorded over glacier, snow and ice was contained within a LST grid cell which contains high percentages of talus, is the probable reason for the large y-intercepts. As expected the MMM LST product had a lower RMS error likely due to the moderating effect of the minimum temperature. However, the daytime MODIS LST product has a larger confidence in identifying daytime cloud cover cloud mask compared to the night-time mask [36], which should act to minimize cloud contamination in the ICM product compared to the MMM product.

Error could include the moderating effect on LST caused by surface water in tundra environments [13] or soil heat or soil water storage [14,22]. Neither the validation sites, nor the Haines Junction Environment Canada monitoring site contained large bodies of water within the MODIS grid cells. Nevertheless, the effects of small and ephemeral surface water bodies on soil heat or water storage remain unaccounted for in both models. Another error could be the larger error in day time LST caused by angular anisotropy [37] when compared with night time LST. The larger RMS error values reported for the ICM product (compared to the MMM which included pre-dawn minimum LST values) likely accentuates angular effects, where the MMM minimizes them.

The differences in slope and y-intercept values between the daily average and the 8-day aggregation results from cloud cover influence of LST observations used in the 8-day period. The daily MMM values are calculated from the same day that has both clear sky minimum and maximum values, which suggests there is little cloud influence compared to the 8-day aggregation. Because the majority of values used in the 8-day aggregation did not have matching daily maximum or minimum values, we can deduce that cloud cover is influencing the LST value to a greater extent than the MMM method. This interpretation is supported by the 8-day values being warmer than the daily values when the temperature was below 0 °C and colder than the daily values when the temperature was above 0 °C. This effect should be considered when aggregating LST values for comparison to air temperature.

Application of the two models presented here can provide a relatively high resolution (1 km) spatial/climatic product with which to evaluate the glacier surface processes and tundra dynamics in the sub-arctic, either on their own or in conjunction with other modelled products in a gap filling capacity. For low statured arctic vegetation, LST is likely a better indicator of plant temperature than air temperature [38] and there are well documented relationships between leaf temperature and photosynthesis of arctic tundra plants [39].

Although improvements in cloud masking will reduce the negative effects of cloud contamination in LST data as a whole, the cost in terms of the reduction of data availability and quality could be significant in areas of the cryosphere that experience large amounts of cloud cover. Future work should be devoted to the identification of cloud contaminated LST data and the degree to which contaminated LST values are affected. By quantifying the degree of contamination, data volume will be preserved and quality improved. Furthermore, the effect that different LST aggregation schemes have on the aggregated products should also be investigated with consideration being given to land cover and regional climate. Statistical techniques such should be investigated to make the most use of limited ground station data. The implementation of the ICM method for every LST observation could provide a powerful technique to fill gaps and aggregate LST data. Integration with interpolated air temperature products or other forms of surface temperature, such as passive microwave, could provide other sources of data. However, mismatches in spatial and temporal scale will pose challenges, especially when trying to mitigate the influence of persistent cloud cover.

6. Conclusions

Mean Daily Surface Temperature modelled from (1) MODIS LST tile data and interpolated daily air temperature curve shape (ICM) and (2) MODIS maximum and minimum data were compared to mean air temperatures from seven alpine locations covering tundra, talus and glacier land cover and elevations between 1,408 m and 2,319 m. The correlations at the seven locations for the ICM method were strong (R² values between 0.72 and 0.85) and all the correlations displayed positive slopes and y-intercepts. RMS errors ranged between 4.09 and 4.90 °C. The average LST calculated from the MMM algorithm had an even stronger correlation with air temperature (R² = 0.90) and smaller variability (RMSE = 2.67 K) compared with the ICM. However, very few daily averages were produced with the MMM method because of almost continuous cloud cover over the study area. The aggregation of MMM to 8-day LST averages produced similar results as the daily averages (R² = 0.84; RMSE = 1.54 K), but had a warmer y-intercept and flattened slope, which is most likely the result of the insolating effect of clouds. Our results indicate that the ICM is a viable gap filling technique for situations where only single observations of day-time LST are available, but the daily average of maximum and minimum LST values produce a better result compared to air temperature. Finally, the common practice of aggregation of daily LST values to 8-day averages requires further research to determine the effect of cloud cover and land cover on this product and how it is applied.