An algorithm to retrieve Land Surface Temperature using Landsat-8 Dataset

Soil moisture, surface temperature, and vegetation are variables that play an important role in our environment which in turn increases the demand for accurate estimation of certain geophysical parameters such as weather, flooding, and land classification. However, for accurate Land Surface Temperature (LST) estimation, remotely sensed data of key environmental forms were considered and applied in this research. The goal of this study was to apply a suitable algorithm for LST estimation from the Landsat-8 dataset that gives a great accuracy when compared with in-situ observations. Spatial and temporal Landsat-8 data were acquired which provided the analytical structure for linking specific data successfully due to fine resolutions. The data were then applied to determine brightness temperatures, vegetation cover, and surface emissivity which demonstrated the effectiveness of the Split-Window Algorithm as an optimum method for LST retrieval from satellite. The results show temperature variation over a long period of time can be used in observing varying temperature values based on terrain i.e. High temperatures in fully built up areas and low temperatures in the well-vegetated regions. Finally, accurate LST estimation is important for land classification, energy budget estimations as well as agricultural production.


Introduction
Remote Sensing which is the science of using instrument-based techniques in obtaining, analysis and interpreting data about objects or areas from a considerable distance without being in direct contact as against on-site observations (Lillesand et al., 2014).These techniques can vary from manual interpretation of aerial photos to the computer-based analysis of non-visible radiation collected by satellite-based sensors.Optical sensors can be broken down into three basic types: panchromatic, multispectral and hyperspectral.According to (Li et al., 2013, Qin andKarnieli, 1999), the extensive requirement of temperature data on a large scale for environmental studies and management activities of the Earth's resources by Environmentalists, Geologists and Health Practitioners has made the remote sensing of LST an important issue in recent decades and as such, many efforts have been devoted to the establishment of methodology for its retrieval.
Remote sensing of LST using infrared radiation gives the average surface temperature of the scene covered by the radiometer's Field Of View (FOV).The spatial variation in LST acquisition causes a decrease in a measurement accuracy which in turn makes it difficult to validate LST from space using in-situ measurements (Minnett and Barton, 2010).To, therefore, retrieve LST from satellite Thermal Infrared (TIR) data, factors such as radiometric and atmospheric calibrations, surface emissivity correction as well as the characterization of spatial variability over land cover for vegetation has to be quantified.These factors are dependent on atmospheric conditions and emissivity of the land surface materials.Therefore, the effects of both atmosphere and Emissivity must be corrected through an algorithm for accurate LST estimation (Sobrino and Raissouni, 2000).
According to Ustin et al., (2004), there is a growing awareness among environmental scientists that remote sensing can and must play a role in providing the data needed to assess ecosystems conditions and to monitor change at all spatial scales.Thus, acquiring LST from remotely sensed data becomes one of the significant factors in this study.As a key parameter of the surface energy budget, LST is directly related to surface energy fluxes and to the latent heat flux, evapotranspiration and water stress Torrion et al., (2014), surface longwave emission and computing soil moisture Cammalleri and Vogt, (2015) and for understanding meteorological and hydrological processes in a changing climate (Duan et al., 2014).This paper attempts to apply a Split Window (SW) algorithm, which requires two essential parameters (emissivity and transmittance) by using the Landsat-8 dataset.Keeping the accuracy of the LST estimate, we avoid the complicated expression of the algorithm and the difficult calculation of its parameters.A complete and detailed description is given to the derivation of the algorithm, which includes the theoretical basis for remote sensing of LST.

Review of Landsat for LST Retrieval
The Landsat Programme, launched in the early 1970s had its capacity innovatively expanded over time which in turn increased the quality and size of captured data by onboard sensors.So far, there have been eight Landsat satellites launched successfully with the commencement of the Landsat-1 in 1972.The Landsat satellites have a near-earth orbit operation, scanning a ground area with 705km altitude.Landsat 1-3 with both the infrared and visible bands had an 80m resolution for the 185x185km imagery.Landsat 4-5 were equipped with double multispectral sensors; Thematic Mapper and Multispectral Scanner, with seven high-resolution spectral bands.Landsat 6 encountered a failure due to its immediate loss after launching while the Landsat 7 ETM + was simply an advanced variant of the TM with a panchromatic band with high resolutions.Landsat-8 has two major instruments (TIRS and OLI), at the 30m spatial resolution for the SWIR, visible and NIR; 100m for the TIR and 15m for the panchromatic bands (Markham et al., 2013).The presence of double TIR bands in the Landsat-8 OLI/TIRS satellite in the atmospheric window (10-12µm) is the main advantage over previous Landsat series.
In the research by Mallick et al., (2008), the potential of remote sensing to study the urban morphology is presented by estimating the spatial distribution and intensities of geophysical parameters using Landsat 7 ETM+.Their results show that the satellite derived emissivity values are in the acceptable range and the NDVI and fractional vegetation cover are effective in deriving surface emissivity.Sobrino et al., (2004) used Landsat TM 5 data to compare the proposed algorithms of Jiménez-Muñoz and Sobrino, (2003) and Qin et al., (2001) for LST estimation.Their results showed RMSD values of 2K and 0.9K respectively when 9x9 pixel samples are considered.Using Landsat TM/ETM + data with a radiative transfer model, Li et al., (2004) suggested that the average difference is 0.98℃ for Landsat 7 and 1.47℃ for Landsat 5 when compared.However, all these studies used data providing daily to half-hourly coverage.They are appropriate for operational monitoring, but their low spatial resolution limits their use for detailed spatial analysis of land cover changes and intensively fragmented areas.

Split Window Algorithm
LST estimation from satellite TIR radiometers has been achievable through different systems; Single-Channel (SC), SW and Dual-Angle (DA) algorithms.Various studies concentrate on the utilization of sensors onboard the polar-orbiting satellites on account of their high spatial resolutions (Pandya et al., 2014).Accurate LST estimation by utilizing the SC algorithm require high-quality air transmittance codes in estimating air quantities, Land Surface Emissivity data, precise air profile and topographic impacts Sobrino et al., (2004) while involving the inversion of the Radiative Transfer Equation as shown in Equation [1].A substitute methodology utilized over ocean uses the differential climatic absorption of both adjoining channels at 10-12μm in the SW algorithm initially developed by McMillin, (1975) which does not need air profile data at acquisition time.The SW algorithm predicts that LSE in the TIR regions of the 10-12μm region is established (Atitar and Sobrino, 2009, Sobrino andRaissouni, 2000).The DA algorithm as proposed by Sobrino et al., (2004) consolidates both the emissivity at forward and nadir perspective as it incorporates just emissivity reliance with no unequivocal WV reliance.This method requires that estimations are made for a fundamentally lengthy path; if not, the calculation becomes uncertain.
where Lsensor is the radiance registered by the sensor, also referred to as Top of Atmosphere (TOA) radiance,    is the black body radiance related to the surface temperature by Planck's law and   is the LST,   and   are the upwelling and down-welling atmospheric radiances respectively in Wsr - 1 m -2 μm -1 ,  is the atmospheric transmissivity and  is the LSE.The selection of the LST algorithm for Landsat-8 data is based on a literature review from Markham et al., (2013) and Tang and Li, (2014) and the characteristics of the Landsat-8 satellite data.The single-window method requires a highquality atmospheric profile and is sensitive to uncertainties in the atmospheric corrections.Considering the simultaneous retrieval of the LST and the LSE, the Temperature and Emissivity Separation (TES) method may be a candidate.However, significant errors in the LST and LSE for the surfaces with low spectral contrast emissivity such as water, snow, vegetation can be caused by the TES (Gillespie et al., 2011).Since the Landsat-8 satellite observes the land almost at-nadir, the dual-angle algorithm was discarded.Among the various methods proposed for LST determination, the quadratic SW algorithm as shown in Equation [2] from Coll et al., (1994) and Sobrino et al., (1993) has received considerable attention because of its simplicity.
Where   is the LST,   and   are the at-sensor brightness temperatures for two TIR channels,  and  are coefficients and  is a constant.To make Equation [2] applicable, various research has been carried out.Such include Sobrino and Raissouni, (2000) that proposed modifying constant  as the linear combination of , ∆ and atmospheric water vapour content.François and Ottlé, (1996) presented different coefficients for different  values.Sun and Pinker, (2003) addressed the SW coefficients according to different surface types to account for LSE effect.However, it is observed that the constant  is a combination of the Mean: and difference; from the LSE of the TIRs channels, thereby keeping the other coefficients independent of LSE.It should be noted that the land surface is complex and that the LSE may be quite different from unity and depends on the channel as stated in (Coll and Caselles, 1997).Considering the effect of the emissivity, the quadratic SW algorithm may not work well.As François and Ottlé, (1996) noted in their work, when ε is greater than 0.95, good accuracy can be obtained using the quadratic method with the emissivity-dependent coefficients.That means, when ε is low, the quadratic relationship no longer performs well or is changed.Therefore, the quadratic method has to be re-examined closely, especially for low emissivity.
An efficient and convenient method for LST retrieval must be explored and investigated as it is influenced by numerous elements such as surface layout, land cover, and aerosol.Qin and Karnieli, (1999) which shows the detailed processes.
This result is further verified from Jin et al., (2015) which shows that SWA coefficients for Sea Surface Temperature estimation can be utilized for LST retrieval if there is an estimation of the ε effect.A modified C is therefore combined with the results of A and B in developing the LST algorithm.
A major factor to be considered is, however, the ambient temperature of the study area through the Atmospheric Transmittance (AT).This is important as the saturated water vapour content affects the satellite imagery in the form of cloud cover.Hence, an estimation of the AT values becomes a necessity.To derive the AT, two main calculations are done namely Water Vapour Content and Mean Atmospheric Temperature.Relative Humidity (RH), as well as the temperature of upper-air estimation, are two fundamental estimations utilized as a part of the numerical climate model analysis for operational climatic forecasts.Equations [3] -[6] as described in (Jarraud, 2008, McRae, 1980, Wagner and Pruß, 2002) were applied in showing the relationship between saturation water vapour pressure and water vapour at a specific temperature.Where   (, ) is the saturation vapour pressure of moist air; () is the pressure function;   () is the saturation vapour pressure in pure phase;  is the ambient temperature;  is the atmospheric pressure;  is the water vapour content and   is the standard atmospheric pressure of 1013.25ℎ .

Application of algorithm
For the application and validation of the SWA, Landsat-8 data were downloaded from EarthExplorer as well as observational data (atmospheric pressure, relative humidity, and temperature) were obtained from the South African Weather Station (SAWS) for the corresponding dates for which the satellite imagery was acquired.Data from the satellite imagery were radiometrically normalized, converted to surface reflectance, and quantitatively analysed to enable the establishment of generic and consistent classification rules after which the main processing began accurate LST retrieval through the SWA.
As shown in Figure [1], the following processing was then done in deriving the LST:

TOA spectral transmittance
The purpose of this was to calibrate the produced noise by the sensors which measure reflectance from the Earth surface in the form of Digital Numbers (DN) representing each pixel.The calibration as used by (Vicente-Serrano et al., 2008) was implemented for the TIR bands as shown in Equation [7]: Where   and   are spectral radiances,   is the maximum DN and   is the TOA.

Brightness Temperature
After the conversion of DNs to reflection, the TIR data were then converted to Brightness Temperatures (  ) which "is the microwave radiance traveling upward from the top of Earth's atmosphere" Yang et al., (2014) using the thermal constants provided in the metadata file.In order to achieve accuracy for the   conversion, equation [8] was implemented through the ATCOR module of ERDAS Imagine.The result is converted to Celsius from Kelvin by adding the absolute zero.
Where  1 and  2 are band-specific constants

Normalized Difference Vegetation Index
The NDVI is an indicator used to analyse the biomass or greenness of the observed area.As Weng et al., (2014) stated, estimating NDVI is essential since the amount of vegetation present is a factor for LST retrieval.Equation [9] was used in deriving the NDVI.

Fractional Vegetation Cover
Spatial variation in the radiometric temperature of surfaces is related to variations of the soil-water concentration vertically and surface greenness detailed by the FVC.The variations of detailed explanations of FVC utilized by remote sensing and modeling, be that as it may, can prompt an error between what is utilized as part of a model and what is estimated utilizing remote sensing.Hence, the formula as described in Wang et al., (2015) was used in deriving the FVC from the NDVI image as shown in Equation [10].

