# Vector soup: high-throughput identification of Neotropical phlebotomine sand flies using metabarcoding

Arthur Kocher, Jean-charles Gantier, Pascal Gaborit, Lucie Zinger, Helene Holota, Sophie Valiere, Isabelle Dusfour, Romain Girod, Anne-Laure Bañuls, Jérôme Murienne

► To cite this version:

Arthur Kocher, Jean-charles Gantier, Pascal Gaborit, Lucie Zinger, Helene Holota, et al. Vector soup: high-throughput identification of Neotropical phlebotomine sand flies using metabarcoding. Molecular Ecology Resources, 2016, 17 (2), pp.172-182. 10.1111/1755-0998.12556. hal-04917754

HAL Id: hal-04917754

https://hal.science/hal-04917754v1

Submitted on 28 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
1 Vector soup: high-throughput identification of neotropical phlebotomine

2 sand flies using metabarcoding

3 ARTHUR KOCHER*§, JEAN-CHARLES GANTIER†, PASCAL GABORIT¶, LUCIE ZINGER*, HELENE

4 HOLOTA*, SOPHIE VALIERE‡, ISABELLE DUSFOUR¶, ROMAIN GIROD¶, ANNE-LAURE BAÑULS§ AND

5 JEROME MURIENNE*.

6 *CNRS, Université Toulouse III Paul Sabatier, ENFA; UMR5174 EDB (Laboratoire Evolution et Diversité Biologique), Toulouse,

7 France; †Laboratoire des Identifications Fongiques et Entomo-parasitologiques, Mennecy, France; ‡GeT-PlaGe, Genotoul,

8 INRA Auzeville, 31326 Castanet-Tolosan, France; §UMR MIVEGEC (IRD 224 – CNRS 5290 – Université de Montpellier), 911

9 Avenue Agropolis, F34394 Montpellier, France; ¶Medical Entomology Unit, Institut Pasteur de la Guyane, 23 Avenue

10 Pasteur, BP 6010, 97306 Cayenne Cedex, French Guiana.

11
12 Keywords: Leishmaniasis, mitochondrial 16S rRNA gene, DNA barcode, vector monitoring, French

13 Guiana, Psychodidae

14 Running title: Sand fly identifications using metabarcoding

15
16 Abstract

17 Phlebotomine sand flies are hematophagous dipterans of primary medical importance. They
18 represent the only proven vectors of leishmaniasis worldwide and are involved in the transmission of
19 various other pathogens. Studying the ecology of sand flies is crucial to understand the epidemiology
20 of leishmaniasis and further control this disease. A major limitation in this regard is that traditional
21 morphological-based methods for sand fly species identifications are time-consuming and require
22 taxonomic expertise. DNA metabarcoding holds great promise in overcoming this issue by allowing
23 the identification of multiple species from a single bulk sample. Here, we assessed the reliability of a
24 short insect metabarcode located in the mitochondrial 16S rRNA for the identification of neotropical
25 sand flies, and constructed a reference database for 40 species found in French Guiana. Then, we
26 conducted a metabarcoding experiment on sand flies mixtures of known content and showed that
27 the method allows an accurate identification of specimens in pools. Finally, we applied
28 metabarcoding to field samples caught in a 1 ha forest plot in French Guiana. Besides providing
29 reliable molecular data for species-level assignations of phlebotomine sand flies, our study proves
30 the efficiency of metabarcoding based on the mitochondrial 16S rRNA for studying sand fly diversity
31 from bulk samples. The application of this high-throughput identification procedure to field samples
32 can provide great opportunities for vector monitoring and eco-epidemiological studies.
33

33
34 **Introduction**

35 Sand flies (Diptera, Psychodidae, Phlebotominae) form a subfamily of hematophagous dipterans that currently comprises close to 1,000 species distributed worldwide. They are of high 36 medical importance as they include all the proven vectors of leishmaniasis and are involved in the 37 transmission of other viral and bacterial diseases (Depaquit et al. 2010; Ready 2013). Three sand fly 38 genera are present in the New World, but only species of the *Lutzomyia* genus França, 1924 have 39 been described as vectors of leishmaniasis [the classification proposed by Young & Duncan (1994) 40 will be used throughout the text, but see Supporting Information for the alternative classification of 41 Galati (1995, 2014)]. The *Lutzomyia* genus currently includes 477 described species, among which 56 42 are proven or suspected vectors of the disease (Killick-Kendrick 1990; Maroli et al. 2013). 43

44 Investigating the ecology of sand flies is essential to better understand the epidemiology of 45 Leishmaniasis and further control this disease, particularly in regions undergoing important 46 environmental changes (Desjeux 2001; Rotureau 2006; Ready 2008; de Oliveira Miranda et al. 2015). 47 A major difficulty in this regard lies in the identification of sand fly specimens, especially in areas such 48 as the neotropics that harbour a high diversity of species. Traditional morphological identifications of 49 these minute insects are mainly based on digestive and reproductive organs (Young & Duncan 1994; 50 Munstermann 2004; El-Hossary 2006) and represent meticulous and time-consuming tasks 51 (dissection, mounting, examination), even for taxonomic experts.

52 Molecular methods based on various genetic markers have been developed to overcome this issue 53 (Depaquit 2014). So far, two techniques have been mainly used for molecular identifications of sand 54 flies: the polymerase chain reaction-restriction fragment-length polymorphism (PCR-RFLP; Aransay et 55 al. 1999; Terayama et al. 2008; Latrofa et al. 2012; Tiwary et al. 2012; Minter et al. 2013) and DNA 56 barcoding (Krüger et al. 2011; Kumar et al. 2012; Contreras Gutiérrez et al. 2014; Al-Ajmi et al. 2015; 57 Nzelu et al. 2015; Maia et al. 2015; Romero-Ricardo et al. 2016). DNA barcoding consists in the
58 sequencing of a short standard DNA marker for the taxonomic assignment of a sample by looking for
59 the closest match in a reference database (Hebert et al. 2003a). In the last decade, it has emerged as
60 the main alternative for species identification in a wide range of taxa, including medically important
61 arthropods (Ondrejicka et al. 2014). More recently, an extended version of DNA barcoding, referred
62 to as DNA metabarcoding, has enabled the identification of multiple species from bulk samples by
63 mass-PCR-amplification and high-throughput sequencing (Taberlet et al. 2012). This method has
64 been successfully applied to arthropod samples (Yu et al. 2012; Ji et al. 2013; Gibson et al. 2014),
65 opening new avenues for the study of insect biodiversity.

66 Here, the reliability of DNA metabarcoding to identify sand flies fly species contained in bulk samples
67 is assessed. The current standard barcode for animals is a c. 650-bp portion of the mitochondrial
68 cytochrome c oxidase subunit 1 gene (COI) which was already proved useful for the identification of
69 neotropical sand flies (Contreras Gutiérrez et al. 2014; Nzelu et al. 2015; Romero-Ricardo et al.
70 2016). The standard COI barcode is however not well suited for metabarcoding applications. First, it
71 is too long to be sequenced on the current high-throughput sequencing platforms such as the
72 Illumina Miseq. Moreover, it is virtually impossible to find highly conserved primer binding sites in
73 this coding gene, because of high mutation rates at the third codon position. This may not be a
74 problem when sequencing single specimens, but can lead to significant amplification bias or
75 detection failure when it comes to metabarcoding (Taberlet 2012, Deagle et al. 2014). Recently, a
76 short fragment of the mitochondrial 16S ribosomal RNA gene was shown to have good properties for
77 insect metabarcoding (Ins16S_1; Clarke et al. 2014). We therefore evaluated the reliability of the
78 Ins16_1 insect metabarcode for the identification of sand flies and constituted the first DNA
79 reference library for species found in French Guiana. We then conducted a metabarcoding assay on
80 artificially constituted sand fly mixtures to assess the accuracy of the method. Finally, we tested the
81 metabarcoding approach on pools of specimens freshly caught in the field.

82 **Material and method**
83 Reference library of DNA metabarcodes

84 Specimens of sand flies were collected with Center for Disease Control (CDC) light traps in various
85 locations in French Guiana between 2007 and 2011 (see Supporting Information for details).

86 Specimens of the same species collected in Trinidad and Tobago in 2012 were also included to assess
87 the validity of the marker over wide geographic ranges. When possible, specimens from both sexes
88 were included in each species, in order to provide molecular evidence for the current acceptations of
89 male-female pairings based on morphological criteria. The samples were processed at Institut
90 Pasteur de la Guyane (Cayenne, French Guiana) and stored at the LIFE laboratory (Mennecy, France).

91 After dissection, thoraxes were kept for molecular study while heads, wings and genitalia were used
92 for morphological identifications as follows: the pieces were treated according to Abonnenc (1972)
93 and mounted in Euparal medium (Réf 7356.1 Carl Roth GmbH). The morphological and
94 morphometrical data was obtained with the DM5500B LeicaTM microscope associated with the
95 calibration software Leica Application Suite (Version LAS 3.8.0.). The morphological keys of Young &
96 Duncan (1994) were then used for species identifications and specimens were stored as vouchers.

97 DNA was extracted using the DNeasy Blood and Tissue kit (Qiagen. Valencia. CA. USA) and amplified
98 using the Ins16S_1 primer pair (Clarke et al. 2014; Ins16S_1-F. 5'- TRRGACGAGAAGACCCTATA-3';
99 Ins16S_1-R. 5'- TCTTAATCCAACATCGAGGTC-3'). PCR amplification was performed in 25 µL mixtures
100 containing 2 µL of DNA template, 0.2 µL of AmpliTaq Gold® (5U/ µL; Applied Biosystems, Foster City,
101 CA, USA), 2.5 µL 10X PCR buffer (provided with Ampitaq Gold®), 0.5 µL dNTPs (2.5 mM each,
102 Promega, Madison, WI, USA), 1 µL of each primer (10 µM), 0.25 bovine serum albumin (10mg/mL,
103 Promega), 2.5 µL MgCl2 (25 mM, Applied Biosystems) and nuclease-free water (Promega).

104 Thermocycling conditions consisted in a first denaturation step at 95°C (10 min) followed by 40 cycles
105 of 30s at 95°C, 30s at 50°C and 30s at 72°C, and a final elongation step at 72°C (10 min). To enable
106 the sequencing of multiple PCR products in a single sequencing run, different combinations of tags
added to the 5' end of each primer were used (Binladen et al. 2007). Tags were eight base pairs long with at least five differences between them to minimize ambiguities in the downstream analyses.

PCR products were pooled and sent to the GeT-PlaGe core facilities of GenoToul (Toulouse, France) for paired-end sequencing on a MiSeq platform (Illumina, San Diego, CA, USA). Prior to sequencing, DNA was purified using the HighPrep PCR system (MAGBIO GENOMICS, Gaithersburg, MD, USA) and used for library construction with the Illumina NEXTflex PCR-Free DNA sequencing kit following the supplier's instructions (Bioo Scientific corp., Austin, TX, USA). Sequencing data was pre-processed and stored on the NG6 platform (Mariette et al. 2012) and all downstream bioinformatics analyses were performed on computing facilities of the GenoToul bioinformatics platform (Toulouse, France).

Sequence reads were analysed using the OBITools package (Boyer et al. 2016) following Riaz et al. (2011). Pair-end reads were aligned and merged, taking into account the phred quality scores for consensus construction and alignment score computation. Reads were then assigned to their corresponding sample, based on the tagged primer sequences. Low quality reads (with alignment scores below 50, containing Ns or shorter than 50 bp) were removed. Identical reads were dereplicated while keeping the coverage information (i.e. number of reads merged). For each sample, only the most abundant sequence was kept. Those that were represented by less than 50 reads were not considered. The resulting DNA reference database was then formatted to be used for taxonomic assignment of sequences. The bash script used for these bioinformatics steps and the DNA reference database in OBITools-extended fasta format are available in the Supporting Information.

The sequences were then aligned with ClustalW2 (Larkin et al. 2007) using default parameters. To allow comparisons with other DNA barcoding studies, pairwise genetic distances were computed using the Kimura-2-parameter (K2P) distance model (Kimura 1980). A neighbour-joining tree was generated (Saitou & Nei 1987) to analyse the clustering patterns among species. The R package
131 'spider' (R Development Core Team 2014, Brown et al. 2012) was then used to compute distance
132 statistics and evaluate the reliability of the metabarcode in assigning samples to the correct species
133 using the 'best match' criteria (Meier et al. 2006). In order to assess the robustness of the results to
134 various analytical procedures, we performed the same analyses with a different alignment algorithm,
135 muscle (Edgar 2004), or by taking into account the secondary structure of the 16S rRNA using
136 RNAsalsa (Stocsits et al. 2009), as well as using different genetic distance models (JC89 and TN93;
137 Jukes & Cantor 1969; Tamura & Nei 1993).

138 *Metabarcoding validation*

139 Next, a metabarcoding assay was conducted on DNA mixtures of known content, hereafter referred
to as "mock communities" (Fig. 1). The first type of mock community consisted of specimens from
140 two distinct species (L. (Trichophoromyia) ininii and L. (Trichopygomyia) trichopyga) that were pooled
141 at ratios of one to 10 or one to 50 prior to bulk DNA extraction. We tested two DNA extraction
142 protocols, adapted from Qiagen extraction kits and Chelex-based protocols (Casquet et al. 2012; Ji et
143 al. 2013, see Supporting Information for more details). The first one is designed to maximize DNA
144 quality whereas the second is cheap and requires little manipulation. Each extraction protocol was
145 performed twice for both abundance ratios. Because the concentration of PCR inhibitors and
146 template DNA can influence amplification success and the formation of PCR artefacts (Qiu et al.
147 2001), each resulting DNA solution was diluted 1, 10 and 100 fold in nuclease-free water (Promega)
148 prior to PCR amplifications. The second type of mock community consisted of a pool of DNA from 10
149 morphologically identified specimens belonging to distinct species. In this case the DNA was
150 extracted individually prior to pooling because the specimens were also included in the reference
151 library.
152

153 Three independent PCR replicates were made for each sample following the protocol described
154 above. For PCR multiplexing, we used the Latin Square Design as described by Esling et al. (2015).
155 This design, in which one over two tagged primer combinations are used, has been shown to
represent the optimal trade-off between multiplexing capacities and minimizing tag switching events.

Sequencing of amplicons as well as sequence read quality filtering and sample assignment were conducted as described above for the individual specimens. Sequences represented by only one read (singletons) were discarded because they are most likely to arise from artefacts. Sequences were then clustered within each sample using *Sumaclust* (Mercier *et al.* 2013) with a similarity threshold of 99%, as an initial denoising step to remove PCR or sequencing errors (Huse *et al.* 2010).

Taxonomic assignment was then performed on our reference DNA database to define operational taxonomic units (OTUs). We further curated the data in order to remove obvious contaminants as well as artefacts such as PCR chimeras, remaining sequencing errors, and missassigned sequence due to tag-switching events (Esling *et al.* 2015). The filtering process was implemented in R and consisted of three steps: (i) a reference based filtering where we kept only the OTUs that had at least 97% identity with its best match in the reference database (ii) filtering of mistagging errors based on OTUs assigned to non-used tag combos (non-critical mistags, see Esling *et al.* 2015 for details) (iii) removal of sequences that were not recovered in at least 2 PCR replicates. The results were then compared to the genuine content of samples to identify false negative and false positive OTUs. For bulk-extracted samples, read abundance ratios were computed to assess whether they were consistent with the specimen ratios. The scripts used for these steps, as well as the sequencing data, are available in the Supporting Information so that the analyses can be reproduced.

Ultimately, metabarcoding was applied to unidentified specimens freshly caught in the field. Sand flies were collected using ten CDC light traps over three nights in a 1 ha forest plot in the region of Iracoubo in French Guyana in October 2015 (WGS84 coordinates: lat=5.25156, long=-53.32059). After each night, the sand flies caught in the traps were killed by freezing, counted and put in 2 mL microcentrifuge tubes (one tube per trap night). Bulk-DNA extraction was performed directly in these tubes after ethanol evaporation using the Chelex protocols. DNA was diluted 10 fold in nuclease-free
water prior to PCR amplification. The metabarcoding protocol was then applied as previously described. In order to assess the comprehensiveness of the sampling, an incidence-based species accumulation curve (SAC) was constructed using the 'exact' method as implemented in the R package 'Vegan' (Oksanen et al. 2016). The mean Chao2 estimator of species richness was then computed from 100 permutations along the SAC. Assuming that PCR amplification bias should be consistent across samples, the relative read abundance was weighted by the number of specimens collected in each sample and summed to represent the composition of the sand fly fauna in the whole site.

## 188 Results

### 189 Reference library of DNA metabarcodes

About 80 species of sand flies have been found in French Guiana, most of them belonging to the genus *Lutzomyia* (Leger et al. 1977; Young & Duncan 1994; Gantier et al. 2006; a complete checklist of French Guiana species is presented in the Supporting Information). Sequences were obtained for 178 morphologically identified specimens representing 38 species of the genus *Lutzomyia* and two species of the genus *Brumptomyia*. These included *L. (Nyssomyia) umbratilis*, *L. (Nyssomyia) flaviscutellata* and *L. (Psychodopygus) s. maripaensis*, the three proven vectors of *Leishmania* species infecting humans locally (Le Pont & Pajot 1980; Pajot et al. 1982; Killick-Kendrick 1990; Fouque et al. 2007). All the sequences were deposited in GenBank (accessions KU761608-KU761817). On average, 4.45 specimens per species were sequenced, with more than half of the species being represented by at least 3 specimens and 11 species represented by only one. The sequence length varied between 216 and 222 with an average of 218.5, close to the estimated mean length of the marker in insects (216 bp; Clarke et al. 2014). The average base composition was 37.4% adenine (A), 4.5% cytosine (C), 11.6% guanine (G) and 46.5% thymine (T) (GC content 14.8%). 47.3% of sites were identical in all sequences and the mean pairwise identity was 83.3%. Mean intraspecific and interspecific distances by species ranged from 0 to 1.7% and from 0.3 to 22.2% respectively. The
distributions of intra and interspecific K2P distances slightly overlapped (Fig. 2a). Nevertheless, for each species taken independently, the maximum intraspecific distance was smaller than the minimum interspecific distance, except for *L. (Nyssomyia) antunesi* and *L. (Nyssomyia) yulli pajoti* that were not well distinguished (Fig. 2b, Fig. 3). 98% of the specimens would have been correctly identified if assigned to the species of the closest sequence in the database (without accounting for the species represented by only one specimen). In the neighbour-joining tree generated with the 178 sequences, every species formed a monophyletic clade except *L. (N.) antunesi* and *L. (N.) yulli pajoti* (Fig 3.). All male-female pairings were confirmed and species represented by specimens from both French Guiana and Trinidad and Tobago were also grouped unambiguously. Slight differences were found in the topologies of the NJ trees constructed using the different alignment algorithms and genetic distance models. However, the clustering patterns at the species level and the results of the taxonomic assignation test were strictly identical among these methods.

## *Metabarcoding validation*

PCR amplifications were successful for all samples except for three samples extracted with the Chelex protocol that had not been diluted. After sequencing and removal of low quality reads and singletons, 2,626,783 sequencing reads were obtained, representing 59,354 distinct sequences. All expected species were detected in each sample and PCR replicate. Unexpected OTUs were also retrieved but represented only 1.6% of the reads, the majority of which was not identified as any sand fly species. Filtering of the sequences that (i) did not have a close match in the reference database, (ii) resulted from tag switching events and (iii) were not detected in at least 2 PCR replicates, removed 2.6%, 2.4% and 0.5% of the reads respectively. All unexpected OTUs were removed after these steps except in three over 22 bulk-extracted samples where one unexpected species was found. Conversely, all expected species remained detectable in every sample except *L. (Evandromyia) infraspinosa* in the mixture of ten species and *L. (T.) ininii* in one of the bulk-extracted sample whatever the DNA dilution prior to PCR. In bulk-extracted samples, read abundance ratios
revealed a mean abundance bias of 2.3 toward *L. (T.) trichopyga* ($\sigma$=1.5). The ratios were fairly consistent among samples and highly similar across PCR replicates. The results of the metabarcoding experiment are synthesized in Table 1.

A total of 415 sand flies were collected in 28 trap nights. From these 28 samples, 24 were successfully amplified and sequenced. 26 OTUs were identified of which 25 were assigned to species and one was assigned to the *Nyssomyia* subgenus. The species accumulation curve constructed from the cumulative incidence of species across trap nights reaches an asymptote (Fig. 4). The mean Chao2 estimator oscillates around 25 species along the SAC and converges to a value of 25.08 (SE=0.33). After weighting the read counts by the number of specimens in each sample, the more abundant species was *L. infraspinosa* (25.4% of reads). *L. umbratilis* and *L. flaviscutellata*, two locally proven vectors of Leishmaniasis, were also found (1.6 and 3.2% of reads respectively).

## Discussion

DNA metabarcoding represents a promising tool to overcome taxonomic impediments in assessing the sand fly fauna diversity and distribution. Because COI is not the best option for metabarcoding (Deagles et al. 2014), we assessed the efficiency of an insect metabarcode located in the mitochondrial 16S ribosomal RNA gene to discriminate between species of this medically important group. Since the distributions of intra and interspecific genetic distances slightly overlapped, the "barcoding gap" obtained is not as clear as expected for the standard COI barcoding marker (Hebert et al. 2003b). The latter has been claimed to be illusive in dipterans (Meier et al. 2006), but several evaluations of the COI marker in phlebotomine sand flies did report a barcoding gap between intra and interspecific distances (Kumar et al. 2012; Contreras Gutiérrez et al. 2014; Nzelu et al. 2015; Romero-Ricardo et al. 2016). Nevertheless, the assignment of each sequence to the species of its closest match in the database resulted in almost 98% correct identifications. This is remarkable because the choice of a metabarcode implies a compromise between its versatility and its resolution
255 (taxonomic coverage and specificity; Riaz et al. 2011). Except for *L. (N.) antunesi* and *L. (N.) yulli* 256 *pajoti*, the distance tree obtained with the mt 16S rRNA marker clustered every species studied in 257 monophyletic clades. It also unambiguously grouped conspecific individuals sampled in very distant 258 locations (French Guiana and Trinidad and Tobago). This is noteworthy assuming that intraspecific 259 genetic variation correlates with the geographic scale of the sampling (Bergsten et al. 2012). Finally, 260 the marker confirmed the associations of all conspecific males and females, highlighting the 261 reliability of traditional morphological methods in this regard. Most of the subgenera represented in 262 the dataset were not retrieved as monophyletic clades in the NJ tree. It is important to keep in mind 263 that the latter only aims at visualising the clustering of individuals belonging to the same species in 264 order to assess the reliability of the marker for taxonomic assignations based on genetic distances. It 265 should not be confused with a phylogenetic analysis and cannot be used to infer deeper taxonomic 266 relationships. 267

267 Our metabarcoding experiments showed that it was possible to accurately identify the content of 268 sand fly mixtures even in samples exhibiting a relatively high species diversity or abundance 269 imbalance. A meticulous data curation is essential to avoid significant overestimations of diversity in 270 metabarcoding studies. It represents the most critical step of analysis pipelines (Kunin et al. 2010; Yu 271 et al. 2012; Esling et al. 2015). Here, the sequence filtering process was very efficient in removing 272 false positive OTUs. Only one unexpected species was detected in three of the 22 bulk-extracted 273 samples. Because the bulk-samples were constituted of entire specimens that were morphologically 274 identified without dissection and mounting, it is difficult to exclude that these false positive OTUs 275 were due to the inclusion of misidentified specimens rather than resulting from molecular biases. 276 These results highlight the advantage of working with a relatively comprehensive reference 277 database. First, it allowed a reference-based filtering that eliminated the majority of unexpected 278 sequences. Second, defining OTUs in a reference-supervised way ensures OTU taxonomic 279 consistency. De novo OTU picking is far from being a trivial matter (Bonder et al. 2012; Sun et al. 280)
280 2012). We tested several combinations of preclustering and clustering algorithms with various parameters on our data (Sumaclust, Swarm and MCL, Van Dongen 2000; Mahé *et al.* 2014). All led to important discrepancies between the resulting de novo OTUs and taxonomic assignments of the sequences (data not shown). On the other hand, the obvious disadvantage of such a stringent reference-based analysis pipeline is that it is blind to species that are not represented in the reference database. Nevertheless, the same sequencing data can be reanalysed later on with an updated reference database. Moreover, it is possible to check for the presence of unknown OTUs in the data by relaxing the stringency of the filtering process. However, great care should be taken by using this approach, as it will inevitably increase the number of artifactual OTUs (false positives). We have looked in the field samples for sequences that had between 90% and 97% identity with their best match in the reference database and a read abundance of at least 200 in at least 2 PCR replicates. This way, we identified two OTUs that may correspond to *Lutzomyia* species not represented in the reference database: one having 94.1% identity with *L. hirsuta* and another having 94.6% identity with *L. sylvicola*.

294 Both extraction protocols led to very satisfying results. Since the Chelex-based protocol is cheap and time efficient, it represents a good candidate for routine use. Nevertheless, it lacks a purification step, which can lead to PCR amplification failure and quick degradation of the DNA in the samples, making them not suitable for long term storage. This likely explains the amplification failures in three out of the four Chelex-extracted samples that were not diluted prior to PCR. Dilution of the samples in distilled water at a minimum of 1:10 is therefore advisable for samples extracted with the Chelex protocol.

301 Several metabarcoding studies have shown significant incongruences between the genuine content of mock communities and amplicon relative abundances (Amend *et al.* 2010; Yu *et al.* 2012; Piñol *et al.* 2015). Thus, the conclusion is usually that metabarcoding is inappropriate for quantification purposes. In pool-extracted samples, read abundance ratios tended to overestimate the proportion
of *L. (T.) trichopyga* (the more abundant species in the mixture), with a mean bias of 2.3. This is may be due to an actual difference in biomass, even though *L. (T.) trichopyga* and *L. (T.) ininii* have a relatively close typical size. Variability in primer-template mismatches have been shown to result in amplification bias of up to five orders of magnitude (Bru *et al.* 2008; Piñol *et al.* 2015). However, the complete mitochondrial genome of several sand fly species, including *L. (T.) ininii* and *L. (T.) trichopyga*, were sequenced, and all exhibit 100% similarity at the *Ins16S_1* primers binding sites (Kocher *et al.* 2015, Ye *et al.* 2015 and further unpublished data). Various others technical and biological factors can impact the quantitative signal and lead to discrepancies (Pompanon *et al.* 2012), but this is beyond the scope of this paper. However, the bias found here was less than one order of magnitude and was fairly consistent across samples and PCR replicates. Relative read abundance may therefore be indicative when studying sand fly mixtures, at least for inter-sample comparisons.

The results obtained on our field samples indicate that metabarcoding was efficient for characterizing the local sand fly community. The species accumulation curve reaches an asymptote and the *Chao2* estimator converges to a value extremely close to the observed number of species, giving good confidence in the comprehensiveness of the sampling. In total, 25 distinct species were identified in the plot, which is comparable to the results of other studies conducted in the region (Rotureau *et al.* 2006; Fouque *et al.* 2007). As mentioned above, great caution has to be taken for quantitative interpretation of metabarcoding. Assuming that PCR amplification bias is consistent across samples, the results presented in Fig. 4b can be regarded as a proxy of species composition for between site comparisons, but are not to be confused with real proportions of specimens.

The use of metabarcoding allowed us to obtain these results with little laboratory work (28 samples processed instead of 415 if we had treated each specimen individually). However, the number of sand flies caught here was relatively low. The added value of metabarcoding in terms of cost and time-effectiveness depends on the number of specimens pooled prior to bulk DNA extraction and on
the number of samples multiplexed for high-throughput sequencing. Using the present protocol, we estimate that 10,000 specimens can be analysed by a single person in 20 days and for a total cost of less than 5,000 euros. As a mater of comparison, just mounting those specimens for morphological identification would take more than 200 days (considering 10 min per specimen and 8 hours of work per day).

In conclusion, this study highlight that the insect metabarcode developed by Clarke et al. (2014) is reliable for species-level assignations of phlebotomine sand flies and that metabarcoding using this marker can accurately identify the composition of sand fly mixtures. The application of the present protocol to field samples can significantly improve the cost and time effectiveness of data generation on sand fly fauna. Moreover, the primers used here have been estimated to amplify over 90% of insect taxa among the publicly available data at the time of their design and can therefore be employed widely. With adequate reference databases, this approach could be extended to concomitantly analyse other insect vectors communities such as mosquitoes, opening new perspectives for eco-epidemiological studies.

## Acknowledgments

This work was supported by "Investissement d'Avenir" grants managed by Agence Nationale de la Recherche (CEBA: ANR-10-LABX-25-01; TULIP: ANR-10-LABX-41, ANR-11-IDEX-0002-02) as well as project METABAR (ANR-11-BSV7-0020). The collection of sand flies was performed in the framework of the French National Project IAEL: ANR-06-SEST-20. We would like to thank Romuald Carinci, Jean Issaly, Virginie Rougeron and Benoit de Thoisy for their involvement in the collection and processing of the samples, Pierre Taberlet and Eric Coissac for fruitful discussions about metabarcoding and Pieternel Verschuren for her kind revision of the English text.

## References

353 Abonnenc E (1972) Les phlébotomes de la région éthiopienne (Diptera, Psychodidae). *Mémoires ORSTOM*, **55**.

354 Al-Ajmi RA, Ayaad TH, Al-Anazi M, Al-Ahdal MN, Al-Qahtani AA (2015) Molecular identification of natural sand fly
355 species populations inferred from ITS2 rRNA gene in Saudi Arabia. *Journal of Environmental Biology*, **36**.

356 Alvar J, Vélez ID, Bern C et al. (2012) Leishmaniasis Worldwide and Global Estimates of Its Incidence. *PLoS ONE*, **7**, e35671.

358 Amend AS, Seifert KA, Bruns TD (2010) Quantifying microbial communities with 454 pyrosequencing: does read abundance count? *Molecular Ecology*, **19**, 5555–5565.

360 Aransay AM, Scoulica E, Chaniotis B, Tselentis Y (1999) Typing of sandflies from Greece and Cyprus by DNA polymorphism of 18S rRNA gene. *Insect Molecular Biology*, **8**, 179–184.

362 Bergsten J, Bilton DT, Fujisawa T et al. (2012) The Effect of Geographical Scale of Sampling on DNA Barcoding. *Systematic Biology*, sys037.

364 Binladen J, Gilbert MTP, Bollback JP et al. (2007) The Use of Coded PCR Primers Enables High-Throughput Sequencing of Multiple Homolog Amplification Products by 454 Parallel Sequencing. *PLoS ONE*, **2**.

366 Bonder MJ, Abeln S, Zaura E, Brandt BW (2012) Comparing clustering and pre-processing in taxonomy analysis. *Bioinformatics*, **28**, 2891–2897.

368 Boyer F, Mercier C, Bonin A et al. (2016) obitools: a unix-inspired software package for DNA metabarcoding. *Molecular Ecology Resources*, **16**, 176–182.

370 Brown SD, Collins RA, Boyer S et al. (2012) Spider: An R package for the analysis of species identity and evolution, with particular reference to DNA barcoding. *Molecular Ecology Resources*, **12**, 562–565.

372 Bru D, Martin-Laurent F, Philippot L (2008) Quantification of the detrimental effect of a single primer-template mismatch by real-time PCR using the 16S rRNA gene as an example. *Applied and environmental microbiology*, **74**, 1660–1663.

375 Casquet J, Thebaud C, Gillespie RG (2012) Chelex without boiling, a rapid and easy technique to obtain stable amplifiable DNA from small amounts of ethanol-stored spiders. *Molecular Ecology Resources*, **12**, 136–141.

377 Clarke LJ, Soubrier J, Weyrich LS, Cooper A (2014) Environmental metabarcodes for insects: *in silico* PCR reveals potential for taxonomic bias. *Molecular Ecology Resources*, **14**, 1160–1170.

379 Contreras Gutiérrez MA, Vivero RJ, Vélez ID, Porter CH, Uribe S (2014) DNA Barcoding for the Identification of Sand Fly Species (Diptera, Psychodidae, Phlebotominae) in Colombia. *PLoS ONE*, **9**, e85496.

381 Deagle BE, Jarman SN, Coissac E, Pompanon F, Taberlet P (2014) DNA metabarcoding and the cytochrome c oxidase subunit I marker: not a perfect match. *Biology Letters*, **10**, 20140562.

383 Depaquit J (2014) Molecular systematics applied to Phlebotomine sandflies: Review and perspectives. *Infection*, **28**, 744–756.

385 Depaquit J, Grandadam M, Fouque F, Andry PE, Peyrefitte C (2010) Arthropod-borne viruses transmitted by Phlebotomine sandflies in Europe: a review. *Euro Surveill*, **15**, 19507.

387 Desjeux P (2001) The increase in risk factors for leishmaniasis worldwide. *Transactions of The Royal Society of Tropical Medicine and Hygiene*, **95**, 239–243.

389 Edgar RC (2004) MUSCLE: multiple sequence alignment with high accuracy and high throughput. *Nucleic acids research*, **32**, 1792–1797.

391 El-Hossary S (2006) Morphological characteristics for sand fly taxonomy. *Research and Training Center on Vector of Diseases, Ain Shams University, Ain Shams*, 1–25.

393 Esling P, Lejzerowicz F, Pawlowski J (2015) Accurate multiplexing and filtering for high-throughput amplicon-sequencing. *Nucleic Acids Research*, **43**, 2513–2524.

395 Gantier JC, Gaborit P, Rabarison P (2006) Les phlébotomes (Diptera: Psychodidae) de la Guyane Française. Description du mâle de *Lutzomyia (Trichopygomyia) depaquiti* n. sp. *Parasite*, **13**, 11–15.

397 Galati, E. A. B. (1995). Phylogenetic systematics of Phlebotominae (Diptera, Psychodidae) with emphasis on American groups. *Bol Dir Malariol San Amb*, **35**(Supl 1), 133-142.

399 Galati EAB (2014) Phlebotominae (Diptera, Psychodidae) Classificação, Morfologia, Terminologia e Identificação de Adultos. Apostila. Departamento de Epidemiologia Faculdade de Saúde Pública Universidade de São Paulo, Brasil, 133 pp.

402 Gibson J, Shokralla S, Porter TM et al. (2014) Simultaneous assessment of the macrobiome and microbiome in a bulk sample of tropical arthropods through DNA metasystematics. *Proceedings of the National Academy of Sciences*, **111**, 8007–8012.

405 Hebert PD, Cywinska A, Ball SL, others (2003a) Biological identifications through DNA barcodes. *Proceedings of the Royal Society of London B: Biological Sciences*, **270**, 313–321.

407 Hebert PDN, Ratnasingham S, Waard JR de (2003b) Barcoding animal life: cytochrome c oxidase subunit 1 divergences among closely related species. *Proceedings of the Royal Society of London B: Biological Sciences*, **270**, S96–S99.
410 Huse SM, Welch DM, Morrison HG, Sogin ML (2010) Ironing out the wrinkles in the rare biosphere through
411 improved OTU clustering. *Environmental Microbiology*, **12**, 1889–1898.

412 Ji Y, Ashton L, Pedley SM *et al*. (2013) Reliable, verifiable and efficient monitoring of biodiversity via
413 metabarcoding. *Ecology Letters*, **16**, 1245–1257.

414 Jukes TH, Cantor CR (1969) Evolution of protein molecules. *Mammalian protein metabolism*, **3**, 132.

415 Killick-Kendrick R (1990) Phlebotomine vectors of the leishmaniases: a review. *Medical and Veterinary
416 Entomology*, **4**, 1–24.

417 Kimura M (1980) A simple method for estimating evolutionary rates of base substitutions through comparative
418 studies of nucleotide sequences. *Journal of molecular evolution*, **16**, 111–120.

419 Kocher A, Gantier J-C, Holota H *et al*. (2015) Complete mitochondrial genome of *Lutzomyia* (*Nyssomyia*)
420 *umbratilis* (Diptera: Psychodidae), the main vector of *Leishmania guyanensis*. *Mitochondrial DNA*, 1–3.

421 Krüger A, Strüven L, Post RJ, Faulde M (2011) The sandflies (Diptera: Psychodidae, Phlebotominae) in military
422 camps in northern Afghanistan (2007–2009), as identified by morphology and DNA "barcoding". *Annals of
423 Tropical Medicine & Parasitology*, **105**, 163–176.

424 Kumar NP, Srinivasan R, Jambulingam P (2012) DNA barcoding for identification of sand flies (Diptera:
425 Psychodidae) in India. *Molecular Ecology Resources*, **12**, 414–420.

426 Kunin V, Engelbrektson A, Ochman H, Hugenholtz P (2010) Wrinkles in the rare biosphere: pyrosequencing errors
427 can lead to artificial inflation of diversity estimates. *Environmental Microbiology*, **12**, 118–123.

428 Larkin MA, Blackshields G, Brown NP *et al*. (2007) Clustal W and Clustal X version 2.0. *Bioinformatics*, **23**, 2947–
429 2948.

430 Latrofa MS, Annoscia G, Dantas-Torres F, Traversa D, Otranto D (2012) Towards a rapid molecular identification
431 of the common phlebotomine sand flies in the Mediterranean region. *Veterinary parasitology*, **184**, 267–270.

432 Leger N, Abonnenc E, Pajot F-X, Kramer R, Claustre J (1977) Liste commentée des phlébotomes de la Guyane
433 Française. *Cahiers ORSTOM. Série Entomologie Médicale et Parasitologie*, **15**, 217–232.

434 Lightburn E, Meynard JB, Morand J et al. (2001) Surveillance épidémiologique des leishmanioses tégumentaires
435 en Guyane Synthèse de données militaires sur 10 ans. *Méd Trop*, **62**, 545–553.

436 Mahé F, Rognes T, Quince C, De Vargas C, Dunthorn M (2014) Swarm: robust and fast clustering method for
437 amplicon-based studies. *PeerJ*, **2**, e593.

438 Maia C, Parreirab R, Cristóvão JM, Afonso MO, Campino L (2015) Exploring the utility of phylogenetic analysis of
439 cytochrome oxidase gene subunit I as a complementary tool to classical taxonomical identification of
440 phlebotomine sand fly species (Diptera, Psychodidae) from southern Europe. *Acta Tropica*, **144**, 1–8.

441 Mariette J, Escudié F, Allias N *et al*. (2012) NG6: Integrated next generation sequencing storage and processing
442 environment. *BMC genomics*, **13**, 462.

443 Maroli M, Feliciangeli MD, Bichaud L, Charrel RN, Gradoni L (2013) Phlebotomine sandflies and the spreading of
444 leishmaniases and other diseases of public health concern. *Medical and Veterinary Entomology*, **27**, 123–147.

445 Meier R, Shiyang K, Vaidya G, Ng PK (2006) DNA barcoding and taxonomy in Diptera: a tale of high intraspecific
446 variability and low identification success. *Systematic biology*, **55**, 715–728.

447 Mercier C, Boyer F, Bonin A, Coissac E (2013) SUMATRA and SUMACLUST: fast and exact comparison and
448 clustering of sequences. In: *Programs and Abstracts of the SeqBio 2013 workshop. Abstract*, pp. 27–29.

449 Minter LM, Yu T, Florin DA *et al*. (2013) Molecular Identification of Sand Flies (Diptera: Psychodidae) in Eastern
450 North America by Using PCR-RFLP. *Journal of Medical Entomology*, **50**, 920–924.

451 Munstermann LE (2004) Phlebotomine sand flies, the Psychodidae. *Biology of disease vectors*, 141–151.

452 Nzelu CO, Cáceres AG, Arrunátegui-Jiménez MJ *et al*. (2015) DNA barcoding for identification of sand fly species
453 (Diptera: Psychodidae) from leishmaniasis-endemic areas of Peru. *Acta tropica*, **145**, 45–51.

454 de Oliveira Miranda DE, da Silva Sales KG, da Gloria Faustino MA *et al*. (2015) Ecology of sand flies in a low-
455 density residential rural area, with mixed forest/agricultural exploitation, in north-eastern Brazil. *Acta tropica*,
456 **146**, 89–94.

457 Oksanen J, Blanchet FG, Kindt R *et al*. (2016). vegan: Community Ecology Package. R package version 2.3-3.

458 Ondrejicka DA, Locke SA, Morey K, Borisenko AV, Hanner RH (2014) Status and prospects of DNA barcoding in
459 medically important parasites and vectors. *Trends in Parasitology*, **30**, 582–591.

460 Peçanha Brazil R, Fuzzari Rodrigues AA, Andrade Filho JD (2015) Sand Fly Vectors of *Leishmania* in the Americas -
461 A Mini Review. *Entomology, Ornithology & Herpetology: Current Research*, **4**.

462 Piñol J, Mir G, Gomez-Polo P, Agustí N (2015) Universal and blocking primer mismatches limit the use of high-
463 throughput DNA sequencing for the quantitative metabarcoding of arthropods. *Molecular Ecology Resources*,
464 **15**, 819–830.
465 Pompanon F, Deagle BE, Symondson WOC et al. (2012) Who is eating what: diet assessment using next generation sequencing. *Molecular Ecology*, **21**, 1931–1950.

467 Qiu X, Wu L, Huang H et al. (2001) Evaluation of PCR-generated chimeras, mutations, and heteroduplexes with 16S rRNA gene-based cloning. *Applied and Environmental Microbiology*, **67**, 880–887.

469 R Development Core Team (2014) R: A language and environment for statistical computing. *R Foundation for Statistical Computing*, Vienna, Austria, 2012. ISBN 3-900051-07-0.

471 Ready PD (2008) Leishmaniasis emergence and climate change. *Revue Scientifique Et Technique-Office International Des Epizooties*, **27**.

473 Ready PD (2013) Biology of Phlebotomine Sand Flies as Vectors of Disease Agents. *Annual Review of Entomology*, **58**, 227–250.

475 Riaz T, Shehzad W, Viari A et al. (2011) ecoPrimers: inference of new DNA barcode markers from whole genome sequence analysis. *Nucleic Acids Research*, **39**, e145–e145.

477 Romero-Ricardo L, Lastre-Meza N, Pérez-Doria A, Bejarano EE (2016) DNA barcoding to identify species of phlebotomine sand fly (Diptera: Psychodidae) in the mixed leishmaniasis focus of the Colombian Caribbean.

479 *Acta Tropica*, **159**, 125–131.

480 Rotureau B, Gaborit P, Issaly J et al. (2006) Diversity and ecology of sand flies (Diptera: Psychodidae: Phlebotominae) in coastal French Guiana. *American Journal of Tropical Medicine and Hygiene*, **75**, 62.

482 Saitou N, Nei M (1987) The neighbor-joining method: a new method for reconstructing phylogenetic trees. *Molecular biology and evolution*, **4**, 406–425.

484 Stocsits RR, Letsch H, Hertel J, Misof B, Stadler PF (2009) Accurate and efficient reconstruction of deep phylogenies from structured RNAs. *Nucleic Acids Research*, **37**(18), 6184–6193.

486 Sun Y, Cai Y, Huse SM et al. (2012) A large-scale benchmark study of existing algorithms for taxonomy-independent microbial community analysis. *Briefings in Bioinformatics*, **13**, 107–121.

488 Taberlet P, Coissac E, Pompanon F, Brochmann C, Willerslev E (2012) Towards next-generation biodiversity assessment using DNA metabarcoding. *Molecular ecology*, **21**, 2045–2050.

490 Tamura K, Nei M (1993) Estimation of the number of nucleotide substitutions in the control region of mitochondrial DNA in humans and chimpanzees. *Molecular biology and evolution*, **10**, 512–526.

492 Terayama Y, Kato H, Gomez EA et al. (2008) Molecular Typing of Sand Fly Species (Diptera, Psychodidae, Phlebotominae) from Areas Endemic for Leishmaniasis in Ecuador by PCR-RFLP of 18S Ribosomal RNA Gene. *Journal of Veterinary Medical Science*, **70**, 907–913.

495 Tiwary P, Kumar D, Rai M, Sundar S (2012) PCR-RFLP Based Method for Molecular Differentiation of Sand Fly Species *Phlebotomus argentipes*, *Phlebotomus papatasi* and *Sergentomyia babu* Found in India. *Journal of Medical Entomology*, **49**, 1515–1518.

498 Van Dongen S (2000) A cluster algorithm for graphs. *Report-Information systems*, 1–40.

499 Ye F, Liu T, King SD, You, P (2015). Mitochondrial genomes of two phlebotomine sand flies, *Phlebotomus chinensis* and *Phlebotomus papatasi* (Diptera: Nematocera), the first representatives from the family Psychodidae.
501 *Parasites & vectors*, **8**(1), 1–13.

502 Young DG, Duncan MA (1994) *Guide to the identification and geographic distribution of Lutzomyia sand flies in Mexico, the West Indies, Central and South America (Diptera: Psychodidae)*. Gainesville, FL: Associated Publishers.

505 Yu DW, Ji Y, Emerson BC et al. (2012) Biodiversity soup: metabarcoding of arthropods for rapid biodiversity assessment and biomonitoring. *Methods in Ecology and Evolution*, **3**, 613–623.

507 **Data accessibility**

508 *Sand flies specimens*: the specimens used for the DNA reference database are kept as voucher at the

509 LIFE laboratory (Mennecy, France). Pictures can be sent on demand.
510 **DNA reference database:** All sequences were published in GenBank (accessions KU761608-KU761817). The sequences are also available in the OBITools extended fasta format together with the corresponding taxonomic database in the Supporting information (S4).

513 **Metabarcoding data:** sequencing data of mock communities and field samples is available in the Supporting information (S5). We provide the data in fasta format after dereplication and removal of singletons and in tabular format after taxonomic assignment of the sequences.

516 **Analysis pipeline:** bash and R script used for the bioinformatic treatment of the data are available in the Supporting Information (S6 and S7).

518 **Author contributions**

519 A.K. and J.M. designed the study. JC.G., P.G., I.D., R.G. and AL.B. were involved in the collection of samples on the field. JC.G. and P.G. conducted the morphological identifications. A.K., H.H. and S.V.
520
521 performed the laboratory work. A.K., J.M. and L.Z. participated in the analysis of the data. A.K. wrote
522 the manuscript and all authors contributed in its improvement.

523 **Supporting Information**

524 Additional Supporting Information may be found in the online version of this article.

525 **S1:** Checklist of the phlebotomine sand fly species reported in French Guiana. The classifications of Young & Duncan (1994) and Galati (1995, 2014) are presented.

527 **S2:** List of the samples used for the metabarcode reference database along with the GenBank accession numbers of the corresponding sequences

529 **S3:** DNA extraction protocols for pools of sand flies

530 **S4:** DNA reference database in OBITools extended fasta format and modified NCBI taxonomic database in ecopcr format

532 **S5:** Metabarcoding data for mock communities and field samples. We included the sequencing data in fasta format after dereplication and removal of singletons, as well as in tabular format after taxonomic assignment of the sequences. The sample description files used for the assignment of the sequences to their corresponding samples with the OBITools are also included

536 **S6:** Bash script used for the bioinformatic treatment of the sequencing data

537 **S7:** R script used for the filtering of metabarcoding data and output as community matrices
539 Tables

540 Table 1: Results of the metabarcoding experiment. Mean numbers of false negative (FN) and false

541 positive (FP) OTUs are presented for each type of sample. For bulk-extracted samples, L. (T.)

542 trichopyga vs L. (T.) ininii read abundance ratios obtained after complete filtering process are

543 presented along with the resulting bias toward L. (T.) trichopyga.

<table><thead><tr><th rowspan="2"></th><th colspan="2">No filtering</th><th colspan="2">Reference-based filtering</th><th colspan="2">Mistags-based filtering</th><th colspan="2">Intersection of PCR replicates</th><th rowspan="2">Reads abundance bias<sup>4</sup></th></tr><tr><th>FN</th><th>FP</th><th>FN</th><th>FP</th><th>FN</th><th>FP</th><th>FN</th><th>FP</th></tr></thead><tbody><tr><td>Mix 10<sup>1</sup></td><td>0.0</td><td>7.0</td><td>0.0</td><td>2.0</td><td>0.0</td><td>0.0</td><td>1.0</td><td>0.0</td><td>NI<sup>5</sup></td></tr><tr><td>Bulk 10<sup>2</sup></td><td>0.0</td><td>5.7</td><td>0.0</td><td>2.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>2.2</td></tr><tr><td>Bulk 50<sup>3</sup></td><td>0.0</td><td>6.2</td><td>0.0</td><td>3.0</td><td>0.0</td><td>0.6</td><td>0.3</td><td>0.3</td><td>2.5</td></tr></tbody></table>

544

545 ¹mixture of DNA extracted from ten distinct sand fly species.

546 ²bulk-extracted samples including 2 species with a specimen ratio of 10

547 ³bulk-extracted samples including 2 species with a specimen ratio of 50

548 ⁴abundance bias: bias found between the relative abundance of L. (T.) trichopyga and L. (T.) ininii included in the mock

549 mixtures and corresponding read abundance ratio

550 ⁵non-informative
<figure><img src="image_2.png" /><figcaption>Figure 1: Schematization of the metabarcoding experiment workflow.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>Figure 2: (A) Distribution of intraspecific (in blue) and interspecific (in red) genetic distances among sand fly samples. (B) Maximum intraspecific genetic distance against minimum interspecific distance for each sample. The identity line is plotted.</figcaption></figure>

554

555

556

557

558
<figure><img src="image_4.png" /></figure>

559
560
561
562

Figure 3: Neighbour-joining tree based on Kimura-2-parameter pairwise distances of the 16S metabarcode from sand fly specimens. Numbers in brackets indicate the number of specimens sequenced per species. Bootstrap support values higher than 50% are indicated on the nodes
<figure><img src="image_5.png" /><figcaption>Figure 4: Results of metabarcoding on field samples. (A) Species accumulation curve (black solid line) based on the incidence of species trap nights constructed with the 'exact' method along with its confidence interval (black dotted lines). The solid red line indicates the mean Chao2 estimator computed from 100 permutations. (B) Overall OTU composition after weighting the read counts by the number of specimens in each sample.</figcaption></figure>

563

564

565

566

567

568