Surface Temperature of the Planet Earth from Satellite Data

: The Intergovernmental Panel on Climate Change regular scientiﬁc assessments of global warming is based on measurements of air temperature from weather stations, buoys or ships. More speciﬁcally, air temperature annual means are estimated from their integration into climate models, with some areas (Africa, Antarctica, seas) being clearly underrepresented. Present satellites allow estimation of surface temperature for a full coverage of our planet with a sub-daily revisit frequency and kilometric resolution. In this work, a simple methodology is developed that allows estimating the surface temperature of Planet Earth with MODIS Terra and Aqua land and sea surface temperature products, as if the whole planet was reduced to a single pixel. The results, through a completely independent methodology, corroborate the temperature anomalies retrieved from climate models and show a linear warming trend of 0.018 ± 0.007 ◦ C / yr.


Introduction
The Intergovernmental Panel on Climate Change (IPCC) publishes Assessment Reports periodically, with the aim of analyzing the current state of the climate. The IPCC's Fifth Assessment Report [1] was the last of these reports and provided the scientific input into the Paris Agreement, which aims to strengthen the global response to the threat of climate change by holding the increase in the global average temperature to well below 2 • C above pre-industrial levels and to pursue efforts to limit the temperature increase to 1.5 • C above pre-industrial levels. The reference global temperatures datasets are derived from air temperature in situ data, with various sources, such as ships, drifting buoys or Argo floats [2], in the case of sea and ocean extensions, and permanent ground-based meteorological stations, which typically measure air temperature at 1.5-2 m, [3] on land. These in situ stations do not cover the entire Earth surface, meaning a lack of representative measures from certain areas of the planet [4].
On the other hand, satellite data provide now temporal series of an adequate span for climate analysis independent of the in situ data. After the launch of the Terra platform, on the 18th of December of 1999, and then the Aqua satellite, on the 4th of May of 2002, four daily measures became available from the MODIS instrument for each 1 km by 1 km pixel (at 10:30, 13:30, 22:30, 01:30 nominal time) [5]. Moreover, the international scientific community has developed surface temperature algorithms since the 80s, which are now capable of estimating the surface Earth temperature from satellite-derived data. Studies show global surface temperature (including land and ocean surfaces) estimations with an associated uncertainty between 1 and 2 K [6].
The aim of this study was to determine if remote sensing is able to retrieve reliably our planet surface temperature changes. We used the MODIS Terra and Aqua Sea Surface Temperature (SST) and Land Surface Temperature (LST) products to estimate the Earth Surface Temperature (EST-the surface temperature of Planet Earth including LST and SST observations), which we then compared with the global annual air temperatures estimated by the NOAA's National Climatic Data Center (NOAA NCDC) [7]. Finally, the different behavior and trend of the surface temperature were analyzed regionally and locally.

Materials and Methods
The MODIS Terra and Aqua SST MOD25 and MYD25 products [8] and the LST MOD11C1 and MYD11C1 products [9] have been downloaded for this study. SST products derived from MODIS are developed by the NASA Ocean Biology Processing Group (OBPG) and the MODIS Ocean Team Computing Facility (MOTCF). Calibrating coefficients, which vary in time, are continuously estimated and verified by the Rosenstiel School of Marine and Atmospheric Science (RSMAS) at the University of Miami. The data used in this study (product MOD25 Version 2.0) have been processed by applying the long wave algorithm [10], which considers MODIS 31 and 32 bands.
MODIS SST spatial resolution is 4.63 km, resulting in a global image of 8640 × 4320 pixels. Pixels have been labelled according to their quality. In this study, only SST pixels with an associated QC (Quality Control) value of 0 and 1, meaning good (clear data with associated scan angles less than 55 • ) and acceptable (good/acceptable data in glint or high scan angle) quality [10], respectively, have been included in estimations. Those pixels in which the atmospheric correction has failed, where cloud or ice have been detected or have not been processed for other reasons, were not included in our computations, and neither was its associated area proportion.
LST products MOD11C1 (Terra) and MYD11C1 (Aqua) at 0.05 • resolutions have been used, resulting in a 7200 × 3600 pixel image size. As with the MOD25 product, LST data has been processed using MODIS 31 and 32 bands as input for a long wave split window algorithm. The average LST error product is 1 • C for temperatures in the range from −10 to 50 • C [11]. LST pixels have also been filtered according to quality control tags explained on the product ATBD (Algorithm Theoretical Basis Document), meaning that only the pixels with QC values of 00 and 01 have been considered, ignoring those not produced due to cloud effects or other reasons [12].The accuracy of the MOD11 version 6 LST product has been tested and was found to mostly fall within 1 K [12][13][14], while MOD25 SST product is provided with an error lower than 0.45 K in most cases [10,15].
The NOAA NCDC uses the Merged Land-Ocean Surface Temperature Analysis (MLOST) in its operational and climate assessment activities, integrated by a land and sea air temperature data. The land dataset contains historical records for over 7000 locations worldwide, requiring adjustments in station locations, temperature instrumentation and land use. The source of sea surface temperature of MLOST is the International Comprehensive Ocean-Atmosphere Data Set (ICOADS), consistent on integrated marine meteorological observations, primarily, from ships and buoys. Both LST and SST grids are used to construct a coherent picture of historical temperature variations across the globe since the late nineteenth century [16].
SST and LST global means have been computed from MODIS Terra at 1:30, 13:30 and Aqua at 10:30 and 22:30, using 8-day composite images, according to SST ij and LST ij , which the values of SST and LST, respectively, for each pixel (i,j) on a given day of the year (t). Moreover, m is the column pixel dimension and n is, its line dimension; A ij , is the area of each pixel and A total is the sum of the area of the pixels considered in the calculation. Therefore, we obtained an LST and an SST value for each 8-day composite image for the whole span of the data (2003-2016).

