Kévin Hidalgo

C Montazeau

D Siaussat

V Braman

Marie Trabalon

F Simard

D Renault

R K Dabiré

K Mouline

Kevin Hidalgo
email: kevin.hidalgo@gmx.com

Distinct Physiological, biochemical and morphometric adjustments in the malaria vectors Anopheles gambiae and An. coluzzii as means to survive to dry season conditions in Burkina Faso

Keywords: Aestivation, migration, adipokinetic hormone (AKH), cuticular hydrocarbons, metabolic rate

come

Distinct physiological, biochemical and morphometric adjustments in the malaria vectors Anopheles gambiae and A. coluzzii as means to survive dry season conditions in Burkina Faso

Introduction

In Western Africa, malaria is mainly transmitted by the infectious bite of two anthropophilic female mosquito species, Anopheles coluzzii and An. gambiae. Both species can occupy a large range of ecological niches, but their densities follow the pace of the availability of their larval breeding-sites, i.e. standing water. In dry savannah areas, these larval habitats dry-up at the onset of the dry season, leading to the fall of mosquito populations [START_REF] Dao | Signatures of aestivation and migration in Sahelian malaria mosquito populations[END_REF][START_REF] Mamai | Monitoring Dry Season Persistence of Anopheles gambiae s.l. Populations in a Contained Semi-Field System in Southwestern Burkina Faso, West Africa[END_REF]. Only the mosquito An. coluzzii can locally persist all year long, as it can also exploit permanent larval breeding-sites like dams, rice fields or river edges, which do not entirely disappear during the dry season [START_REF] Baldet | Malaria transmission in 1999 in the rice field area of the Kou Valley (Bama), (Burkina Faso)[END_REF][START_REF] Costantini | Living at the edge: biogeographic patterns of habitat segregation conform to speciation by niche expansion in Anopheles gambiae[END_REF][START_REF] Diabate | The role of agricultural use of insecticides in resistance to pyrethroids in Anopheles gambiae sl in Burkina Faso[END_REF][START_REF] Diabate | The spread of the Leu-Phe kdr mutation through Anopheles gambiae complex in Burkina Faso: genetic introgression and de novo phenomena[END_REF]. The wider range of habitats exploited by An. coluzzii for larval breeding (temporary drying-up at the onset of the dry season, or permanent) may be explained by its aptitude to express different adaptive phenotypes related to the local breeding-site characteristics (Gimonneau et al., 2012a, 2012b, Hidalgo et al., 2015a, 2016). Recent works suggested that An. coluzzii mosquitoes collected from temporary larval breeding-sites might be "strong aestivators" programmed to enter a dormant state following specific (yet unknown) stimuli at the onset of the dry season [START_REF] Hidalgo | Comparative physiological plasticity to desiccation in distinct populations of the malarial mosquito Anopheles coluzzii[END_REF][START_REF] Wilps | Hormonal regulation of carbohydrate metabolism in the blowfly Phormia terraenovae[END_REF]. However, when An. coluzzii mosquitoes exploit permanent breeding-sites, they are considered as "weak aestivators" entering a facultative dormant state depending, at least in part, on the breeding abilities of the areas they exploit. By contrast, the mosquito An. gambiae exploits temporary larval breeding-sites only. Using population dynamics studies, it has been proposed that these mosquitoes did not enter a dormant state, but rather display high dispersal abilities to rapidly colonize new favourable habitats at the onset of the dry season [START_REF] Dao | Signatures of aestivation and migration in Sahelian malaria mosquito populations[END_REF][START_REF] Huestis | Ecophysiology of Anopheles gambiae s.l.: Persistence in the Sahel[END_REF][START_REF] Mamai | Monitoring Dry Season Persistence of Anopheles gambiae s.l. Populations in a Contained Semi-Field System in Southwestern Burkina Faso, West Africa[END_REF]. Consequently, physiological traits, in particular those associated with desiccation resistance or dispersal abilities, may differ between the two species at the onset of the dry season [START_REF] Arcaz | Desiccation tolerance in Anopheles coluzzii: the effects of spiracle size and cuticular hydrocarbons[END_REF][START_REF] Hidalgo | Novel insights into the metabolic and biochemical underpinnings assisting dryseason survival in female malaria mosquitoes of the Anopheles gambiae complex[END_REF][START_REF] Mamai | Morphological changes in the spiracles of Anopheles gambiae s.l (Diptera) as a response to the dry season conditions in Burkina Faso (West Africa)[END_REF].

Aestivation is a well-known survival process in insects inhabiting arid areas, but the physiological mechanisms and regulative metabolism controlling this process remain controversial [START_REF] Denlinger | Mosquito Diapause[END_REF]. As suggested for insect overwintering diapause, metabolic rate depression should be one of the most significant physiological traits characterising aestivation [START_REF] Benoit | Water management by dormant insects: comparisons between dehydration resistance during summer aestivation and winter diapause[END_REF][START_REF] Storey | Metabolic rate depression and biochemical adaptation in anaerobiosis, hibernation and estivation[END_REF]. Decrease of the metabolic rate has already been observed to prevent water losses and increase fitness of several insect species inhabiting arid environments [START_REF] Chown | Respiratory water loss in insects[END_REF][START_REF] Hadley | Water relations of terrestrial arthropods[END_REF][START_REF] Rourke | Geographic and altitudinal variation in water balance and metabolic rate in a California grasshopper, Melanoplus sanguinipes[END_REF]. Yet, this pattern has not been observed in field-sampled An. coluzzii [START_REF] Huestis | Variation in metabolic rate of Anopheles gambiae and A. arabiensis in a Sahelian village[END_REF][START_REF] Huestis | Seasonal variation in metabolic rate, flight activity and body size of Anopheles gambiae in the Sahel[END_REF], even in A c c e p t e d m a n u s c r i p t populations assumed to be "strong aestivators". When field-sampled during the rainy season and at the onset of the dry season, the high inter-individual variability of mosquitoes' phenotypes, resulting from uncontrolled age, trophic and reproductive states, may have prevented the drawing of a clear seasonal phenotypic pattern [START_REF] Nespolo | Intrapopulational variation in the standard metabolic rate of insects: repeatability, thermal dependence and sensitivity (Q10) of oxygen consumption in a cricket[END_REF][START_REF] Rogowitz | Energy metabolism of eucalyptus-boring beetles at rest and during locomotion: gender makes a difference[END_REF].

In addition to changes in metabolic rate, decreased permeability of the cuticle via changes in hydrocarbon composition for instance, have already been linked to desiccation resistance in several arthropod species [START_REF] Blomquist | Introduction: History and overview of insect hydrocarbons. Insect Hydrocarb[END_REF][START_REF] Gibbs | Cuticular lipids and water balance[END_REF][START_REF] Stinziano | Rapid desiccation hardening changes the cuticular hydrocarbon profile of Drosophila melanogaster[END_REF]. In particular, amounts of long-chain hydrocarbons correlate to increased cuticle hydrophobicity, and we expect that such adjustments would occur in anopheline species at the onset of the dry season. In parallel, the catabolism of specific organic compounds, like stored glycogen which generates significantly higher amounts of bounded water than any other energetic substrates [START_REF] Schmidt-Nielsen | Animal Physiology: Adaptation and Environment[END_REF], can contribute to reduce body water losses. The catabolism of glycogen by the glycogen phosphorylase is under the control of the adipokinetic hormones (AKH) [START_REF] Wilps | Hormonal regulation of carbohydrate metabolism in the blowfly Phormia terraenovae[END_REF][START_REF] Gäde | Regulation of Intermediary Metabolism and Water Balance of Insects by Neuropeptides[END_REF][START_REF] Ziegler | The putative AKH receptor of the tobacco hornworm, Manduca sexta, and its expression[END_REF], whose variations may signal changes in metabolism activity at the entry of the winter for dormant insects [START_REF] Hahn | Energetics of insect diapause[END_REF]. Similarly, variations of AKH levels may depict changes in desiccation resistance over the year in anopheline species, and increased amounts of this peptide can be expected in An. coluzzii at the onset of the dry season [START_REF] Hidalgo | Comparative physiological plasticity to desiccation in distinct populations of the malarial mosquito Anopheles coluzzii[END_REF].

AKH is also particularly important for flight maintenance, as it supplies the energetic substrates needed for striated muscle contraction in insects [START_REF] Arrese | Insect fat body: energy, metabolism, and regulation[END_REF][START_REF] Van Der Horst | Insect adipokinetic hormones: release and integration of flight energy metabolism[END_REF]. Thus, an increased expression of AKH gene expression levels in female An. gambiae at the onset of the dry season could be part of the "dispersive" phenotype.

In the present study, we aimed to experimentally trigger and characterize the specific phenotypes associated with current assumptions about survival strategies of An. coluzzii and An. gambiae: (1) females from An. gambiae species are expressing a dispersive behaviour to escape the harsh conditions of the dry season, whereas An. coluzzii females only virtually disappear and aestivate, and (2) the intensity of the aestivation phenotype in An. coluzzii females depends on the larval breeding-site characteristics (permanent to temporary pattern).

To address these assumptions, female mosquitoes were experimentally exposed to the environmental conditions they experience during the rainy (RS) and the onset of the dry season (ODS) in Burkina-Faso (West Africa). Variations in their metabolic rates, water contents, allometric measures (i.e. body size and dry mass), AKH peptides mRNA expression levels A c c e p t e d m a n u s c r i p t (Anoga-AKH-I, Anoga-AKH-II and their putative receptor Anoga-AKH-R), and cuticular hydrocarbon fingerprints were measured and compared (i) between RS and ODS conditions for each mosquito population, and (ii) among populations. We assumed that (i) An. coluzzii, females would have enhanced abilities to prevent body water losses at ODS, (ii) a decreased metabolic rate (which would also help preventing body water losses) should be observed at ODS in the aestivating females of An. coluzzii, but not in An. gambiae, where we expected to measure an increase in CO2 release, (iii) differential AKH gene expression levels between populations and experimental conditions that would reveal different metabolic mobilisation and, hence, different survival strategies (aestivating vs. dispersive), and (iv) seasonal changes in hydrocarbon composition of the cuticle should be observed in both species, notably in aestivating An. coluzzii specimens. Finally, we assumed that distinct ecophysiological phenotypes between populations of An. coluzzii should be observed according to their site of collection (i.e. temporary or permanent breeding-sites areas).

Materials and methods

Mosquito populations and rearing conditions

Mosquito populations

