Diversity and ecological niche model of malaria vector and non-vector mosquito species in Covè, Ouinhi, and Zangnanado, Southern Benin

The present study aimed to assess mosquito species diversity, distribution, and ecological preferences in the Covè, Ouinhi, and Zangnanado communes, Southern Benin. Such information is critical to understand mosquito bio-ecology and to focus control efforts in high-risk areas for vector-borne diseases. Mosquito collections occurred quarterly in 60 clusters between June 2020 and April 2021, using human landing catches. In addition to the seasonal mosquito abundance, Shannon's diversity, Simpson, and Pielou's equitability indices were also evaluated to assess mosquito diversity. Ecological niche models were developed with MaxEnt using environmental variables to assess species distribution. Overall, mosquito density was higher in the wet season than in the dry season in all communes. A significantly higher Shannon's diversity index was also observed in the wet season than in the dry seasons in all communes (p < 0.05). Habitat suitability of An. gambiae s.s., An. coluzzii, Cx. quinquefasciatus and Ma. africana was highly influenced by slope, isothermality, site aspect, elevation, and precipitation seasonality in both wet and dry seasons. Overall, depending on the season, the ecological preferences of the four main mosquito species were variable across study communes. This emphasizes the impact of environmental conditions on mosquito species distribution. Moreover, mosquito populations were found to be more diverse in the wet season compared to the dry season.


Collection, processing and occurrences of mosquito species
Collections of adult mosquitoes were carried out every three months across study clusters between June 2020 and April 2021.In each cluster, four households located 15-20 m apart were selected to ease the overnight supervision of mosquito collections that occurred using the human landing catches (HLC) technique.Four collectors were assigned to each household, with two collectors (one indoor and one outdoor) collecting mosquitoes from 19:00 to 01:00, and the second group collecting from 01:00 to 07:00 in the morning.Mosquito specimens collected through HLCs were morphologically identified and separated as Culicinae and Anophelinae using a binocular loupe.Anophelinae mosquitoes were identified at the species level using the keys of Gillies and De Meillon 23 and Gillies and Coetzee 24 .Culicinae mosquitoes were also identified at the species level using the identification keys of Doby (1955) and Edwards (1941).Specimens of An. gambiae s.l. were preserved using silica gel and molecularly identified using the protocol of Santolamazza et al. 25 .
Two databases, one with data collected during the wet season (June to November 2020), and the second with data collected over the dry season (December 2020 to April 2021) were generated.Each database included five mosquito species, of which three (An.gambiae sensu lacto, Cx. quinquefasciatus, and Ma.africana) were morphologically identified, and two (An.gambiae s.s., An. coluzzii) were molecularly identified.The databases contained both occurrence and mosquito abundance data.The geographical coordinates of mosquito collection households were used as occurrence data points for the different mosquito species.Households where a species was absent were excluded from the occurrence data for that species.Species occurrence points were rarefied at a resolution of 500 m to avoid pseudo-replication of data points and reduce clustering.

Climatic and environmental data
Overall, a total of 23 topographical, climatic and landscape variables were obtained to build the models (Supplementary Table S1).The 19 bioclimatic variables (bio1 to bio19), slope, Shuttle Radar Topography Mission (SRTM) elevation, soil type and aspect datasets were downloaded from the WorldClim database at a spatial resolution of 30 arc seconds (~ 1 km) 26 .Landsat 8 and 9 satellite data were downloaded from the USGS website 27 and the normalized differential vegetation index (NDVI) was calculated using ArcGIS software.All raster layers were cropped to the extent of the Covè-Zangnanado-Ouinhi communes, resampled to 1000 m resolution using bilinear interpolation, and transformed to the projected WGS 1984 World Mercator coordinate system using the Extract by Mask function in ArcGIS 10.8.The correlation between bioclimatic variables was checked, and highly correlated variables (Pearson's correlation coefficients r > 0.7) were excluded from the model.To maximise model performance, minimise model overfitting (minimise the number of variables in the model) and multicollinearity (minimise the correlation between variables in the models), we used Pearson's correlation coefficient (r) to select the final predictor variables from the initial pool of 19 climate variables described above.Annual mean temperature (bio1), isothermality (bio3), maximum temperature of warmest month (bio5), annual precipitation (bio12), seasonality of precipitation (bio15), Precipitation of Wettest Quarter (bio16) were the less correlated bioclimatic variables (r < 0. 7) that were retained to build the ecological niche models for An.gambiae s.s., An. coluzzii, Cx. quinquefasciatus and Ma.africana during the dry and wet seasons.

Ecological niche models
Due to its robust predictive accuracy 28 , predictive performance using presence-only data 29 and insensitivity to sample size 30 , the maximum entropy (MaxEnt) modelling approach is frequently used to develop ecological niche models.
To generate habitat suitability maps for An.gambiae s.s., An. coluzzii, Cx. quinquefasciatus and Ma.africana, 12 topographic, landscape and climatic variables were applied as predictor environmental variables alongside species occurrence data using MaxEnt version 3.4.4 31.The selected output option was "cloglog" which produces habitat suitability estimates in a range between 0 and 1 for each pixel.The modelling process for each species used k-fold cross-validation replicated over 20 iterations to partition the occurrence data and use all presence points for training and testing 32 .MaxEnt parameters were configured with a maximum of 10,000 iterations, keeping other parameters as default.Linear, quadratic and product transformations of the environmental variables were used with a regularisation multiplier of 1.
Several methods including the jackknife test and permutation techniques, were used to assess the importance of the variables in the model.The jackknife procedure was used to assess the regularized training gain of models built using each variable individually, and when each variable was excluded from the model in turn.Variables for which the model gain decreased the most when omitted or which increased the model gain the most when considered in isolation were identified as the most informative variables in the model that contain information not found in the other variables.Measuring importance by permutation involves randomly permuting the values of each variable between the presence and background data, and then re-evaluating the model.The resulting variation in fit performance is calculated for each variable, where a more significant variation indicates an increased importance of that variable.
Finally, the assessment of model performance for each species used the area under the receiver operating characteristics curve (AUC) based on the testing data, calculated as an average of the AUCs obtained over the 20 iterations of the model.An AUC of 1 indicates perfect prediction, while an AUC of 0.5 suggests prediction equivalent to chance.

Data analysis
To scrutinize the structure of the mosquito population within each commune, three indicators were determined: • The Shannon's diversity index (H′) assesses the different types of mosquito species present.It ranges between 0 (homogeneous population consisting of only one species) and Hʹmax equals to log S (heterogeneous population for which all individuals of all species are equally distributed), with S being the species richness (total number of species identified in a given area).This was calculated as follows: where "ni" denotes the number of individuals belonging to a given species and "N" represents the overall number of individuals of all mosquito species in a commune.• The Simpson index (1 − D) evaluates the predominance of a given species within the population.It ranges from 0 (lower diversity with relatively few dominant species in the population) to 1 (higher diversity with several equally distributed species in the population).The calculation of this index was done as follows: where "ni" is the number of individuals of the given species and "N" is the total number of individuals.• The Pielou's equitability index (Jʹ) measures the evenness of species distribution in a population and derives from the Shannon diversity index (Hʹ).The Pielou's equitability index ranges from 0 (minimum evenness, suggesting a highly skewed distribution where one or a few species dominate the population) to 1 (maximum evenness which means that all species in the population have equal abundance).The formula to calculate this index is the following: where "Hʹ" is Shannon's diversity index and "Hʹmax" is the maximum possible value of Hʹ if each species has the same probability.
The diversity indices were calculated using Past version 3.14 software 33 .They were compared between communes using the Tukey test (package Agricolae).This test was performed using R software version 4.3.0,with the statistical significance threshold set at 0.05.
The mean density of mosquitoes per season was calculated for An.gambiae s.l, Cx. quinquefasciatus and Ma.africana at the cluster level.This was compared between seasons (dry and wet) using a mixed effect generalized linear model with a negative binomial distribution.Collection rounds were included in the model as random effects.The seasons were included as a fixed effect.Chi-square test was used to analyse the proportions of sibling species identified in An. gambiae s.l.Stata 15.0 software (Stata Corp., College Station, TX, USA) was used for these analyses.
Ecological niche models were produced using MaxEnt version 3.4.4software 31 .The maps were produced using ArcGIS version 10.8 software 34 .

Ethics approval
The study protocol received ethical approval from both Benin's National Ethics Committee for Health Research (Reference N°30/MS/DC/SGM/DRFMT/CNERS/SA, Approval n°6 of 04/03/2019) and the Ethics Committee of the London School of Hygiene and Tropical Medicine (Approval 16237-1).All participants provided informed consent prior to their involvement.All collectors underwent mosquito capture training that will allow them to avoid bites during the collection.Health measures taken included yellow fever vaccination before the study and access to local healthcare facilities for treatment of confirmed malaria cases or other ailments with similar symptoms.
In Mansonia spp., Ma. africana was the most predominant, accounting for 90% (74,945 out of 83,241) of the total number of individuals of this genus.In Culex spp., Cx. quinquefasciatus and Cx.annulioris accounted for 85.9% (40,058 out of 46,641) and 11.1% (5185 out of 46,641) respectively.In Anopheles spp. the most frequent species mosquito species was An. gambiae s.l.(96.9%).

Seasonal abundance and diversity of mosquito species
The species diversity within mosquito populations was assessed in the 3 communes.The highest species richness was observed in Zangnanado with 24 mosquito species collected, compared to 20 species for both Covè and Ouinhi (Fig. 3).A total of 18 mosquito species were common to all three study communes (Fig. 3).Overall, lower species richness was observed during the dry season compared to the wet season (21 vs 25 respectively).
Overall, there is no evidence for a higher density of An. gambiae s.l in the dry season compared to the wet season in Covè ( 31 1).
The Shannon-Weaver index assessing mosquito diversity was 0.90-1.18,with the highest value recorded in Covè.Analysis of this index indicated a significant difference in mosquito population diversity between wet and dry seasons in each of the three study communes (p < 0.05) (Table 1).

Others species
Others species The Piélou equitability index ranged between 0.50 (Ouinhi) and 0.69 (Zangnanado) and was similar (p > 0.05) between dry and wet seasons in Covè and Ouinhi.However, in Zangnanado, the Piélou equitability index was significantly higher in the dry season as compared to the wet season (p = 0.0061) (Table 1).
The Simpson index varied between 0.49 (Ouinhi) and 0.59 (Covè), and was similar (p > 0.05) between the wet and dry seasons in each of the three communes (Table 1).

Ecological niche models for mosquito species
During the wet season, An. gambiae s.s., An. coluzzii, Cx. quinquefasciatus, and Ma.africana showed very high habitat suitability in the south of Covè, south-west and east of Zangnanado, and centre of Ouinhi (Fig. 4).However, the north of Covè and north-western of Zangnanado showed moderate to high suitability for An.gambiae  s.s.(Fig. 4A).The models also revealed high habitat suitability in the north of Zangnanado and south-west of Ouinhi during the same period for An.coluzzii, Cx. quinquefasciatus and Ma.africana (Fig. 4C,E,G, respectively).
In the dry season, very high habitat suitability was observed throughout the commune of Covè for An.gambiae s.s., An. coluzzii and Cx.quinquefasciatus (Fig. 4B,D,F).The centre of Ouinhi was a highly favourable habitat for An.coluzzii, Cx. quinquefasciatus and Ma.africana (Fig. 4D,F,H).The west and south-east of Zangnanado showed high habitat suitability for all four mosquito species (Fig. 4B,D,F,H).

Contribution of variables to mosquito habitat suitability in the Covè-Zangnanado-Ouinhi health zone
Figure 5 shows the impact of the different evaluated bioclimatic and environmental data on the habitat suitability of An. gambiae s.s, An. coluzzii, Cx. quinquefasciatus and Ma.africana.
For An. gambiae s.s, slope and isothermality were the most important predictors of habitat suitability (Fig. 5).Thus, habitat suitability for this mosquito species increases when slope decreases (wet and dry seasons), and isothermality increases (wet season).
Regarding An. coluzzii, habitat suitability was positively correlated with isothermality (wet season), and negatively correlated with precipitation seasonality (dry season).
For Cx. quinquefasciatus, the habitat was most suitable when the isothermality increased (wet season), and the precipitation seasonality decreased (dry season).
Ma. africana habitat suitability increased as the precipitation seasonality decreased in both wet and dry seasons (Fig. 5).

Discussion
The present study focussed on the diversity, spatio-temporal distribution and ecological preferences of mosquito populations from the communes of Covè, Zangnanado and Ouinhi, areas endemic for malaria and lymphatic filariasis 16,35 .Overall, all studied mosquito species were more abundant in the wet season, except for An.gambiae s.l.whose density was higher in the dry season.Regardless of the commune, An. coluzzii predominated over An. gambiae s.s. in all seasons (wet and dry).Mosquito species diversity was higher in the wet season.The habitat suitability of An. gambiae s.s., An. coluzzii, Ma. africana and Cx.quinquefasciatus were mainly dependent on slope and precipitations, with some variations between species and seasons.
Overall, a total of 25 mosquito species belonging to five different genera were collected.The Culicinae composed mainly of Mansonia spp., and Culex spp., were significantly more abundant than the Anophelinae made up largely of An. gambiae s.l.This trend is similar to that observed in several entomological monitoring trials carried out both in Benin 36,37 and in other countries 38,39 .Mosquitoes of the Aedes genus were collected in low proportions, which could be due to the collection methods used, as well as the diurnal activity of this mosquito genus.Indeed, previous studies have shown that the collection of Aedes mosquitoes was greater with BG Sentinel Traps [40][41][42] , with peaks of activity usually observed early in the morning and in the afternoon 43 .
Our results showed a seasonal variation in the density and diversity of mosquito species in each study commune.The high densities of Cx. quinquefasciatus and Ma.africana during the wet season, a period conducive to the development of mosquito breeding sites, are consistent with previous observations made in India, Kenya, and Cambodia [44][45][46] .The unusual observation of An. gambiae s.l. that revealed to be more abundant in the dry season, could be because of the presence of flowing rivers in the low altitude areas of the study communes.These rivers are suitable for irrigated crops such as rice growing and market gardening, which create suitable habitats www.nature.com/scientificreports/for the development of larvae of An. gambiae s.l. 47In addition, the withdrawal of the Ouémé and Zou rivers from their beds during the dry season created numerous breeding sites favourable to this mosquito species.A similar observation was made around the Sanaga River in Southern Cameroon 48 .The diversity of mosquito populations is a paramount indicator as it provides information on the risk of mosquito-borne diseases 49 .In the three study communes, this index was higher in the wet season than in the dry season.This suggests that the mosquito populations of the study communes had few dominant species.Indeed, Ma. africana, Cx. quinquefasciatus and An.gambiae s.l. were the three most abundant mosquito species as they accounted for 85% (98,169/115,550) of the total number of mosquitoes sampled during the wet season.This suggests that the breeding sites in these communes were more suitable for the proliferation of these three mosquito species.This high diversity could increase the risk of transmission of mosquito-borne diseases such as malaria, and lymphatic filariasis 49 .
The study of ecological niche patterns revealed crucial information about the seasonal distribution of mosquito species.During the wet season, An. gambiae s.s., An. coluzzii, Cx. quinquefasciatus and Ma.africana were mainly present in the south of Covè, the southwest and east of Zangnanado, as well as the centre of Ouinhi.This shows that these species share a similar ecological niche in the wet season.Indeed, these regions are low-altitude floodplains where water and vegetation coexist, thus forming semi-permanent breeding sites favourable to the development of An. coluzzii 50 , Cx. quinquefasciatus and Ma.africana.This concurs with the higher proportion of An. coluzzii compared to An. gambiae s.s.recorded in the area.The colonisation of An. gambiae s.s. in the north of Covè, and the northwest of Zangnanado (high altitude areas) during this period, could be due to the numerous temporary breeding sites created by rainwater, as previously reported by Minakawa et al. 51 in Kenya.
During the dry season, the entire commune of Covè, the centre of Ouinhi, as well as the west and southeast of Zangnanado (areas consisted of plains previously flooded and where the recession started) offered highly suitable habitats for all four species of mosquitoes due to the recession of the flood.The same observation was made around the Gambia 52 and Ghibe 53 rivers.From the wet season to the dry season, the observation made is that suitable habitats for mosquitoes seem to move from the northeast to the southwest.This could be due to the direction of the flow of the Ouémé River which is influenced by the slope.
According to the work of Accrombessi et al. 54 conducted in the study area, sites with a high density of An. gambiae s.l. were associated with a high prevalence of malaria, so it would be judicious to focus malaria vector control interventions in areas where ecological niche models revealed suitable habitats for this mosquito species.
Past studies have highlighted several environmental factors that influenced the distribution of different mosquito species 55,56 .In the present study, our models revealed that precipitation seasonality, slope, isothermality and elevation were particularly relevant variables in the distribution of An. gambiae s.s., An. coluzzii, Cx. quinquefasciatus and Ma.africana.Indeed, the slope was important in determining in the distribution of An. gambiae s.s., while the elevation significantly influenced the distribution of An. coluzzii, Cx. quinquefasciatus, and Ma.africana.These results corroborate those obtained with Cx. quinquefasciatus in the Galápagos 57 .
Although the present study provides useful information on the abundance and diversity of mosquito species, as well as the ecological factors that influence their distribution, it has some limitations such as: the lack of consideration of socio-economic factors (demography, level of urbanization, agricultural activity, etc.), and the availability of hosts in the analyses; the difficulty in generalizing the results obtained to other regions due to the variability in the ecology, the short study period (one year); the failure to take into account complex interactions between species (predation, parasitism, competition) in the calculation of diversity indices; and the reliance on one mosquito collection technique.Another major gap for the present study is also the lack of Plasmodium sporozoite infection data.

Conclusion
The present study is a significant contribution to the understanding of the ecology of the main mosquito species encountered in the communes of Covè, Ouinhi and Zangnanado in Southern Benin.Indeed, it shows that the abundance and diversity of these mosquito species were much more pronounced in the wet season than in the dry season, with ecological preferences varying depending on the season.These findings will enable the implementation of targeted control strategies that take into account both the ecological preferences of each mosquito species as well as their seasonality.

Figure 1 .
Figure 1.Map of the study area.The map was drawn by C.J.A. using ArcGIS Desktop 10.8.1 software, study and online data was provided by Noce et al. 22 .

Figure 2 .
Figure 2. Composition and relative abundance of mosquito species.

Figure 4 .
Figure 4. Habitat suitability models in wet and dry seasons for An.gambiae s.s.(A and B), An. coluzzii (C and D), Cx. quinquefasciatus (E and F) and Ma.africana (G and H) under current climatic conditions.

Figure 5 .
Figure 5. Contribution of variables in mosquito species models during the wet and dry seasons.

Table 1 .
Seasonal abundance of mosquito species per commune.Indicator values with different superscripts within the same commune for the same mosquito species are significantly different (p < 0.05).The mean density was expressed in terms of number of bites/person/night (b/p/n).

Table 2 .
Seasonal abundance of Anopheles gambiae complex species per commune per season.An, Anopheles; N, number of; %, proportion; CI, confidence interval.