Spatial and Temporal Heterogeneities of Aedes albopictus Density in La Reunion Island: Rise and Weakness of Entomological Indices

Following the 2006 Chikungunya disease in La Reunion, questions were raised concerning the monitoring survey of Aedes albopictus populations and the entomological indexes used to evaluate population abundance. The objectives of the present study were to determine reliable productivity indexes using a quantitative method to improve entomological surveys and mosquito control measures on Aedes albopictus. Between 2007 and 2011, 4 intervention districts, 24 cities, 990 areas and over 850,000 houses were used to fulfil those objectives. Four indexes including the classical Stegomyia index (House Index, Container Index, Breteau Index) plus an Infested Receptacle Index were studied in order to determine whether temporal (year, month, week) and/or spatial (districts, cities, areas) heterogeneities existed. Temporal variations have been observed with an increase of Ae. albopictus population density over the years, and a seasonality effect with a highest population during the hot and wet season. Spatial clustering was observed at several scales with an important autocorrelation at the area scale. Moreover, the combination among these results and the breeding site productivity obtained during these 5 years allowed us to propose recommendations to monitor Aedes albopictus by eliminating not the most finding sites but the most productive ones. As the other strategies failed in La Reunion, this new approach should should work better.


Introduction
The emergence of arboviruses and their vectors raised health and economic problems during the recent past years. The latest major arbovirus pandemic was due to the Chikungunya arbovirus which was not listed as a major health problem in Africa, Asia and to a lesser extent Europe. Chikungunya expansion started in Kenya in 2005, reached the South Western region of the Indian Ocean [1,2] and in parallel Central Africa in 2005 [3], and then the European [4] and Asian continents [5]. The worldwide dispersal of this arbovirus was due to the lack of Human population immunity, and the efficient competence [6], ecology [7], dynamic [8] and density of the vectors [9]. The major vector incriminated during this pandemy was Aedes albopictus (Skuse) [10,11]. The Chikungunya virus outbreak in [2005][2006] in La Reunion Island, located in the South West of the Indian Ocean, has initiated the monitoring of Ae. albopictus in urban areas of the Island by the vector control service of ARS (Agence Regionale de Sante), local representation of the French Ministry of Health.
Ae. albopictus is the major Aedes species in La Reunion since 1953 and its distribution can be observed in urban, suburban and natural areas [12]. The island landscape is divided by more than 300 ravines crossing urban environments, often composed of houses with gardens. Female Ae. albopictus can lay eggs in a high variety of natural and artificial breeding sites. In La Reunion, bamboo stump, rock holes, flower plates and small recipients are particularly attractive breeding sites [13]. In La Reunion, environmental conditions control the abundance of Aedes albopictus and adults generally emerge between 11 and 14 days after a rainy event. (Boyer et al. In prep). The meteorological data, the availability and the number of breeding sites, the dimension and the urbanization of the cities and the presence of green areas could affect the density of Ae. albopictus in the urban environment in La Reunion. Distribution, abundance, and individual fitness of the mosquitoes depend on biotic and abiotic factors [14,15,16,17], i.e. the presence of unforested areas, water depth, and temperature. In addition, predators and competitors have also significant effects on mosquito larvae [18,19]. Population abundance is markedly cyclical according to the season: the population of Aedes albopictus in New Jersey appeared to be very seasonal correlated with high daily temperature [20]. Moreover the productivity indexes in different urban breeding sites have been investigated on Aedes aegypti but few studies were carried on Ae. albopictus. A lack of information still exist on containers productivity and their relative contribution to the adult population size in the Island.
The indices used to evaluate Stegomyia population densities, such as the house index (HI: percentage of houses with at least one active breeding site), container index (CI: percentage of containers with larvae) and Breteau index (BI: number of active breeding sites per 100 houses) are widely used as empirical standard parameters [21]. These indices approximate the population size and determine monitoring threshold for control programs [11]. The Pan American Health Organization described 3 levels of risk for Dengue transmission with Ae. aegypti based on the HI (low (HI,0.1%), medium (0.1%,HI,5%), and high (HI.5%) [22]) and 2 risk levels for Yellow Fever transmission, with HI.1% and BI.5% [23]. These entomological indexes developed on Ae. aegypti can also allow to identify Dengue transmission areas [23]. The infested receptacle index (IRI: mean number of positive containers per house) was also used to measure public motivation in cleaning up their yards. Based on the strong correlation among the number of pupae and adults, pupae per person index [24], pupae per premises index (PPI) [21] and pupae per hectare index (PHI) [21] were used. In Italy, the number of laid eggs, adults and the pupal demographic monitoring were correlated [25]. But according to the vector control agency, neither pupal indexes or Adult Productivity Index (API) [26], can be daily or weekly monitored in La Reunion. However Stegomyia indexes failed to describe the risk of arbovirus transmission. For example, the Breteau index puts equal weight on all positive containers, however less than 20% of the total containers can be responsible for more than 80% of the Aedes adults [27,28].
The biology of Ae. albopictus is related to water-filled containers distributed throughout suburban landscapes. Production levels of immatures were correlated with the abundance of mosquitopositive containers [29] and the amount of discarded containers in the habitat [30]. The distribution of the containers can create «hot spots» of mosquito production that can be point sources for reinfestation of neighborhoods [31]. The monitoring need to be made at the community level as a single non-complying household can infest an entire neighborhood [20]. In Peru, a clusterisation of the pupae was also related to Ae. aegypti dispersal [32]. Moreover, the mosquito production varied with containers' type in urban and suburban areas [29]. Based on the Index of Container Importance (ICI), a 60 increase-fold of the mosquito production was described due to the presence of positive containers in the neighborhood [29]. The occurrence of these hot spots may spatially vary over time, but the existence or absence of neighborhood clustering of Ae. albopictus can indicate locations where control efforts should be focused [29]. Indeed, the dispersal of Ae. albopictus is driven by the availability of oviposition sites [29]. More, the proximity of containers and other larval habitats increases the probability of pupae abundance [33]. The operational challenge is to adapt the better strategy thorough space and time [33].
We developed an approach with a geographic information system (GIS) combined with spatial statistics with temporal aspects. Indeed, vector control programs could be more efficient if the spatial locations of highly productive areas of Ae. albopictus can be predictable. The purpose of the present study was to identify Ae. albopictus breeding sites, to study the different productivity indexes and to improve the entomological survey monitoring. Without any a priori knowledge on the significance of each of these indexes on the Ae. albopictus population density, we will understand whether the density of mosquitoes is homogeneous or heterogeneous at a spatial and temporal scales.

Ethics statement
During one day randomly chosen per week, ARS agents asked for the permission to examine the water-holding containers in houses and checked every container. Residents provided informed consent to the agents to have their residence used in the study. In private properties, the owner gave permission to conduct the study on its site.The field studies did not involve endangered species. No study on vertebrate animals was carried out.
Additionally owners/residents gave permission to conduct a study on their land in order to determine the productivity of Ae. albopictus breeding sites. These visits occurred in a mission of Public Health and the agents from the Vector Control Service possessed annual District Ordinance allowing interventions in private places. In practice, the agents asked people to visit and it goes well in 98% of the cases with the people agreement. These visits represented around 70,000 visits per year. If people did not agree, the agents continued on to the next house unless the situation justified intervention with the municipal police (which was never the case for this research study). There is therefore no ethics committee since the Vector Control Service is not in research area but in a national mission of Public Health. The inspection of natural and artificial breeding sites were carried out in urban and rural areas to identify and destroy breeding sites, and also to inform people about prophylaxis and vector diseases. Traditional Stegomyia indexes (House Index HI, Container Index CI, Breteau Index BI) and Infested Receptacle Index (IRI) have been selected from this entomological study (Table 1). Public spaces such as schools and parks were also monitored regularly.

Entomological survey and Productivity indexes
The pupae per premises index (PPI), pupae per hectare index (PHI) and the adult density index (ADI) were not used during the monitoring. From a practical point of view, Culex and Aedes pupae are present in the same breeding sites in La Reunion, and it is not possible for the agents to distinguish easily and rapidly the different species. The PHI was also not used because it existed an important heterogeneity inside the same city or zone in term of human population density and habitat density. Moreover, the following of the habitat study was not possible in term of agent complement.

Monitoring of breeding sites
During the routine service request inspection process, breeding sites were destroyed and larvicide was applied only rarely. Yard sanitation and education of the residents are also part of the service request routine: agents demonstrated how to find and definitively eliminate mosquitoes' breeding sites and explained the part of each owner in the integrated community-based control strategy for preventing La Réunion Island from outbreaks. The agents also give some advices in case of active arbovirus transmission for the protection against mosquito bites (use of repellents, trousers and long sleeves, and see quickly a doctor in a case of fever).
Between 2007 and 2011, 60 houses per zone were controlled. Every container containing water was recorded in both the ravines and the urban areas every week. First, immature stages were collected with a pipette and the number of larvae (3 rd and 4 th instars), pupae were counted, and productivity (pupae+larvae) was calculated. Every single container with water was described (physical characteristic: size, depth, amount of water, water quality and appearance (clear, tinted, turbid, polluted), and shade (shadow, mi-shadow, sunny).
Four meteorological stations in four cities (Saint-Pierre, Saint-Denis, Saint-Paul, Saint-Benoît) in the 4 different districts were chosen to determine the effect of precipitations and temperature on the index variation. The zones were chosen near to the meteorological stations (less than 2 kms) with an altitude reflecting the altitude of the meteorological station (less than 300 m) and with a continuous coherent landscape ( Figure 2).

Weekly monitoring
A weekly monitoring was decided to estimate the natural variation within homogeneous breeding sites. One site was followed during 32 weeks in order to estimate productivity fluctuations. Nine different bamboo holes, as a natural habitat, were followed at the same place with the same environmental and meteorological conditions (same temperature, precipitations, relative humidity, wind). Bamboo holes represented one of the most frequent Ae. albopictus natural breeding sites encountered in  La Reunion: identical bamboo stumps (same size, same depth) spaced out 50 cms and located in the ravine ''Chemin Severe'' in Saint-Benoit were chosen. Each bamboo stump was emptied every week and larvae and pupae were counted and then returned to the hole. Natural water was supplied every week.

Data analysis
All data were exported and carried out with JMP 8.0 (Statsoft inc., Paris, France), R 2.14.1 freeware (The R Development Core Team, 2008) and ArcGIS 10.0 (GIS software: ESRI ArcGIS TM , Redlands, CA, United States).
Spatial analysis. The spatial analyses were carried out on the same dataset of 6548 zones used during the temporal analyses. Analyses were carried out on the four indexes per year. We used the global Moran index to test the spatial autocorrelation choice distance (ArcGis Analyst Spatial Tools; 500, 800 and 1000 meters were tested). The Moran's coefficient of autocorrelation test determines whether the distribution of these indexes is spatially clustered, dispersed or random in space [34]. The location of clusters and outliers was detected with the local Moran test [35] (ArcGIS, Spatial Analyst Tools). LISA (Local Indicators of Spatial Association) identified the core of a cluster analysing the values of sites per pair. A local Moran index was determined for each site. The sum of all the local indices was observed proportionality with the global Moran index, which take into account the average value. Then, local Moran indices were compared to the global Moran index in order to determine the nature of the spatial liaison (High-High, High-Low, Low-High or Low-Low). The significance was calculated after a normalization of the local Moran indices (P = 0.05). Analyses were parameterized with the neighbour sites in fixed distance band of 1000 meters. Map representations were specifically developed using the built-in GIS software functions (Arc-Info GIS 9.3 grid module, ESRI Inc).
Meteorological analysis. The meteorological analyses were carried out on a dataset of 957 data in the North (Saint-Denis), 457 in the South (Saint-Pierre), 329 in the West (Saint-Paul) and 501 in the East (Saint-Benoît). The step function results allowed to create nested models. We chose the best model with effecting likelihood-ratio test. After, a Pearson test was performed to determine whether the model was correctly adjusted. The residuals were also analyzed to make sure that normality and homoscedasticity of residuals were followed to carry out an ANOVA.
Temporal analysis. Meteorological data were used to make a time series by day-decades for Saint-Denis zona (The 3.88% missing data were extrapolated from the average of the previous and next indice's values). The stationary (Dickey Fuller test) and the autocorrelation (Ljung-Box test) were tested for each index. The time series were divided into three parts (general, seasonal, residuals). Then a global linear model was carried out with the meteorological data to explain the seasonal results the residuals. Productivity analysis. Analyses were performed on a dataset coming from 1323 breeding sites to study the effect of environmental factors on the productivity. The differences among the types of the container were indicated by an ANOVA. The correlations between the numbers of immature Ae. albopictus and 9 factors (width, depth, type of the breeding site, the town, the organic matter, the effect of the sunshine, the presence of Anopheles and Culex and the quality of the water) were tested with Pearson tests. The influence of the significant factors on the presence of immature stages has been tested using an ANOVA.
Results 850,804 houses, equally distributed among the four districts and representing 25% of the total houses, were visited (Table 2). However, 50.5% of these houses were not inspected due to the absence of the inhabitants or the decline of the owners. 16.6% of houses were positive, representing 69,958 positive houses with 148,741 positive containers ( Table 2).

Temporal analysis
Temporal analysis divided the real values into three different parts, i.e. the global trend of the considered index, the seasonality of the index and the residuals (Figure 3). In general the indexes increased from 2007 to 2010, reaching a plateau after 2010. In details, House Index (HI) values oscillated between 0.14 and 0.18 before 2009, from 0.14 to 0. 22  Ci and IRI values increased continuously from 0.18 up to 0.36 and from 1.8 up to 2.5, respectively. A significant seasonality was observed with the four indices: An annual cycle was apparent in all the zones with low index values during the winter (dry and cold) and high index values during the summer (humid and hot), but with substantial heterogeneity in the relative strength of the multiannual cycle, especially the timing of annual picks. In details, HI and CI values were high from December to February, and low from April to November for HI and between July and November with BI, with a only high pick in the middle of August. BI values were high in January and February and low from May to November. The seasonality of IRI was significant but the seasonality was less marked than the other indices. The residuals were tested and correlated with climatic factors including the precipitations and the temperatures (Table 3).
One site containing identical bamboo stumps breeding sites was followed up every week during 32 weeks in order to estimate the weekly productivity variation. 6

Meteorological analysis
The BI cannot be explained by the meteorological factors: the model cannot be adjusted, meaning that other factor(s) should explain this index. The interactions between the temperatures (minimum, maximum or average) and the precipitations had significant influences on HI, CI and IRI ( Table 4). The geographical differences into temperature and humidity between  (Table 5). In 2008, CI is the only clustered index whereas in 2009 all the indices were clustered except HI ( Table 3). Variations of the clusterisation also existed for each index depending of the year (Figure 4). For example, we can observe a majority of 'Low-Low' clusters in Saint-Paul in 2007 indicating a low density of mosquito, but a majority of 'High-High' clusters in 2008 indicating an important variation in the same city depending of the year. This pattern of index value over time can be observed in several zones: only Cilaos, Salazie and Plaine-des-Palmistes conserved clusters with low levels or no clusterisation.
Breeding site productivity 1,948 breeding sites have been assessed with 1,323 positive breeding sites (67.9%) in which 33,332 larvae and 4,554 pupae were counted. The number of pupae, larvae (instar 3 & 4) and the productivity (pupae + larvae) were correlated with length, diameter, sun exposition, types of breeding sites, and depth of the breeding sites. Moreover, the number of pupae was correlated with the quality of the water. These 3 variables were not correlated with the area, cities or with the presence/absence of Anopheles and Culex individuals. In details, the number of Aedes was significantly higher in the deepest and widest shaded breeding sites. In the same way, more larvae and pupae were found in the largest surface breeding-sites (observed maximum area = 8 m64 m). The water quality was correlated only with the number of pupae: more pupae were counted in troubled and polluted water than in clear water.
The productivity was significantly higher in the breeding sites constantly under shadow (60.10) than in the breeding sites always exposed to the sun (5.43). The breeding sites alternatively under the shadow or exposed to the sun presented an intermediate productivity (22.54). The same pattern was observed with larvae and the pupae.
Ae. albopictus breeding sites were classified into 7 classes: natural breeding sites such as bamboo stumps, tree and rock holes generally found in ravines, flower plates, small and big containers, tires, basins and reservoirs, and vegetal ones. The number of immature Ae. albopictus per breeding sites ranged from 5.3 up to 26.8 in natural breeding sites, vegetal breeding sites and flower plates. Intermediate productivity values were obtained in small and big containers (71.5-72.5) ( Figure 5). The highest numbers of immature were found in basins and reservoirs (in average 122.2) and particularly in tires (in average 109.8).

Aedes albopictus and Stegomyia indices
The four Stegomyia indicators were able to detect seasonal variations of Ae. albopictus, as already demonstrated for Ae. aegypti, and highlighted the summer as the period at major risk of high mosquito density. However these indexes applied to Ae. albopictus showed important limitations as predictors of possible outbreaks in La Reunion, likely due to the high spatial and temporal heterogeneity values which overtake the standards for other vector-borne diseases. Since 1996, it is classically accepted that epidemics of Yellow Fever, Dengue and other arboviruses would not occur when HI,1% and BI,5 [31]. In Cuba, BI.1 is threshold for action [23]. In Trinidad the indexes for Ae. aegypti during a Dengue Fever crisis between 2002 and 2004 [38], ranged from 20 to 34 during the dry season and from 30 to 46 during the wet season. In La Reunion BI for Ae. albopictus were always superior to 20. Similarly, HI was always superior to 12% and CI to 20%. These thresholds based on Ae. aegypti biology, may not be applied for a highly invasive, competitive and opportunistic species such as Ae. albopictus. Indeed the thresholds for defining the alert would be disease-specific depending on the transmission potential of the considered disease and of the involved mosquito species. The classical Stegomyia indexes revealed little use to orientate the choices of the control agency. In addition, socio-economic and environmental conditions of La Reunion can affect the ability of these indexes to characterize the variation of Ae. albopictus densities. Particularly, BI which was defined with Aedes agypti cannot be explained by the season or the meteorological data, indicating that this index is not an entomological index but more a socioentomological index, representative of the population involvement. IRI was the less discriminant index suggesting not a great difference in human behaviour. The CI and the HI were the most discriminant index from this study, and CI is also the one (such as IRI) highlighting the HH cluster where the Chikungunya human cases appeared in 2010 ( Figure S1). However an alert threshold corresponding to the biological characteristic of La Reunion remains to be defined. But the final problem concerning the prophylaxis indicators described or developed within this study remains that they are based on larval indicators, and not on adult indicators, which are the real indicators of a potential transmitted disease risk.

Prophylaxis indicators
A better involvement of more qualified technicians, a higher number of houses visited, and a new and more efficient protocol of  [36], leading to a continuous increase of Ae. albopictus, despite a good knowledge of the vectorial risk for 78.9% of them. This lack of interest impacts the BI. Some widespread behaviors (water in potted plants, regular watering) favor the persistence of breeding sites in over 22% of the controlled houses in 2011. Human behavior more than the weather can allow the maintenance of the Ae. albopictus population.
The most productive breeding sites are the largest ones, i.e. basins, reservoirs, and tires. Key containers are defined as individual containers of any type that contribute significantly to an urban Aedes problem [22]. Often, the mean number of pupae per container was considered most representative of the production of adults because of the low pupal mortality and a high correlation was reported for the number of pupae and the number of adults for Ae. aegypti and Ae. albopictus [39] [40]. In our study case, we observed no difference between pupal and larval productivity.
The productivity of breeding sites needed to introduce indexes which could take into account the real productivity depending of the type of breeding sites. These results led to the determination of a breeding site productivity ranking, defining the most productive sites needing more intervention. For the control of Ae. albopictus  populations, we attributed a productivity coefficient for every type of breeding sites class as defined by the control vector unit. The monitoring of entomological indexes need to be simple, reproducible, reliable and pertinent. It requires the use of several indexes and a routine monitoring of adults must be developed. We propose a monitoring based on the abundance of containers and their characterization (7 classes) without counting of the number of larvae and pupae in urban breeding site. A quantifiable productivity index can help the mosquito control units to be more efficient by targeting the most productive breeding sites. This is in line with the proposition of Yebakima to define a weighted index for Aedes species which consider the specific productivity of each breeding site [41]. Yebakima developed a weighted Breteau Index to specifically characterize the mosquito breeding sites according to their characteristic productivity.

Spatiotemporal heterogeneities
The mosquito population varied monthly, increasing during the summer months, according to the meteorological data. Interesting variations over the week between breeding-sites few centimeters apart were also observed not related to meteorological data, human intervention or population behaviour. The intrinsic productivity value of a breeding site is not constant despite the same environmental conditions. This random production added complexity to the analysis of the breeding site productivity and highlighted the actual limits of the monitoring. Important spatial variations of the productivity were observed from macro-scale (districts, city) to micro-scale (zones). This spatial discontinuity could be explained by the absence of ecological continuity in La Reunion, because of the alternation ravine/city and a lack of an overall/global land use. Among the 5 most productive cities, 4 are most inhabited cities and are located along the coast, at low altitudes. Properties often have small gardens with potted flowers, and high number of inhabitants increases the environmental damages with more garbage, tires and wastes leading to more mosquito breeding sites. In contrast, the 3 cities with the lowest index values are located inland at high altitudes and are isolated by high altitude geographical barriers.
In addition to this variations, an important clusterisation was found, which shows the complexity of the distribution of Ae. albopictus within the different localities and the difficulty of vector control. A consequence of this highly clustered, local spatial pattern is that missing some houses during vector control operations can leave intact mosquito clusters that could repopulate the area [37]. This clusterisation is also temporal as some zones can change from HH intensity to LL. This has a direct impact on vector control, preventing the localization of highly productive premises that could be targeted [37]. Because of all these variations, the data of this study can only be interpreted as a snapshot but cannot be predictive of outbreaks: the highly productive locations could not be localized in advance with this spatial approach in La Reunion.

Rise and Weakness
We can conclude on advisedly use and misuse of these indicators. Among the good surprises, the 4 Stegomyia index described well the seasonal dynamic and density of the vector, even if the Breteau Index should not be considered anymore as an entomological index. These indices were fluctuant, seasonal and discriminant. We were also able to highlight some limits including that these indicators are not epidemiological risks' indicators. Moreover, these 4 indicators are not suitable to the prediction of the epidemiological risk because the thresholds were not well defined, and the WHO thresholds were not suitable. Moreover, with the actual sampling of the control agency, the cryptic and public larvae breeding ground, tree holes or ravines' breeding sites were not detected. And, due to the specificity of La Reunion landscape, an exhaustive sampling would imply a too high number of agents. So far, the systematic destruction of Ae. albopictus habitats has failed to control the species population and a continuous recolonization of areas from which individuals have been eliminated was observed. In La Reunion, the insecticide spraying and the larval breeding sites suppression were ineffective to control Ae. albopictus population. A specific suppression of the most productive breeding sites is likely to achieve better results, at the same time reducing the cost of the intervention. Figure S1 Variation of the clusterisation for each index depending of the considered yea (with Dengue and Chilungunya human cases. This clusterisation is the Anselin representation of the clusterisation of local Moran index. HH represents High-High clusters which is a cluster with a high index surrounded with clusters with high index; LL Low-Low, LH Low-High, and HL High-Low. (TIF)