Experiments were conducted using two mosquito populations of An. coluzzii, and one population of An. gambiae. The two populations of An. coluzzii were established in September-October 2012 from 50 gravid females collected from two localities of Burkina Faso (Bama and Soumousso). At Bama [11°23'N, 04°24'W], An. coluzzii is the prevailing mosquito species (Gimonneau et al., 2012a), and it is breeding all year long in rice fields (i.e. permanent population). By contrast, at Soumousso [11°01N, 04°02'W], population dynamics of An. coluzzii match the seasonal variations of local climatic conditions [START_REF] Dabiré | Anopheles funestus (Diptera: Culicidae) in a Humid Savannah Area of Western Burkina Faso: Bionomics, Insecticide Resistance Status, and Role in Malaria Transmission[END_REF]. The population of An. gambiae used in this study was established at the IRSS insectary in 2009 from females of temporary populations collected in the village of Soumousso. This population thrives in sympatry with populations of An. coluzzii.

The three populations were maintained under the same insectary conditions (27 ºC, 70% humidity, light/dark cycles of 12/12h) prior to being used for the analysis (approximatively 6 generations). The three populations were maintained in separate rearing rooms to avoid any crossbreeding.

Experimental rearing conditions

A c c e p t e d m a n u s c r i p t

The three mosquito populations were reared from eggs to adults into two climatic chambers (Binder KBF720 VWR international S.A.S), programmed to reproduce the natural daily fluctuations of rearing conditions experienced by mosquitoes during the rainy season (RS; N=1 climatic chamber) and at the onset of the dry season (ODS; N=1 climatic chamber).

Programmed daily temperature and relative humidity conditions mimicked field conditions, as earlier used by (Hidalgo et al., 2015b[START_REF] Hidalgo | Novel insights into the metabolic and biochemical underpinnings assisting dryseason survival in female malaria mosquitoes of the Anopheles gambiae complex[END_REF][START_REF] Hidalgo | Comparative physiological plasticity to desiccation in distinct populations of the malarial mosquito Anopheles coluzzii[END_REF][START_REF] Mamai | Metabolomic and Ecdysteroid Variations in Anopheles gambiae s.l. Mosquitoes Exposed to the Stressful Conditions of the Dry Season in Burkina Faso, West Africa[END_REF]. Moreover, variation of day length was also accounted by using hourly-recorded light intensity data during the RS and ODS conditions with a weather monitoring station (Weatherlink; Davis Instruments, Hayward, CA, U.S.A.) (Supplementary data 1). Water temperature in rearing trays was similar to air temperature [START_REF] Hidalgo | Novel insights into the metabolic and biochemical underpinnings assisting dryseason survival in female malaria mosquitoes of the Anopheles gambiae complex[END_REF][START_REF] Hidalgo | Comparative physiological plasticity to desiccation in distinct populations of the malarial mosquito Anopheles coluzzii[END_REF].

For each mosquito population, three independent batches of eggs were synchronously collected from more than 50 different caged females, and reared in RS or ODS conditions. The climatic conditions inside climatic chambers were switched between each batch of eggs to account for any potential "chamber effect". Upon emergence, only females were maintained within the climatic chambers in small cages at a density of 30 females per cage, until they reached four to six-day old. Emergent females were fed ad libitum with glucose (10% w/v) and water. Twenty-four hours before the experiments (i.e. when females were three to five days old), females were food-and water-deprived in order to avoid any bias in carbon dioxide (VC02) and water loss rate (VH2O) measurements. The first batch corresponded to mosquitoes used for VCO2, VH2O and allometric measurements. The second batch was used for analysing the mRNA expression of the AKH peptides. The third batch corresponded to mosquitoes used for cuticular hydrocarbons analyses.

Gas exchange, flight and allometric measurements

2.2.1.

Metabolic (VCO2) and water loss (VH2O) rates

The VCO2 and VH2O rates of females from the three populations were measured using a flow through respirometry and a CO2-H2O gas analyzer (Sable System, Li-7000 CO2/H2O infrared gas analyzer, Li-Cor-Biosciences, Lincoln, NE). For each population and experimental condition (i.e. RS; ODS), six to nine runs were performed on pools of four randomized females (so that the minimum CO2 identification threshold within the gas analyzer system was reached).

Four to six-day old females were placed into a 200 x 400 mm (L x l) glass chamber. The chamber was connected to the gas analyzer, and flushed with a constant flow of air at a rate of 200 mL/min. The temperature and relative humidity of the air entering the chamber were modulated using a Peltier effect temperature controller and a Dew point generator/relative A c c e p t e d m a n u s c r i p t humidity controller system, respectively. The programmed conditions corresponded to the hottest and driest period of the day for each season during which females are resting (from noon to four p.m.). Hence, the relative humidity and temperature of the flowing air were respectively 18% and 34 ºC for the ODS conditions, and 70% and 29 ºC for the RS conditions.

Each run lasted 100 min during which both the rates of CO2 and H2O releases (VCO2 and VH2O) were recorded. Following an acclimation period of 25 min, measurements of VCO2 and VH2O were performed for 15 min. The measurement sequence was repeated three times.

Control measures were performed in an identical parallel system without mosquitoes for 10 min before each measurement sequence (i.e. 3 x 10 min). Data were analyzed using the Expedata software (Sable System International V.1.0.1). Amounts of mosquitoes' VCO2 and VH2O were calculated by subtracting the control values to those measured with mosquitoes in the chamber.

After each run, viability of the females was checked before they were snap-frozen and dried for 3 days at 60 ºC, and weighed (Sartorius, 0.1mg accuracy). For each run, VCO2 and VH2O were averaged per 15 min period and expressed in µl/min/mg of mosquito dry mass to assess for mass variability between species and conditions [START_REF] Gray | Evidence from mosquitoes suggests that cyclic gas exchange and discontinuous gas exchange are two manifestations of a single respiratory pattern[END_REF][START_REF] Huestis | Variation in metabolic rate of Anopheles gambiae and A. arabiensis in a Sahelian village[END_REF][START_REF] Huestis | Seasonal variation in metabolic rate, flight activity and body size of Anopheles gambiae in the Sahel[END_REF].

Flight activity

The flight activity, which might influence the rate of CO2 and H2O released by females, was recorded during the VCO2 and VH2O analyzes. Females were recorded using a HD Camera (Sanyo, model). A total of 44 batches of four pooled female mosquitoes (i.e. 176 femalesrepresenting from six to nine batches by experimental condition and population) were analyzed using the ImageJ1.41.0 software.

Allometric measures

The females' body size was estimated using their wing length (WL). Immediately after VCO2 and VH2O running sequences, females were snap-frozen. The right wing was removed before females were dried and weighed. Wings were mounted on a microscope slide, and pictures were taken under a stereomicroscope (x 20, Leica DFC425). Wing length was measured from the allula to the wing tip as described by [START_REF] Charlwood | Biological variation in Anopheles darlingi root[END_REF]. The measures were performed using ImageJ1.41.0 software (Wayne Rasband, National Institute of Health, Rockville, MD, U.S.A.) to an accuracy of ±0.001 mm. Damaged wings were discarded from the analysis. WL was divided by the dry mass of mosquitoes to estimate the surface to volume A c c e p t e d m a n u s c r i p t ratio of the specimens. A total of 158 wings was measured (i.e. from 24 to 36 per experimental condition and mosquito population).

AKH-related genes expression variations

RNA extraction and cDNA synthesis

For each population and experimental condition (RS; ODS), three to four samples of ten randomly pooled four day-old female mosquitoes were analysed. Collected mosquitoes were immediately snap-frozen in liquid nitrogen, and conserved at -80 °C with RNAlater solution (Ambion, USA) until further processing.

Before RNA extraction, the RNAlater solution was removed, and the total RNAs of each sample were extracted using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) coupled with the RNeasy Kit (QIAGEN). Samples were thereafter treated with DNase I (Ambion, USA)

following the manufacturer's instructions. RNA amounts were quantified by spectrophotometry at 260 nm (Nanodrop2000, Thermo Scientific). Single-stranded cDNAs were then synthesized with Superscript II reverse transcriptase (Gibco BRL, Invitrogen) as described in Bigot et al.

(2012) and according to the manufacturer's instructions.

2.3.2.

Real-time quantitative PCRs

All real-time quantitative PCRs were conducted as described in [START_REF] Bigot | Peripheral regulation by ecdysteroids of olfactory responsiveness in male Egyptian cotton leaf worms, Spodoptera littoralis[END_REF]. Ten genes (Actine, Rps13, Rps7, Rpl5, h3a, Cytp450, Tubulin, hsp83, EGFR, and 18s) were tested as putative housekeeping genes, following a BestKeeper analysis [START_REF] Pfaffl | Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper -Excel-based tool using pair-wise correlations[END_REF]. The Rps13 was selected as the reference gene since its expression was stable in all samples, whatever the experimental condition or the anopheline population tested. Further, mRNA of the two AKHs genes Anoga-AKH-I and Anoga-AKH-II, and the AKH receptor gene Anoga-AKH-Receptor (Kaufmann and Brown, 2006) were examined. Specific primers (reverse and forward) for both housekeeping genes and target genes were designed using the Eprimer3 software (http://emboss.bioinformatics.nl/cgi-bin/emboss/eprimer3; Table 1).

Each PCR reaction was technically triplicated and consisted of 6 µL of absolute Blue SYBR Green Fluor (ROCHE Molecular Systems Inc., USA), 2 µL of cDNA (25 ng.µl -1 ), 0.5 µL of each reverse and forward primers (10 µM), and 3 µL of RNA free water. The qPCR program and conditions have been described in [START_REF] Bigot | Peripheral regulation by ecdysteroids of olfactory responsiveness in male Egyptian cotton leaf worms, Spodoptera littoralis[END_REF]. normalise candidate gene expression levels to the geometric mean of the reference gene level using the Q-Gene software [START_REF] Simon | Q-Gene: processing quantitative real-time RT-PCR data[END_REF].

Cuticular hydrocarbon fingerprints

For each mosquito population exposed to RS and ODS conditions, the external cuticle (i.e. epicuticular) hydrocarbons were extracted from six to eight samples of four randomly pooled four day-old females, so that the minimum samples' mass reached 1 mg. Samples were washed for 10 min three times in 100 µL of dichloromethane. The three elution solutions were pooled and dried up under nitrogen. Dried extracts were conserved at -20 °C until further processing. Upon analysis, samples were re-suspended in 50 µL of dichloromethane containing 60 ng of n-tetradecane (internal standard). Samples were analyzed by gas-chromatography (GC) on a HP 6850A chromatograph equipped with a HP-5MS fused silica column (30m x 0.25mm; film thickness 0.25mm; Chrompack, Les Ulis, France) and an ionization flame detector (FID).

We used the procedure described by [START_REF] Trabalon | Changes in cuticular compounds composition during the gregarious period and after dispersal of the young in Tegenaria atrica (Araneae, Agelenidae)[END_REF]. Briefly, the oven temperature increased from 100 °C to 200 °C at a pace of 5 °C/min, and from 200 °C to 320 °C at a pace of 3 °C/min. The composition of the total hydrocarbon extract of each sample was identified by GC-MS on a HP 1890A equipped with a HP-5MS column. Only clearly defined peaks were used for profile characterization and quantitative analyzes (surface greater than 0.1% of the total amount of chemical compounds detected for each sample). The compounds were identified by comparing their mass spectra and GC retention times with those of pure standards. The results were expressed as mass/µL extract injected (ng/µL). Quantification was done by GC-FID analysis using internal and external standards:fatty acids (C16:0, C18:0, C18:1 and C18:2), methyl esters (palmitate, stearate) and n-alkanes (C13 to C34), Sigma ® .

Data analysis

All statistical procedures were conducted using the R 3.1.1 statistical software (R Development Core Team, 2008). Before analysis, the normal data distribution and homoscedasticity of variables were verified using Shapiro-Wilk and Bartlett's tests, respectively. Accordingly, analyzes were further processed using parametric or non-parametric tests. For all analyzes, main effects and all relevant first and second order interactions were tested in full models. Model simplification used stepwise removal of terms, where the significance of terms was estimated using the difference in Akaike's information criterion (AIC). When needed, a Tukey HSD procedure was used to perform post-hoc comparisons among the levels of significant factors.

A c c e p t e d m a n u s c r i p t

A two-way ANOVA was first performed to test whether allometric measures (dry mass, wing length and surface to volume ratio) and the flight activity of females varied according to the rearing conditions (RS, ODS), or among their population origins.

Generalised linear models (GLM, quasi-poisson error and logit function) were used to investigate the variations of metabolic rate (VCO2; expressed as µL/min/mg of dry mass), including the rearing conditions, the anopheline populations, the flight activity level, the VH2O, the wing size and the surface to volume ratio as explanatory variables. A second GLM was used to investigate the amount of water loss in females (VH2O, expressed as µL/min/mg of dry mass) with the same above-mentioned variables (except that the VH2O was replaced by VCO2).

Because of a significant effect of the surface to volume ratio on VH2O (ddl=1, χ 2 =91.04, P<0.001), a linear regression model was plotted to observe the direction and amplitude of this correlation. This relationship significantly differed according to the rearing conditions (RS or ODS) the mosquitoes were exposed to (ddl=2, χ 2 =5.91, P<0.05), but was similar among their population origins (ddl=1, χ 2 =1.0, P=0.59). Values of the three populations were then concatenated, and linear regressions between surface to volume ratio and VH2O were plotted distinctly for RS and ODS conditions.

The expression levels of the genes Anoga-AKH-I, Anoga-AKH-II and Anoga-AKH-R were analyzed using two way ANOVA procedures, where the experimental condition and the mosquito population were used as explanatory variables.

The cuticular hydrocarbon concentrations were first log-transformed (x = log10 [X + 1])

to meet the assumption of normally distributed residuals. An ANOVA was performed to test for the effect of mosquito population and experimental conditions on cuticular hydrocarbon fingerprints. In addition, a MANOVA was performed to examine qualitative differences in cuticular hydrocarbon profiles among anopheline populations and experimental conditions.

Two-way ANOVAs (where the experimental conditions and population of origin were considered as explanatory variables) were conducted for each individual hydrocarbon. A linear discriminant analysis (LDA) was performed to further assess qualitative and quantitative differences among the mosquito populations and rearing conditions. The between-and withingroup(s) degrees of freedom together with the F-value were reported for each LDA axis. The distribution's significance of the LDA was assessed using a Monte Carlo test with 10,000 permutations (P < 0.001).

Results

A c c e p t e d m a n u s c r i p t

3.1. Gas exchange, flight and allometric measurements 3.1.1. Metabolic (VCO2) and water loss (VH2O) rates

Overall, the three populations displayed equivalent VCO2 (GLM, ddl=2, χ 2 = 0.15, P=0.46). However, this parameter was significantly influenced by the rearing conditions (GLM, ddl=1, χ 2 = 7.45, P<0.01; Fig. 1). The difference was particularly striking for An. coluzzii mosquitoes from the temporary population exposed to ODS conditions: females were characterized by an important increase of their VCO2 (Fig. 1). A similar trend, yet nonsignificant, was observed for the permanent population of An. coluzzii (GLM, ddl=2, χ 2 = 5.11

P=0.07).

Rates of VH2O significantly differed among the three populations, and according to the environmental conditions (Table 2). VH2O of all populations significantly decreased when mosquitoes were exposed to ODS conditions (Table 2; Fig. 1). The decrease of VH2O from RS to ODS was higher in females of An. gambiae (Fig. 1).

Flight activity

Flight activity was similar among the three populations (ANOVA, ddl=2, F=0.33, P=0.72). An increase of the flight activity was observed for the three populations subjected to ODS conditions, but this trend was not significant (ANOVA, ddl=1, F=2.79, P=0.11, Fig. 2).

Allometric measures

Analyses showed that both wing size and female's dry mass significantly differed among the three populations (ANOVA, ddl=2, Fsize=12.91, Fmass=16.02, P<0.001; Fig. 3A-B).

Only dry mass was significantly influenced by rearing conditions (ANOVA, ddl=1, F=5.4, P<0.05; Fig. 3B), and this effect was mainly due to the significant increase of the dry mass of An. gambiae females reared in ODS conditions compared to those reared in RS ones (Fig. 3B).

Overall, the surface to volume ratio of females differed among populations (ANOVA, ddl=1, F=13.0, P<0.001, Fig. 3C) but not between rearing conditions (ANOVA, ddl=1, F=0.26,

P=0.61).

A significant interaction between environmental conditions and populations was observed (ANOVA, ddl=2, F=4.8, P<0.05), suggesting that the impact of the environmental conditions on the surface to volume ratio differed among the populations. Indeed, only females of An. gambiae exhibited significant phenotypic modifications, with the ones reared under ODS conditions exhibiting a 1.8 fold lower surface to volume ratio (Fig. 3C). This decrease is due to the significant increase of dry mass in the corresponding females (see above, Fig. 3B). VCO2 was not influenced by wing size, flight activity, nor the surface to volume ratio (GLM, ddl=1, χ 2 = 1.62, P=0.20; χ 2 = 0.47, P=0.49; χ 2 = 0.08, P=0.49, respectively).

VH2O was influenced by the variation of the surface to volume ratio of the mosquitoes, but not by the variation of their wing size and flight activity (Table 2). For all populations and for both rearing conditions, VH2O increased linearly with the surface to volume ratio of the females (Fig. 4). This linear relationship was different in mosquitoes exposed to RS than in those exposed to ODS (Fig. 4), as supported by the significant interaction terms between surface to volume ratio and rearing conditions (Table 2).

AKH-related genes expression variations

Expression levels of the Anoga-AKH-R and Anoga-AKH-I genes differed among the populations (ANOVA, ddl=2, FAKH-R=5.95, FAKH-I=2.80, P<0.05), but not those of Anoga-AKH-II (ANOVA, ddl=2, F=2.11, P=0.16). Expression levels of the three genes were not significantly influenced by the rearing conditions, while significant interaction term between rearing conditions and population was found (ANOVA, ddl=2, Anoga-AKH-Receptor:

F=15.76, P<0.001; Anoga-AKH-I: F=8.33, P<0.01; Anoga-AKH-II: F=5.53, P<0.05).

Hence, Anoga-AKH-R and Anoga-AKH-I genes displayed similar patterns of expression, with a significant increase of their mRNA amounts between RS and ODS conditions for the permanent population of An. coluzzii. Conversely, a decrease was observed for An. coluzzii from temporary sites (Fig. 5A-B). A distinct pattern was reported for the Anoga-AKH-II gene in which mRNA levels significantly increased for both populations from temporary sites, whereas a decrease (although marginally significant) was observed for An. coluzzii from permanent breeding sites (Fig. 5C).

Cuticular hydrocarbon fingerprints

Thirteen compounds were identified as mosquito cuticular hydrocarbons. Total amount of cuticular hydrocarbons significantly differed among the three mosquito populations (ANOVA, ddl=2, F=31.81, P<0.001), but there were no significant effects of the rearing conditions (ANOVA, ddl=1, F=0.78, P=0.38; Fig. 6). In particular, females An. gambiae expressed 8.75to 15.30-fold more cuticular hydrocarbon amounts than An. coluzzii populations under RS and ODS conditions, respectively. No difference between the two An. coluzzii populations was observed (Fig. 6).

A c c e p t e d m a n u s c r i p t

Distinct cuticular hydrocarbon profiles were measured among the three populations (MANOVA, F2,26=6.65, P<0.001), the two rearing conditions (MANOVA, F1,13=17.43, P<0.01), and the interaction between these two terms (MANOVA, F2,26=3.80, P<0.01). In the LDA analysis, the first axis (LD1) accounted for 25.99% of the total inertia, and the variation among groups was 19.31 times higher than the variation within groups (Fig. 7). LD1 mainly separated female An. coluzzii from permanent breeding sites reared under ODS conditions from all other groups. This clear cut-off was mainly characterized by the increased amounts of 11tetracosane and 12-tetracosane in ODS-reared females of this permanent An. coluzzii population, and 11-methyltricosane, 13-methylhexacosane, 13-methylpentacosane, 3methyltricosane, 3-methylpentacosane, cholesterol, n-hentriacontane, n-pentacosane, ntetracosane and n-tricosane in the other groups (Fig. 7 and Supplementary data 2). The second axis (LD2) accounted for 22.43% of the total inertia, and the variation among groups was 5.06 times higher than the variation within groups. LD2 mainly separated female An. coluzzii from temporary breeding sites reared under ODS conditions from all other groups (Fig. 9). Increased amounts of n-nonacosane were measured in ODS-reared temporary population of An. coluzzii, and n-heptacosane in RS-reared temporary and permanent An. coluzzii (Fig. 7 and Supplementary data 2).

Discussion

The present work aimed at characterizing and comparing the influence of environmental conditions found in sub-Saharan Africa at the onset of the dry season (ODS), and at the rainy season (RS), on phenotypic adjustments displayed by females of An. gambiae and An. coluzzii sampled from localities where larval breeding-sites are permanent or temporary. Our results

showed that all measured phenotypic traits displayed significant changes depending on the environmental climatic conditions (i.e. temperature, RH and photoperiod) experienced by the mosquitoes during their development. Phenotypic adjustments also varied according to the considered population, and the larval breeding-site dynamics (temporary vs. permanent).

Seasonal metabolic rate variations are not consistent with the idea of an aestivation strategy in An. coluzzii

According to studies lead on dormant insect species, and more particularly on those that overwinter as diapausing organisms [START_REF] Hahn | Meeting the energetic demands of insect diapause: Nutrient storage and utilization[END_REF][START_REF] Kambule | Metabolic rate in diapause and nondiapause brown locust eggs correlated with embryonic development[END_REF], we predicted that metabolic rate would have been lowered in An. coluzzii reared under ODS Conversely, metabolic rate should have been increased in female An. gambiae, as they are supposed to exhibit a dispersive flight strategy. Surprisingly, VCO2 emitted by An. gambiae did not show any significant variation from RS to ODS conditions, and significantly increased, or tended to do so, in temporary and permanent populations of An. coluzzii, respectively. The physiological coercions that insects have to face during aestivation may partly differ from those faced during overwintering diapause, as revealed by the unexpected VCO2 fingerprints we are reporting with this work. [START_REF] Denlinger | Mosquito Diapause[END_REF] pointed out that resting metabolic rate of malarial mosquitoes is not lower, but rather exhibits a very complex pattern of change during the course of the dry season, depending on the temperature and on the female's reproductive status. The high temperatures encountered during aestivation may impose unique physiological constraints [START_REF] Denlinger | Mosquito Diapause[END_REF], and the sensitivity of aerobic metabolism to temperature changes has been evidenced in earlier studies [START_REF] Clarke | Seasonal Acclimatization and Latitudinal Compensation in Metabolism: Do They Exist?[END_REF][START_REF] Huestis | Seasonal variation in metabolic rate, flight activity and body size of Anopheles gambiae in the Sahel[END_REF][START_REF] Terblanche | The effects of acclimation on thermal tolerance, desiccation resistance and metabolic rate in< i> Chirodica chalcoptera</i>(Coleoptera: Chrysomelidae)[END_REF], including in field-sampled mosquitoes [START_REF] Huestis | Variation in metabolic rate of Anopheles gambiae and A. arabiensis in a Sahelian village[END_REF][START_REF] Huestis | Seasonal variation in metabolic rate, flight activity and body size of Anopheles gambiae in the Sahel[END_REF]. The extreme temperatures and low humidity levels encountered during the dry season seem to promote an increased metabolic activity of female mosquitoes, maybe underlining the increased demand of energetic substrates (i.e. ATP) to flee the desiccating conditions (migration) or resist locally (aestivation). The increase of flight activity -only marginally significantof the An. coluzzii females at ODS may have also contributed to increase VCO2 values.

By contrast, female An. gambiae did not show significant variation of their metabolic rate from RS to ODS conditions, underlying that they express lower plasticity levels than An.

coluzzii and put in place distinct physiological mechanisms to cope with hot and dry conditions, as we initially suggested.

Mosquitoes exhibit different ways to prevent body water loss at ODS

The increase in VCO2 observed at ODS in An. coluzzii females should participate to increase its desiccation level, as gas exchanges are supposed to increase water losses in insects.

Yet, our results did not support this assumption, as VCO2 variations observed in all mosquitoes at both RS and ODS conditions did not influence their rates of water losses (VH2O). Previous studies showed that insect's metabolic rate fluctuations do not necessarily lead to significant changes in water losses [START_REF] Bradley | Physiological Responses to Selection for Desiccation Resistance in Drosophila melanogaster[END_REF][START_REF] Chown | Respiratory water loss in insects[END_REF][START_REF] Djawdan | Does selection for stress resistance lower metabolic rate?[END_REF][START_REF] Rourke | Geographic and altitudinal variation in water balance and metabolic rate in a California grasshopper, Melanoplus sanguinipes[END_REF]Williams et al., 1998;Williams and Bradley, 1998). In particular, the control of water losses is mainly achieved by adjustments of cuticular permeability, as cuticular transpiration is

A c c e p t e d m a n u s c r i p t the principle way through which insects experience desiccation [START_REF] Benoit | Water management by dormant insects: comparisons between dehydration resistance during summer aestivation and winter diapause[END_REF]). Yet, mosquito's water losses significantly decreased in all populations at ODS, suggesting that mechanisms other than adjustments of the metabolic rate participate to reduce anopheline transpiration under such climatic conditions.

The reduction of the surface to volume ratio is described in several species as one of the mechanisms preventing body water losses [START_REF] Benoit | Suppression of water loss during adult diapause in the northern house mosquito, Culex pipiens[END_REF][START_REF] Hadley | Water relations of terrestrial arthropods[END_REF].

Overall, our data showed a clear linear relationship between the rate of VH2O in female An.

gambiae and their surface to volume ratio. Such seasonal allometric change in An. gambiae is due to an increase of females' dry mass at ODS, highlighting potential differences in their abilities to seasonally harvest, store and/or use their metabolic resources, as shown for Aedes albopictus [START_REF] Reiskind | Is bigger really bigger? Differential responses to temperature in measures of body size of the mosquito, Aedes albopictus[END_REF].

Variation levels of adipokinetic peptides (AKH) may also help to enhance desiccation resistance in anopheline mosquitoes, as shown in overwintering dormant insects [START_REF] Hahn | Energetics of insect diapause[END_REF]). Here, we measured the expression levels of two AKH peptide genes (Anoga-AKH-I and Anoga-AKH-II) and their putative receptor (Anoga-AKH-Receptor). Overall, we highlighted that both Anoga-AKH-I and Anoga-AKH-Receptor genes showed very similar patterns of expression, suggesting that their expression levels are regulated in parallel, as expected for receptor/peptide couples. By contrast, Anoga-AKH-Receptor expression levels did not match with those of the Anoga-AKH-II. Few years ago, Anoga-AKH-II was described to code for an intermediary AKH and corazonin protein named AKH/corazonin-related (ACP) (Kaufmann and Brown, 2006). Functional roles of corazonin remain obscure in insects, even if they seem to have cardio-acceleratory effects, and should thus increase metabolic rate, as observed in the American cockroach [START_REF] Sláma | Effect of corazonin and crustacean cardioactive peptide on heartbeat in the adult American cockroach (Periplaneta americana)[END_REF]. The protein corazonin is involved in migratory processes of gregarious locusts [START_REF] Tawfik | Identification of the gregarization-associated darkpigmentotropin in locusts through an albino mutant[END_REF], and is released under nutritional stress [START_REF] Veenstra | Does corazonin signal nutritional stress in insects? Insect[END_REF], thus corroborating the hypothesis of a dispersal strategy at ODS in An.

gambiae. In addition, variation of ACP mRNA levels correlates to those observed for glycogen phosphorylase mRNA in Anopheles species [START_REF] Hidalgo | Comparative physiological plasticity to desiccation in distinct populations of the malarial mosquito Anopheles coluzzii[END_REF], suggesting a possible physiological interaction between these two actors in An. gambiae. In particular, the increase in ACP expression in An. gambiae at ODS could help degrading large amounts of stored carbohydrates, thus providing non-negligible amounts of metabolic water and energetic substrates to sustain an increased dispersal [START_REF] Huestis | Ecophysiology of Anopheles gambiae s.l.: Persistence in the Sahel[END_REF][START_REF] Mamai | Monitoring Dry Season Persistence of Anopheles gambiae s.l. Populations in a Contained Semi-Field System in Southwestern Burkina Faso, West Africa[END_REF].

A c c e p t e d m a n u s c r i p t

By contrast, under ODS conditions, the AKH-I/AKH-Receptor couple was overexpressed in the permanent population of An. coluzzii, whereas it was under-expressed in the temporary one, and did not show any seasonal variation in An. gambiae. Expression levels of AKH-I mRNA characterized the distinct molecular response of permanent An. coluzzii populations, which continue to be active and reproduce at ODS. Oddly, these results are in contrast to those we previously obtained using 7-day-old females, in which no seasonal differences were observed with the same An. coluzzii populations [START_REF] Hidalgo | Comparative physiological plasticity to desiccation in distinct populations of the malarial mosquito Anopheles coluzzii[END_REF]. Such a difference may be explained by the younger age of our mosquitoes, in addition to taking in account photoperiod fluctuations in our experimental rearing. AKH peptides are suspected to mediate physiological responses to desiccation, as they allow the release and transport, from the fat body to the hemolymph, of organic compounds with osmoprotectant functions (i.e.

trehalose, proline), notably in overwintering insects [START_REF] Gäde | Regulation of Intermediary Metabolism and Water Balance of Insects by Neuropeptides[END_REF][START_REF] Hahn | Energetics of insect diapause[END_REF][START_REF] Isabel | AKH-producing neuroendocrine cell ablation decreases trehalose and induces behavioral changes in Drosophila[END_REF][START_REF] Wilps | Hormonal regulation of carbohydrate metabolism in the blowfly Phormia terraenovae[END_REF][START_REF] Ziegler | The putative AKH receptor of the tobacco hornworm, Manduca sexta, and its expression[END_REF]. Further studies are required to confirm this hypothesis, but An. coluzzii populations may require increased levels of osmoprotectants to deal with the desiccating conditions of the dry season while remaining reproductively active. Overall, seasonal variations of AKH-I and ACP expressions may represent valuable markers of the survival strategies and/or stress response mechanisms elicited in An. coluzzii and An. gambiae populations.

Seasonal changes in cuticular hydrocarbon fingerprints suggest distinct desiccation resistance strategies in anopheline mosquitoes

Changes in cuticular hydrocarbon composition are expected in anopheline mosquitoes reared under ODS conditions as cuticular transpiration is the main way of water loss in insects [START_REF] Chown | Respiratory water loss in insects[END_REF][START_REF] Chown | Insect Physiological Ecology: Mechanisms and Patterns[END_REF][START_REF] Hadley | Water relations of terrestrial arthropods[END_REF][START_REF] Johnson | Effect of mating stage on water balance, cuticular hydrocarbons and metabolism in the desert harvester ant, Pogonomyrmex barbatus[END_REF]. Using gas chromatography techniques, we showed quantitative differences in the total amount of cuticular hydrocarbons between the two sibling species, with An. gambiae females displaying the highest amounts, whatever the environmental conditions, while no difference was evidenced between the two An. coluzzii populations. Interestingly, total amount of cuticular hydrocarbons did not vary in mosquitoes according to the rearing conditions. In insects, desiccation resistance is not related to the total amount of cuticular hydrocarbons a species exhibits [START_REF] Arcaz | Desiccation tolerance in Anopheles coluzzii: the effects of spiracle size and cuticular hydrocarbons[END_REF][START_REF] Gibbs | Physiological mechanisms of evolved desiccation resistance in Drosophila melanogaster[END_REF][START_REF] Kwan | Adaptation to Desiccation Fails to Generate Pre-and Postmating Isolation in Replicate Drosophila Melanogaster Laboratory Populations[END_REF], but rather to their relative composition [START_REF] Arcaz | Desiccation tolerance in Anopheles coluzzii: the effects of spiracle size and cuticular hydrocarbons[END_REF][START_REF] Gibbs | Physiological mechanisms of evolved desiccation resistance in Drosophila melanogaster[END_REF][START_REF] Gibbs | Cuticular lipids and water balance[END_REF][START_REF] Nelson | Cuticular lipids and desiccation resistance in overwintering larvae of the goldenrod gall fly, Eurosta solidaginis (Diptera: Tephritidae)[END_REF].

We highlighted qualitative hydrocarbon differences between the two An. coluzzii populations on one hand, and An. gambiae on the other, whatever the season we considered.

A c c e p t e d m a n u s c r i p t

The profile of An. gambiae females is mainly characterized by higher amounts of (1) cholesterol, a precursor of insects' steroids, and (2) methyl-branched alkanes, known as important semiochemical cues for social and sexual insect recognition [START_REF] Howard | Ecological, Behavioral, and Biochemical Aspects of Insect Hydrocarbons[END_REF]. Such qualitative differences might be part of the biological basis leading to assortative mating between the species and, to a lesser extent, to distinct innate biology and ecology [START_REF] Niang | Does extreme asymmetric dominance promote hybridization between Anopheles coluzzii and Anopheles gambiae s.s. in seasonal malaria mosquito communities of West Africa?[END_REF].

Our data underlined distinct cuticular hydrocarbon adjustments between the two species in response to seasonal changes. While cuticular hydrocarbons did not change with environmental conditions in An. gambiae females, seasonal differences were observed in An.

coluzzii. This finding gives credit to our earlier observations, which suggested the occurrence of seasonal changes in cuticular composition (i.e. aromatic amino acids and RR2-proteins) at ODS in these females but not in An. gambiae [START_REF] Hidalgo | Novel insights into the metabolic and biochemical underpinnings assisting dryseason survival in female malaria mosquitoes of the Anopheles gambiae complex[END_REF].

Lastly, distinct seasonal adjustments of cuticular hydrocarbons were observed between

An. coluzzii exploiting temporary or permanent larval breeding-sites. Temporary population of An. coluzzii displayed an increased amount of n-nonacosane (C29H60) at ODS, a compound already observed to increase under desiccating winter in the scorpion Centruroides sculpturatus [START_REF] Toolson | Seasonal effects on cuticular permeability and epicuticular lipid composition inCentruroides sculpturatus ewing 1928 (Scorpiones: Buthidae)[END_REF]. Conversely, permanent population of this species showed increased amounts of the 11-and 12-tetracosane (C24H50). The three compounds melt at hot temperatures, but n-nonacosane has a higher melting point (62 to 66 °C) than n-tetracosane (48 to 54 °C). Insect literature suggests that diapausing species over-express long chain hydrocarbons [START_REF] Benoit | Water management by dormant insects: comparisons between dehydration resistance during summer aestivation and winter diapause[END_REF][START_REF] Kankare | Transcriptional Differences between Diapausing and Non-Diapausing D. montana Females Reared under the Same Photoperiod and Temperature[END_REF], because it improves the cuticle impermeability avoiding desiccation [START_REF] Chung | Wax, sex and the origin of species: Dual roles of insect cuticular hydrocarbons in adaptation and mating[END_REF][START_REF] Gibbs | Thermodynamics of cuticular transpiration[END_REF]. Thus, the increased amount of longer hydrocarbon chains in An. coluzzii females collected from temporary breeding-sites (n-nonacosane) could support the "strong" aestivation strategy suggested for these mosquitoes during the dry season. Interestingly, n-tetracosane accumulated in temporary population of An. coluzzii at ODS is known to increase during the reproductive period of Dipterans [START_REF] Jurenka | Hydrocarbon profiles of diapausing and reproductive adult face flies (Musca autumnalis)[END_REF] supporting the expectation of a maintained reproductive activity in these mosquitoes at ODS [START_REF] Baldet | Malaria transmission in 1999 in the rice field area of the Kou Valley (Bama), (Burkina Faso)[END_REF][START_REF] Costantini | Living at the edge: biogeographic patterns of habitat segregation conform to speciation by niche expansion in Anopheles gambiae[END_REF][START_REF] Diabate | The role of agricultural use of insecticides in resistance to pyrethroids in Anopheles gambiae sl in Burkina Faso[END_REF][START_REF] Diabate | The spread of the Leu-Phe kdr mutation through Anopheles gambiae complex in Burkina Faso: genetic introgression and de novo phenomena[END_REF]. Further studies are still needed to tear apart the exact roles of these candidate hydrocarbons, but, at least, their distinct profiles between conditions and populations consolidate the expectation of distinct survival strategies in An. coluzzii as a function of the ecology of their larval breeding-sites.

Conclusion

A c c e p t e d m a n u s c r i p t

Understanding the survival strategies put in place by malaria mosquitoes and other tropical insects to cope with the dry season conditions is a major challenge [START_REF] Denlinger | Mosquito Diapause[END_REF][START_REF] Wilps | Hormonal regulation of carbohydrate metabolism in the blowfly Phormia terraenovae[END_REF]. Our study underlined distinct ecophysiological adjustments in An. coluzzii and An. gambiae females to prevent water losses under ODS conditions. These adjustments were also driven by the type of breeding site exploited by the mosquito population (i.e. permanent or temporary larval ecotype). In addition, specific changes in epicuticle hydrocarbons and AKH mRNA expression levels between mosquitoes highlights the high phenotypic diversity in An. gambiae s.l. species complex. These results are in line with the hypothesis of distinct survival strategies among these species [START_REF] Wilps | Hormonal regulation of carbohydrate metabolism in the blowfly Phormia terraenovae[END_REF]. Further investigations are required to examine the exact roles of these markers in the mosquito metabolism and biology, but we believe they could constitute relevant markers to highlight specific dry season strategies developed by "permanent" and "temporary" anopheline populations during the dry season. Altogether, our results bring better understanding of mosquito's biology during the dry season, opening the way for alternative method of malaria control in sub-Saharan Africa. Different letters represent significant differences between the experimental modalities at P value < 0.05 (ANOVA). n= 32, 24 and 36 samples (one female per sample) for RS conditions, and n= 28, 32 and 24 samples (one female) for ODS conditions in, respectively, the permanent and temporary populations of An. coluzzii and An. gambiae.

Figure captions

  The cycle threshold values (Ct-values) for both reference and target genes were determined using the Light-Cycler® 480 software (Roche, FRANCE). The average Ct value of each technical triplicate was used to A c c e p t e d m a n u s c r i p t

  Can flight activity and allometry influence females' VCO2 and VH2O?

  these mosquitoes are assumed to enter a dormant state named aestivation.

Fig. 1

 1 Fig.1Mean ± s.e.m. rate of CO2 release (in µl/min/mg of mosquito dry mass) (A) and H2O loss (in µl/min/mg of mosquito dry mass) (B) by the three anopheline populations under the rainy season (RS, white bar) and the onset of the dry season (ODS, grey bar) environmental conditions. Different letters represent significant differences between experimental modalities at P value<0.05 (GLM). n= 8, 6 and 9 samples of four pooled females for RS, and n= 7, 8 and 6 samples of four pooled females for ODS in, respectively, the permanent and temporary populations of An. coluzzii, and An. gambiae.

Fig. 2

 2 Fig.2 Mean ± s.e.m. rate of flight activity frequency (in second/minute) recorded in the three anopheline populations under the environmental conditions of the rainy season (RS, white bar) and the onset of the dry season (ODS, grey bar). Different letters represent significant differences between experimental modalities at P value < 0.05 (ANOVA). Samples are the same as those in figure 1.

Fig. 3

 3 Fig.3 Mean ± s.e.m. of A: wing size (in mm), B: Dry mass (in mg), and C: surface to volume ratio (in mm/mg of dry mass) of the three anopheline populations under the environmental conditions of the rainy season (RS, white bar) and the onset of the dry season (ODS, grey bar).

Fig. 4

 4 Fig.4 Linear regression between the females' surface to volume ratio (in mm/mg of dry mass) and H2O losses (in µl/min/mg of mosquito dry mass) when mosquitoes are reared under the environmental conditions of the rainy season (RS; filled line), or the onset of the dry season (ODS; dashed line). Samples are the same as those described in Fig. 1 and 3.

Fig. 5

 5 Fig.5Relative normalised mRNA expression level ± s.e. of A: Anoga-AKH-receptor; B: Anoga-AKH-I; and C: Anoga-AKH-II in the three anopheline populations reared under the environmental conditions of the rainy season (RS, white bar) and the onset of the dry season (ODS, grey bar). Different letters represent significant differences between experimental modalities at P value < 0.05 (ANOVA). n= 3, 3 and 4 samples of ten pooled females for RS; n= 3, 3 and 3 samples of ten pooled females for ODS in, respectively, the permanent and temporary populations of An. coluzzii and An. gambiae.

Fig. 6 Fig. 7

 67 Fig.6 Mean ± s.e.m. of the total amount of cuticular hydrocarbon (in ng/mosquito) in the three anopheline populations reared under environmental conditions of the rainy season (RS, white bar) and the onset of the dry season (ODS, grey bar). Different letters represent significant differences between experimental modalities at P value < 0.05 (ANOVA). n= 4, 4 and 3 samples of four pooled females for RS; n= 5, 7 and 5 samples of four pooled females for ODS

Table 1 .

 1 Nucleotide sequences of the primers used in qRT-PCR reactions for the amplification of Actin, Rps13, Rps7, Rpl5, h3a, CytP450, Tubulin, hsp83, EGFR, 18s, Anoga-AKH-I, Anoga-AKH-II, and Anoga-AKH-R in An. gambiae.

Table 2 .

 2 Results of the GLM computed on the VH2O with the experimental conditions, the anopheline populations, VCO2, flight activity and surface to volume ratio as explicative variables. Signif. code: < 0.001 "***"; <0.05 "*".

Acknowledgements

The authors would like to thank Mamai Wadaka and Sougrinoma Zoungrana (IRSS, Bobo-Dioulasso) for their help in mosquito rearing and sample preparation.

Funding

This research was funded by the Agence Nationale de la Recherche through grant ANR-08-MIEN-006 to FS.


Philippe Bich

Nash equilibrium existence for some discontinuous games

Answering to an open question of Herings et al. (see [3]), one extends their fixed point theorem to mappings defined on convex compact subset of R n , and not only polytopes. Such extension is important in non-cooperative game theory, where typical strategy sets are convex and compact. An application in game theory is given.

Introduction

In [START_REF] Herings | A fixed point theorem for discontinuous functions[END_REF], Herings et al. prove the following new fixed point theorem for possibly discontinuous mappings: Theorem 1.1 Let P a non empty polytope, i.e. the convex hull of a finite subset of R n ; let f : P → P which is "locally gross direction preserving" in the following sense: for every x ∈ P such that f (x) = x, there exists V x , an open neighborhood of x in P such that for every u and v in V x , f (u) -u, f (v) -v ≥ 0.

Then f admits a fixed-point, i.e. there exists x ∈ P such that f (x) = x.

Here, and throughout this paper, for every x ∈ R n and y ∈ R n , x, y denotes the eulidean scalar product of x and y.

This theorem is a generalization of Brouwer fixed point theorem (see [START_REF] Brouwer | Über Abbildung von Mannigfaltigkeiten[END_REF]) which says that every continuous mapping from the unit closed ball of R n to itself admits a fixed point.

The authors provide an application to the existence of equilibria for games with possibly discontinuous payoffs. There is a huge literature on discontinuous games (see, for example, [START_REF] Reny | On the existence of pure and mixed strategy Nash equilibria in discontinuous games[END_REF]); a usual issue in such games is that classical fixed point theorems as Brouwer's one cannot be directly applied to yield Nash equilibria existence. Thus, Theorem 1.1 could be a answer to this problem for some classes of discontinuous games.

Yet, there is a restriction in Theorem 1.1: the set P must be a polytope. But typical strategy sets in game theory are rather compact and convex sets. Thus, an important question in practice is to know if Theorem 1.1 holds true for such subset of R n . In [START_REF] Herings | A fixed point theorem for discontinuous functions[END_REF], one can read: "whether locally gross direction preserving is sufficient to guarantee the existence of a fixed point on an arbitrary non empy convex and compact set is still an open question".

The aim ot this paper is to answer to this question. A first natural idea to extend Theorem 1.1 to the case where P is any compact and convex subset of R n is the following: first, approximate P by a sequence of polyhedra P n ⊂ P ; then, Theorem 1.1 applied to the mappings proj Pn • f |Pn provides a sequence of fixed points x n of these mappings; finally, one could hope that the sequence (x n ) converges (up to an extraction) to a fixed point of f . But, as the authors say: "... the discontinuities of f on the boundary of P prevent us from taking the limit of the sequence of polyhedra. So to resolve this problem, a different approach is needed". In this paper, we provide such a different approach.

As a matter of fact, our approach drives us to generalize strictly the "locally gross direction preserving" property. Thus, we extend Herings et al. in two ways: we replace polytopes by any non empty and convex subet of R n , and extend their continuity assumption. Finally, we prove that our fixed point theorem can be applied to yield the existence of Nash equilibria in possibly discontinuous games.

The main theorem

The following theorem extend Herings et al. result to compact and convex subsets of R n . Theorem 2.1 Let C be a non empty convex and compact subset of R n , and let f : C → C. Assume that f is "locally gross direction preserving" in the following sense: for every x ∈ C such that f (x) = x, there exists V x , an open neighborhood of x in C such that for every u and

v in V x , f (u) -u, f (v) -v ≥ 0. Then, there exists x ∈ C such that f (x) = x.
Proof. Throughout this paper, for every x ∈ R n and r ∈ R + , B(x, r) denotes the open ball centered at x of radius r, for the euclidean norm of R n .

• Step 1: one proves that any locally gross direction preserving mapping satisfy the following property that we call half-continuity:

Definition 2.2 A mapping f : C → C is said to be half-continuous if: ∀x ∈ C, x = f (x) ⇒ ∃p ∈ R n , ∃ > 0 such that: ∀x ∈ B(x, ) ∩ C, x = f (x ) ⇒ p, f (x ) -x > 0.
To prove Step 1, suppose that f : C → C is locally gross direction preserving, and prove that it is half-continuous. Let x ∈ C such that f (x) = x, and let V x be an open neighborhood of x such that for every u and

v in V x , f (u) -u, f (v) -v ≥ 0. Let {f (x 1 )-x 1 , ..., f (x k )-x k } be a basis of the vector space F :=span{f (y)- y, y ∈ V x }, where k ∈ N * , and x 1 , ...,x k are in V x . Then define p = k i=1 (f (x i ) -x i ). Let x ∈ V x such that f (x ) = x . One clearly have p, f (x ) -x ≥ 0 (1)
from "locally gross direction preserving" property and from the definition of p.

Besides, since one has f (x i ) -x i , f (x ) -x ≥ 0 for every i = 1, .., k, Inequation 1 is an equality if and only if one has

∀i = 1, ..., k, f (x ) -x , f (x i ) -x i = 0.
This last property would imply f (x ) -x ∈ F ⊥ ∩ F = {0}, a contradiction with the assumption that f (x ) = x . Thus, Inequality 1 is strict, and Step 1 if proved.

• Step 2: supposing that f : C → C has no fixed point and is halfcontinuous, one builds a continuous mapping p : C → R n such that for every x ∈ C, p(x), f (x) -x > 0.

Suppose that f admits no fixed point and is half-continuous. For every x ∈ C, and from half-continuity, there exists p x ∈ R n and x > 0 such that

for all x ∈ B(x, x ) ∩ C, p x , f (x ) -x > 0. One has C ⊂ ∪ x∈C B(x, ). Since C is compact, there exists x 1 , ..., x n in C such that C ⊂ ∪ n i=1 B(x i , x i ). Consider λ 1 , ...,
λ i (x) = 1 and such that ∀i = 1, ..., n, ∀x ∈ C, x / ∈ B(x i , x i ) ⇒ λ i (x ) = 0.
Define the mapping p :

C → R n by ∀x ∈ C, p(x) = n i=1 λ i (x)p x i .
Clearly, from the the properties of the λ i , one has

∀x ∈ C, p(x), f (x) -x > 0 and
Step 2 is proved.

• Step 3: Proof of Theorem 2.1
Suppose that f : C → C, a locally gross direction preserving mapping, has no fixed point. From Step one, f is half-continuous. From Step two, there exists a continuous mapping p :

C → R n such that for every x ∈ C, p(x), f (x) -x > 0. Now, define g : C → C by ∀x ∈ C, g(x) = proj C (x + p(x)),
where for every y ∈ R n , proj C (y) denotes the orthogonal projection of y on the convex set C. Clearly, g is a continuous mapping. Thus, from Brouwer fixed point theorem, there exists x ∈ C such that g(x) = x, or equivalently proj C (x + p(x)) = x. Besides, from a standard characterization of projection of y ∈ R n on a convex set C, one has ∀c ∈ C, y -proj C (y), c -proj C (y) ≤ 0 Applying this inequality to y = x + p(x) and c = f (x), one obtains

x + p(x) -proj C (x + p(x)), f (x) -proj C (x + p(x)) ≤ 0, or equivalently, since proj C (x + p(x)) = x, p(x), f (x) -x ≤ 0.
a contradiction with the definition of p. This ends the proof of the Theorem 2.1.

Extension of Theorem 2.1

In the proof given above, one can notice that one only needs that the mapping f is half-continuous. Thus, from Step two and three of the previous proof, one obtains the following result:

Theorem 3.1 Let C a non empty convex compact subset of R n , let f : C → C half-continuous, which means: ∀x ∈ C, x = f (x) ⇒ ∃p ∈ R n , ∃ > 0 such that: ∀x ∈ B(x, ) ∩ C, x = f (x ) ⇒ p, f (x ) -x > 0.
Then f admits a fixed point.

In fact, half-continuity is strictly weaker than locally gross direction preserving assumption. Indeed, geometrically, Herings et al. assumption requires that for every x ∈ C which is not a fixed point, then f (y) -y and f (z) -z make a sharp angle or are orthogonal for y and z on some neighborhood of x. Our half-continuity assumption requires that for every x ∈ C which is not a fixed point, then all the vectors f (y) -y are in a same strict half-space for y on some neighborhood of x. This is clearly a weaker assumption. In the following example, one exhibits a large class of half-continuous mappings that may be not locally gross direction preserving:

Example Let C a closed subset of R n and F a multivalued mapping from C to R n . Suppose that F has a closed graph, non empty convex values and that F (C) is bounded. Suppose that for every x ∈ C one has x / ∈ F (x), i.e. F has no fixed point. Then we assert that any selection f of F (i.e. for every x ∈ C, f (x) ∈ F (x)) is half-continuous. To prove this result, let x ∈ C. Since F (x) is convex and compact (because F has a closed graph and F (C) is bounded) and since, by assumption, x / ∈ F (x), a separation theorem implies that there exists p ∈ R n such that ∀y ∈ F (x), p, y -x > 0.

(

) 2
and in particular, since f is a selection of F ,

p, f (x) -x > 0. (3)
Now, if f is not half-continuous at x, then from Equation 3, there exists a sequence (x n ) n∈N of C converging to x ∈ C, and such that

∀n ∈ N, p, f (x n ) -x n ≤ 0. (4)
Let us define, for every integer n ∈ N, y n = f (x n ). Since f is a selection of F , one has y n ∈ F (x n ) for every integer n. Since F (C) is bounded, and since F has a closed graph, the sequence (y n ) n∈N converges (up to an extraction) to ȳ ∈ F (x).

Passing to the limit in Equation 4, one obtains p, ȳ -x ≤ 0 (5) which contradicts Equation 2, because ȳ ∈ F (x). Now, one constructs a multivalued mapping satisfying the assumptions of Example above, and a particular selection f of F which is not locally gross direction preserving. Let C be the closed unit ball of R 2 , and define

f : C → R 2 by f (x, y) = (x, y) + (1, 0) if x ≤ 0 and f (x, y) = (x, y) + (-1, 1) if x > 0. Define a multivalued mapping F from C to R 2 by ∀(x, y) ∈ C, x = 0 ⇒ F (x) = {f (x)} and ∀(0, y) ∈ C, F (0, y) = {(0, y) + (-1 + 2t, 1 -t), t ∈ [0, 1]}.
Clearly, f is a selection of F , F has non empty, convex values, and F (C) is compact. We let the reader prove that F has a closed graph, and clearly F has no fixed points. So, f , which is a selection of F , is half-continuous (which can be checked directly). Now, f is not locally gross direction preserving at every (0, ȳ) for ȳ ∈ [-1, 1]. Indeed, one has f (x, y) -(x, y), f (x , y ) -(x , y ) = -1 for x ≤ 0 and x > 0. Thus one can not have f (x, y) -(x, y), f (x , y ) -(x , y ) ≥ 0 for (x, y) and (x , y ) in a neighborhood of (0, ȳ), and so f is not locally gross direction preserving at (0, ȳ).

An extension of Nash equilibrium existence theorem

In this subsection, we apply our main fixed point theorem to Game theory. Consider n players, and for every i = 1, ..., n, let X i be the strategy space of player i. Suppose each X i is included in R N where N ∈ N * is fixed. For every i = 1, ..., n, let b i : Π n j=1 X j → X i be a best reply function of player i and b = (b 1 , ..., b n ). For every x = (x 1 , ..., x n ) ∈ Π n j=1 X j , b i (x) is a best strategy of player i, given the strategy x j of all the others players (j = i). A game G is defined by the couple G = ((X i ) n i=1 , (b i ) n i=1 ). A Nash equilibrium of G is x = (x 1 , ..., x n ) ∈ Π n j=1 X j such that for every i = 1, ..., n, b i (x) = x i .

We now define half-continuous games:

Definition 4.1 A game G = ((X i ) n i=1 , (b i ) n i=1
) is said to be half-continuous if for every x ∈ Π n i=1 X i which is not a Nash equilibrium, there exists p ∈ (R N ) n and a neighborhood V x of x in (R N ) n such that for every x ∈ V x which is not a Nash equilibrium, one has p, b(x ) -x > 0.

A possible interpretation of this definition is the following: suppose there are N goods, and suppose a strategy x i of each player i = 1, ..., n is to specify the quantity x i (j) he wants to buy for every good j = 1, ..., N . Suppose that the players do not know the real prices of the N goods, but each player i = 1, ..., N is ready to pay p i for having one unit of each good. Consider a strategy profil x = (x 1 , ..., x n ). For every i = 1, ..., n, p i (b i (x) -x i ) would be the best profit that player i can do if he replaces his strategy x i by b i (x), given the strategy of the others. Call this the best profit of player i reachable from x.

Then, a game G is half-continuous if each time x is not a Nash equilibrium, then there exists a profil of prices p = (p 1 , ..., p n ) such that the sum of best profits of all payers reachable from x ∈ Π n i=1 X i is strictly positive, for every non Nash equilibrium x in a neighborhood of x.

In the following, a game G = ((X i ) n i=1 , (b i ) n i=1 ) is compact (resp. convex) if all the X i are compact (resp. convex): Theorem 4.2 Every compact, convex and half-continuous game G = ((X i ) n i=1 , (b i ) n i=1 ) admits a Nash-equilibrium.

Proof. Let X = Π n i=1 X i . Clearly, since G is half-continuous, the mapping b : X → X is half-continuous. From Theorem 3.1, b admits a fixed point which is a Nash equilibrium.

  λ n a partition of unity subordinate to this open covering. That means that each λ i is a continuous mapping from C to [0, 1] such that

	∀x ∈ C,

n i=1

Centre d'Economie de la Sorbonne, 106/112 Boulevard de l'Hôpital, 75013 Paris. Email: bich@univ-paris1.fr.


Thomas Bonald

Bertrand Charpentier

Alexis Galland

Alexandre Hollocou

Ant, Arthropod, Spider, Bee Animal

Bat, Vampire George Dog

Byron

Baron Byron

Parrot, Gull, Surtsey Herring Chordate

Gull

Hierarchical Graph Clustering using Node Pair Sampling

Keywords: Scientific classification, Animal, Chordate, Binomial nomenclature, Bird 231 Iron, Oxygen, Electron, Hydrogen, Phase (matter) Cretaceous, Jurassic 51 Binomial nomenclature, Bird, Carolus Linnaeus, Insect, Bird migration 24 Paleocene

We present a novel hierarchical graph clustering algorithm inspired by modularity-based clustering techniques. The algorithm is agglomerative and based on a simple distance between clusters induced by the probability of sampling node pairs. We prove that this distance is reducible, which enables the use of the nearest-neighbor chain to speed up the agglomeration. The output of the algorithm is a regular dendrogram, which reveals the multi-scale structure of the graph. The results are illustrated on both synthetic and real datasets. * Part of this work has been done while Bertrand Charpentier was a Master student at KTH, Sweden.

Introduction

Many datasets can be represented as graphs, being the graph explicitely embedded in data (e.g., the friendship relation of a social network) or built through some suitable similarity measure between data items (e.g., the number of papers co-authored by two researchers). Such graphs often exhibit a complex, multi-scale community structure where each node is invoved in many groups of nodes, so-called communities, of different sizes.

One of the most popular graph clustering algorithm is known as Louvain in name of the university of its inventors [START_REF] Blondel | Fast unfolding of communities in large networks[END_REF]. It is based on the greedy maximization of the modularity, a classical objective function introduced in [START_REF] Newman | Finding and evaluating community structure in networks[END_REF]. The Louvain algorithm is fast, memory-efficient, and provides meaningful clusters in practice. It does not enable an analysis of the graph at different scales, however [Fortunato andBarthelemy, 2007, Lancichinetti andFortunato, 2011]. The resolution parameter available in the current version of the algorithm 1 is not directly related to the number of clusters or the cluster sizes and thus is hard to adjust in practice.

In this paper, we present a novel algorithm for hierarchical clustering that captures the multi-scale nature of real graphs. The algorithm is fast, memory-efficient and parameter-free. It relies on a novel notion of distance between clusters induced by the probability of sampling node pairs. We prove that this distance is reducible, which guarantees that the resulting hierarchical clustering can be represented by regular dendrograms and enables a fast implementation of our algorithm through the nearest-neighbor chain scheme, a classical technique for agglomerative algorithms [START_REF] Murtagh | Algorithms for hierarchical clustering: an overview[END_REF].

The rest of the paper is organized as follows. We present the related work in Section 2. The notation used in the paper, the distance between clusters used to aggregate nodes and the clustering algorithm are presented in Sections 3, 4 and 5. The link with modularity and the Louvain algorithm is explained in Section 6. Section 7 shows the experimental results and Section 8 concludes the paper.

Related work

Most graph clustering algorithms are not hierarchical and rely on some resolution parameter that allows one to adapt the clustering to the dataset and to the intended purpose [START_REF] Reichardt | Statistical mechanics of community detection[END_REF], Arenas et al., 2008, Lambiotte et al., 2014, Newman, 2016]. This parameter is hard to adjust in practice, which motivates the present work.

Usual hierarchical clustering techniques apply to vector data [Ward, 1963, Murtagh andContreras, 2012]. They do not directly apply to graphs, unless the graph is embedded in some metric space, through spectral techniques for instance [Luxburg, 2007, Donetti and[START_REF] Donetti | [END_REF].

A number of hierarchical clustering algorithms have been developped specifically for graphs. The most popular algorithms are agglomerative and characterized by some distance between clusters, see [Newman, 2004, Pons and Latapy, 2005, Huang et al., 2010, Chang et al., 2011]. None of these distances has been proved to be reducible, a key property of our algorithm. Among non-agglomerative algorithms, the divisive approach of [START_REF] Newman | Finding and evaluating community structure in networks[END_REF] is based on the notion of edge betweenness while the iterative approach of [START_REF] Sales-Pardo | Extracting the hierarchical organization of complex systems[END_REF] and [START_REF] Lancichinetti | Detecting the overlapping and hierarchical community structure in complex networks[END_REF] look for local maxima of modularity or some fitness function; other approaches rely on statistical interence [START_REF] Clauset | Hierarchical structure and the prediction of missing links in networks[END_REF], replica correlations [START_REF] Ronhovde | Multiresolution community detection for megascale networks by information-based replica correlations[END_REF] and graph wavelets [START_REF] Tremblay | Graph wavelets for multiscale community mining[END_REF]. To our knowledge, non of these algorithms has been proved to lead to regular dendrograms (that is, without inversion).

Finally, the Louvain algorithm also provides a hierarchy, induced by the successive aggregation steps of the algorithm [START_REF] Blondel | Fast unfolding of communities in large networks[END_REF]. This is not a full hierarchy, however, as there are typically a few aggregation steps. Moreover, the same resolution is used in the optimization of modularity across all levels of the hierarchy, while the numbers of clusters decrease rapidly after a few aggregation steps. We shall see that our algorithm may be seen as a modified version of Louvain using a sliding resolution, that is adapted to the current agglomeration step of the algorithm.

Notation

Consider a weighted, undirected graph G = (V, E) of n nodes, with V = {1, . . . , n}. Let A be the corresponding weighted adjacency matrix. This is a symmetric, non-negative matrix such that for each i, j ∈ V , A ij > 0 if and only if there is an edge between i and j, in which case A ij is the weight of edge {i, j} ∈ E. We refer to the weight of node i as the sum of the weights of its incident edges,

w i = j∈V A ij .
Observe that for unit weights, w i is the degree of node i. The total weight of the nodes is:

w = i∈V w i = i,j∈V A ij .
We refer to a clustering C as any partition of V . In particular, each element of C is a subset of V , we refer to as a cluster.

Node pair sampling

The weights induce a probability distribution on node pairs,

∀i, j ∈ V, p(i, j) = A ij w ,
and a probability distribution on nodes,

∀i ∈ V, p(i) = j∈V p(i, j) = w i w .
Observe that the joint distribution p(i, j) depends on the graph (in particular, only neighbors i, j are sampled with positive probability), while the marginal distribution p(i) depends on the graph through the node weights only.

We define the distance between two distinct nodes i, j as the node pair sampling ratio 2 :

d(i, j) = p(i)p(j) p(i, j) , (1)
with d(i, j) = +∞ if p(i, j) = 0 (i.e., i and j are not neighbors). Nodes i, j are close with respect to this distance if the pair i, j is sampled much more frequently through the joint distribution p(i, j) than through the product distribution p(i)p(j). For unit weights, the joint distribution is uniform over the edges, so that the closest node pair is the pair of neighbors having the lowest degree product.

Another interpretation of the node distance d follows from the conditional probability,

∀i, j ∈ V, p(i|j) = p(i, j) p(j) = A ij w j .
This is the conditional probability of sampling i given that j is sampled (from the joint distribution). The distance between i and j can then be written

d(i, j) = p(i) p(i|j) = p(j) p(j|i) .
Nodes i, j are close with respect to this distance if i (respectively, j) is sampled much more frequently given that j is sampled (respectively, i).

Similarly, consider a clustering C of the graph (that is, a partition of V ). The weights induce a probability distribution on cluster pairs,

∀a, b ∈ C, p(a, b) = i∈a,j∈b p(i, j),
and a probability distribution on clusters,

∀a ∈ C, p(a) = i∈a p(i) = b∈C p(a, b).
We define the distance between two distinct clusters a, b as the cluster pair sampling ratio:

d(a, b) = p(a)p(b) p(a, b) , (2)
with d(a, b) = +∞ if p(a, b) = 0 (i.e., there is no edge between clusters a and b). Defining the conditional probability

∀a, b ∈ C, p(a|b) = p(a, b) p(b) ,
which is the conditional probability of sampling a given that b is sampled, we get

d(a, b) = p(a) p(a|b) = p(b) p(b|a) .
This distance will be used in the agglomerative algorithm to merge the closest clusters. We have the following key results.

Proposition 1 (Update formula) For any distinct clusters a, b, c ∈ C, d(a ∪ b, c) = p(a) p(a ∪ b) 1 d(a, c) + p(b) p(a ∪ b) 1 d(b, c) -1 .
2 The distance d is not a metric in general. We only require symmetry and non-negativity.

Proof. We have:

p(a ∪ b)p(c)d(a ∪ b, c) -1 = p(a ∪ b, c), = p(a, c) + p(b, c), = p(a)p(c)d(a, c) -1 + p(b)p(c)d(b, c) -1 ,
from which the formula follows. By the reducibility property, merging clusters a and b cannot decrease their minimum distance to any other cluster c.

Clustering algorithm

The agglomerative approach consists in starting from individual clusters (i.e., each node is in its own cluster) and merging clusters recursively. At each step of the algorithm, the two closest clusters are merged. We obtain the following algorithm:

1. Initialization C ← {{1}, . . . , {n}} L ← ∅ 2. Agglomeration For t = 1, . . . , n -1, • a, b ← arg min a ,b ∈C,a =b d(a , b ) • C ← C \ {a, b} • C ← C ∪ {a ∪ b} • L ← L ∪ {{a, b}} 3. Return L
The successive clusterings C 0 , C 1 , . . . , C n-1 produced by the algorithm, with C 0 = {{1}, . . . , {n}}, can be recovered from the list L of successive merges. Observe that clustering C t consists of n -t clusters, for t = 0, 1, . . . , n -1. By the reducibility property, the corresponding sequence of distances d 0 , d 1 , . . . , d n-1 between merged clusters, with d 0 = 0, is non-decreasing, resulting in a regular dendrogram (that is, without inversions) [START_REF] Murtagh | Algorithms for hierarchical clustering: an overview[END_REF].

It is worth noting that the graph G does not need to be connected. If the graph consists of k connected components, then the clustering C n-k gives these k connected components, whose respective distances are infinite; the k -1 last merges can then be done in an arbitrary order. Moreover, the hierarchies associated with these connected components are independent of one another (i.e., the algorithm successively applied to the corresponding subgraphs would produce exactly the same clustering). Similarly, we expect the clustering of weakly connected subgraphs to be approximately independent of one another. This is not the case of the Louvain algorithm, whose clustering depends on the whole graph through the total weight w, a shortcoming related to the resolution limit of modularity (see Section 6).

Aggregate graph. In view of (2), for any clustering C of V , the distance d(a, b) between two clusters a, b ∈ C is the distance between two nodes a, b of the following aggregate graph: nodes are the elements of C and the weight between a, b ∈ C (including the case a = b, corresponding to a self-loop) is i∈a,j∈b A ij . Thus the agglomerative algorithm can be implemented by merging nodes and updating the weights (and thus the distances between nodes) at each step of the algorithm. Since the initial nodes of the graph are indexed from 0 to n -1, we index the cluster created at step t of the algorithm by n + t. We obtain the following version of the above algorithm, where the clusters are replaced by their respective indices:

1. Initialization V ← {1, . . . , n} L ← ∅ 2. Agglomeration For t = 1, . . . , n -1, • i, j ← arg min i ,j ∈V,i =j d(i , j ) • L ← L ∪ {{i, j}} • V ← V \ {i, j} • V ← V ∪ {n + t} • p(n + t) ← p(i) + p(j) • p(n + t, u) ← p(i, u) + p(j, u) for u ∈ V \ {n ∪ t} 3. Return L
Nearest-neighbor chain. By the reducibility property of the distance, the algorithm can be implemented through the nearest-neighbor chain scheme [START_REF] Murtagh | Algorithms for hierarchical clustering: an overview[END_REF]. Starting from an arbitrary node, a chain a nearest neighbors is formed. Whenever two nodes of the chain are mutual nearest neighbors, these two nodes are merged and the chain is updated recursively, until the initial node is eventually merged. This scheme reduces the search of a global minimum (the pair of nodes i, j that minimizes d(i, j)) to that of a local minimum (any pair of nodes i, j such that d(i, j) = min j d(i, j ) = min i d(i , j)), which speeds up the algorithm while returning exactly the same hierarchy. It only requires a consistent tie-breaking rule for equal distances (e.g., any node at equal distance of i and j is considered as closer to i if and only if i < j). Observe that the space complexity of the algorithm is in O(m), where m is the number of edges of G (i.e., the graph size).

Link with modularity

The modularity is a standard metric to assess the quality of a clustering C (any partition of V ). Let δ C (i, j) = 1 if i, j are in the same cluster under clustering C, and δ C (i, j) = 0 otherwise. The modularity of clustering C is defined by [START_REF] Newman | Finding and evaluating community structure in networks[END_REF]:

Q(C) = 1 w i,j∈V (A ij - w i w j w )δ C (i, j), (3)
which can be written in terms of probability distributions,

Q(C) = i,j∈V (p(i, j) -p(i)p(j))δ C (i, j).
Thus the modularity is the difference between the probabilities of sampling two nodes of the same cluster under the joint distribution p(i, j) and under the product distribution p(i)p(j). It can also be expressed from the probability distributions at the cluster level,

Q(C) = a∈C (p(a, a) -p(a) 2 ).
It is clear from (3) that any clustering C maximizing modularity has some resolution limit, as pointed out in [START_REF] Fortunato | Resolution limit in community detection[END_REF], because the second term is normalized by the total weight w and thus becomes negligible for too small clusters. To go beyond this resolution limit, it is necessary to introduce a multiplicative factor γ, called the resolution. The modularity becomes:

Q γ (C) = i,j∈V (p(i, j) -γp(i)p(j))δ C (i, j), (4)
or equivalently,

Q γ (C) = a∈C (p(a, a) -γp(a) 2 ).
This resolution parameter can be interpreted through the Potts model of statistical physics [START_REF] Reichardt | Statistical mechanics of community detection[END_REF], random walks [START_REF] Lambiotte | Random walks, Markov processes and the multiscale modular organization of complex networks[END_REF], or statistical inference of a stochastic block model [Newman, 2016].

For γ = 0, the resolution is minimum and there is a single cluster, that is C = {{1, . . . , n}}; for γ → +∞, the resolution is maximum and each node has its own cluster, that is C = {{1}, . . . , {n}}.

The Louvain algorithm consists, for any fixed resolution parameter γ, of the following steps:

1. Initialization C ← {{1}, . . . , {n}}

Iteration

While modularity Q γ (C) increases, update C by moving one node from one cluster to another.

Aggregation

Merge all nodes belonging to the same cluster, update the weights and apply step 2 to the resulting aggregate graph while modularity is increased.

Return C

The result of step 2 depends on the order in which nodes and clusters are considered; typically, nodes are considered in a cyclic way and the target cluster of each node is that maximizing the modularity increase.

Our algorithm can be viewed as a modularity-maximizing scheme with a sliding resolution. Starting from the maximum resolution where each node has its own cluster, we look for the first value of the resolution parameter γ, say γ 1 , that triggers a single merge between two nodes, resulting in clustering C 1 . In view of (4), we have:

γ 1 = max i,j∈V p(i, j) p(i)p(j)
.

These two nodes are merged (corresponding to the aggregation phase of the Louvain algorithm) and we look for the next value of the resolution parameter, say γ 2 , that triggers a single merge between two nodes, resulting in clustering C 2 , and so on. By construction, the resolution at time t (that triggers the t-th merge) is γ t = 1/d t and the corresponding clustering C t is that of our algorithm. In particular, the sequence of resolutions γ 1 , . . . , γ n-1 is non-increasing.

To summarize, our algorithm consists of a simple but deep modification of the Louvain algorithm, where the iterative step (step 2) is replaced by a single merge, at the best current resolution (that resulting in a single merge). In particular, unlike the Louvain algorithm, our algorithm provides a full hierarchy. Moreover, the sequence of resolutions γ 1 , . . . , γ n-1 can be used as an input to the Louvain algorithm. Specifically, the resolution γ t provides exactly n -t clusters in our case, and the Louvain algorithm is expected to provide approximately the same number of clusters at this resolution.

Experiments

We have coded our hierarchical clustering algorithm, we refer to as Paris3 , in Python. All material necessary to reproduce the experiments presented below is available online4 .

Qualitative results. We start with a hierarchical stochastic block model, as described in [START_REF] Lyzinski | Community detection and classification in hierarchical stochastic blockmodels[END_REF]. There are n = 160 nodes structured in 2 levels, with 4 blocks of 40 nodes at level 1, each block of 40 nodes being divided into 4 blocks of 10 nodes at level 2 (see Figure 1). We also show in Figure 3 the number of clusters with respect to the resolution parameter γ for Paris (top) and Louvain (bottom). The results are very close, and clearly show the hierarchical structure of the model (vertical lines correspond to changes in the number of clusters). The key difference between both algorithms is that, while Louvain needs to be run for each resolution parameter γ (here 100 values ranging from 0.01 to 20), Paris is run only once, the relevant resolutions being direct outputs of the algorithm, embedded in the dendrogram (see Section 6).

We now consider four real datasets, whose characteristics are summarized in Table 1. The first dataset, extracted from OpenStreetMap5 , is the graph formed by the streets of the center of Paris. To illustrate the quality of the hierarchical clustering returned by our algorithm, we have extracted the two "best" clusterings, in terms of ratio between successive distance merges in the corresponding dendrogram; the results are shown in Figure 4. The best clustering gives two clusters, Rive Droite (with Ile de la Cité) and Rive Gauche, the two banks separated by the river Seine; the second best clustering divides these two clusters into sub-clusters. The second dataset, extracted from OpenFlights6 , is the graph of airports with the weight between two airports equal to the number of daily flights between them. We run Paris and extract the best clusterings from the largest component of the graph, as for the OpenStreet graph. The first two best clusterings isolate the Island/Groenland area and the Alaska from the rest of the world, the corresponding airports forming dense clusters, lightly connected with the other airports. The following two best clusterings are shown in Figure 5, with respectively 5 and 10 clusters corresponding to meaningful continental regions of the world. The third dataset is the graph formed by links between pages of Wikipedia for Schools7 , see [START_REF] West | Wikispeedia: An online game for inferring semantic distances between concepts[END_REF]Leskovec, 2012]. The graph is considered as undirected. Observe that the ability of selecting the clustering associated with some target number of clusters is one of the key advantage of Paris over Louvain. Moreover, Paris gives a full hiearchy of the pages, meaning that each of these clusters is divided into sub-clusters in the output of the algorithm. Table 2 (bottom) gives for instance, among the 500 clusters found by Paris (that is, in C n-500 ), the 10 largest clusters that are subclusters of the first cluster of Table 2 (top), related to animals / taxonomy. The subclusters tend to give meaningful groups of animals, revealing the multi-scale structure of the dataset.

The fourth dataset is the subgraph of Wikipedia restricted to pages related to humans. We have done the same experiment as for Wikipedia for Schools. The results are shown in Table 3. Again, we observe that clusters form relevant groups of people, with cluster #1 corresponding to political figures for instance, this cluster consisting of meaningful subgroups as shown in the bottom of Table 3. All this information is embedded in the dendrogram returned by Paris. Quantitative results. To assess the quality of the hierarchical clustering, we use the cost function proposed in [Dasgupta, 2016] and given by:

where the sum is over all left and right clusters a, b attached to any internal node of the tree representing the hierarchy. This is the expected size of the smallest subtree containing two random nodes i, j, sampled from the joint distribution p(i, j) introduced in Section 4. If the tree indeed reflects the underlying hierarchical structure of the graph, we expect most edges to link nodes that are close in the tree, i.e., whose common ancestor is relatively far from the root. The size of the corresponding subtree (whose root is this common ancestor) is expected to be small, meaning that ( 5) is a relevant cost function. Moreover, it was proved in [START_REF] Cohen-Addad | Hierarchical clustering: Objective functions and algorithms[END_REF] that if the graph is perfectly hierarchical, the underlying tree is optimal with respect to this cost function.

The results are presented in Table 5 for the graphs considered so far and the graphs of Table 4, selected from the SNAP datasets [START_REF] Leskovec | SNAP Datasets: Stanford Large Network Dataset Collection[END_REF]. The cost function is normalized by the number of nodes n so as to get a value between 0 and 1. We compare the performance of Paris to that of a spectral algorithm where the nodes are embedded in a space of dimension 20 by using the 20 leading eigenvectors of the Laplacian matrix L = D -A (D is the diagonal matrix of node weights) and applying the Ward method in the embedding space. The spectral decomposition of the Laplacian is based on standard functions on sparse matrices available in the Python package scipy.

Observe that we do not include Louvain in these experiments as this algorithm does not provide a full hierarchy of the graph, so that the cost function ( 5) is not applicable. The results are shown when the algorithm runs within some time limit (less than 1 hour on a computer equipped with a 2.8GHz Intel Core i7 CPU with 16GB of RAM). The best performance is displayed in bold characters. Observe that both algorithms have similar performance on those graphs where the results are available. However, Paris is much faster than the spectral algorithm, as shown by Table 6 (for each algorithm, the initial load or copy of the graph is not included in the running time; running times exceeding 1 hour are not shown). Paris is even faster than Louvain in most cases, while providing a much richer information on the graph.

Graph

Graph

Conclusion

We have proposed an agglomerative clustering algorithm for graphs based on a reducible distance between clusters. The algorithm captures the multi-scale structure of real graphs. It is parameter-free, fast and memory-efficient. In the future, we plan to work on the automatic extraction of clusterings from the dendrogram, at the most relevant resolutions.

Proposition 2 (

 2 Reducibility) For any distinct clusters a, b, c ∈ C, d(a ∪ b, c) ≥ min(d(a, c), d(b, c)). Proof. By Proposition 1, d(a ∪ b, c) is a weighted harmonic mean of d(a, c) and d(b, c), from which the inequality follows.

Figure 1 :

 1 Figure 1: A hierachical stochastic block model with 2 levels of hierarchy.

Figure 2 :

 2 Figure 2: Dendrogram associated with the clustering of Paris on a hierachical stochastic block model of 16 blocks.

Figure 3 :

 3 Figure 3: Number of clusters with respect to the resolution parameter γ for Paris (top) and Louvain (bottom) on the hierachical stochastic block model of Figure 1.

Graph

Figure 4 :

 4 Figure 4: Clusterings of the OpenStreet graph by Paris.

Figure 5 :

 5 Figure 5: Clusterings of the OpenFlights graph by Paris.

  b)(|a| + |b|),

Table 1 :

 1 Summary of the datasets.

		# nodes # edges Avg. degree
	OpenStreet	5,993	6,957	2.3
	OpenFlights	3,097	18,193	12
	Wikipedia Schools	4,589	106,644	46
	Wikipedia Humans 702,782 3,247,884	9.2

Table 2 (

 2 top) shows the 10 largest clusters of C n-100 , the 100 last (and thus largest) clusters found by Paris. Only pages of highest degrees are shown for each cluster.

Table 2 :

 2 Hierarchical clustering of the graph Wikipedia for Schools by Paris.

Table 3 :

 3 Hierarchical clustering of the graph Wikipedia Humans by Paris.

Table 4 :

 4 Summary of the considered graphs from SNAP.

		# nodes	# edges Avg. degree
	Facebook	4,039	88,234	44
	Amazon	334,863	925,872	5.5
	DBLP	317,080	1,049,866	6.6
	Twitter	81,306	1,342,310	33
	Youtube	1,134,890 2,987,624	5.2
	Google	855,802	4,291,352	10

Table 5 :

 5 Performance comparison of a spectral algorithm and Paris in terms of normalized Dasgupta's cost.

		Spectral		Paris
	OpenStreet		0.0103	0.0102
	OpenFlights		0.125		0.130
	Facebook		0.0479	0.0469
	Wikipedia Schools	0.452	0.402
	Amazon		-	0.0297
	DBLP		-	0.110
	Twitter		-	0.0908
	Youtube		-	0.185
	Wikipedia Humans	-		131
	Google		-	0.0121
	Graph	Spectral	Louvain	Paris
	OpenStreet	0.86s	0.19s		0.17s
	OpenFlight	0.51s	0.31s		0.33s
	Facebook	5.9s	1s		0.71s
	Wikipedia Schools	16s	2.2s		1.5s
	Amazon	-	45s		43s
	DBLP	-	52s		31s
	Twitter	-	35s		21s
	Youtube	-	8 min	16 min 30s
	Wikipedia Humans	-	2 min 30s 2 min 10s
	Google	-	3 min 50s 1 min 50s

Table 6 :

 6 Comparison of running times.

Paris = Pairwise AgglomeRation Induced by Sampling.

See https://github.com/tbonald/paris

https://openstreetmap.org

https://openflights.org

https://schools-wikipedia.org