Next Article in Journal
Combined Effects of Soil Silicon and Host Plant Resistance on Planthoppers, Blast and Bacterial Blight in Tropical Rice
Next Article in Special Issue
Vast Gene Flow among the Spanish Populations of the Pest Bactrocera oleae (Diptera, Tephritidae), Phylogeography of a Metapopulation to Be Controlled and Its Mediterranean Genetic Context
Previous Article in Journal
Ultrastructure of Antennal Sensory Organs in Nine Flesh Flies (Diptera: Sarcophagidae): New Insight into the Definition of Family Sarcophagidae
Previous Article in Special Issue
Fruit Fly Larval Survival in Picked and Unpicked Tomato Fruit of Differing Ripeness and Associated Gene Expression Patterns
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Insects
Volume 13
Issue 7
10.3390/insects13070603
insects-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

The Impact of Fast Radiation on the Phylogeny of Bactrocera Fruit Flies as Revealed by Multiple Evolutionary Models and Mutation Rate-Calibrated Clock

by
Federica Valerio
1,
Nicola Zadra
2,3,
Omar Rota-Stabelli
2,3,4 and
Lino Ometto
1,*
1
Department of Biology and Biotechnology, University of Pavia, 27100 Pavia, Italy
2
Research and Innovation Centre, Fondazione Edmund Mach, 38010 San Michele all’Adige, Italy
3
Department of Cellular, Computational and Integrative Biology—CIBIO, University of Trento, 38123 Trento, Italy
4
Center Agriculture Food Environment—C3A, University of Trento, 38010 San Michele all’Adige, Italy
*
Author to whom correspondence should be addressed.
Insects 2022, 13(7), 603; https://doi.org/10.3390/insects13070603
Submission received: 29 April 2022 / Revised: 22 June 2022 / Accepted: 28 June 2022 / Published: 30 June 2022
(This article belongs to the Special Issue Genetics and Ecological Evolution of Dipteran Pest Species)
Browse Figures
Review Reports Versions Notes

Abstract

:

Simple Summary

Comparative approaches are widely used to investigate the traits that underlie particular biological and ecological traits. They are effective, however, only if we know the correct relationships among the species under consideration. Here, we aimed at reconstructing a robust phylogeny for selected true fruit flies of the genus Bactrocera, which are well known as important agricultural pests worldwide. Existing phylogenetic inferences are still ambiguous, especially concerning the relationship between the two major pests B. dorsalis and B. tryoni. In this study, we employed a genome-scaled dataset and different models of molecular evolution to reconstruct the phylogenetic relationships of ten Bactrocera species and two outgroups and further date their divergence times. The resulting phylogeny fully supports B. dorsalis as more closely related to B. latifrons than to B. tryoni, opposite to what was supported by previous works. This incongruence likely derives from a fast divergence of these lineages, as revealed by our clock analysis, which can lead to conflicting results when using few genetic markers. Our results thus highlight the utility of using large datasets and of exploring different evolutionary models to study the evolutionary history of species of economic importance.

Abstract

Several true fruit flies (Tephritidae) cause major damage to agriculture worldwide. Among them, species of the genus Bactrocera are extensively studied to understand the traits associated with their invasiveness and ecology. Comparative approaches based on a reliable phylogenetic framework are particularly effective, but several nodes of the Bactrocera phylogeny are still controversial, especially concerning the reciprocal affinities of the two major pests B. dorsalis and B. tryoni. Here, we analyzed a newly assembled genomic-scaled dataset using different models of evolution to infer a phylogenomic backbone of ten representative Bactrocera species and two outgroups. We further provide the first genome-scaled inference of their divergence by calibrating the clock using fossil records and the spontaneous mutation rate. The results reveal a closer relationship of B. dorsalis with B. latifrons than to B. tryoni, contrary to what was previously supported by mitochondrial-based phylogenies. By employing coalescent-aware and heterogeneous evolutionary models, we show that this incongruence likely derives from a hitherto undetected systematic error, exacerbated by incomplete lineage sorting and possibly hybridization. This agrees with our clock analysis, which supports a rapid and recent radiation of the clade to which B. dorsalis, B. latifrons and B. tryoni belong. These results provide a new picture of Bactrocera phylogeny that can serve as the basis for future comparative analyses.

1. Introduction

Tephritidae, commonly known as true fruit flies, are an incredibly diverse group of phytophagous insects. This family includes more than 5000 species in over 500 genera [1], making it one of the largest among Diptera. Several of these species cause extensive damage to agriculture worldwide, since females lay eggs within the fruit flesh (i.e., pericarp) and the larvae rapidly develop in the fruit, eating it and inducing bacterial and fungal decay. The accompanying economic impacts are considerable: for instance, the olive fruit fly, Bactrocera oleae (Rossi), has been estimated to cause an annual economic loss of approximately $800 million in the Mediterranean basin alone [2]. Moreover, such impact is exacerbated by these species’ ability to invade new areas and establish viable populations. For example, Bactrocera dorsalis and Bactrocera latifrons, both native to South-East Asia, have recently expanded their distribution range in the Hawaiian Islands [3,4], in several sub-Saharan African countries [5,6,7,8,9,10,11,12], and also in Europe [9,13,14], where they are posing an immediate or potential threat to local agriculture.
Most relevant tephritid pests belong to five genera: Anastrepha Schiner, Ceratitis MacLeay, Rhagoletis Loew, Zeugodacus Hendel and Bactrocera Macquart. The genus Bactrocera is the most economically important and comprises over 450 species [15], at least 50 of which are considered important pests [16,17]. Identifying the key ecological and biological traits relevant for their invasiveness and host preference is extremely important to drive effective control strategies. To this aim, comparative approaches based on a reliable phylogenetic framework are particularly effective, as they allow to trace the evolution of species-specific and shared traits and evaluate whether and to what extent control measures may be applied in related species [18]. The Bactrocera genus is subdivided into subgenera, species groups and species complexes, and the relationships between these groups have been repeatedly revised (see, e.g., [19,20,21,22,23,24,25]). Different molecular phylogenies of Bactrocera resulted in different relationships depending on the type (nuclear and/or mitochondrial) and number of markers, and on the number of taxa analyzed [23,24,25,26,27,28,29]. For instance, based on mtDNA data, B. dorsalis is usually considered to be more closely related to B. tryoni than to B. latifrons (e.g., [30,31]), but recent studies suggested a poorly supported closer relationship between B. dorsalis and B. latifrons [24,26,32] (Figure 1a). This uncertainty may also be important for studies relevant for their management, since different phylogenetic reconstructions may affect the interpretation of comparative morphological and genomic studies that focus on pest-related biological traits.
Here, we employed a rigorous phylogenomic approach to infer the phylogeny of ten representative Bactrocera species, including seven of the nine major agricultural pests of worldwide importance (see Table 1 in [16]), and we calibrated a relaxed clock to estimate their divergences. Rather than resolving the full Bactrocera phylogeny, our aim is to provide a robust phylogenetic and chronological backbone that can clarify the yet unresolved (B. dorsalis, B. tryoni, B. latifrons) clade, and which can then serve as basis for powerful comparative analyses. Our dataset also includes species characterized by distinct diets, which range from polyphagy (e.g., B. dorsalis) to monophagy (e.g., B. oleae), and, therefore, understanding the evolutionary history of this group may also provide precious information on the mode and tempo of the evolution of their biology and ecology.

2. Materials and Methods

2.1. Datasets

In our analysis, we included Bactrocera and Dacini species for which genomic and/or transcriptomic resources were available: Bactrocera jarvisi (Tryon), Bactrocera oleae (Rossi), Bactrocera minax (Enderlein), Bactrocera bryoniae (Tryon), Bactrocera correcta (Bezzi), Bactrocera dorsalis (Hendel), Bactrocera latifrons (Hendel), Bactrocera musae (Tryon), Bactrocera tryoni (Froggatt), Bactrocera zonata (Saunders) and Zeugodacus cucurbitae (Coquillett). In particular, we analyzed datasets of coding sequences (CDS) and of RNA (transcript) sequences for Ceratitis capitata [33], Z. cucurbitae [34], B. dorsalis [35], B. latifrons (NCBI accession: MIMC00000000.1), B. minax [36] and B. oleae [37]. For the six remaining Bactrocera species (B. bryoniae, B. correcta, B. jarvisi, B. musae, B. tryoni and B. zonata), we downloaded the available RNA-Seq raw reads (see Table S1 for SRA accession numbers) and assembled the corresponding transcriptomes using default parameters with Trinity v. 2.7.0 [38].

2.2. Orthologous Gene Set Identification and Alignment

Orthologs across the ten Bactrocera species and the outgroups Z. cucurbitae and C. capitata were identified using a reciprocal-best-hit approach using pairwise BLASTn searches [39] between the C. capitata CDS sequences and each of the other datasets. Putative 1:1 orthologs were first aligned using MAFFT [40] and any incomplete codon (based on the C. capitata sequence) was removed. We then re-aligned the ortholog sets using the PRANK algorithm [41] implemented in the tool TranslatorX [42]. We minimized bias in our datasets by 1) removing alignments containing sequences with internal stop codons and 2) using a custom perl script to remove problematic and ambiguous alignment regions [43]. Using this pipeline, we ultimately identified 110 orthologous gene sets across all twelve species (single alignments are available as Supplementary File S1).
We concatenated orthologs to generate an alignment of 189,891 nucleotides (nts) and then translated it using the standard genetic code to obtain an alignment of 63,297 amino acids (aas). We further generated an alignment of 24,885 nts containing only 4-fold degenerate sites retrieved using MEGA7 [44].

2.3. Phylogenetic Analyses

We inferred phylogenetic relationships using both a maximum likelihood (ML) and a Bayesian approach. To explore possible sources of systematic errors, we employed homogenous, heterogeneous and coalescent-aware models of evolution.
We ran ML analyses on both the concatenated aa alignment using the PROTGAMMAGTR model, and on the nt alignment partitioned in first, second and third (1 + 2 + 3) position using the GAMMAGTR model implemented in RAxML [45]. We also ran an ML analysis based on the 4-fold degenerate site alignment using the GAMMAGTR model. In all cases, node support was calculated by the rapid bootstrap feature of RAxML (100 replicates). We also estimated bootstrap supports using the coalescent-aware analysis of ASTRAL [46], which was based on all single ML gene trees obtained by RAxML using the same models of the concatenated analyses for either the nt or the aa sequences. Bootstrap values were estimated by performing either 100 multi-locus bootstrap replicates or gene + site resampling (using the -g option).
The same three datasets, aas, nts (1 + 2 + 3), and 4-fold degenerate sites, were used to run Bayesian analyses in BEAST v. 2.5.1 [47]. The aminoacidic dataset was analyzed with a LG + G4 substitution model. The 4-fold degenerate site dataset and the codon dataset were analyzed with a GTR + G4 substitution model. The codon dataset was split into three partitions, corresponding to the codon positions, setting linked trees across them. We performed additional Bayesian analysis on the aminoacidic dataset using a CAT + GTR model and on the 4-fold degenerate site dataset using the among-site heterogeneous CAT model with gamma distribution with PhyloBayes [48]. In both cases, we ran two independent MCMC chains and checked for convergence using the associated tracecomp and bpcomp commands. For the aminoacidic dataset we let both chains run until parameters were stabilized with maxdiff = 0, reldiff < 0.2 (except for loglik, alpha, stat and rrent, with reldiff between 0.62 and 0.77) and the effective sample size was effsize > 170 (except for loglik, alpha, stat, rrent and allocent, with effsize between 16 and 119). For the 4-fold degenerate site dataset, both chains ran until parameters were stabilized with maxdiff = 0, reldiff < 0.25 (except for stat, with reldiff < 0.4) and the effective sample size was effsize > 170 (except for nmode, statalpha, kappa and allocent, with effsize between 17 and 85). The summary statistics reached stability at the end of the runs and were periodically visualized with the script graphphylo (https://github.com/wrf/graphphylo accessed on 15 April 2022).
Bayesian analyses were also ran using StarBEAST2 [49], which employs a multispecies coalescent method to estimate species trees from multiple sequence alignments, as implemented in the BEAST2 package and according to the tutorial provided by Taming the BEAST2 [50]. For this analysis we used either the nucleotide or the aminoacidic alignments, linking the site models across the gene sets.
BEAST2 and StarBEAST2 analyses had chains run for 108 generations (or until convergence for at least 4 × 107 generations), sampling trees and parameters every 1000 generations and inspecting convergence and likelihood plateauing in Tracer. Posterior consensus trees were generated after discarding the first 10% of generations as burn-in. For the StarBEAST2 analysis, single-gene trees were loaded and visualized by DensiTree [51].

2.4. Dating Analysis

Because of the numerous incongruences between the species tree obtained by the multi-locus analyses and the single-gene trees (see Section 3), which could bias the dating analysis, we produced a conservative dataset by: (i) limiting the species sample to 10 representative species (C. capitata, Z. cucurbitae, B. bryoniae, B. dorsalis, B. jarvisi, B. latifrons, B. minax, B. musae, B. oleae, B. tryoni) and (ii) considering only those 37 genes that produced a ML tree supporting the consensus species tree with minimum ML bootstrap values of 50 at each node.
Divergence times were then estimated by BEAST v. 2.5.1 using the 4-fold degenerate sites of the resulting concatenated dataset (11,768 nt). We calibrated the crown of the Dacinae (Ceratitis-Bactrocera split) by setting a root prior uniformly distributed between 6 and 65 million years ago, which correspond to the age of a Ceratitis fly fossil [52] and of the Schizophora radiation [53,54], respectively. The exact placement of fossils is often disputed and, in general, fossils are also misused because phylogenetic models generally use them to calibrate nodes even though they represent extinct lineages on the stem and not the exact ancestors of the extant ones. In light of this, the mutation rate represents an alternative valid calibration. Indeed, the 4-fold degenerate site dataset also allowed us to use a spontaneous (neutral) mutation rate as an additional prior. Since mutation rate in Tephritidae is not known yet, we assumed it to be similar to that of Drosophila (another Diptera) and used the estimate of 0.0346 (SD = 0.0028) substitutions per base pair per million years, provided in [55,56]. This assumption is reasonable, as the mutation rate of insects is quite conserved, with the spontaneous mutation rate of Heliconius melpomene (Lepidoptera) being very close to that of D. melanogaster (2.9 × 10−9 and 3 × 10−9, respectively [57,58]). In addition, the large SD derives from the error associated with the mutation rate estimate and, thus, also accounts for the between-branch variation in substitution rates. Because in Bactrocera we assumed eight generations per year (in nature they range from 3 to 5 of B. oleae and sub-tropical B. dorsalis populations, to >12 for the tropical species; [59,60,61,62]) and to account for uncertainty, we finally set as prior a normally distributed mean of 0.028 (SD = 0.03). In a second approach, we set a mutation rate log-normally distributed with “mean in real space” M = 0.028 and S = 0.82 (to produce the same 95% quantile, 0.077, as the normal distribution). For both approaches we performed a model selection to choose the most fitting clock and demographic prior based on the marginal likelihood values with the nested sampling approach implemented in the NS package [63]. We tested the strict and the LOGN relaxed clock and the Yule and the birth–death models, for a total of eight different combinations. Following the recommendations provided by the dedicated Taming the Beast tutorial [50], the sub-chain length was set at 50,000, which corresponds to the length of the MCMC run (i.e., 5 × 107) divided by the smallest ESS value observed across the eight model runs (i.e., ~1000), and the number of particles was set at 10. A model was considered favored over another if the difference between the two marginal likelihoods (i.e., the Bayes factor (BF) in log space) was more than twice the sum of the corresponding standard deviations. We ran eight different analyses that used different combinations of priors and model settings and performed a model selection to identify the most appropriate for our data. In particular, the nested sampling approach allowed us to estimate the marginal likelihoods of the different models and make pairwise comparisons using the associated Bayes factor. All models had marginal likelihoods with a standard deviation ranging from 2.5 to 2.9, which was small enough to assess whether a model was favored over another one. In all cases, we employed a GTR + G replacement model. Because the Bayesian phylogenetic analysis on the concatenated 4-fold degenerate sites resulted in a topology incongruent to the one supported by all other ML and Bayesian analyses (see Section 3), the species tree was fixed according to the consensus topology.
All analyses were performed twice, with chains ran for 5 × 107 generations, sampling trees and parameters every 1000 generations and inspecting convergence and likelihood plateauing in Tracer. Both chains resulted well mixed, with average effective sample size (ESS) values across posterior values being well above 200. The consensus trees (maximum clade credibility trees) were generated after discarding the first 20% of generations as burn-in.

3. Results

3.1. Phylogenetic Analysis

The results of our analyses strongly indicate that B. dorsalis is more closely related to B. latifrons than to B. tryoni (Figure 1). This relationship is highly supported according to both the ML bootstrap values (≥98, Figure 1a and Figure S1) and the Bayesian posterior probabilities obtained using both BEAST and PhyloBayes analyses (PP = 1, Figure 1a and Figures S2–S4), no matter whether they were based on the codon or aminoacidic alignments. The only support for the relationship ((B. dorsalis, B. tryoni), B. latifrons) comes from the Bayesian analysis run in BEAST2 and was based on the concatenated alignment of the 4-fold degenerate sites (Figure S5). This latter dataset contains sites under (nearly) neutral evolution and, therefore, is suitable for divergence time estimates because it allows a calibration in which a neutral spontaneous mutation rate can be considered equal to the substitution rate [18,56]. This type of dataset is, however, more prone to saturation, and hence, artefacts due to systematic errors. When the 4-fold degenerate site dataset was analyzed in PhyloBayes (which is also a Bayesian implementation) using the among-site heterogeneous CAT model, instead of the among-site GTR homogenous model in BEAST2, we retrieved the same topology obtained by the aforementioned Bayesian and ML analyses on the nucleotide and aminoacidic datasets (Figure S6; see also Figure S7 for the results of the ML analysis on the same dataset). Heterogeneous models such as CAT can accommodate for systematic errors related to site-specific variation of the replacement pattern: this suggests that the discordant topology obtained by BEAST2, and also supported by mitochondrial data, is likely an artefact due to model inadequacy.

3.2. Dating Analysis

Our dating analysis is based on the best combination of priors according to a model selection whose results (Table S2) indicated as the favored model the one where we set the mutation prior with a log-normal distribution, a strict clock, and a birth–death model. The Bayes factor values, even after correcting for uncertainty by subtracting the corresponding standard deviations, were well above two, which provides overwhelming support for that model [64]. The fact that a strict clock is favored over a relaxed clock is consistent with the low mean value of the coefficient of variation parameter (i.e., the standard deviation of branch rates divided by the mean rate), which equals 0.24. Therefore, in the main text we report and discuss the results obtained by this analysis.
The clade including Zeugodacus and Bactrocera split from Ceratitis between 6 and 45 million years ago (mean 19 million years ago, mya; 95% highest posterior density, 95%HPD, of 6–45 mya; Figure 2), with Zeugodacus then diverging ca. 7.6 mya (95%HPD: 2.3–18 mya). The diversification of Bactrocera followed immediately after, ca. 5.9 mya (95%HPD: 1.8–14 mya), with a rapid radiation subtending the species-rich clade [24] that includes the polyphagous species B. dorsalis, B. latifrons and B. tryoni happening between ca. 2.1–1.9 mya (95%HPDs 0.6–5 mya).

4. Discussion

4.1. Phylogenetic Analyses Reveal a Closer Affinity of B. dorsalis to B. latifrons Than to B. tryoni

The results of the phylogenetic analyses reveal relationships that are in contrast with all phylogenies inferred from mitochondrial sequences (e.g., [30,31]), which support B. dorsalis as being more closely related to B. tryoni than to B. latifrons. The results are instead consistent with a closer relationship between B. dorsalis and B. latifrons, as only partially supported by three previous studies. The first of such phylogenetic analyses, based on 167 Dacini species including Bactrocera [24], however, was not conclusive in determining the relationships between these three species, since the phylogeny had many unresolved nodes, including those relative to the most common ancestors of B. dorsalis, B. latifrons and B. tryoni. The low power to disentangle such relationships likely derives from the small dataset—seven nuclear genes—used to produce their phylogeny. A second study by Choo and colleagues [32] analyzed 116 orthologous genes across 11 Bactrocera species: despite the larger dataset, their results were also not adequately supported, as the split between (B. dorsalis, B. latifrons) and B. tryoni in a ML analysis had a bootstrap value of 70, indicating lack of statistical confidence. We also note that the incongruence does not extend to the species that in our analyses are identified as the closest relatives to B. dorsalis (i.e., B. musae), B. latifrons (i.e., B. bryoniae) and B. tryoni (i.e., B. jarvisi), and which are characterized by a longer divergence history (see below). Similarly, a third study [26] that produced an even larger phylogenomic dataset using highly multiplexed amplicon sequencing was not conclusive regarding the relationships between these three species, with a polymorphism-aware phylogenetic model [65] analysis favoring B. dorsalis being more closely related to B. tryoni than to B. latifrons, and an ASTRAL analysis in which this relationship was unresolved.
The discordance between the results obtained by previous studies and our study points to a complex evolutionary history of the group, as suggested by coalescent-aware analyses. Combining single-gene analyses into a coalescent framework also supports B. dorsalis as being more closely related to B. latifrons than to B. tryoni, both when using ASTRAL and StarBEAST2 (Figures S8–S13). It is evident, however, that many genes have phylogenies not consistent with the inferred species phylogeny. For instance, the StarBEAST2 approach reveals a high number of genes having an alternative phylogeny within the ((B. dorsalis, B. latifrons), B. tryoni), suggesting incomplete lineage sorting due to fast radiation or (ongoing) hybridization. For example, in the nucleotide datasets, 58 genes have a ((B. dorsalis, B. latifrons), B. tryoni) topology, 38 a ((B. dorsalis, B. tryoni), B. latifrons) topology and 14 a ((B. tryoni, B. latifrons), B. dorsalis) topology. This uncertainty is also apparent in the ASTRAL results, whereby the support for such a clade falls to <92 when bootstrapping by gene resampling (compare Figures S8 and S9). In the light of these results, we cannot exclude the possibility that these species experienced and still experience hybridization events, which could then result in widespread introgression events. Indeed, hybrids have been reported for several closely related Bactrocera species [66,67,68,69,70,71,72,73,74,75,76], and although none of the published studies involved the pair of species analyzed in our analyses, possible introgression can occur via direct hybridization or via intermediate hybridization events involving other closely related species. Finally, our results are also consistent with biogeography data: B. dorsalis and B. latifrons are originally from South Asia, whereas B. tryoni is native of Australia, suggesting that their ancestral distributional range also reflects their degree of evolutionary divergence.

4.2. Dating Analysis Suggests Fast and Recent Radiation in Bactrocera

Incomplete lineage sorting is expected for rapid radiations (e.g., [77]): this is exactly what is revealed by our molecular clock analyses, which support more recent divergences for the Bactrocera radiation compared with previous estimates [23,32,78,79]. Consistent with a rapid radiation of the (B. dorsalis, B. latifrons, B. tryoni) clade, the results of the clock analysis place its origin in the mid-Pliocene, with a mean at ~2.08 mya (95%HPD: 0.6–5 mya), and a subsequent, very close cladogenesis centered at ~1.87 mya (95%HPD: 0.6–4.5 mya) separating B. dorsalis and B. latifrons (Figure 2; see also Figures S14–S21 for the time trees obtained using the different models reported in Table S2). Interestingly, during this period the sea rose at peak levels [80] and, thus, increased distances between islands and island groups, possibly facilitating allopatric speciation (the three species have native ranges in southeastern Asia and Australia). The proximity of the two cladogenetic events and the large overlap of their 95% confidence intervals agree with a rapid radiation, which could have resulted not only in frequent incomplete lineage sorting but also in possible hybridization events as discussed above. This would also explain the discordant results between the nuclear and the mitochondrial phylogenies, a finding that is reported in many organisms, including insects [81,82,83,84].
Our results also revealed that the split between Zeugodacus + Bactrocera and C. capitata could have occurred between 6 and 45 million years ago (mean of ca. 19 mya). Previously published works proposed more ancient C. capitata and Bactrocera splits at around 24.9 mya (no confidence intervals reported, [78]), 83 mya (95% highest posterior density: 64–103 mya [79]), and 110.9 mya (95% confidence interval: 91.2–131.4 mya [23]). The difference in the estimates could be due to several factors such as the type of molecular markers (nuclear or mitochondrial), the choice of clock and demographic models, and more importantly, the type of calibrations employed. Here we employed a strategy based on a calibration that combined fossil data and mutation rate, in which we assumed a defined number of generations per year on the basis of available demographic evidence [59,60,61,62]. Generation time can, however, greatly influence the divergence time estimation: for example, if in our analysis we had assumed five generations per year (instead of eight) without changing the other parameters, the mean divergence time would increase from ~20 to ~32 mya for the Bactrocera-Ceratitis split. The latter is close to what was estimated by Choo and colleagues [32], who performed phylogenetic and dating analyses concatenating nuclear (n = 116) and mitochondrial (n = 13) genes of a similar species dataset and whose divergence time confidence intervals for this node partially overlap with ours (31.21 mya, 95%HPD: 41.27–21.61; this study: 18.86 mya, 95%HPD: 6–44.98). We also employed a mutation rate prior [55] for the 4-fold degenerate sites that is different and, in our opinion, more conservative than the calibration used by Choo and colleagues. The molecular clock is known to change among sites and lineages [85] and through time [86,87]. Therefore, when choosing a rate prior it is important that it refers to a similar timescale, and especially to sites that evolve similarly. Our approach used a molecular rate estimated from Drosophila, which has a generation time and a phylogeny timescale similar to Bactrocera; hence, the 4-fold degenerate site rate prior may indeed be a reliable assumption. In contrast, Choo and colleagues used a different approach, whereby they specified as prior only the divergence between Rhagoletis (Diptera: Tephritidae) and Drosophila, previously inferred in [88].
Finally, we would also like to point out that the mutation rate prior distribution can be a strong determinant of the divergence time estimates and, therefore, need to be carefully tested against alternative models (Table S2): for example, mean divergence time between C. capitata and Bactrocera is estimated between 18.3 and 18.9 mya when using a log-normal distribution (Figure 2 and Figures S13–S16), whereas it is estimated between 13 and 13.2 mya using a normal distribution (Figures S17–S21), and the same reduction of ~30% holds for all other divergence times across the phylogeny. Nevertheless, model selection indicates that a log-normal is the most fitting prior distribution of the rate; therefore, the divergence times depicted in Figure 2 should be regarded as more likely. Irrespectively of the divergence estimate for the deep splits, our results clearly indicate that the clade to which B. dorsalis, B. latifrons and B. tryoni belong underwent a rapid and recent diversification.

5. Conclusions

Overall, our clock analyses reveal a recent and fast radiation scenario of Bactrocera evolution and our phylogeny provides a useful framework for future comparative genomics and comparative biology studies in the major Bactrocera pest species. The possibility that hybridization can still occur between closely related species also warns about the possibility that selective events in one species (for instance, resistance to insecticides) may be readily transferred to other species by introgression. More generally, our results once again highlight the importance of comparing different evolutionary models to understand complex phylogenies.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/insects13070603/s1, Figure S1: Phylogeny of Bactrocera obtained by RAxML on the concatenated codon alignments; Figure S2: Phylogeny of Bactrocera obtained by BEAST2 using a Bayesian analysis on the concatenated aminoacidic alignments; Figure S3: Phylogeny of Bactrocera obtained by BEAST2 using a Bayesian analysis on the concatenated codon alignments; Figure S4: Phylogeny of Bactrocera obtained by PhyloBayes using a Bayesian analysis on the concatenated aminoacidic alignments; Figure S5: Phylogeny of Bactrocera obtained by BEAST2 using a Bayesian analysis on the 4-fold degenerate sites of the concatenated alignments; Figure S6: Phylogeny of Bactrocera obtained by PhyloBayes using a Bayesian analysis on the 4-fold degenerate sites of the concatenated alignments; Figure S7: Phylogeny of Bactrocera obtained by RAxML using a maximum likelihood analysis on the 4-fold degenerate sites of the concatenated alignments; Figure S8: Phylogeny of Bactrocera inferred using ASTRAL and based on single-gene aminoacidic alignments obtained by RAxML and with bootstrap values estimated by performing multi-locus bootstrap replicates; Figure S9: Phylogeny of Bactrocera inferred using ASTRAL and based on single-gene aminoacidic alignments obtained by RAxML and with bootstrap values estimated by performing gene + site resampling; Figure S10: Phylogeny of Bactrocera inferred using ASTRAL and based on single-gene codon alignments obtained by RAxML and with bootstrap values estimated by performing multi-locus bootstrap replicates; Figure S11: Phylogeny of Bactrocera inferred using ASTRAL and based on single-gene codon alignments obtained by RAxML and with bootstrap values estimated by performing gene + site resampling; Figure S12: Multi-locus phylogenesis of Bactrocera. Bayesian analyses obtained by StarBEAST2 for each of the aminoacidic orthologous genes alignments; Figure S13: Multi-locus phylogenesis of Bactrocera. Bayesian analyses were obtained by StarBEAST2 for each of the codon orthologous genes alignments; Figure S14: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate log-normally distributed as prior, a log-normal clock and a birth–death model; Figure S15: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate log-normally distributed as prior, a strict clock and a birth–death model; Figure S16: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate log-normally distributed as prior, a log-normal clock and a Yule model; Figure S17: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate log-normally distributed as prior, a strict clock and a Yule model; Figure S18: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate normally distributed as prior, a log-normal clock and a birth–death model; Figure S19: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate normally distributed as prior, a strict clock and a birth-–death model; Figure S20: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate normally distributed as prior, a log-normal clock and a Yule model; Figure S21: Molecular time tree of Bactrocera obtained by BEAST2 using the 4-fold degenerate sites using a mutation rate normally distributed as prior, a strict clock and a Yule model; Table S1: SRA databases used to reconstruct the transcriptome of five Bactrocera species; Table S2: Results of the model selection of the different BEAST analyses used to estimate divergence times; File S1: alignments of the 110 orthologous gene sets.

Author Contributions

Conceptualization, L.O.; methodology, F.V., N.Z., O.R.-S. and L.O.; formal analysis, F.V., N.Z. and L.O.; investigation, F.V., N.Z., O.R.-S. and L.O.; writing—original draft preparation, F.V. and L.O.; writing—review and editing, F.V., L.O. and O.R.-S.; supervision, L.O. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Data Availability Statement

Publicly available datasets were analyzed in this study. Raw data can be found at https://www.ncbi.nlm.nih.gov/sra accessed on 24 April 2020, with accession numbers indicated in the additional Table S2, whereas available CDS data can be found at https://www.ncbi.nlm.nih.gov/genome/ accessed on 24 April 2020. The datasets generated and analyzed during the current study, when not available as supplementary information, and the associated scripts are available from the corresponding author on request.

Acknowledgments

We acknowledge the support of the Italian Ministry of Education, University and Research (MIUR) Dipartimenti di Eccellenza Program (2018–2022), Department of Biology and Biotechnology “L. Spallanzani”, University of Pavia.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Bickel, D.; Pape, T.; Meier, R. Diptera Diversity: Status, Challenges and Tools; Brill: Leiden, The Netherlands, 2009; ISBN 978-90-04-18100-7. [Google Scholar]
  2. Daane, K.M.; Johnson, M.W. Olive Fruit Fly: Managing an Ancient Pest in Modern Times. Annu. Rev. Entomol. 2010, 55, 151–169. [Google Scholar] [CrossRef] [PubMed]
  3. Barr, N.B.; Ledezma, L.A.; Leblanc, L.; San Jose, M.; Rubinoff, D.; Geib, S.M.; Fujita, B.; Bartels, D.W.; Garza, D.; Kerr, P.; et al. Genetic Diversity of Bactrocera dorsalis (Diptera: Tephritidae) on the Hawaiian Islands: Implications for an Introduction Pathway Into California. J. Econ. Entomol. 2014, 107, 1946–1958. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  4. Vargas, R.I.; Nishida, T. Survey for Dacus latifrons (Diptera: Tephritidae). J. Econ. Entomol. 1985, 78, 1311–1314. [Google Scholar] [CrossRef]
  5. Lux, S.A.; Copeland, R.S.; White, I.M.; Manrakhan, A.; Billah, M.K. A New Invasive Fruit Fly Species from the Bactrocera dorsalis (Hendel) Group Detected in East Africa. Int. J. Trop. Insect Sci. 2003, 23, 355–361. [Google Scholar] [CrossRef]
  6. Drew, R.A.I.; Tsuruta, K.; White, I. A New Species of Pest Fruit Fly (Diptera: Tephritidae: Dacinae) from Sri Lanka and Africa. Afr. Entomol. 2004, 13, 149–154. [Google Scholar]
  7. Vayssières, J.-F.; Goergen, G.; Lokossou, O.; Dossa, P.; Akponon, C. A New Bactrocera Species in Benin among Mango Fruit Fly (Diptera: Tephritidae) Species. Fruits 2005, 60, 371–377. [Google Scholar] [CrossRef] [Green Version]
  8. Goergen, G.; Vayssières, J.-F.; Gnanvossou, D.; Tindo, M. Bactrocera invadens (Diptera: Tephritidae), a New Invasive Fruit Fly Pest for the Afrotropical Region: Host Plant Range and Distribution in West and Central Africa. Env. Entomol. 2011, 40, 844–854. [Google Scholar] [CrossRef] [PubMed]
  9. EPPO Global Database. Available online: https://gd.eppo.int (accessed on 21 June 2022).
  10. Crop Protection Compendium. Available online: https://www.cabi.org/cpc (accessed on 21 June 2022).
  11. Mwatawala, M.; Makundi, R.; Maerere, A.P.; De Meyer, M. Occurrence of the Solanum Fruit Fly Bactrocera latifrons (Hendel) (Diptera: Tephritidae) in Tanzania. J. Afrotrop. Zool. 2010, 6, 83–89. [Google Scholar]
  12. Mziray, H.A.; Makundi, R.H.; Mwatawala, M.; Maerere, A.; De Meyer, M. Host Use of Bactrocera latifrons, a New Invasive Tephritid Species in Tanzania. J. Econ. Entomol. 2010, 103, 70–76. [Google Scholar] [CrossRef]
  13. Nugnes, F.; Russo, E.; Viggiani, G.; Bernardo, U. First Record of an Invasive Fruit Fly Belonging to Bactrocera dorsalis Complex (Diptera: Tephritidae) in Europe. Insects 2018, 9, 182. [Google Scholar] [CrossRef] [Green Version]
  14. Gargiulo, S.; Nugnes, F.; Benedetta, F.; Bernardo, U. Bactrocera latifrons in Europe: The Importance of the Right Attractant for Detection. Bull. Insectol. 2021, 74, 311–320. [Google Scholar]
  15. Doorenweerd, C.; Leblanc, L.; Norrbom, A.L.; San Jose, M.; Rubinoff, D. A Global Checklist of the 932 Fruit Fly Species in the Tribe Dacini (Diptera, Tephritidae). ZooKeys 2018, 730, 19–56. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  16. Vargas, R.I.; Piñero, J.C.; Leblanc, L. An Overview of Pest Species of Bactrocera Fruit Flies (Diptera: Tephritidae) and the Integration of Biopesticides with Other Biological Approaches for Their Management with a Focus on the Pacific Region. Insects 2015, 6, 297–318. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  17. White, I.M.; Elson-Harris, M.M. Fruit Flies of Economic Significance: Their Identification and Bionomics; CAB International: Wallingford, UK, 1992; ISBN 0-85198-790-7. [Google Scholar]
  18. Ometto, L.; Cestaro, A.; Ramasamy, S.; Grassi, A.; Revadi, S.; Siozios, S.; Moretto, M.; Fontana, P.; Varotto, C.; Pisani, D.; et al. Linking Genomics and Ecology to Investigate the Complex Evolution of an Invasive Drosophila Pest. Genome Biol. Evol. 2013, 5, 745–757. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  19. Clarke, A.R.; Armstrong, K.F.; Carmichael, A.E.; Milne, J.R.; Raghu, S.; Roderick, G.K.; Yeates, D.K. Invasive Phytophagous Pests Arising through a Recent Tropical Evolutionary Radiation: The Bactrocera dorsalis Complex of Fruit Flies. Annu. Rev. Entomol. 2004, 50, 293–319. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  20. Drew, R.A.I. The Tropical Fruit Flies (Diptera: Tephritidae: Dacinae) of the Australasian Region and Oceanic Regions. Mem. Qld. Mus. 1989, 26, 1–521. [Google Scholar]
  21. Drew, R.A.I.; Hancock, D.L. Phylogeny of the Tribe Dacini (Dacinae) Based on Morphological, Distributional, and Biological Data; CRC Press: Boca Raton, FL, USA, 1999; pp. 509–522. ISBN 978-0-429-12467-9. [Google Scholar]
  22. Drew, R.A.I.; Romig, M.C. Keys to the Tropical Fruit Flies of South-East Asia; CABI: Wallingford, UK, 2016. [Google Scholar]
  23. Krosch, M.N.; Schutze, M.K.; Armstrong, K.F.; Graham, G.C.; Yeates, D.K.; Clarke, A.R. A Molecular Phylogeny for the Tribe Dacini (Diptera: Tephritidae): Systematic and Biogeographic Implications. Mol. Phylogenet. Evol. 2012, 64, 513–523. [Google Scholar] [CrossRef] [Green Version]
  24. San Jose, M.; Doorenweerd, C.; Leblanc, L.; Barr, N.; Geib, S.; Rubinoff, D. Incongruence between Molecules and Morphology: A Seven-Gene Phylogeny of Dacini Fruit Flies Paves the Way for Reclassification (Diptera: Tephritidae). Mol. Phylogenet. Evol. 2018, 121, 139–149. [Google Scholar] [CrossRef]
  25. Virgilio, M.; Jordaens, K.; Verwimp, C.; White, I.M.; De Meyer, M. Higher Phylogeny of Frugivorous Flies (Diptera, Tephritidae, Dacini): Localised Partition Conflicts and a Novel Generic Classification. Mol. Phylogenet. Evol. 2015, 85, 171–179. [Google Scholar] [CrossRef]
  26. Dupuis, J.R.; Bremer, F.T.; Kauwe, A.; San Jose, M.; Leblanc, L.; Rubinoff, D.; Geib, S.M. HiMAP: Robust Phylogenomics from Highly Multiplexed Amplicon Sequencing. Mol. Ecol. Resour. 2018, 18, 1000–1019. [Google Scholar] [CrossRef]
  27. Muraji, M.; Nakahara, S. Phylogenetic Relationships among Fruit Flies, Bactrocera (Diptera, Tephritidae), Based on the Mitochondrial rDNA Sequences. Insect Mol. Biol. 2001, 10, 549–559. [Google Scholar] [CrossRef] [PubMed]
  28. Nakahara, S.; Muraji, M. Phylogenetic Analyses of Bactrocera Fruit Flies (Diptera: Tephritidae) Based on Nucleotide Sequences of the Mitochondrial COI and COII Genes. Res. Bull. Plant Prot. Jpn. 2007, 44, 1–12. [Google Scholar]
  29. Smith, P.T.; Kambhampati, S.; Armstrong, K.A. Phylogenetic Relationships among Bactrocera Species (Diptera: Tephritidae) Inferred from Mitochondrial DNA Sequences. Mol. Phylogenet. Evol. 2003, 26, 8–17. [Google Scholar] [CrossRef]
  30. Yong, H.-S.; Song, S.-L.; Lim, P.-E.; Eamsobhana, P.; Suana, I.W. Complete Mitochondrial Genome of Three Bactrocera Fruit Flies of Subgenus Bactrocera (Diptera: Tephritidae) and Their Phylogenetic Implications. PLoS ONE 2016, 11, e0148201. [Google Scholar] [CrossRef]
  31. Zhang, B.; Liu, Y.H.; Wu, W.X.; Wang, Z.L. Molecular Phylogeny of Bactrocera Species (Diptera: Tephritidae: Dacini) Inferred from Mitochondrial Sequences of 16S RDNA and COI Sequences. Fla. Entomol. 2010, 93, 369–377. [Google Scholar] [CrossRef]
  32. Choo, A.; Nguyen, T.N.M.; Ward, C.M.; Chen, I.Y.; Sved, J.; Shearman, D.; Gilchrist, A.S.; Crisp, P.; Baxter, S.W. Identification of Y-Chromosome Scaffolds of the Queensland Fruit Fly Reveals a Duplicated Gyf Gene Paralogue Common to Many Bactrocera Pest Species. Insect Mol. Biol. 2019, 28, 873–886. [Google Scholar] [CrossRef]
  33. Papanicolaou, A.; Schetelig, M.F.; Arensburger, P.; Atkinson, P.W.; Benoit, J.B.; Bourtzis, K.; Castañera, P.; Cavanaugh, J.P.; Chao, H.; Childers, C.; et al. The Whole Genome Sequence of the Mediterranean Fruit Fly, Ceratitis capitata (Wiedemann), Reveals Insights into the Biology and Adaptive Evolution of a Highly Invasive Pest Species. Genome Biol. 2016, 17, 192. [Google Scholar] [CrossRef] [Green Version]
  34. Sim, S.B.; Geib, S.M. A Chromosome-Scale Assembly of the Bactrocera cucurbitae Genome Provides Insight to the Genetic Basis of White Pupae. G3 Genes Genomes Genet. 2017, 7, 1927–1940. [Google Scholar] [CrossRef] [Green Version]
  35. Geib, S.M.; Calla, B.; Hall, B.; Hou, S.; Manoukis, N.C. Characterizing the Developmental Transcriptome of the Oriental Fruit Fly, Bactrocera dorsalis (Diptera: Tephritidae) through Comparative Genomic Analysis with Drosophila melanogaster Utilizing ModENCODE Datasets. BMC Genom. 2014, 15, 924. [Google Scholar] [CrossRef] [Green Version]
  36. Wang, J.; Xiong, K.-C.; Liu, Y.-H. De Novo Transcriptome Analysis of Chinese Citrus Fly, Bactrocera minax (Diptera: Tephritidae), by High-Throughput Illumina Sequencing. PLoS ONE 2016, 11, e0157656. [Google Scholar] [CrossRef]
  37. Bayega, A.; Djambazian, H.; Tsoumani, K.T.; Gregoriou, M.E.; Sagri, E.; Drosopoulou, E.; Mavragani-Tsipidou, P.; Giorda, K.; Tsiamis, G.; Bourtzis, K.; et al. De Novo Assembly of the Olive Fruit Fly (Bactrocera oleae) Genome with Linked-Reads and Long-Read Technologies Minimizes Gaps and Provides Exceptional y Chromosome Assembly. BMC Genom. 2020, 21, 259. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  38. Grabherr, M.G.; Haas, B.J.; Yassour, M.; Levin, J.Z.; Thompson, D.A.; Amit, I.; Adiconis, X.; Fan, L.; Raychowdhury, R.; Zeng, Q.; et al. Full-Length Transcriptome Assembly from RNA-Seq Data without a Reference Genome. Nat. Biotechnol. 2011, 29, 644–652. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  39. Camacho, C.; Coulouris, G.; Avagyan, V.; Ma, N.; Papadopoulos, J.; Bealer, K.; Madden, T.L. BLAST+: Architecture and Applications. BMC Bioinform. 2009, 10, 421. [Google Scholar] [CrossRef] [Green Version]
  40. Katoh, K.; Standley, D.M. MAFFT Multiple Sequence Alignment Software Version 7: Improvements in Performance and Usability. Mol. Biol. Evol. 2013, 30, 772–780. [Google Scholar] [CrossRef] [Green Version]
  41. Löytynoja, A.; Goldman, N. A Model of Evolution and Structure for Multiple Sequence Alignment. Philos. Trans. R. Soc. B Biol. Sci. 2008, 363, 3913–3919. [Google Scholar] [CrossRef]
  42. Abascal, F.; Zardoya, R.; Telford, M.J. TranslatorX: Multiple Alignment of Nucleotide Sequences Guided by Amino Acid Translations. Nucleic Acids Res. 2010, 38, W7–W13. [Google Scholar] [CrossRef] [Green Version]
  43. Ramasamy, S.; Ometto, L.; Crava, C.M.; Revadi, S.; Kaur, R.; Horner, D.S.; Pisani, D.; Dekker, T.; Anfora, G.; Rota-Stabelli, O. The Evolution of Olfactory Gene Families in Drosophila and the Genomic Basis of Chemical-Ecological Adaptation in Drosophila suzukii. Genome Biol. Evol. 2016, 8, 2297–2311. [Google Scholar] [CrossRef] [Green Version]
  44. Kumar, S.; Stecher, G.; Tamura, K. MEGA7: Molecular Evolutionary Genetics Analysis Version 7.0 for Bigger Datasets. Mol. Biol. Evol. 2016, 33, 1870–1874. [Google Scholar] [CrossRef] [Green Version]
  45. Stamatakis, A. RAxML Version 8: A Tool for Phylogenetic Analysis and Post-Analysis of Large Phylogenies. Bioinformatics 2014, 30, 1312–1313. [Google Scholar] [CrossRef]
  46. Zhang, C.; Rabiee, M.; Sayyari, E.; Mirarab, S. ASTRAL-III: Polynomial Time Species Tree Reconstruction from Partially Resolved Gene Trees. BMC Bioinform. 2018, 19, 153. [Google Scholar] [CrossRef] [Green Version]
  47. Bouckaert, R.; Vaughan, T.G.; Barido-Sottani, J.; Duchêne, S.; Fourment, M.; Gavryushkina, A.; Heled, J.; Jones, G.; Kühnert, D.; De Maio, N.; et al. BEAST 2.5: An Advanced Software Platform for Bayesian Evolutionary Analysis. PLoS Comput. Biol. 2019, 15, e1006650. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  48. Lartillot, N.; Philippe, H. A Bayesian Mixture Model for Across-Site Heterogeneities in the Amino-Acid Replacement Process. Mol. Biol. Evol. 2004, 21, 1095–1109. [Google Scholar] [CrossRef] [PubMed]
  49. Ogilvie, H.A.; Bouckaert, R.R.; Drummond, A.J. StarBEAST2 Brings Faster Species Tree Inference and Accurate Estimates of Substitution Rates. Mol. Biol. Evol. 2017, 34, 2101–2114. [Google Scholar] [CrossRef]
  50. Barido-Sottani, J.; Bošková, V.; Du Plessis, L.; Kühnert, D.; Magnus, C.; Mitov, V.; Müller, N.F.; Pečerska, J.; Rasmussen, D.A.; Zhang, C.; et al. Taming the BEAST—A Community Teaching Material Resource for BEAST 2. Syst. Biol. 2018, 67, 170–174. [Google Scholar] [CrossRef] [PubMed]
  51. Bouckaert, R.R.; Heled, J. DensiTree 2: Seeing Trees Through the Forest. bioRxiv 2014, 012401. [Google Scholar] [CrossRef] [Green Version]
  52. Norrbom, A. New Genera of Tephritidae (Diptera) from Brazil and Dominican Amber, with Phylogenetic Analysis of the Tribe Ortalotrypetini. Insecta Mundi 1994, 8, 1–15. [Google Scholar]
  53. Wiegmann, B.M.; Trautwein, M.D.; Winkler, I.S.; Barr, N.B.; Kim, J.-W.; Lambkin, C.; Bertone, M.A.; Cassel, B.K.; Bayless, K.M.; Heimberg, A.M.; et al. Episodic Radiations in the Fly Tree of Life. Proc. Natl. Acad. Sci. USA 2011, 108, 5690–5695. [Google Scholar] [CrossRef] [Green Version]
  54. Junqueira, A.C.M.; Azeredo-Espin, A.M.L.; Paulo, D.F.; Marinho, M.A.T.; Tomsho, L.P.; Drautz-Moses, D.I.; Purbojati, R.W.; Ratan, A.; Schuster, S.C. Large-Scale Mitogenomics Enables Insights into Schizophora (Diptera) Radiation and Population Diversity. Sci. Rep. 2016, 6, 1–13. [Google Scholar] [CrossRef] [Green Version]
  55. Keightley, P.D.; Trivedi, U.; Thomson, M.; Oliver, F.; Kumar, S.; Blaxter, M.L. Analysis of the Genome Sequences of Three Drosophila melanogaster Spontaneous Mutation Accumulation Lines. Genome Res. 2009, 19, 1195–1201. [Google Scholar] [CrossRef] [Green Version]
  56. Obbard, D.J.; Maclennan, J.; Kim, K.-W.; Rambaut, A.; O’Grady, P.M.; Jiggins, F.M. Estimating Divergence Dates and Substitution Rates in the Drosophila Phylogeny. Mol. Biol. Evol. 2012, 29, 3459–3473. [Google Scholar] [CrossRef] [Green Version]
  57. Keightley, P.D.; Pinharanda, A.; Ness, R.W.; Simpson, F.; Dasmahapatra, K.K.; Mallet, J.; Davey, J.W.; Jiggins, C.D. Estimation of the Spontaneous Mutation Rate in Heliconius melpomene. Mol. Biol. Evol. 2015, 32, 239–243. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  58. Keightley, P.D.; Ness, R.W.; Halligan, D.L.; Haddrill, P.R. Estimation of the Spontaneous Mutation Rate per Nucleotide Site in a Drosophila melanogaster Full-Sib Family. Genetics 2014, 196, 313–320. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  59. Vargas, R.I.; Walsh, W.A.; Kanehisa, D.; Jang, E.B.; Armstrong, J.W. Demography of Four Hawaiian Fruit Flies (Diptera: Tephritidae) Reared at Five Constant Temperatures. Ann. Entomol. Soc. Am. 1997, 90, 162–168. [Google Scholar] [CrossRef]
  60. Stephens, A.E.A.; Kriticos, D.J.; Leriche, A. The Current and Future Potential Geographical Distribution of the Oriental Fruit Fly, Bactrocera dorsalis (Diptera: Tephritidae). Bull. Entomol. Res. 2007, 97, 369–378. [Google Scholar] [CrossRef]
  61. Theron, C.D.; Manrakhan, A.; Weldon, C.W. Host Use of the Oriental Fruit Fly, Bactrocera dorsalis (Hendel) (Diptera: Tephritidae), in South Africa. J. Appl. Entomol. 2017, 141, 810–816. [Google Scholar] [CrossRef] [Green Version]
  62. Li, X.; Yang, H.; Wang, T.; Wang, J.; Wei, H. Life History and Adult Dynamics of Bactrocera dorsalis in the Citrus Orchard of Nanchang, a Subtropical Area from China: Implications for a Control Timeline. ScienceAsia 2019, 45, 212–220. [Google Scholar] [CrossRef] [Green Version]
  63. Russel, P.M.; Brewer, B.J.; Klaere, S.; Bouckaert, R.R. Model Selection and Parameter Inference in Phylogenetics Using Nested Sampling. Syst. Biol. 2019, 68, 219–233. [Google Scholar] [CrossRef]
  64. Kass, R.E.; Raftery, A.E. Bayes Factors. J. Am. Stat. Assoc. 1995, 90, 774–795. [Google Scholar] [CrossRef]
  65. Schrempf, D.; Minh, B.Q.; De Maio, N.; von Haeseler, A.; Kosiol, C. Reversible Polymorphism-Aware Phylogenetic Models and Their Application to Tree Inference. J. Theor. Biol. 2016, 407, 362–370. [Google Scholar] [CrossRef] [Green Version]
  66. Drew, R.A.I.; Lambert, D.M. On the Specific Status of Dacus (Bactrocera) aquilonis and D. (Bactrocera) tryoni (Diptera: Tephritidae). Ann. Entomol. Soc. Am. 1986, 79, 870–878. [Google Scholar] [CrossRef]
  67. Cruickshank, L.; Jessup, A.J.; Cruickshank, D.J. Interspecific Crosses of Bactrocera tryoni (Froggatt) and Bactrocera jarvisi (Tryon) (Diptera: Tephritidae) in the Laboratory. Aust. J. Entomol. 2001, 40, 278–280. [Google Scholar] [CrossRef]
  68. Pike, N.; Wang, W.Y.S.; Meats, A. The Likely Fate of Hybrids of Bactrocera tryoni and Bactrocera neohumeralis. Heredity 2003, 90, 365–370. [Google Scholar] [CrossRef] [PubMed]
  69. Wee, S.-L.; Tan, K.-H. Evidence of Natural Hybridization between Two Sympatric Sibling Species of Bactrocera dorsalis Complex Based on Pheromone Analysis. J. Chem. Ecol. 2005, 31, 845–858. [Google Scholar] [CrossRef] [PubMed]
  70. Ebina, T.; Ohto, K. Morphological Characters and PCR-RFLP Markers in the Interspecific Hybrids between Bactrocera carambolae and B. papayae of the B. dorsalis Species Complex (Diptera: Tephritidae). Res. Bull. Plant Prot. Serv. Jpn. 2006, 42, 23–34. [Google Scholar]
  71. Shearman, D.C.A.; Frommer, M.; Morrow, J.L.; Raphael, K.A.; Gilchrist, A.S. Interspecific Hybridization as a Source of Novel Genetic Markers for the Sterile Insect Technique in Bactrocera tryoni (Diptera: Tephritidae). J. Econ. Entomol. 2010, 103, 1071–1079. [Google Scholar] [CrossRef]
  72. Schutze, M.K.; Jessup, A.; Ul-Haq, I.; Vreysen, M.J.B.; Wornoayporn, V.; Vera, M.T.; Clarke, A.R. Mating Compatibility Among Four Pest Members of the Bactrocera dorsalis Fruit Fly Species Complex (Diptera: Tephritidae). J. Econ. Entomol. 2013, 106, 695–707. [Google Scholar] [CrossRef]
  73. Augustinos, A.A.; Drosopoulou, E.; Gariou-Papalexiou, A.; Bourtzis, K.; Mavragani-Tsipidou, P.; Zacharopoulou, A. The Bactrocera dorsalis Species Complex: Comparative Cytogenetic Analysis in Support of Sterile Insect Technique Applications. BMC Genet. 2014, 15, S16. [Google Scholar] [CrossRef] [Green Version]
  74. Bo, W.; Ahmad, S.; Dammalage, T.; Tomas, U.S.; Wornoayporn, V.; Ul Haq, I.; Cáceres, C.; Vreysen, M.J.B.; Schutze, M.K. Mating Compatibility between Bactrocera invadens and Bactrocera dorsalis (Diptera: Tephritidae). J. Econ. Entomol. 2014, 107, 623–629. [Google Scholar] [CrossRef]
  75. Gilchrist, A.S.; Shearman, D.C.; Frommer, M.; Raphael, K.A.; Deshpande, N.P.; Wilkins, M.R.; Sherwin, W.B.; Sved, J.A. The Draft Genome of the Pest Tephritid Fruit Fly Bactrocera tryoni: Resources for the Genomic Analysis of Hybridising Species. BMC Genom. 2014, 15, 1153. [Google Scholar] [CrossRef] [Green Version]
  76. Yeap, H.L.; Lee, S.F.; Robinson, F.; Mourant, R.G.; Sved, J.A.; Frommer, M.; Papanicolaou, A.; Edwards, O.R.; Oakeshott, J.G. Separating Two Tightly Linked Species-Defining Phenotypes in Bactrocera with Hybrid Recombinant Analysis. BMC Genet. 2020, 21, 132. [Google Scholar] [CrossRef]
  77. Pollard, D.A.; Iyer, V.N.; Moses, A.M.; Eisen, M.B. Widespread Discordance of Gene Trees with Species Tree in Drosophila: Evidence for Incomplete Lineage Sorting. PLoS Genet. 2006, 2, e173. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  78. Yaakop, S.; Ibrahim, N.J.; Shariff, S.; Zain, B.M.M. Molecular Clock Analysis on Five Bactrocera Species Flies (Diptera: Tephritidae) Based on Combination of COI and NADH Sequences. Orient. Insects 2015, 49, 150–164. [Google Scholar] [CrossRef]
  79. Nardi, F.; Carapelli, A.; Boore, J.L.; Roderick, G.K.; Dallai, R.; Frati, F. Domestication of Olive Fly through a Multi-Regional Host Shift to Cultivated Olives: Comparative Dating Using Complete Mitochondrial Genomes. Mol. Phylogenet. Evol. 2010, 57, 678–686. [Google Scholar] [CrossRef] [PubMed]
  80. Zhong, G.; Geng, J.; Wong, H.K.; Ma, Z.; Wu, N. A Semi-Quantitative Method for the Reconstruction of Eustatic Sea Level History from Seismic Profiles and Its Application to the Southern South China Sea. Earth Planet. Sci. Lett. 2004, 223, 443–459. [Google Scholar] [CrossRef]
  81. DeSalle, R.; Giddings, L.V. Discordance of Nuclear and Mitochondrial DNA Phylogenies in Hawaiian Drosophila. Proc. Natl. Acad. Sci. USA 1986, 83, 6902–6906. [Google Scholar] [CrossRef] [Green Version]
  82. Beltrán, M.; Jiggins, C.D.; Bull, V.; Linares, M.; Mallet, J.; McMillan, W.O.; Bermingham, E. Phylogenetic Discordance at the Species Boundary: Comparative Gene Genealogies Among Rapidly Radiating Heliconius Butterflies. Mol. Biol. Evol. 2002, 19, 2176–2190. [Google Scholar] [CrossRef] [Green Version]
  83. Putnam, A.S.; Scriber, J.M.; Andolfatto, P. Discordant Divergence Times among Z-Chromosome Regions between Two Ecologically Distinct Swallowtail Butterfly Species. Evolution 2007, 61, 912–927. [Google Scholar] [CrossRef]
  84. Toews, D.P.L.; Brelsford, A. The Biogeography of Mitochondrial and Nuclear Discordance in Animals. Mol. Ecol. 2012, 21, 3907–3930. [Google Scholar] [CrossRef]
  85. Drummond, A.J.; Ho, S.Y.W.; Phillips, M.J.; Rambaut, A. Relaxed Phylogenetics and Dating with Confidence. PLoS Biol. 2006, 4, e88. [Google Scholar] [CrossRef]
  86. Ho, S.Y.W.; Lanfear, R.; Bromham, L.; Phillips, M.J.; Soubrier, J.; Rodrigo, A.G.; Cooper, A. Time-Dependent Rates of Molecular Evolution. Mol. Ecol. 2011, 20, 3087–3101. [Google Scholar] [CrossRef]
  87. Ho, S.Y.W.; Lo, N. The Insect Molecular Clock. Aust. J. Entomol. 2013, 52, 101–105. [Google Scholar] [CrossRef]
  88. Misof, B.; Liu, S.; Meusemann, K.; Peters, R.S.; Donath, A.; Mayer, C.; Frandsen, P.B.; Ware, J.; Flouri, T.; Beutel, R.G.; et al. Phylogenomics Resolves the Timing and Pattern of Insect Evolution. Science 2014, 346, 763–767. [Google Scholar] [CrossRef] [PubMed]
Insects 13 00603 g001 550
Figure 1. Phylogeny of Bactrocera inferred from the aminoacidic alignments of 110 orthologous nuclear genes. (a) Alternative phylogenetic relationships with supported (*) and non-supported (ns) nodes between B. dorsalis (d), B. latifrons (l) and B. tryoni (t) as presented in the literature and based on mitochondrial (Mt) DNA (e.g., [30,31]) or nuclear (Nuc) DNA (e.g., [24,26,32]). (b) Phylogeny obtained using a maximum likelihood phylogenetic analysis on the concatenated aminoacidic alignments (63,297 aa). Support at nodes is given as bootstrap values for the ML analyses (both for the aminoacidic and nucleotide alignments), bootstrap values estimated by performing 100 multi-locus gene + site resampling using a multi-locus coalescent-aware phylogenetic analysis in ASTRAL across all 110 genes (As) and as posterior probabilities for the Bayesian BEAST2 and PhyloBayes analyses on the aminoacidic dataset (BaB and BaP, respectively). Bactrocera dorsalis, B. latifrons and B. tryoni tips are color-coded as in panel A. (c) Bayesian analysis obtained by StarBEAST2, which employs a multispecies coalescent method to estimate species trees from multiple sequence alignments (i.e., one for each of the 110 orthologous gene sets). For this analysis we used the aminoacidic alignments, linking the site models across the gene sets. Note the numerous discordant gene trees, especially within the B. dorsalis–B. latifrons–B. tryoni clade, compared to the species tree (supported by the gene trees in blue).
Figure 1. Phylogeny of Bactrocera inferred from the aminoacidic alignments of 110 orthologous nuclear genes. (a) Alternative phylogenetic relationships with supported (*) and non-supported (ns) nodes between B. dorsalis (d), B. latifrons (l) and B. tryoni (t) as presented in the literature and based on mitochondrial (Mt) DNA (e.g., [30,31]) or nuclear (Nuc) DNA (e.g., [24,26,32]). (b) Phylogeny obtained using a maximum likelihood phylogenetic analysis on the concatenated aminoacidic alignments (63,297 aa). Support at nodes is given as bootstrap values for the ML analyses (both for the aminoacidic and nucleotide alignments), bootstrap values estimated by performing 100 multi-locus gene + site resampling using a multi-locus coalescent-aware phylogenetic analysis in ASTRAL across all 110 genes (As) and as posterior probabilities for the Bayesian BEAST2 and PhyloBayes analyses on the aminoacidic dataset (BaB and BaP, respectively). Bactrocera dorsalis, B. latifrons and B. tryoni tips are color-coded as in panel A. (c) Bayesian analysis obtained by StarBEAST2, which employs a multispecies coalescent method to estimate species trees from multiple sequence alignments (i.e., one for each of the 110 orthologous gene sets). For this analysis we used the aminoacidic alignments, linking the site models across the gene sets. Note the numerous discordant gene trees, especially within the B. dorsalis–B. latifrons–B. tryoni clade, compared to the species tree (supported by the gene trees in blue).
Insects 13 00603 g001
Insects 13 00603 g002 550
Figure 2. Molecular time tree of Bactrocera. Bactrocera (plus Zeugodacus) likely originated during the Miocene optimum (around 19 mya) and experienced recent, fast cladogenetic events less than 5 mya. The analysis was conducted setting the mutation rate log-normally distributed as prior, a strict clock and a birth–death model. Mean and 95% highest posterior density of the inferred age (blue bars) are reported for each node.
Figure 2. Molecular time tree of Bactrocera. Bactrocera (plus Zeugodacus) likely originated during the Miocene optimum (around 19 mya) and experienced recent, fast cladogenetic events less than 5 mya. The analysis was conducted setting the mutation rate log-normally distributed as prior, a strict clock and a birth–death model. Mean and 95% highest posterior density of the inferred age (blue bars) are reported for each node.
Insects 13 00603 g002
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Valerio, F.; Zadra, N.; Rota-Stabelli, O.; Ometto, L. The Impact of Fast Radiation on the Phylogeny of Bactrocera Fruit Flies as Revealed by Multiple Evolutionary Models and Mutation Rate-Calibrated Clock. Insects 2022, 13, 603. https://doi.org/10.3390/insects13070603

AMA Style

Valerio F, Zadra N, Rota-Stabelli O, Ometto L. The Impact of Fast Radiation on the Phylogeny of Bactrocera Fruit Flies as Revealed by Multiple Evolutionary Models and Mutation Rate-Calibrated Clock. Insects. 2022; 13(7):603. https://doi.org/10.3390/insects13070603

Chicago/Turabian Style

Valerio, Federica, Nicola Zadra, Omar Rota-Stabelli, and Lino Ometto. 2022. "The Impact of Fast Radiation on the Phylogeny of Bactrocera Fruit Flies as Revealed by Multiple Evolutionary Models and Mutation Rate-Calibrated Clock" Insects 13, no. 7: 603. https://doi.org/10.3390/insects13070603

APA Style

Valerio, F., Zadra, N., Rota-Stabelli, O., & Ometto, L. (2022). The Impact of Fast Radiation on the Phylogeny of Bactrocera Fruit Flies as Revealed by Multiple Evolutionary Models and Mutation Rate-Calibrated Clock. Insects, 13(7), 603. https://doi.org/10.3390/insects13070603

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop