Diverse novel resident Wolbachia strains in Culicine mosquitoes from Madagascar

Wolbachia endosymbiotic bacteria are widespread throughout insect species and Wolbachia transinfected in Aedes mosquito species has formed the basis for biocontrol programs as Wolbachia strains inhibit arboviral replication and can spread through populations. Resident strains in wild Culicine mosquito populations (the vectors of most arboviruses) requires further investigation given resident strains can also affect arboviral transmission. As Madagascar has a large diversity of both Culicine species and has had recent arboviral outbreaks, an entomology survey was undertaken, in five ecologically diverse sites, to determine the Wolbachia prevalence. We detected diverse novel resident Wolbachia strains within the Aedeomyia, Culex, Ficalbia, Mansonia and Uranotaenia genera. Wolbachia prevalence rates and strain characterisation through Sanger sequencing with multilocus sequence typing (MLST) and phylogenetic analysis revealed significant diversity and we detected co-infections with the environmentally acquired bacteria Asaia. Mosquitoes were screened for major arboviruses to investigate if any evidence could be provided for their potential role in transmission and we report the presence of Rift Valley fever virus in three Culex species: Culex tritaeniorhynchus, Culex antennatus and Culex decens. The implications of the presence of resident Wolbachia strains are discussed and how the discovery of novel strains can be utilized for applications in the development of biocontrol strategies.

ScIentIfIc REPORTS | (2018) 8:17456 | DOI: 10.1038/s41598-018-35658-z been confirmed after Wolbachia has established, despite local transmission events every year for the previous 13 years with a trend of increasing imported cases 24 . Further releases are ongoing in DENV endemic countries and mathematical models suggest the wMel strain of Wolbachia could reduce the basic reproduction number, R0, of DENV transmission by 66-75% 25 . Wolbachia strains also inhibit other medically important arboviruses including chikungunya virus (CHIKV) 14,26 , YFV 27 and recently ZIKV 26,28 .
There is also evidence that resident Wolbachia strains in laboratory mosquito colonies can inhibit arbovirus transmission [29][30][31][32] . Therefore, the prevalence and diversity of natural Wolbachia strains in field populations of mosquitoes requires further investigation to fully understand how widespread this bacterium is in wild Culicine populations. For example, recent studies undertaken to determine if Ae. aegypti contains resident strains report contrasting results with no evidence seen in a large collection from 27 countries 33 but detection of a resident strain in populations from the Southeastern United States 34 . Additional bacterial species that make up the mosquito microbiota may influence the presence of resident Wolbachia strains. In Ae. aegypti laboratory colonies, the acetic acid bacterium Asaia infects the gut and salivary glands and there is competition with transinfected Wolbachia strains for colonisation of mosquito reproductive tissues 35 . Asaia has been detected in field populations of Cx. quinquefasciatus 36 suggesting a complex association between these two bacterial species in wild mosquito populations.
The Culicinae subfamily of mosquitoes (Family Culicidae) contains a large variety of species (>3,000) particularly concentrated in tropical regions. Genera that transmit, or are implicated in transmitting, arboviruses include Aedeomyia (Ad.), Aedes, Culex, Mansonia (Ma.) and Uranotaenia (Ur.). Madagascar, a large island located off the southeast coast of mainland Africa, represents an optimal location to investigate Wolbachia prevalence in Culicines given its geographic isolation and large diversity of both mosquito species and arboviral diseases. Currently there have been 237 mosquito species morphologically identified in Madagascar in 15 different genera 37,38 , with 59% of these species thought to be endemic species in Madagascar. There have also been several arbovirus outbreaks in recent history including Rift Valley fever virus (RVFV) 39 , DENV (2006), CHIKV (2006CHIKV ( , 2007CHIKV ( , 2009CHIKV ( & 2010 and West Nile virus (WNV) 40 . Previous studies have identified 64 species in Madagascar that are implicated in transmission of medical or veterinary pathogens 37 . In this study, we undertook an entomology survey in five ecologically diverse sites in Madagascar to determine the prevalence of resident Wolbachia strains in adult female mosquitoes and discovered multiple novel Wolbachia strains. The characterisation, typing and phylogeny of these novel Wolbachia strains was analysed through Sanger sequencing, with multilocus sequence typing (MLST), to determine the diversity and evolutionary history of resident strains. We also compared Wolbachia and Asaia prevalence rates and screened mosquitoes for major medically important arboviruses to investigate if there was evidence for the involvement of any of these species as potential arbovirus vectors in these sites.

Results
Mosquito abundance and diversity at study sites. Overall, 1174 individuals belonging to seven genera were captured from five collection sites ( Fig. 1; Table 1). Qualitative collection site characteristics, median humidity and temperature and altitude are described in Supplementary Table S1. Cx. antennatus was collected in large numbers from Tsaramandroso and Ivato Aeroport and Cx. univitattus was the dominant species in Ankazobe. 820 (69.85%) of the mosquitoes were caught in CDC light traps, while 354 (30.15%) were caught with Zebubaited traps. Tsaramandroso had the largest number of individuals caught (n = 644) and highest diversity, which may have been influenced by both a high median relative humidity (79%) and temperature (24 °C) for the collection period. In contrast, Bemokotra had the smallest number of individuals collected (n = 19) and this could have been influenced by a lower median relative humidity (50%) despite a median temperature of 26 °C. Ivato Aeroport had the lowest diversity of species collected in the study and this could possibly be due to a low median temperature (11 °C) despite a median relative humidity of 88%. Aedes, Aedeomyia and Ficalbia (Fi.) species were collected only at lower elevation sites (Tsaramandroso, Bemokotra) and Culex species were collected predominantly at low and high elevation sites (Tsaramandroso, Ankazobe and Ivato Aeroport). However, it should be noted that there are many more factors that influence mosquito species abundance and diversity.
Detection of resident Wolbachia strains and mosquito species identification. We screened individual mosquitoes for the presence of resident Wolbachia strains using three conserved Wolbachia genes (16 rRNA, wsp and ftsZ) and detected resident strains in six mosquito genera ( Table 2). The prevalence of Wolbachia infection within mosquito species ranged from 3% (1/32) of Cx. antennatus, to 100% for Ad. madagascarica (9/9). PCR analysis resulted in congruency with the same individual samples amplifying fragments for all three genes, with the exception of Ma. uniformis, Fi. circumtestacea, Cx. antennatus and Cx. duttoni, in which the wsp gene was not amplified in any individuals. As Madagascar contains diverse mosquito species, we amplified a fragment of the mitochondrial cytochrome c oxidase subunit I (CO1) gene 41 of Wolbachia-infected individuals to provide as much molecular confirmation of species as possible, given the lack of available CO1 sequences in certain regions of the gene for some species. We were unable to amplify any CO1 gene fragments of Fi. circumtestacea but the remaining Wolbachia-infected species were confirmed by successful CO1 sequencing (Fig. 2a). CO1 sequences were deposited in Genbank (accession numbers MK033247-MK033269, Supplementary Table S1) including the first CO1 sequences for Ad. madagascarica, and two Uranotaenia species from Madagascar, as well as the first Cx. decens sequences covering this region of the CO1 gene.
Surprisingly, we found no evidence of resident Wolbachia strains in members of the Cx. pipiens complex indicating the absence of the wPip strain which is usually ubiquitous in Cx. pipiens pipiens and Cx. pipiens quinquefasciatus sibling species 42 . As species within the Cx. pipiens complex are morphologically indistinguishable, we selected a sub-sample of individuals from different collection sites to determine the species using both species-specific PCR analysis and sequencing targeting a CO1 gene fragment that discriminates species within the Cx. pipiens complex 43 . PCR-based species identification targeting the ace-2 gene 44  Cx. pipiens pipiens within the Cx. pipiens complex (Fig. 2b). We also analysed three females that were morphologically identified as Cx. pipiens complex individuals but did not result in any species confirmation through PCR (Fig. 2b). CO1 sequencing indicated these individuals were genetically diverse and, with the sequences currently available for comparison, were most closely related to Cx. rima and sequences named Cx. sp.GP-A from Kenya -Genbank accession number KU380413.1 (94% nucleotide identities, but only 89-94% coverage to the gene fragment amplified), suggesting the possibility of unknown Culex species present in Madagascar. We also screened individuals from all mosquito species for Asaia and detected the presence of this competing bacteria in 13 of 17 species (Table 2). Variable prevalence rates were observed and we detected co-infections of Wolbachia and Asaia in Ad. madagascarica (n = 8), Cx. antennatus (n = 1), Cx. decens (n = 2), Cx. duttoni (n = 1), Ma. uniformis (n = 5) and Uranotaenia species (n = 4). Fisher's exact post hoc tests revealed no association between Wolbachia and Asaia prevalence rates in species that had at least one individual infected with either bacteria ( Table 2).
Wolbachia strain characterisations. Sanger sequencing of 16S rRNA, wsp and multilocus sequence typing (MLST) gene fragments was undertaken with Wolbachia-infected individuals to characterise, type and determine strain phylogenies. As Madagascar contains some mosquito species that are only present on the island, we predicted the presence of divergent Wolbachia strains and this was reflected by variation in the success of amplifying gene target fragments using standard primers and alternative protocols that use M13 sequencing adaptors and/or primers comprising increased degeneracies 45 . Despite using a variety of primers, we were unable to amplify wsp or any of the five MLST gene loci for Cx. antennatus and Cx. duttoni. For the remaining resident Wolbachia strains, we amplified wsp from all strains except wUnif-Mad (Table 3). Phylogenetic analysis based on the 16S rRNA gene indicates most of these novel Wolbachia strains are clustering with Supergroup B Wolbachia strains, such as wPip and wDei (Fig. 3). Interestingly, only wUra1 clusters with Supergroup A Wolbachia strains, such as wMel and wRi. As the wsp gene has been evolving at a faster rate and provides more informative strain phylogenies, with improved phylogenetic resolution 46 , we compared the phylogenetic relationships for strains that amplified a wsp fragment (Fig. 4). Similar phylogenetic relationships are shown with wsp, with wUra1 sequences clustering with wMel (Supergroup A), but the remaining strains (wUra2, wDec and wMad) clustered with Supergroup B Wolbachia strains. MLST resulted in successful sequencing of the majority of gene loci, although we were unable to amplify hcpA from wUnif-Mad despite using alternative protocols with degenerate primers (Table 3). ftsZ universal primers 47 were also needed to obtain a sequence for wMad. Each of the wMad gene sequences produced an exact match to alleles already present in the database, with the resulting allelic profile producing an exact match with strain type 125. This allelic profile matches four isolates currently present in the MLST database, detected in three Lepidopteran species, from French Polynesia, India, Japan and Tanzania (Isolate numbers 40, 247, 270 and 322 respectively), but no isolates from mosquito or other Dipteran hosts. The remaining strains had new alleles for at least one of the MLST gene loci (sequences differed from those currently present in the database by at least one nucleotide difference per locus), resulting in new allelic profiles highlighting the diversity and novelty of these resident strains. Wolbachia gene sequences were deposited in GenBank (accession numbers MK026554-MK026563 and MK033270-MK033320, Supplementary Table S2).  The phylogeny of these novel strains based on concatenated sequences of all five MLST gene loci (where possible) (Fig. 5), or four MLST gene loci for analysis of wUnif-Mad (hcpA sequence not available) alongside the other strains (Fig. 6), also reveals significant strain divergence. As expected, due to the general incongruence between host phylogenetic relationships and their resident Wolbachia strains 48,49 , the majority of the novel Wolbachia strains in our study provide no evidence for clustering with other strains that infect mosquitoes. One exception is the wMad strain which, although distinct, appears closely related to wPip (the resident strain normally present in the Cx. pipiens complex). The other exception is the wUnif-Mad strain (based on 4 MLST loci, without hcpA), which is closely related to the Wolbachia strain previously detected in Ma. uniformis in Kenya (Isolate number 500), and where the pattern of amplification success was also matched as there was absence of amplification for both hcpA and wsp. However, wUnif-Mad differs by two nucleotides within the fbpA sequence, producing a novel fbpA allele and therefore a new allelic profile ( Fig. 6 and Table 3).
Site locations for Wolbachia-infected mosquitoes. The novel resident Wolbachia strain in Cx. decens, named wDec, was only found in individuals collected from Tsaramandroso with no detectable infections in individuals from Anivorano Nord or Bemokotra. Resident Wolbachia strains in both Cx. antennatus (wAnt) and Cx. duttoni (wDut) were from individuals collected in Tsaramandroso. The novel wMad Wolbachia strain was detected in all Ad. madagascarica females collected in Tsaramandroso. We detected the wUnif-Mad strain in Ma. uniformis from both Tsaramandroso and Anivorano Nord. Our single Wolbachia-infected Fi. circumtestacea (wCir) was collected from Tsaramandroso and novel Wolbachia strains in Uranotaenia species (wUra1 and wUra2) were detected in individuals from Anivorano Nord, Tsaramandroso and Bemokotra.

Arbovirus detection in mosquitoes.
In total, 704 mosquitoes of species considered potential vectors of medically important arboviruses in Madagascar were screened using arbovirus-specific PCR assays. Mosquitoes were screened as either individuals or pools of 3-5 individuals from the same location and trapping method ( Table 1). RVFV RNA was detected in two individually screened female adults of Cx. decens: one collected from Tsaramandroso and one from Bemokotra. We also detected RVFV in ten individual Cx. tritaeniorhynchus females and five individual Cx. antennatus females from Tsaramandroso. We sequenced RVFV PCR products from all positive individuals to confirm correct target amplification and also used CO1 sequencing to confirm mosquito species (representative samples are shown in Fig. 2a). No evidence of infection for the remaining arboviruses was present from any mosquitoes collected from the five sites in our study.

Discussion
Madagascar contains some of the most diverse and unique flora and fauna given its geographical isolation. There are currently 237 known species of mosquitoes present on the island, but recent newly described species 37 would suggest the possibility of more species present only in Madagascar. In our study, we assessed the Wolbachia prevalence of mosquito species that have the potential to be arbovirus vectors from five diverse ecological sites. To our knowledge, we are reporting for the first time resident Wolbachia strains in seven mosquito species: Ad. The novel resident Wolbachia strain in Cx. decens, named wDec, was only found in individuals collected from one location (Tsaramandroso) and an overall prevalence of 20% indicates that prevalence rates are likely to be variable across the island. Cx. decens is present in all biogeographic domains of Madagascar and is particularly  50 . In Madagascar, this species was previously shown to be infected with WNV and Babahoyo virus 50 and in mainland Africa has been shown to be infected with additional arboviruses including Moussa virus 51 . The wMad strain was found in all Ad. madagascarica individuals we collected from Tsaramandroso. Ad. madagascarica was previously found naturally infected with WNV and is zoophilic, with a preference for feeding on avian blood 52 , demonstrating its likely involvement in zoonotic arbovirus transmission in Madagascar. In previous entomological surveys in north-western Madagascar where WNV is endemic, Ad. madagascarica was the most abundant species and WNV RNA was detected in one pool suggesting this species may play a role in WNV maintenance or transmission 53 . Novel resident Wolbachia strains were also detected in two Uranotaenia species, of which there is very limited knowledge of these species in Madagascar. The lack of available CO1 sequences prevented molecular species confirmation but Ur. anopheloides is endemic to Madagascar and the Comoros archipelago and is most abundant is warmer regions of the western domain 38 . Ur. alboabdominalis is thought to occur mainly in the eastern and western domains of Madagascar and Ur. neireti in the central and eastern domains between 900-2000 m above sea level 37 .

Mosquito species Number
There is also limited knowledge of Fi. circumtestacea in Madagascar, with reports of its presence in the eastern and western domains 37 , but this species is not currently confirmed as being involved in arbovirus transmission. As we were only able to amplify a fragment of the 16S rRNA gene for resident strains in Cx. antennatus (wAnt)  Table 3. Novel resident Wolbachia strain wsp and MLST gene allelic profiles. Exact matches to existing alleles present in the database are shown in bold, novel alleles are shown in italics and denoted by the allele number of the closest match (number of single nucleotide differences to the closest match). *Alternative degenerate primers used to generate sequence. and Cx. duttoni (wDut), further studies are needed to determine if the lack of additional Wolbachia gene amplification is due to strain variability or low-density infections preventing successful amplification. An alternative explanation could be that amplification resulted from Wolbachia genes integrated into the mosquito host genomes resulting from ancient horizontal gene transfers 54,55 . Although most genes transferred between Wolbachia and their hosts are non-functional in the recipient genome 56 , expression of genes associated with Wolbachia prophage regions was previously observed in Ae. aegypti 55 . As we extracted RNA from field caught mosquitoes, we were able to demonstrate expression of multiple Wolbachia genes (including the ftsZ cell cycle gene that plays a central role during bacterial cytokinesis and the Wolbachia surface protein gene for most strains) which would indicate amplification is more likely of bacterial gene origin rather than integrated into the host genome.
Ma. uniformis is present in all biogeographical domains of Madagascar, is considered an anthropophilic vector of numerous arboviruses such as RVFV 57 and WNV 52 , and is a known vector of Wuchereria bancrofti filarial nematodes in numerous countries. Resident Wolbachia strains in Ma. uniformis have previously been shown in southeast Asia 58 , Kenya 59 and more recently in Sri Lanka 60 . Interestingly in Sri Lanka, 3/3 (100%) of Ma. uniformis individuals amplified the wsp gene 60 but in our study all seven individuals failed to amplify wsp. Ma. uniformis was also one of the four Wolbachia-infected mosquito species from populations in Kenya 59 and this provides the only comparative MLST and reports the requirement for nested PCR to amplify hcpA. Our Wolbachia prevalence rate of 29% (7/24) for Ma. uniformis is similar to the 29% (5/19) reported in Kenyan populations 59 suggesting this species has variable prevalence rates in wild populations.
The absence of Wolbachia infections in the Cx. pipiens complex (particularly Cx. pipiens pipiens which were confirmed by PCR and CO1 sequencing) was surprising given high infection rates of the wPip strain are often  [61][62][63][64][65] . Furthermore, our results contrast with a study undertaken in Madagascar as part of a wider Incompatible Insect Technique (IIT) from southwestern Indian Ocean islands in which the wPip strain was detected in Cx. pipiens quinquefasciatus mosquitoes from Antananarivo 66 . This suggests the possibility of variable infection rates between Cx. pipiens pipiens and its sibling species Cx pipiens quinquefasciatus in Madagascar populations (which we did not collect in our study). We also undertook molecular identification of three Culex individuals (also Wolbachia-uninfected) that were phylogenetically similar to Cx. rima (Fig. 2b) which is consistent with a previous entomological survey that describes an unknown species morphologically similar to Cx. rima 67 . Further phylogenetic studies are warranted to determine the diversity of Culex species and the variability of Wolbachia prevalence within the Cx. pipiens complex. The presence of resident Wolbachia strains in mosquito vector species which can transmit human arboviruses could be influencing arboviral transmission dynamics in field populations of Madagascar. In laboratory studies, resident Wolbachia strains have been shown to impact arboviral transmission. For example, Wolbachia was shown to reduce DENV infection of salivary glands and limit transmission in Ae. albopictus 29 . Outbreaks of RVFV in Madagascar are thought to result from infected domestic animals imported from east Africa 68 and Cx. antennatus has previously been identified as an important RVFV vector in Madagascar 39 . Our results would indicate Cx. tritaeniorhynchus is also likely contributing to transmission in Tsaramandroso given we also detected RVFV RNA in multiple non-blood fed females of this species. RVFV has not previously been detected in this species in Madagascar, however, Cx. tritaeniorhynchus was implicated as a major vector during a large outbreak of RVFV in Saudi Arabia 69 so this species contribution to transmission in Madagascar could potentially be currently underestimated. The detection of RVFV in three Culex species at Tsaramandroso and Cx. decens at Bemokotra would correlate with previous studies that have shown that RVFV is one of the most abundant and widely distributed arboviruses across the island and there have been several recent RVFV outbreaks 39,70 . In our study, we detected RVFV in Wolbachia-uninfected Cx. decens and Cx. tritaeniorhynchus individuals but our single Wolbachia-infected Cx. antennatus was co-infected with RVFV. This female was non-blood fed indicating potential dissemination of RVFV beyond the bloodmeal. We were also unable to amplify any Wolbachia genes other than 16S rRNA from this individual which could be explained by a very low-density resident strain present in Cx. antennatus potentially resulting in the low prevalence rate of 3.1% (1/32). Wolbachia tissue tropism influences the potential to inhibit arboviruses 15,71 as some resident Wolbachia strains are present predominantly in reproductive tissue and have little effects on arboviruses 72 but transinfected strains that are present in tissues such as salivary glands have strong inhibitory effects 15,16 . However, as arbovirus infection rates are normally low in mosquito populations, correlating the prevalence of Wolbachia and RVFV would require a significantly larger number of mosquitoes.
The co-detection of resident Wolbachia strains and Asaia in seven species also indicates the potential for complex tripartite interaction with arboviruses. In our study, Asaia was also detected in the single Wolbachia-infected Cx. antennatus individual so further studies should be undertaken to determine the dynamics of these microbes within wild mosquito populations. The prevalence of Asaia is highly dependent on the environment given this bacterium can be acquired through various routes including from males to females during mating and from both larval and adult mosquito sugar feeding 73 . For example, prevalence rates of 46% were seen in field-caught Ae. albopictus from Madagascar 74 but the overall composition and diversity of the mosquito microbiota in field populations of Ae. aegypti and Ae. albopictus in Madagascar is also highly dependent on environmental factors 75,76 . It has also recently been shown that there is no evidence for Wolbachia-Asaia co-infections in a wide range of Anopheles species and the environment likely influences Asaia prevalence and density in wild mosquito populations 77 . The discovery of novel resident Wolbachia strains in mosquito species from Madagascar may also impact future attempts to extend Wolbachia biocontrol strategies by using these strains for applied use through transinfection. The naturally occurring wAlbA and wAlbB strains of Wolbachia that infect Ae. albopictus have been successfully transferred to Ae. aegypti 16,18,20 and significantly inhibit DENV 16,20 . This suggests the introduction of novel resident Wolbachia strains from other Culicine species may generate inhibitory effects on arboviruses although this will be dependent on several factors and there is variability between Wolbachia strains on the strength of inhibition 15,16,20,21 . Further experiments should elucidate both the density and tissue tropism of these novel strains to determine candidate strains for transfer to species such as Cx. quinquefasciatus that are major arbovirus vectors. As Cx. quinquefasciatus contain a resident Wolbachia strain (wPip) 78 , introduction of transinfected strains would require a stable association, with introduced strains growing to higher densities in specific tissues which result in inhibition of pathogen transmission 71 . Our discovery of novel resident strains in Culex species closely related to Cx. quinquefasciatus may improve transinfection success and ultimately lead to Wolbachia biocontrol strategies targeting Culex species. Mosquito sampling. Mosquito sampling was undertaken utilizing both CDC light traps and a net trap baited with Zebu (local species of cattle) to attract zoophilic species at each of the five study sites 70 . Trapping was performed over two consecutive nights per site, with CDC light traps set up at dusk to be operational during sunset and overnight, before being collected at dawn. Eight CDC light traps were placed around each site nearby dwellings, near poultry coops, inside areas of denser vegetation (typically forested areas or crop fields) and near potential breeding locations which varied depending on each site. At each site, one Zebu trap was placed in a way which simulated local animal husbandry practices (i.e. inside of a shed, nearby to a large corral, or on its own tied to a post). Zebu traps were set up following sunset and mosquitoes were collected by mouth aspiration before dawn. Temperature and relative humidity were recorded at each site using Lascar Electronics data loggers and the median was calculated over the collection period.

Methods
Morphological identification and RNA extraction. Mosquitoes were anesthetized with chloroform vapor and morphological identification was performed using a variety of mosquito genera keys 37 . Following identification, mosquitoes were placed into 96 well storage plates or Eppendorf tubes and RNAlater (Sigma, UK) was added prior to storage at 4 °C or lower to prevent viral RNA degradation. RNA was extracted from individual whole mosquitoes (or additional pools for further detection of human pathogens) using QIAGEN RNeasy 96 kits according to manufacturer's instructions. RNA was eluted in 40 μl of RNase-free water and stored at −80 °C. A QIAGEN QuantiTect Reverse Transcription kit was used to reverse transcribe RNA, generating cDNA from all RNA extracts, according to manufacturer's instructions.

Molecular species identification.
Sanger sequencing targeting a fragment of the mitochondrial gene cytochrome oxidase 1 (CO1) gene was undertaken on mosquito cDNA extracts to confirm species identification 41 . Selected samples morphologically identified as part of the Cx. pipiens complex were further analyzed using multiplex PCR targeting the ace-2 locus 44 . Sanger sequencing was also undertaken targeting a different fragment of the CO1 gene specifically shown to provide the optimal discrimination of species within the Culex pipiens complex 43 . Wolbachia MLST. Multilocus sequence typing (MLST) was undertaken to characterize Wolbachia strains using the sequences of five conserved genes as molecular markers to genotype each strain. In brief, 450-500 base pair fragments of the coxA, fbpA, hcpA, gatB and ftsZ Wolbachia genes were amplified from one individual from each Wolbachia-infected mosquito species using previously optimised protocols 45 . A Culex pipiens gDNA extraction (previously shown to be infected with the wPip strain of Wolbachia) was used a positive control for each PCR run, in addition to no template controls (NTCs). If no amplification was detected using standard primers, further PCR analysis was undertaken using alternative protocols which included utilisation of M13 sequencing adaptors and/or primers with increased degeneracies 45 .
Sanger sequencing. PCR products were separated and visualised using 2% E-gel EX agarose gels (Invitrogen) with SYBR safe and an Invitrogen E-gel iBase Real-Time Transilluminator. PCR products were submitted to Source BioScience (Source BioScience Plc, Nottingham, UK) for PCR reaction clean-up, followed by Sanger sequencing to generate both forward and reverse reads. Sequencing analysis was carried out in MEGA7 82 as follows. Both chromatograms (forward and reverse traces) from each sample was manually checked, analysed, and edited as required, followed by alignment by ClustalW and checking to produce consensus sequences. Consensus sequences were used to perform nucleotide BLAST (NCBI) database queries, and Wolbachia gene loci sequences were used for searches against the Wolbachia MLST database (http://pubmlst.org/wolbachia) 83 . If a sequence produced an exact match in the MLST database, we assigned the appropriate allele number.
Phylogenetic analysis. Maximum Likelihood phylogenetic trees were constructed from Sanger sequences as follows. The evolutionary history was inferred by using the Maximum Likelihood method based on the Tamura-Nei model 84 . The tree with the highest log likelihood in each case is shown. The percentage of trees in which the associated taxa clustered together is shown next to the branches. Initial tree(s) for the heuristic search were obtained automatically by applying Neighbor-Joining and BioNJ algorithms to a matrix of pairwise distances estimated using the Maximum Composite Likelihood (MCL) approach, and then selecting the topology with superior log likelihood value. The trees are drawn to scale, with branch lengths measured in the number of substitutions per site. Codon positions included were 1st + 2nd + 3rd + Noncoding. All positions containing gaps and missing data were eliminated. The phylogeny test was by Bootstrap method with 1000 replications. Evolutionary analyses were conducted in MEGA7 82 .
Arbovirus screening. Arbovirus screening using published PCR assays included the major arboviruses of public health importance, suspected or having the potential of being transmitted in Madagascar (Supplementary  Table S4). PCR reactions for all assays except ZIKV were prepared using 5 µl of Qiagen SYBR Green Master mix, a final concentration of 1 µM of each primer, 1 µl of PCR grade water and 2 µl template cDNA, to a final reaction volume of 10 µl. Prepared reactions were run on a Roche LightCycler ® 96 System and PCR cycling conditions are described in Supplementary Table S1. Amplification was followed by a dissociation curve (95 °C for 10 seconds, 65 °C for 60 seconds and 97 °C for 1 second) to ensure the correct target sequence was being amplified. ZIKV screening was undertaken using a Taqman probe-based assay using 5 µl of Qiagen QuantiTect probe master mix, a final concentration of 1 µM of each primer, 1 µl of PCR grade water and 2 µl template cDNA, to a final reaction volume of 10 µl. PCR results were analysed using the LightCycler ® 96 software (Roche Diagnostics). Synthetic long oligonucleotide standards (Integrated DNA Technologies) of the amplified PCR product were generated in the absence of biological virus cDNA positive controls and each assay included negative (no template) controls.

Statistical analysis.
Fisher's exact post hoc test in Graphpad Prism 6 was used to compare Wolbachia and Asaia prevalence rates in species that had at least one individual infected with Wolbachia.
Ethics approval and consent to participate. Ethical clearance for the use of Zebu cattle in baited traps was obtained from the ethical committee of the Livestock ministry of Madagascar which is the sole relevant authority for animal care in Madagascar. The ethical committee number is: 2012/WN/Minel/3. We followed the European guidelines (European directives EU 86/609-STE123 and 2010/63/EU) for animal handling. To minimize the risk of infection of mosquito borne diseases, an internal net was used to protect the cattle against mosquito bites.