SST MODIS = SST
where SST t mean and LST t mean are given by: In addition, surface temperature has been estimated for the northern (HN) and southern (HS) hemispheres, and for six latitudinal zones: the low latitudes from the Equator to the Tropic of Cancer in the HN and Tropic of Capricorn in the HS (0-23.5 • ), mid latitudes from the Tropics to the Arctic Circle in the HN and Antarctic Circle in the HS (23.5 • -66.5 • ) and high latitudes from the Arctic and Antarctic Circles to the Poles (66.5 • -90 • ).
Once SST and LST have been estimated, EST has been calculated by weighting both parameters according to the percentage of area covered by sea and land respectively and taking into consideration the different LST and SST image resolution: No gap-filling was carried out for cloud contaminated data. Manual inspection of retrieved SST and LST 8-day averages showed the presence of outliers, which we removed by a classical Z-score approach (all estimates further away than five standard deviations of the mean, calculated from the whole time series-EST, LST and SST-were not considered). We then estimated the annual temperature averages for LST and SST for the period 2003-2016. Trends were retrieved by linear regression. Estimates of Sen slope results led to similar results. Trend significance was estimated with the Mann-Kendall trend test.

MODIS Earth Surface Temperature Matches NOAA NCDC Global Air Temperature Estimations
The comparison between the surface temperature values of the Planet Earth obtained from MODIS products (EST) with the methodology developed in this work and the NOAA-NCDC data is shown in Figure 1. We obtain a correlation coefficient of 0.93, with linear trends of 0.018 ±0.007 • C/yr for MODIS EST and 0.021 • C/yr for NOAA NCDC anomalies. The NOAA NCDC anomalies have been calculated by using a reference temperature of 13.9 • C, which is the average combined land and ocean air temperature of the 1901-2000 period [7].The given correlation is statically significant at the 99.9% confidence level, also been estimated by the Sen's slope method, obtaining 0.018 • C/yr. and a confidence level of 97.5 while the range of differences between the MODIS derived EST and the NOAA NCDC air temperature anomalies is of −0.07 to 0.1 • C and the mean difference for the total period is 0.0043 • C. The median evaluated uncertainty in the period 2003-2016 is 0.10 • C for the SST, 0.20 • C for LST and 0.11 • C for EST. Moreover, trends have also been estimated by the Sen's slope method, obtaining 0.018 ± 0.007 • C/yr. and a confidence level of 97.5%, strengthening our results.
The high correlation found between the estimations of EST obtained from MODIS products in the present work and the in situ NOAA NCDC database emphasizes the high potential of thermal infrared satellite data as a source of accurate input data of great value for global surface temperature estimates and, therefore, for being applied to climatological studies [17].

