CAN RECONSTRUCTED LAND SURFACE TEMPERATURE DATA FROM SPACE PREDICT A WEST NILE VIRUS OUTBREAK ?

Temperature is one of the main drivers of ecological processes. The availability of temporally and spatially continuous temperature time series is crucial in different research and application fields, such as epidemiology and control of zoonotic diseases. In 2010, several West Nile virus (WNV) outbreaks in humans were observed in Europe, with the largest number of cases recorded in Greece. Human cases continued to occur for four more years. The occurrence of the 2010’s outbreak in Greece has been related to positive anomalies in temperature. Currently available remote sensing time series might provide the temporal and spatial coverage needed to assess this kind of hypothesis. However, the main problem with remotely sensed temperature are the gaps caused by cloud cover. With the objective of testing the former hypothesis, we reconstructed daily MODIS Land Surface Temperature (LST) data and derived several indices that are known or hypothesized to be related to mosquito populations, WNV transmission or risk of disease since they might constitute proxies for favoring or limiting conditions. We present the first results of the comparisons of time series of LST-derived indices among locations with WNV human cases and municipalities with and without reported WNV infection in Greece between 2010 and 2014.


INTRODUCTION
West Nile virus (WNV) is one of the mosquito-borne flavivirus most widely distributed in the world.The virus causes a variety of symptoms to humans: from an asymptomatic infection to severe and even fatal encephalitis (Petersen et al., 2013).Wild birds are the major reservoirs of the virus and the main transmission route is through Culex mosquito-vectors.Most patients with WNV infections do not present any symptoms or only display mild flu-like symptoms, but in some cases people develop severe illness involving the central nervous system (Marka et al., 2013).Human cases have been reported from several countries since the 1960s.However, the frequency of reported outbreaks with severe symptoms has increased over the last 15-20 years (Calistri et al., 2010), especially in Europe and the United States.
In 2010, several WNV outbreaks in humans were observed all over Europe, however the largest number of cases was recorded in Greece (Danis et al., 2011).After that, human cases were reported for four more years, during a transmission cycle that goes from June to September (Hadjichristodoulou et al., 2015).The first cases in 2010 were reported in northern Greece (central Macedonia) near the city of Thessaloniki.In subsequent years, the disease further spread both southwards and eastwards, and cases were reported even in the highly-populated Greek capital city of Athens.More than 600 laboratory confirmed cases and 73 deaths were reported between the years 2010 and 2014, making this outbreak the largest of WNV disease in Europe (Pervanidou et al., 2014).
Several factors may be involved in the determination of a WNV outbreak onset.Environmental conditions affecting both avian reservoir hosts and the mosquito vector populations, as well as human behavior and habitat modifications, may regulate WNV amplification and subsequently trigger an outbreak (Gibbs et al., 2006).The environmental conditions associated with the onset of WNV human cases in Europe, despite being investigated by multiple studies, are partly unknown (Paz and Semenza, 2013, Tran et al., 2014, Marcantonio et al., 2015, Semenza et al., 2016).The deficiency of complete and homogeneous continental-wide datasets is one of the underlying reasons for this lack of knowledge.On the contrary, datasets collected at a more local scale, even though representative of only a fraction of the European continent, may provide more robust results.These results may, in certain instance, be transferred to the continental scale.Some of the variables that have been found to be correlated with WNV epidemiology and might be considered risk factors are: positive temperature anomalies in July, high summer temperatures, summer drought, positive anomalies in water index in June, number of districts with WNV infections the previous year, presence of wetlands, irrigated croplands and fragmented forest, type of bird migration route crossing the district under study, among others (Gibbs et al., 2006, Pradier et al., 2008, Paz and Semenza, 2013, Tran et al., 2014, Marcantonio et al., 2015).These variables, though related to WNV human cases, might be acting indirectly through their effects over birds or mosquito populations, and over virus amplification.In the case of mosquito-borne pathogens such as WNV, higher than usual temperatures are known to influence vector competence, to accelerate virus replication within mosquitoes (extrinsic incubation period, EIP), to boost mosquitoes reproduction rates, and to prolong their breeding season (Gibbs et al., 2006, Paz andAlbersheim, 2008).As such, temperature is one of the main environmental factors addressed when studying vector-borne viruses carried by mosquitoes (Paz and Albersheim, 2008, Bisanzio et al., 2011, Ros et al., 2014, Chaskopoulou et al., 2016).
The occurrence of the 2010's emergence of WNV in Greece has been indeed related to positive anomalies in temperature that year (Danis et al., 2011).Even though some attempts have been made to test this hypothesis, they have been based on either sparse data from meteorological stations (Paz et al., 2013) or on data derived from coarse spatial and temporal temperature models that do not match the time frame of Greece WNV outbreaks (Valiakos et al., 2014).Currently available remote sensing time series, however, might provide the temporal and spatial coverage needed to better assess the previous hypothesis.
The MODIS LST sensor, on board of Terra and Aqua satellites, registers temperature four times a day with a spatial resolution of approximately 1 km.Recently, NASA made publicly available a new re-processed version 6 of their products, providing several improvements with regards to previous versions.However, the main problem with remotely sensed temperature is the occurrence of gaps, mostly caused by cloud cover.Hence, data must be reconstructed, i.e. gap filled in order to be useful (Neteler, 2010, Metz et al., 2014).In this sense, long standing FOSS4G such as GRASS GIS (Neteler et al., 2012, GRASS Development Team, 2017) provides at present more tools than ever to easily handle and process hundreds of thousands of maps and fill gaps in remote sensing time series (Gebbert and Pebesma, 2014, Gebbert and Pebesma, 2017, Metz and GRASS Development Team, 2017a, Metz and GRASS Development Team, 2017b, Metz and GRASS Development Team, 2017d).
Given that temperature is one of the main drivers of ecological processes related to vector-borne viruses carried by mosquitoes and that we have the data and tools available to get temporally and spatially continuous temperature time series, our question is then: Can MODIS LST (and derived LST indices) inform us about spatial and temporal trends or changes in temperature that might have contributed to the occurrence of WNV outbreaks in Greece?Here we present a first approach to the problem, however our ultimate goal is to study and model the relationship of different environmental variables and the incidence of human WNV infections in Greece.In this first paper, we present a yearly analysis of LST and derived related variables in sites with reported cases, with the intention of identifying temperature related differences in the periods before and after the first outbreak and municipalities with and without reported cases.

MAIN BODY
2.1 Materials and Methods 2.1.1Data used: Our study area is Greece which has a surface of 131.957 km 2 .For the reconstruction of satellite data however, we used a larger area which boundaries are: North: 42°1'12"N, East: 15°24'6"E, South: 34°41'47"N, West: 34°23'53"E.We used the daily LST Day and Night from MOD11A1 and MYD11A1 version 6 products (1 km of spatial resolution) from January 1 st , 2003 to December 31 st , 2016.The emissivity band in the MODIS LST products and the 30 arc-seconds Global Multi-resolution Terrain Elevation Data 2010 (GMTED2010) were used as covariates.Greek meteorological stations from the public database Global Summary of the Month (GSOM1 ) were used as validation set for the obtained reconstructed LST.
The epidemiological dataset was obtained by the European Centre for Disease Prevention and Control (ECDC).The data consisted of the geo-location (latitude and longitude coordinates) of all WNV fever and encephalitis human cases reported in Greece from July 2010 to September 2014.

LST data processing:
The raw data downloaded from MODIS LP DAAC site (4 images per day) were mosaicked and reprojected to Lambert azimuthal equal area (EU LAEA, EPSG code 3035) with nearest neighbor as the resampling method.Nearest neighbor sampling was chosen in order to preserve the bit pattern in the quality assessment (QA) band.The reprojected MODIS LST values were then filtered with the mandatory and LST QA band allowing for an LST error lower than 3 K.Our LST reconstruction method consisted of temporal interpolation to partially reconstruct Day and Night daily LST (fill short gaps in time) followed by spatial interpolation with covariates and outlier filtering to fill remaining gaps in space.This was done for each time series separately (MOD11A1 Day and Night, MYD11A1 Day and Night).Finally, we aggregated the four reconstructed time series into a daily time series in order to obtain average, minimum and maximum temperature estimates.To asses the quality of the reconstruction, we compared our reconstructed monthy aggregated time series with data from Greek meteorological stations from the GSOM database.

Data analysis:
We aggregated the daily reconstructed LST time series into months, seasons and years by means of arithmetic average.This yielded the average temperature, the average minimum temperature and the average maximum temperature for each period (i.e.: month, season, year).Seasons were defined as follows: winter (December, January, February), spring (March, April, May), summer (June, July, August) and autumn (September, October, November).We then estimated a series of biologically meaningful variables that are known or hypothesized to affect mosquito populations, WNV transmission or disease risk.The thresholds used for deriving LST-related variables were drawn from bibliography on both mosquitoes survival and development and virus amplification (Rueda et al., 1990, Vinogradova, 2000, Reisen et al., 2006, Tachiiri et al., 2006, Gong et al., 2011, Ciota et al., 2014).The following list of annual indices were derived from reconstructed daily LST: • Annual averages of minimum, maximum and average temperature.• Annual standardized anomalies for average, maximum and minimum temperature: estimated as the difference between the yearly value and the long-term average, divided by the long-term standard deviation of the corresponding variable.For example: (Tavi -Average(Tav)) / sd(Tav), where i is 2003, 2004... 2016 and Tav is the average temperature.
• Annual spring warming: estimated as the slope of the linear regression of daily average temperature from February 1 st to June 30 th .• Annual autumnal cooling: estimated as the slope of the linear regression of daily average temperature from October 1 st to December 31 st .• Annual number of days with minimum temperature ≤ 0°C.
• Annual absolute minimum temperature.
• Annual first day with minimum temperature ≤ 0°C.
• Annual number of summer days: estimated as the count of days with average temperature ≥ 25°C.following literature we chose a threshold of 10°C to identify annually the potential mosquito growing season and obtain the corresponding length, start and end dates (Metz and GRASS Development Team, 2017c).• Length, start and end dates of WNV transmision season: following literature we chose a threshold of 14.3°C to annually identify the potential WNV transmission season and obtain the corresponding length, start and end dates (Reisen et al., 2006, Metz andGRASS Development Team, 2017c).• Number of potential mosquito cycles per year and start and end date of each cycle: These indices were estimated by accumulating biologically effective degree days (BEDD) from the 4 maps per day reconstructed LST time series.Thresholds used for accumulation were 10°C and a cutoff of 34°C (Rueda et al., 1990, Vinogradova, 2000).After the accumulation was done, the next step was to identify the potential occurrence of mosquito cycles (from egg hatching to adult emergence), following Culex mosquitoes degree days (DD) requirements (Tachiiri et al., 2006), obtain length, start and end date.Finally, we aggregated the time series containing indicators of each cycle, to get the maximum number of potential cycles per pixel.• Number of potential EIP per year and start and end date of each cycle: These indices were estimated as described above, but using different base and cutoff temperatures: 14.3°C and 32°C, respectively (Reisen et al., 2006, Kilpatrick et al., 2008).The identification of potential EIPs was done following (Reisen et al., 2006), who propose a median EIP of 109 DD.
In the case of annual number of days with minimum temperature ≤ 0°C, the annual absolute minimum temperature and the first day of the year with minimum temperature below 0°C, the aggregation is performed considering a year starting in July and finishing in June of the following calendar year, since we were interested in catching the first colds in autumn that might affect the end of growing season in mosquitoes and their survival.
All the described variables were extracted for the geo-locations with reported WNV human cases, to obtain annual time series of each variable and compare the environment before and after the first WNV outbreak.We spatially aggregated the environmental LST-derived indices (average, median, standard deviation, minimum and maximum) by municipalities to compare conditions in areas with and without reported WNV infection.
2.1.4Software: MODIS data was dowloaded by means of pyModis library2 .All the remote sensing processing and analysis was done in GRASS GIS (GRASS Development Team, 2017).

Comparison of LST with GSOM:
We used a linear model to address the relation between monthly reconstructed average LST and and average monthly temperature from GSOM meteorological stations in the period 2003-2016.The global R 2 was 0.97 (Figure 1).Since the average temperature in GSOM data is estimated as: , GSOM TAVG should be lower in the summer months and higher in winter because in summer, days are longer than nights, and in winter, days are shorter than nights.Monthly average MODIS LST therefore, based on 4 values should thus be more accurate than GSOM based on 2 values.Besides, for summer months average reconstructed LST shows higher values than GSOM data and for winter months, the opposite, as expected (Figure 2).

WNV human cases:
In 2010, Greece underwent an epidemic of WNV infection, with 262 clinical human cases.The virus overwintered and spread across the country in the following years (2011)(2012)(2013)(2014) resulting in more than 600 confirmed human infections (Chaskopoulou et al., 2016).Figure 3 shows the the temporal distribution of the total number of cases per year and per month (years aggregated) between 2010 and 2014.

LST-derived indices:
Time series of LST-derived indexes in localities with reported WNV human infections showed (in most cases) high variability among localities when split per  year of onset of the infection.Especially in the cases of 2012 and 2013 infections, a bimodal distribution was observed, coinciding with the fact that during those years there were two main different areas in which the outbreaks took place (Central-East Macedonia and Attica).When aggregated by municipality, no clear pattern was observed either.For example, Figure 4 shows average LST in municipalities with more than 25 human cases.
Yearly LST indexes (average, minimum average, maximum average) did not show any notorious difference in localities with reported human cases, nor when spatially aggregated in municipalities.However, when considering the standardized anomalies, we observed that 2010 showed positive anomalies in the average minimum temperature during 2010, i.e.: minimum temperatures were higher that year (Figure 5).This trend was consistent and continuous from 2012 onwards and it is also evident when reviewing individual points, i.e: the locality with the highest number of cases (Figure 6).
The slope of spring warming was steeper in 2010, compared to that in 2009, being highest (steepest) in 2012, year of the second most important number of reported WNV human cases (Figure 7).Autumnal cooling in 2010 (slope of September, October and November average temperatures, after the main outbreak) was the steepest cooling observed, both in time and space (not shown).
The number of days with average temperature higher than 20 and lower than 30°C was higher in 2010 in Pella municipality (the municipality with the highest number of reported cases).This index also showed higher values for Pella municipality the year before the main outbreak.On the other hand, the number of tropical nights (days with minimum temperature higher than 20°C) appeared to be higher in 2010 in localities in which the first outbreak took place.This was also visible at municipality level both for 2010 and 2012.Moreover, the number of freezing days and the maximum number of consecutive freezing days in 2009-2010 (July 2009 -June 2010) appeared to be lower in the localities and municipalities where the main outbreak was observed later in 2010 (Figure 8).
Both the number of potential mosquito cycles and WNV EIPs, as well as the length of the mosquito growing season and the WNV transmission season were higher in 2012 and 2013.The differences in the length of seasons are visible in both the localities with reported human cases and when aggregated by municipalities.For the localities with reported cases in 2010, however, there seems to be three groups that depict different season lengths (Figure 9).

Discussion
This work constitutes a first approach to the study of WNV outbreaks in Greece as a function of remotely sensed environmental variables.Since temperature is one of the main drivers of insects life cycles and the WNV is a mosquito-borne pathogen, the occurrence of the disease (given that hosts, vectors and pathogens   coexist in the same area) should be related to temperature, too.
It is well-known that higher than usual temperatures might in-Figure 7. Spring warming slope for locations with WNV reported cases split by year of onset.
Remotely sensed temperatures come handy when resources are limited and meteorological stations scarce (most of the times).
Sensors such as MODIS on board of Terra and Aqua satellites provide temporally and spatially continuous data.However, remotely sensed LST is not free of problems, the main issue is the occurrence of gaps in the data due to cloud cover.Hence, LST needs to be gap-filled first to become fully useful (Neteler, 2010, Metz et al., 2014).A general evaluation of our reconstruction showed good global relationship between air temperature from meteorological stations in Greece and reconstructed MODIS LST.Since it was previously suggested that climatic conditions (higher temperatures, especially) might have favored the occurrence of the 2010 WNV outbreak in Greece (Danis et al., 2011) and Europe in general (Paz and Semenza, 2013), our main objective was to identify temperature-related variables that might show temporal differences in localities with reported human cases and spatio-temporal differences in areas with and without WNV infection.To this aim we estimated diverse indices that are known or hypothesized to represent limiting or favoring conditions for either mosquito or pathogen development and survival based on thresholds reported in the literature (Rueda et al., 1990, Vinogradova, 2000, Reisen et al., 2006, Tachiiri et al., 2006, Gong et al., 2011, Ciota et al., 2014).
Annual LST summaries, such as average, minimum and maximum did not show important spatio-temporal variations, however a trend towards higher temperatures was observed when considering the annual standardized anomalies in average and minimum LST.These differences, especially in minimum temperature anomalies, were particularly evident in the year 2010, when the first outbreak was recorded.These findings, though preliminar, The  In the same line, we observed that the slope of spring warming was steeper in 2010 compared to the previous years.Furthermore, indices that depict counts of days with conditions supposed to be limiting (number of consecutive freezing days) or favorable (number of tropical nights and number of days with temperature between 20 and 30°C) for mosquitos and WNV, also showed differences consistent with the former hypothesis.Moreover, the number of days with average temperature between 20 and 30°C showed higher values in Pella municipality the year before the main outbreak (i.e.: 2009), when bird surveys in hunting areas already detected WNV antibodies (Valiakos et al., 2014).The number of potential mosquito cycles and EIPs as well as their median length however, only showed marked differences for years 2012 and 2013.Meanwhile 2012 showed more than 150 human cases,  in 2013 the number of human cases was reduced to approximately one half of the previous year.This reduction in the number of cases, while conditions seem favorable might be highly likely related to human behaviour, control measures over mosquito populations (Papa et al., 2013, Chaskopoulou et al., 2016) and inmunization of the bird and human population (Papa et al., 2013, Hadjichristodoulou et al., 2015, Chaskopoulou et al., 2016).
The inspection and analysis of annual LST-derived indices for localities and municipalities (spatially aggregated) with reported The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLII-4/W2, 2017 FOSS4G-Europe 2017 -Academic Track, 18-22 July 2017, Marne La Vallée, France WNV human cases showed that a higher temporal resolution (i.e.: monthly or weekly) might be needed to better understand the relationship of outbreak onset and number of cases (or incidence) with temperature derived indices.Furthermore, if these indices are intended to be used to trigger early warnings and public health campaigns, higher temporal resolution is necessary to detect more specific and relevant differences.Notwithstanding this, reconstructed LST and the derived indices allowed us to identify likely meaningful differences that support the hypothesis of a favorable environment (even when spatially aggregated by municipalities with the average or the median).
Though some differences in LST-derived indices were observed both when considering localities and municipalities, some others are only (or better) observed in one or the other spatial scale.Not only that but there seems to be also an important level of variability in time series of localities and municipalities with reported cases (Figure 4).Furthermore, it might be possible as well that different variables or combination of variables might have triggered the onset of cases each year.The absence of reported WNV human cases in 2015 and 20163 , even though the environmental conditions seemed to be propicious, remain an open question.
On the other hand, some other variables relevant for mosquito populations such as humidity, amount of sun hours, availability of breeding sites and vegetation, might interact with temperature and should be included in models as well.For example, the area in which the first outbreak started in the region of Central Macedonia is characterized by many wetlands, rivers and lakes that serve as stopping areas for migratory birds during their migration from overwintering areas in Africa to breeding sites in northern Europe and viceversa (Chaskopoulou et al., 2016).The relative importance of different variables has not yet been evaluated.
The lack of differences in some indices at municipality level might be due to the loss of detail caused by the spatial aggregation of already rather coarse data, a phenomenon known as modifiable area unit problem (MAUP) (Atkinson and Graham, 2006).However, since the exact infection sites are unknown and only localities are reported, a resolution of 1 km might compensate for the mismatch between mosquito habitats and uncertainty in the precise identification of the infection sites.On the other hand, it might be as well that the thresholds used to derive some of the variables were not the optimal ones, since most of them (if not all) were obtained under laboratory conditions.This is an aspect that requires further investigation and especially, field surveys or experiments.