Land Surface Emissivity
The LSE must be known in order to estimate LST since it is the relative capacity of a surface to discharge radiative energy and characterized as the proportion of black-body radiative energy (Planck's law) to land surface radiative energy of similar temperatures (Jimenez-Munoz et al., 2014).
Nonetheless, it stays a challenge to absorb satellite radiance within the TIR regions (10-12m) because of the absence of infrared LSE.Fortunately, several methods in remote sensing have been proposed for LSE retrieval.According to (Dutta, 2015), LSE can be calculated using: Where,   and   are the Soil and Vegetation Emissivity constants for TIR bands respectively as shown in Table [1].

Atmospheric Transmittance
Atmospheric transmittance which is characterized as atmospheric capacity in transmitting electromagnetic energy is reliant on air masses by radiation, water vapour, and air dust.It fluctuates for radiation of various wavelengths; smaller light scattering and atmospheric absorption, greater air transmittance.The value of W as shown in Equation ( 6) is incorporated for the AT when using the SW Algorithm needed in the LST retrieval process.

Land Surface Temperature
SWA utilizes both TIR bands normally situated in the 10-12m atmospheric window.According to (Jiménez-Muñoz, 2008, Kamran et al., 2015), "the basis of the SW algorithm is that the radiance attenuation for atmospheric absorption is proportional to the radiance difference of simultaneous measurements at two different wavelengths, each subject to different amounts of atmospheric absorption".Therefore, processes 4.1 -4.7 is combined with the LST retrieval process as described in Li et al., (2013) and Rozenstein et al., (2014)

RESULTS AND DISCUSSION
As described in the processing, the results entailed spatial distribution maps showing the NDVI, FVC and LSE values.LSE was created using the NDVI threshold technique where the NDVI was reclassified into soil and vegetation.Highly elevated regions had more vegetative cover; hence, LSE was high in these regions.temperatures between 13-23℃ while (e) and (f) have relatively low temperatures due to being along the coastline with low inhabitants.The long variation in the temperature values is attributed to taking the mean temperatures over the time stamp covering all seasons.Also, the results were slightly affected due to errors as indicated in the analysis section for the time series was acquired from the satellite.Input parameters when estimating RH and W contributed to slight errors where some regions had either lower or higher temperatures.Figures [3] and [4] which were the validation sites showed higher LST in urban areas when compared to the non-urban areas.
In order to reveal the diversities of LST in different types of Landcover, the temperature characteristics of four Landcover types were analysed.The minimum, maximum and mean temperatures of the land cover patterns were derived by averaging all corresponding pixel values.The preceding results reflect that the spatial variability of LST is basically consistent with Landcover types of the study sites.In summary, the results had a wide range due to the range of temperatures being considered over various seasons.

Analysis
Sensitivity, statistical regression and error analysis on each study site was then carried out to further determine the validity of this SWA and determine error effects from estimations.

Sensitivity Analysis
The sensitivity analysis provides possible error effects of the estimation parameters used for LST retrieval.LSE, at-sensor   and AT were required for the application of this SWA and errors on these parameters subsequently led to LST errors.The analysis was performed using the python scripts as established by Isaya Ndossi and Avdan, (2016) in which various contributions to final LST error is estimated from numerous derivatives according to the classical error theory (Jin et al., 2015, Rozenstein et al., 2014).LST estimation errors are independent of change in temperatures as it changes over 0-60℃ range by less than ± 0.03°C as  10 - 11 = -2.2℃,Table [2] and 0  There is an error increase in LST estimations when the AWVC decreases which in turn causes an increase in the AT.The impact increases as   difference between the increasing TIRS bands.

Linear regression analysis
Statistical regression analysis was used to establish a relationship between the estimated LST values and the measured air temperature values from the SAWS and was achieved by using mean temperature values over the study period (June 2014 to May 2015).R 2 was achieved by utilizing Microsoft Excel and regression model formula of  ≈ (, ); [15] where  is the dependent variable (measured LST) and  is the independent variable (estimated LST).The analysis resulted in an R 2 of 0.952 (Durban), 0.983 (East London) and 0.979 (Polokwane) as shown in Figure [7].

Error analysis
Error from the split window algorithm used for LST retrieval is hard to measure.However, such errors are minimized.For the error with calibration, it showed that the results are within the pre-determined ±0.6K obtained after July 2014 to (http://landsat.usgs.gov/l8handbook_appendixa.php).Error from LSE can originate from the error of estimating NDVI and approximations.However, soil emissivity is highly relative within the study area which is well vegetated for a ~0.979 emissivity applicability.This research had below 0.005 emissivity error which prompted a 0.16℃ maximum LST error (when interpolated with SAWS data).Thus, the overall estimate is 1.04K -1.27K RMSE in the LST estimation for the three study site as shown in Table [3].
Table [3] was then repeated for the other study sites -East London and Polokwane -with RMSE values of 1.27 and 1.04 respectively.

Conclusion
The applied SWA is a dynamic mathematical tool that provides LST information using the   of TIRS bands, LSE derived from FVC of the OLI bands in accurately estimating LST with the addition of Atmospheric Transmittance.The applicability of this SWA from the results of validation and RMSE for the retrieved LSTs is confirmed when validated against observed SAWS data from the three study sites.On the basis of the validation and the sensitivity analysis, it can be concluded that the SWA presented in this paper is able to provide an alternative to ground-truth observations for accurate LST retrieval from the Landsat-8 dataset.

Acknowledgement
Appreciation goes to the University of KwaZulu-Natal for providing an enabling environment for this research.Thanks also go to SAWS and eThekwini Municipality for providing raw data needed for processing and validation.
(, ) = () *  () Figure [2] has been derived using LSE,   and Emissivity difference between band 10 and 11.An analysis of Figure [2] shows region (a) having highest temperatures of 22-28℃, (b) having mixed temperature levels due to the dam presence and low inhabitants, (c) is well populated with average temperatures of 18-23℃, (d) which is the busy Durban central business district having varying

Table 1 :
(Jimenez-Munoz et al., 2014)-Munoz et al., 2014)Over a short distance, there can be a substantial LSE change.The possibility of utilizing the NDVI values for the retrieval of the vegetation and soil proportion; or the minimum and maximum NDVI value; so as to be able to estimate LSE and in turn the mean and difference LSE as required in the SWA application.The mean and difference LSE shows the TIR bands relationship with the applied formula of: to give Equation [14].

Table 2 :
Estimated LSE from NDVI and FVC