Air Temperature Versus Satellite Temperature
We compare two temperature databases, one derived from MODIS Terra and Aqua observations and the other, interpolated from the NOAA NCDC in situ measures. Several differences exist between both types of data sets. First, MODIS satellite observations refer to skin SST while NOAA NCDC measures the air temperature. This fact has been documented to be responsible for a global bias of 1 °C [18].
Cloud coverage affects the information received by thermal sensors, preventing the measure of those pixels situated below clouds and introducing additional biases in respect to air temperature measures. Cloudiness is the main disadvantage that appears when measuring surface temperature with remote sensing sensor: some regions on the planet, such as polar zones are covered by clouds the most part of the year. On the other hand, MODIS allows a continuous, homogenous and global coverage, ensuring four daily representative observations, as previously specified [19]. Despite the more comprehensive data period enabled by AVHRR LST and SST datasets, it has been discarded due to the variability induced by the satellite's orbital drift [20][21][22].
The wide-scale and periodic observations provided by MODIS are a clear point in favor of satellite estimated surface temperature versus classic air temperature datasets. The second one resides in the global measurement, contrasting with the integration of sparse in situ air temperature measures, not representative on a global scale [19] and which measurements accuracy depend on the sensor characteristics, calibration and exposure. Moreover, the daily air temperature averages only take into consideration the maximum and minimum observations per day. This last consideration leads to errors in the 2-3 °C range, while the level of accurate precision should be 1-2 °C [23].
The aim of this study is not to compare absolute skin and air temperature values, as the existence of biases between both data sets is evident: skin temperature has been documented to be 3 °C lower than air temperature at nighttime [24]. This study neither aims to show local results as several previous studies have done [25][26][27][28][29]. The study shows the potential of satellite thermal observations, represented by MODIS-Terra and MODIS-Aqua, to estimate global surface temperature trends with similar accuracy that the traditional air temperature databases and, therefore, its applicability to climate change studies.

Regional Analysis: Northern Hemisphere Land Surface Contributes Most to Temperature Increases
The annual values of LST, SST and EST are presented in Figure 2 for the period 2003-2016 by hemisphere and globally. Stability in surface temperatures prevails until 2011, when temperature

Air Temperature Versus Satellite Temperature
We compare two temperature databases, one derived from MODIS Terra and Aqua observations and the other, interpolated from the NOAA NCDC in situ measures. Several differences exist between both types of data sets. First, MODIS satellite observations refer to skin SST while NOAA NCDC measures the air temperature. This fact has been documented to be responsible for a global bias of 1 • C [18].
Cloud coverage affects the information received by thermal sensors, preventing the measure of those pixels situated below clouds and introducing additional biases in respect to air temperature measures. Cloudiness is the main disadvantage that appears when measuring surface temperature with remote sensing sensor: some regions on the planet, such as polar zones are covered by clouds the most part of the year. On the other hand, MODIS allows a continuous, homogenous and global coverage, ensuring four daily representative observations, as previously specified [19]. Despite the more comprehensive data period enabled by AVHRR LST and SST datasets, it has been discarded due to the variability induced by the satellite's orbital drift [20][21][22].
The wide-scale and periodic observations provided by MODIS are a clear point in favor of satellite estimated surface temperature versus classic air temperature datasets. The second one resides in the global measurement, contrasting with the integration of sparse in situ air temperature measures, not representative on a global scale [19] and which measurements accuracy depend on the sensor characteristics, calibration and exposure. Moreover, the daily air temperature averages only take into consideration the maximum and minimum observations per day. This last consideration leads to errors in the 2-3 • C range, while the level of accurate precision should be 1-2 • C [23].
The aim of this study is not to compare absolute skin and air temperature values, as the existence of biases between both data sets is evident: skin temperature has been documented to be 3 • C lower than air temperature at nighttime [24]. This study neither aims to show local results as several previous studies have done [25][26][27][28][29]. The study shows the potential of satellite thermal observations, represented by MODIS-Terra and MODIS-Aqua, to estimate global surface temperature trends with similar accuracy that the traditional air temperature databases and, therefore, its applicability to climate change studies.

Regional Analysis: Northern Hemisphere Land Surface Contributes Most to Temperature Increases
The annual values of LST, SST and EST are presented in Figure 2

High Latitudes Show the Highest Land Surface Temperature Increases
In Figure 3, surface temperature profiles are shown in more detail for the defined six latitudinal zones (high, mid and low). High latitudinal zones, located between 66.5 and 90 degrees and presenting extreme temperatures, only represent 9% of the total surface taken into account, thus they scarcely contribute to global EST. Low latitudes sea surface zones represent 29% of the total surface analyzed in computations, contributing nearly 1/3 of the proportion to the global result. NH low latitudes show an absolute increase of 0.50 °C from 2011 to 2016 and SH low latitudes an increase of 0.78 °C. Several studies relate positive SST anomalies in the eastern Tropical Pacific Ocean with El Niño Southern Oscillation (ENSO) 2015 event [30,31], presenting an associated ONI above 0.5°. Moreover, these SST anomalies are reported to affect the North Atlantic Ocean and the Indian Ocean [32,33], which spread through mid-latitudes. According to our results, sea surface in these midlatitude zones represents 36% of the data analyzed, contributing more than 1/3 of the proportion of the global results. SST increase during 2011-2016 is quantified as 0.52 °C in the NH mid latitudes, although it is stable in the SH during this period.
Curves associated with SST present a smoother behavior and less variability than LST curves, the sea temperature being less influenced by the summer and winter alternations due to water thermal capacity. The percentage of sea surface is significantly higher than land in the NH low and mid latitudes and the SH low latitudes, thus the final EST curve resembles more the SST behavior than the LST one.
The global EST trend during the period 2003-2016 is 0.018 °C/yr, with high local variability as can be observed in Table 1. By parameter, global trends are of 0.030 °C/yr. for LST and 0.013 °C/yr. for SST. The warming LST trend is not only higher than the SST trend globally, but also for the other

High Latitudes Show the Highest Land Surface Temperature Increases
In Figure 3, surface temperature profiles are shown in more detail for the defined six latitudinal zones (high, mid and low). High latitudinal zones, located between 66.5 and 90 degrees and presenting extreme temperatures, only represent 9% of the total surface taken into account, thus they scarcely contribute to global EST. Low latitudes sea surface zones represent 29% of the total surface analyzed in computations, contributing nearly 1/3 of the proportion to the global result. NH low latitudes show an absolute increase of 0.50 • C from 2011 to 2016 and SH low latitudes an increase of 0.78 • C. Several studies relate positive SST anomalies in the eastern Tropical Pacific Ocean with El Niño Southern Oscillation (ENSO) 2015 event [30,31], presenting an associated ONI above 0.5 • . Moreover, these SST anomalies are reported to affect the North Atlantic Ocean and the Indian Ocean [32,33], which spread through mid-latitudes. According to our results, sea surface in these mid-latitude zones represents 36% of the data analyzed, contributing more than 1/3 of the proportion of the global results. SST increase during 2011-2016 is quantified as 0.52 • C in the NH mid latitudes, although it is stable in the SH during this period.
Curves associated with SST present a smoother behavior and less variability than LST curves, the sea temperature being less influenced by the summer and winter alternations due to water thermal capacity. The percentage of sea surface is significantly higher than land in the NH low and mid latitudes and the SH low latitudes, thus the final EST curve resembles more the SST behavior than the LST one.
The global EST trend during the period 2003-2016 is 0.018 • C/yr, with high local variability as can be observed in Table 1. By parameter, global trends are of 0.030 • C/yr. for LST and 0.013 • C/yr. for SST. The warming LST trend is not only higher than the SST trend globally, but also for the other nine areas defined. The Copernicus Ocean State Report [34] exhibited an SST warming at a rate of 0.016 • C/yr for the period 1993-2015, meaning a total average increase of 0.4 • C over the considered period, by using satellite-derived ESA Climate Change Initiative (ESA-CCI) observational products Hausfather et al. [     Linear trends have been estimated by individual parameters, SST and LST, and for our estimated EST at global, hemispherical and latitudinal scales. In the case of SSTs, negative trends are mainly observed in the North Atlantic Ocean, South of Greenland, in the North Pacific, the Sea of Japan, the Yellow Sea and in the Southern Ocean. On the other hand, positive trends are irregularly distributed along NH low and mid-latitudes and in the Atlantic and Pacific Oceans, next to the North American and Australian Continent and the Arctic [37]. In the SH, comprehensive stable areas are observed, with positive trends in the South Pacific, Indian and Atlantic Oceans, in agreement with the results shown in the regional SST analysis of the Fifth IPCC Assessment Report, which establishes an SST increase during 1950-2009 of 0.027 • C/decade, 0.10 • C/decade and 0.063 • C/decade, respectively. For LST, positive trends are found for vast areas of Central and Eastern Europe, Scandinavia and Siberia. In the American Continent, warming trends appear in large areas of Alaska and Canada, between the peninsulas of California and Florida, and in northeastern Brazil and Patagonia. Negative trends appear mainly in the Indian Peninsula and in eastern lengths of Antarctica ( Figure 4). Retamales et al. [38] found a negative LST trend in eastern Antarctica of about −1 K/decade, using the MODIS LST product between years 2000-2018. They also estimated the trend with the ERA datasets, obtaining a warming trend of 0.7 K/decade, in disagreement with the results offered by the MODIS LST estimation. The differences between both results show the influence of cloudiness in satellite observations. In particular, cloud coverage increases during winter in this area of the planet, a fact which can impact the calculated temperature trends. The global EST trend has also been estimated by the Theil-Sen method, showing a warming rate of 0.017 • C/yr. with a 97.5% Mann-Kendal significance (Figure 4). Results show the same warming/cooling patterns as estimated by the linear method.
Remote Sens. 2020, 12, x FOR PEER REVIEW 7 of 10 Linear trends have been estimated by individual parameters, SST and LST, and for our estimated EST at global, hemispherical and latitudinal scales. In the case of SSTs, negative trends are mainly observed in the North Atlantic Ocean, South of Greenland, in the North Pacific, the Sea of Japan, the Yellow Sea and in the Southern Ocean. On the other hand, positive trends are irregularly distributed along NH low and mid-latitudes and in the Atlantic and Pacific Oceans, next to the North American and Australian Continent and the Arctic [37]. In the SH, comprehensive stable areas are observed, with positive trends in the South Pacific, Indian and Atlantic Oceans, in agreement with the results shown in the regional SST analysis of the Fifth IPCC Assessment Report, which establishes an SST increase during 1950-2009 of 0.027 °C/decade, 0.10 °C/decade and 0.063 °C/decade, respectively. For LST, positive trends are found for vast areas of Central and Eastern Europe, Scandinavia and Siberia. In the American Continent, warming trends appear in large areas of Alaska and Canada, between the peninsulas of California and Florida, and in northeastern Brazil and Patagonia. Negative trends appear mainly in the Indian Peninsula and in eastern lengths of Antarctica (Figure 4). Retamales et al. [38] found a negative LST trend in eastern Antarctica of about −1 K/decade, using the MODIS LST product between years 2000-2018. They also estimated the trend with the ERA datasets, obtaining a warming trend of 0.7 K/decade, in disagreement with the results offered by the MODIS LST estimation. The differences between both results show the influence of cloudiness in satellite observations. In particular, cloud coverage increases during winter in this area of the planet, a fact which can impact the calculated temperature trends. The global EST trend has also been estimated by the Theil-Sen method, showing a warming rate of 0.017 °C/yr. with a 97.5% Mann-Kendal significance ( Figure 4). Results show the same warming/cooling patterns as estimated by the linear method.

Linear method Sen's slope method
Mann-Kendall test

Conclusions
The comparison between the planetary MODIS Earth Surface Temperature (EST) presented in this work and the NOAA-NCDC air temperature data shows a correlation coefficient of 0.93 between

Conclusions
The comparison between the planetary MODIS Earth Surface Temperature (EST) presented in this work and the NOAA-NCDC air temperature data shows a correlation coefficient of 0.93 between the two databases, demonstrating the high potential of thermal infrared satellite data to provide accurate data for climatic and meteorological studies. Satellite data are essential in the monitoring of SST, LST and EST (always considering the presence of cloud coverage) as a solid source of continuous data in space and time. In this way, the randomness and lack of data in certain regions of the planet, when referring to in situ data, is widely solved.
The global EST trend is 0.018 • C/yr. The highest EST trend has been computed for the North Hemisphere high latitudes, 0.032 • C/yr., as well as the highest LST trend, 0.080 • C/yr. On the other hand, the highest SST trend is associated to the South Hemisphere low latitudes: 0.022 • C/yr. Globally, SST trends (0.013 • C/yr.) are lower than LST trends (0.030 • C/yr.). SST shows less variability in time than LST and a higher influence on EST, especially in southern latitudes, in which the sea surface extension dominates over the land surface extension. In summary, this paper pretends to open a reflection in the scientific community, covering not only the trends observed from satellites but also the current estimations of air temperature, as well as the validation of both in situ and satellite data. We propose to use the satellite data to estimate anomalies of the average temperature of the Earth's surface: the thermal data from the satellites could be a proxy to be included in the IPCC assessment reports. Further work should be done in the future, as the extension of MODIS surface temperature temporal series is enlarged with time.