CONCLUSIONS
This study highlights the usefulness of remote sensing, and especially LST, for the characterization and quantification of temperature changes.The major advantage of MODIS LST is that it captures climatic conditions with short revisit time and global coverage.This might reduce the need of meteorological field data, especially when resources are limited.Satellite tracking of LST is indeed useful and might be complementary to mosquito and virus surveillance and monitoring activities.Furthermore, the set of LST indices derived might be extended to other areas and vector-borne pathogens whose ecology and biology is strongly connected to temperature variations.Besides the application and the results per se, which are yet preliminar, the relevance of this study is also related to the use of open data and free and open source software that provide the advantage of automatization and reproducibility.
Further research is planned to actually relate this environmental information with WNV cases occurrence and incidence, in order to identify the temporal windows relevant to enhance the monitoring and early warning systems to trigger vector control programs and public education campaigns.We intend to continue working in this topic by using higher temporal resolution in response and independent variables, other remotely sensed variables and models that allow the inclusion of spatial and temporal dependencies.

Figure 3 .
Figure 3. WNV human cases in the period 2010-2014 in Greece.Source: European Centre for Disease Prevention and Control (ECDC).

The
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLII-4/W2, 2017 FOSS4G-Europe 2017 -Academic Track, 18-22 July 2017, Marne La Vallée, France (a) Standardized anomalies in the average minimum LST for locations with WNV reported cases split by year of onset.(b) Standardized anomalies in the average minimum LST per municipality (2003-2016).Standardized anomalies per pixel were aggregated with the arithmetic average by municipality.

Figure 5 .
Figure 5. Standardized anomalies in the average minimum LST.

Figure 6 .
Figure 6.Yearly average LST [°C] in the locality with the highest number of reported WNV human cases between 2010 and 2014.
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLII-4/W2, 2017 FOSS4G-Europe 2017 -Academic Track, 18-22 July 2017, Marne La Vallée, France (a) Median number of days with average LST between 20 and 30°C.(b) Median number of days with average minimum LST higher than 20°C.(c)Maximum number of consecutive days with average minimum LST below 0°C.

Figure 8 .
Figure 8. Indices showing number of days with average LST between 20 and 30°C, number of days with average minimum LST higher than 20°C and maximum number of consecutive freezing days aggregated by municipality.
(a) Length of mosquito growing season.(b) Length of WNV transmission season.

Figure 9 .
Figure 9. Length (in days) of mosquito and WNV seasons.

•
Annual number of tropical nights: estimated as the count of days with average minimum temperature ≥ 20°C.•Annual number of days with average temperature ≥ 20 and ≤ 30°C, which is assumed to be the optimal development temperature for mosquitoes.• Maximum annual number of consecutive summer days (average temperature ≥ 25°C).• Maximum annual number of consecutive freezing days (average minimum temperature ≤ 0°C).• Average temperature of mosquito growing season: estimated as the average of daily average temperature between April 1 st and October 31 st each year.• Average temperature of WNV transmision season: estimated as the average of daily average temperature between June 1 st and October 31 st each year.• Length, start and end dates of mosquito growing season: