Spatio-temporal modelling of the first Chikungunya epidemic in an intra-urban setting: The role of socioeconomic status, environment and temperature

Three key elements are the drivers of Aedes-borne disease: mosquito infestation, virus circulating, and susceptible human population. However, information on these aspects is not easily available in low- and middle-income countries. We analysed data on factors that influence one or more of those elements to study the first chikungunya epidemic in Rio de Janeiro city in 2016. Using spatio-temporal models, under the Bayesian framework, we estimated the association of those factors with chikungunya reported cases by neighbourhood and week. To estimate the minimum temperature effect in a non-linear fashion, we used a transfer function considering an instantaneous effect and propagation of a proportion of such effect to future times. The sociodevelopment index and the proportion of green areas (areas with agriculture, swamps and shoals, tree and shrub cover, and woody-grass cover) were included in the model with time-varying coefficients, allowing us to explore how their associations with the number of cases change throughout the epidemic. There were 13627 chikungunya cases in the study period. The sociodevelopment index presented the strongest association, inversely related to the risk of cases. Such association was more pronounced in the first weeks, indicating that socioeconomically vulnerable neighbourhoods were affected first and hardest by the epidemic. The proportion of green areas effect was null for most weeks. The temperature was directly associated with the risk of chikungunya for most neighbourhoods, with different decaying patterns. The temperature effect persisted longer where the epidemic was concentrated. In such locations, interventions should be designed to be continuous and to work in the long term. We observed that the role of the covariates changes over time. Therefore, time-varying coefficients should be widely incorporated when modelling Aedes-borne diseases. Our model contributed to the understanding of the spatio-temporal dynamics of an urban Aedes-borne disease introduction in a tropical metropolitan city.


Introduction
The first chikungunya virus (CHIKV) epidemic in Rio de Janeiro city, the second most populated city in Brazil and its leading tourist destination, occurred in 2016 [1].CHIKV is transmitted to humans by the same vectors as dengue viruses (DENV), the Aedes mosquitoes [2].
Vector-control activities have not prevented Rio de Janeiro from being endemic for dengue for years, nor from having experienced large dengue epidemics every three to four years, in general [3][4][5].
For a chikungunya epidemic to occur, three main elements are necessary, represented by the blue area in Fig 1: mosquito population, susceptible human population, and the virus circulating [6][7][8].The Ae. aegypti mosquito is present all over the city of Rio de Janeiro, facilitating a new arbovirus is established.Because CHIKV and DENV belong to different families, previous immunity to DENV does not cross-react with CHIKV, and the population of Rio de Janeiro can be considered equally naïve to CHIKV before 2016.Therefore, given the presence of the mosquito population and susceptible human population, the occurrence of local CHIKV transmission in Rio de Janeiro was conditioned by the entry of the virus.Once the virus is circulating and established, how it will spread within the city depends on multiple factors.As a consequence, some areas of the city will experience the epidemic at different times, and will also have different attack rates.Reliable data regard to the intra-urban level of vector population, the susceptible human population over time, and the entry time of the virus in each location are necessary to understand the spatio-temporal dynamics of an epidemic.However, these data are not usually available.
Alternatively, measured factors -such as temperature, socioeconomic and environmental factorsindirectly associated with the number of chikungunya cases (represented by red arrows in Fig 1) can be used as proxies.These factors have a direct effect (represented by black arrows in Fig 1) on the necessary elements for an epidemic to occur (inside the blue area), which are unmeasured.The temperature varies with time, and also present different effects on the epidemic over time.
Socioeconomic and environmental characteristics take long periods to show important changes and can be considered fixed during the period of an epidemic.However, the way these factors impact the number of cases can change as the epidemic progresses, i.e., with time-varying effects.
The temperature affects the Ae.aegypti population, the virus, and the mosquito-virus interaction [9].Ideal temperature conditions accelerate all stages of the mosquito life cycle, increasing the population in the long term.In the short term, the temperature influences the mosquito activity as well as the length of the virus incubation period, with maximal transmission occurring around 26-29°C [10].Regarding the environment, the Ae.aegypti mosquitoes are highly adapted to urban settings, and the level of urbanisation is inversely correlated with the proportion of green areas [11].The socioeconomic status impacts the mosquito population as disorderly urbanisation and inadequate sanitary conditions favour the presence of the mosquito most common reproduction site: containers filled with water found inside or in the surroundings of domiciles [12,13].Furthermore, high population densities favour the contact between the mosquito and the human, increasing the chance of infection and transmission [14].
In the last decades, models that take into account the spatial dependency structure of the cases have been applied to better estimate covariates' association with Aedes-borne diseases epidemics [15][16][17][18].The application of spatial models for intra-urban settings is growing more recently [19][20][21][22].These model types, as the intrinsic conditional autoregressive (ICAR) models [23], are built under the assumption that adjacent areas share similar characteristics.The inclusion of a latent spatial random effect accounts for both the spatial structure and unmeasurable factors [24].The inclusion of time-varying coefficients allows us to explore how the effect of the covariates change throughout the epidemic.The temperature is usually included in statistical models for Aedes-borne diseases in a linear fashion and with a pre-defined lag.We propose to estimate the temperature effect in a non-linear framework using a transfer function, including an immediate effect and a memory effect that propagates to future times.Another advantage of using a transfer function is that the estimation of the lag of the effect is data-driven [25].
We used this methodological approach (ICAR models with time-varying coefficients and a transfer function) to identify how temperature, socioeconomic and environmental factors are related to the space-time progression of an Aedes-borne disease epidemic in an intra-urban setting.The first chikungunya epidemic in Rio de Janeiro, a large tropical city with environmental and socioeconomic disparities, presents ideal conditions to this end.

Study site
Rio de Janeiro is the second-largest city in Brazil, with 6,3 million inhabitants (2010), and its primary tourist destination.Rio's area is of 1204 km², with 160 neighbourhoods grouped into four large regions (Downtown, South, North and West).These regions are subdivided in 10 health districts called programmatic areas: area 1.0 (Downtown region); areas 2.1 and 2.2 (South region);  With three mountain massifs and 84 km of beaches, Rio has a diverse geography that is directly associated with the history of occupation and with socioeconomic disparities [26].The Downtown region is the historical, commercial and financial centre of the city, with many cultural establishments.The South region is the most popular tourist destination, with famous beaches and wealthy neighbourhoods.In the North region, there are very large slums ("favelas") and nearly 27% of the population, almost 2.4 million people, live in such communities [27].The West region has more heterogeneous characteristics among its neighbourhoods, being the area 5.1 more densely populated, areas 5.2 and 5.3 less urbanised, and area 4.0 wealthier.

Chikungunya cases
Data on chikungunya cases were obtained from the Notifiable Diseases Information System (SINAN) via the Rio de Janeiro Municipal Secretariat of Heath and are publicly available [28].
We analysed notified cases of chikungunya (confirmed by laboratory or clinicalepidemiological criteria) occurring in Rio de Janeiro municipality between January and December 2016, by week and neighbourhood of the patient's residence.
Case definitions follow the Ministry of Health protocols.A suspected case of chikungunya is defined as a patient with sudden fever of over 38.5°C and severe arthralgia or arthritis not explained by other conditions, and who either lives in endemic areas or has visited one up to two weeks before the onset of symptoms or has an epidemiological link with a confirmed case.A confirmed case is a suspected case with at least one positive specific laboratory test for CHIKV or confirmed by clinical-epidemiological criteria [29].

Socioeconomic data
To characterise the socioeconomic status, we obtained the sociodevelopment index data by neighbourhood from the Instituto Pereira Passos [30].This index is based on eight indicators from the 2010 Demographic Census: 1) the percentage of domiciles with adequate water supply; 2) the percentage of domiciles with adequate sewage; 3) the percentage of domiciles with garbage collection; 4) the average number of toilets per resident; 5) the percentage of illiteracy among residents between 10 and 14 years old; 6) per capita income of the domiciles, expressed as minimum wages; 7) the percentage of domiciles with per capita income up to one minimum wage, and 8) the percentage of domiciles with per capita income greater than five minimum wage.The sociodevelopment index is calculated as the arithmetic average of the normalised indicators (each ranging from 0 to 1, being 0 the worst socioeconomic condition and 1, the best) [30].
We also obtained data on the population by neighbourhood from the Instituto Pereira Passos [31].

Environment and temperature data
Land use data for the city of Rio de Janeiro were obtained from the Instituto Pereira Passos as a shapefile [32].We created the category "green areas" by aggregating: agricultural areas, areas with swamps and shoals, areas with tree and shrub cover, and areas with woody-grass cover.After that, we calculated the proportion of green areas for each neighbourhood (Fig 4B).
Temperature information for 2016 was obtained from 38 meteorological weather stations in Rio de Janeiro, from five different meteorological and environmental institutes.The institutes are the Brazilian National Institute of Meteorology [33], the Brazilian Airspace Control Department [34], the Rio de Janeiro State Environmental Institute [35], the Rio de Janeiro Municipal Environmental Secretariat [36] and the Alerta Rio System [37].All measurements are made according to the recommendations of the World Meteorological Organization [38].The institutes make their meteorological data publicly available, being that the frequency of measurements of the first four organisations is hourly, while the Alerta Rio System measurements are made every 15 minutes.
From the temperature measurements, we computed the daily maximum, minimum and mean temperature.We also evaluated the availability of daily data in terms of missing measurements.The daily records that had more than 60% of missing measurements were excluded.We decided to use the minimum temperature as in tropical climates it acts as a limiting factor for the Ae.aegypti activity and population [39,40].We then obtained the minimum temperature for each week and station, and, to obtain the minimum temperature by neighbourhood, we applied universal kriging.
Briefly, kriging is a method that uses a sample of data points to estimate the value of a given variable over a continuous space [41].First, we interpolated the minimum temperature to a grid with each unit measuring 500m x 500m.The grid with the meteorological weather stations locations is displayed in the S1 Fig. Then we obtained the minimum temperature of the neighbourhood by calculating the average of the minimum temperature of the grid units whose centroids were within the boundaries of the neighbourhood.

Statistical analysis
We used the Stan platform [46,47] to fit spatio-temporal models, more specifically ICAR models, to a dataset consisting of neighbourhoods counts of chikungunya cases, exploring the relationship with sociodevelopment index, the proportion of green areas and minimum temperature.
Let Y i,t be the counts of chikungunya cases at neighbourhood i = 1, 2, …, n = 160, and week t = 1, 2, … , T, where Y i ,t ∼Poisson(μ i ,t ) .We explored different structures for μ i,t , presented in the S1 Appendix along with each Watanabe-Akaike information criterion (WAIC) [48].The model selected based on the WAIC has the following structure: The latent spatial effect is represented by ɸ, which, a priori, follows a conditional autoregressive distribution [23]; that is, the conditional distribution of each ɸ i follows a normal distribution whose mean and variance depend on the neighbourhood structure w ij .Assuming a binary neighbourhood structure, where w ij =1 if areas i and j share borders and 0 otherwise, each ɸ i 10 .
follows a conditional normal distribution whose mean is equal to the average of its neighbours, and its variance is inversely proportional to the number of neighbours d i , that is: The expected number of chikungunya cases at neighbourhood i (e i ) represents the number of cases that would have been observed if there were no differences in the incidence of cases across time and space: β 0 is the intercept, X' i represents a vector of k covariates and β k,t is the coefficient of covariate k in week t.We decided to allow for time-varying coefficients for the covariates to explore if their association with the number of cases vary as the epidemic progresses.The covariates included in the X' i vector were sociodevelopment index and proportion of green areas.The proportion of green areas showed a skewed distribution.Therefore, this variable was transformed to the cubic root.We also fitted models including the population density.However, the 90% credible interval (CI) of the population density coefficient included 0 for most weeks, and the inclusion of this variable increased the WAIC .Therefore, the population density was not considered in the final model.
We estimated the temperature effect in a non-linear fashion using a transfer function (U i,t ), considering that the temperature has an immediate effect (ζ i ) and that a proportion (ρ i ) of this effect propagates to future times.This proportion ρ i is called memory effect and we considered it to be any value between 0 and 1 [25].To combine and visualise both effects of the temperature, we obtained the impulse response function of the temperature for each neighbourhood.This function expresses the effect of a 1 unit increase in the temperature of one week propagating in time [25].The temperature was standardised.
The models were fitted under the Bayesian framework using the Stan platform [46,47] to run four chains of 10000 iterations each where the first 5000 were the warmup.Relatively vague, proper prior distributions were used.We visually inspected the chains and used the R-hat statistic to check convergence [46,49].Model validation was performed by comparing the fitted number of cases (mean and 90% CI) with the observed number of cases [48].It is worth mentioning that we also fitted models that considered the reparametrisation of the Besag-York-Mollié (BYM2) as proposed by Riebler et al. [50].However, the random component was over 90% spatial, and the unstructured effect 90% CI included zero for all neighbourhoods.

Ethics Statement
This study was approved by the Research Ethics Committee of Escola Nacional de Saúde Pública Sergio Arouca (ENSP) -Fundação Oswaldo Cruz, approval number 2.879.430.Informed consent was not required as this is a study using secondary data and the data were analysed anonymously.

Results
Between January and December 2016, 13,627 cases of chikungunya were notified in the city of Rio de Janeiro, corresponding to an incidence of 21.6 cases per 10,000 inhabitants.The number of cases peaked at week 17/2016, with 1118 chikungunya cases (Fig 3A).The cumulative number of cases by neighbourhood ranged from 0 (Grumari, area 4.0) to 721 (Realengo, area 5.1).The highest incidence was found in Catumbi (area 1.0), of 211.0 cases per 10,000 inhabitants (Fig 3B).where the instantaneous temperature effect was null, the memory effect is irrelevant.entire study period.The classification of the neighbourhoods in terms of chikungunya relative risk considering the 90% CI is available in the S5 Fig.

Discussion
In this study, we estimated the associaton of socioeconomic status, environment, and temperature with the spatio-temporal distribution of the first chikungunya epidemic in Rio de Janeiro city.The sociodevelopment index and the proportion of green areas were included in the model with time-varying coefficients, allowing us to explore how the effects of these factors changed throughout the epidemic.The temperature was included in the model in a non-linear fashion using a transfer function, considering that the temperature has an immediate effect and that 20 .
a proportion of this effect propagates to future times.To our knowledge, this is the first time a transfer function is applied for temperature when modelling Aedes-borne diseases.
The sociodevelopment index was inversely associated with the risk of chikungunya in all models (Fig 5A,Figs A and B in S2 Fig).
In fact, the sociodevelopment index presented the strongest effect in the models, which strengthens the hypothesis of chikungunya being a disease of social determination.This index is composed of sanitary conditions indicators, among others, and poor sanitary conditions are known to favour the reproduction of the Ae.aegypti mosquitoes.The association of low socioeconomic locations with increased risk of chikungunya was also found in a study in French Guiana [56] and in a study in Barraquilla, a Colombian city [21].Our results indicate that poor neighbourhoods were affected first and hardest by the chikungunya epidemic, highlighting the importance of vector control activities in socioeconomically vulnerable locations.
When the spatial dependency was not included in the model, the proportion of green areas was negatively associated with the number of chikungunya cases (Fig C in S2 Fig).Such association was observed for dengue in São Paulo, where low vegetation cover areas presented higher dengue incidence rates [22].However, with the inclusion of the spatial component, the effect of the proportion of green areas moved towards the null.This is possible due to spatial confounding, which happens when covariates that are spatially smooth are collinear with spatial random effects [57].In Rio de Janeiro, the majority of the green areas are in the mountain massifs (Fig 2), which trespasses different neighbouring borders.
The temperature was associated with an increase in the risk of chikungunya in most neighbourhoods.In our models, we assumed that the temperature presents an effect that is a combination of an instantaneous effect and of a proportion of this that propagates in time (the memory effect).Although the instantaneous association is relatively low across the different neighbourhoods, it can persist for a long period of time for some of the districts.This was only possible because we allowed the parameters present in the transfer function to change across neighbourhoods.Different areas experiencing heterogeneous temperature effects for Aedes-borne diseases were also observed in previous studies [58,59].The inclusion of space-varying coefficients for the temperature improved the fitting of our models.
The instantaneous association represents the effect of the temperature on the activity of the mosquito, human behaviour and the mosquito-virus interaction.The biting rate of the Ae.aegypti increases with the temperature until around 35 °C [10], and people become more exposed to mosquitoes in warm temperatures.The temperature also accelerates the virus extrinsic incubation period in the mosquito, and transmission was estimated to peak around 28.5 °C [10].On the other hand, the temperature affects the population of mosquitoes by increasing fecundity, egg-to-adultsurvival, development rate and lifespan [10].This effect is not only on the same week but also accumulates in time, which is captured by the memory effect.
Interestingly, the CHIKV epidemic in Rio de Janeiro in 2016 did not reach the whole city, with high-risk areas mostly concentrated in the North and Downtown regions.The decrease in the number of cases coincided with the drop in the minimum temperature, around week 17 (Fig 3A and 4D).These two observations combined suggest that the epidemic was interrupted not because of susceptible human population depletion in the city, but because the drop in temperature caused a reduction in the transmission in such a way that the epidemic was not sustained.It is important to note that although the number of cases diminished substantially, there were still chikungunya cases being reported until the end of the year.Rio de Janeiro is a tropical city, and the minimum temperature rarely is below the minimum temperature needed for transmission to occur, of 13.5 °C [10].A previous study conducted in the city showed that the Ae.aegypti population varies seasonally, but the mosquito is endemic all over the year [60].This could explain the long persistence of the temperature effect in time in some neighbourhoods (Fig 9 and S1 Video).
The application of spatial models considering intra-urban scenarios is still growing for Aedes-borne diseases.Our study identified high-risk neighbourhoods for the first chikungunya epidemic in Rio de Janeiro city, which concentrated mostly in the North and Downtown regions (Fig 10 and S2 Video).Such regions were already identified as high-risk locations for dengue [61].
In our previous study, neighbourhoods from these regions were more likely to constitute simultaneous clusters for dengue, Zika and chikungunya [62].These regions have a combination of factors that favours the Ae.aegypti's ecology and transmission of diseases: low vegetation, low socioeconomic status and increased temperature (Fig 4).Vector control activities should be prioritised and intensified in the identified high-risk areas, as they also appear to be the first ones affected by the epidemic.Additionally, the long persistence of the temperature effect (Fig 9 and S1 Video) and a small number of cases even after the decline of the epidemic, indicate that the mosquito continues to circulate and to transmit the disease throughout the year.Therefore, interventions must be designed to be continuous and to work in the long term in these locations.
Our study has some limitations.As for any study using secondary data on Aedes-borne diseases, there is an uncertainty on the diagnosis of the reported cases as well as under-reporting.It is important to consider that, in the same year, the city was also experiencing dengue and Zika epidemics [62].Because of the association between Zika and severe congenital manifestations, the disease awareness around Zika may have improved the reporting rates [63].However, the simultaneous occurrence of three arbovirus epidemics may have impaired the differential diagnosis, as they cause similar symptoms.Another limitation is the spatial unit.We analysed the data aggregated at the neighbourhood level.Although the data are more reliable at this spatial unit, smaller areas inside the same neighbourhood can present different socioeconomic and environmental characteristics.Finer scales such as census tracts should be considered in future studies.Finally, an important limitation is the assumption that the chikungunya risk is related to the neighbourhood of residence, while some people may get infected in other locations.
The model here presented has the potential to be applied to other cities and other urban Aedes-borne diseases.Mosquito population information is expensive to collect and often unreliable.
Therefore not depending on such data is a strength of our model.Another strength is the application of a transfer function to estimate the non-linear effect of the temperature and its duration.By using temperature, socioeconomic status and proportion of green areas data as proxies of the key elements of CHIKV transmission, our model contributed to understanding better the spatio-temporal dynamics of the first chikungunya epidemic in a tropical metropolitan city.Importantly, our results indicate that even considering the environment and the temperature, the socioeconomic status plays a major role in affecting the incidence and distribution of a first Aedes-borne disease epidemic in a large city.This strengthens the importance of improving sanitary conditions and taking measures to diminish social inequality as necessary actions to control and prevent Aedes-borne diseases.

Fig 1 .
Fig 1.A theoretical model for a chikungunya epidemic in a given region.Direct associations are represented by black arrows and indirect associations by red arrows.The blue area includes the necessary elements for the epidemic to occur.

Fig 4 .
Fig 4. Sociodevelopment index in 2010 (A), percentage of green areas in 2015 (B), minimum temperature (ºC) average in 2016 by neighbourhood (C) and boxplot of the minimum temperature (ºC) by neighbourhood and week (D), Rio de Janeiro city, Brazil.

Fig 5 .
Fig 5. Time-varying coefficients (in the log scale, mean and 90% credible interval) for sociodevelopment index (SDI) (A) and proportion of green areas (B) for a spatial model for chikungunya cases from weeks 9 to 52 2016 and controlling for minimum temperature, Rio de Janeiro city, Brazil.

Fig 6 .
Fig 6.Chikungunya cases mean spatial effects (in the log scale) controlling for sociodevelopment index, proportion of green areas, and minimum temperature, weeks 9 to 52 2016, Rio de Janeiro city, Brazil.

Fig 7 .
Fig 7. Minimum temperature instantaneous effect (in the same week) (A) and memory effect (B) on chikungunya cases (in the log scale) by neighbourhood, mean and 90% credible interval, controlling for sociodevelopment index and proportion of green areas, and the latent spatial effect, weeks 9 to 52 2016, Rio de Janeiro city, Brazil.

Fig 8 .
Fig 8. Impulse response function of the minimum temperature effect on chikungunya cases over time, posterior mean and 90% credible interval, controlling for sociodevelopment index and proportion of green areas, and the latent spatial effect, in selected neighbourhoods, Rio de Janeiro city, Brazil.

Fig 9 .
Fig 9. Minimum temperature instantaneous effect on chikungunya cases and its propagation in time by neighbourhood in selected weeks, controlling for sociodevelopment index and green areas proportion, and the latent spatial effect, Rio de Janeiro city, Brazil.

Fig 10 .
Fig 10.Posterior mean chikungunya relative risk by neighbourhood in selected weeks, controlling for sociodevelopment index, proportion of green areas, and minimum temperature, and the latent spatial effect, Rio de Janeiro city, Brazil.