Chloroplast genomic insights into adaptive evolution and rapid radiation in the genus Passiflora (Passifloraceae)

Abstract

Chloroplasts are essential organelles in plants and eukaryotic algae, responsible for photosynthesis, fatty acid synthesis, amino acid production, and stress responses. The genus Passiflora, known for its species diversity and dynamic chloroplast (cp) genome evolution, serves as an excellent model for studying structural variations. This study investigates evolutionary relationships within Passiflora by sequencing 11 new chloroplast genomes, assessing selective pressures on cp genes, and comparing plastid and nuclear phylogenies. Passiflora cp genomes showed significant variations in size, gene content, and structure, ranging from 132,736 to 163,292 base pairs, especially in Decaloba. Structural rearrangements and species-specific repeat patterns were identified. Selective pressure tests revealed significant adaptive evolution in certain lineages, with several genes, including clpP and petL, under positive selection. Phylogenetic analyses confirmed the monophyly of subgenera Astrophea, Passiflora, and Decaloba, while Deidamioides appeared polyphyletic. Nuclear phylogenetic analysis based on 35S rDNA sequences supported the monophyly of Astrophea but showed inconsistencies within subgenus Passiflora compared to cp genome data. This study highlights the evolutionary complexity of Passiflora cp genomes, demonstrating significant structural variations and adaptive evolution. The findings underscore the effectiveness of plastid phylogenomics in resolving phylogenetic relationships and provide insights into adaptive mechanisms shaping cp genome diversity in angiosperms.

Introduction

Chloroplasts are essential organelles found in plants and some algae, playing critical roles in many crucial cellular processes, such as photosynthesis, fatty acid synthesis, amino acid production, and plant stress responses [41, 42, 83, 101]. This important organelle contains its own genetic material, distinct from the nuclear genome, which usually is organized as a circular DNA molecule, but in some cases is can be found as a linear structure [65]. The chloroplast genome (cp) is typically well conserved in angiosperms, ranging from 120 to 160 kilobases in size and encoding approximately 100 genes, including those genes encoding photosynthetic proteins and the ones encoding proteins used in the machinery for the cp gene expression [14, 82, 84].

Chloroplast comparative genomics can shed light on different evolutionary processes and functional adaptations, revealing, for example, the potential of cp genome structure variation in plant adaptation to abiotic stresses or different environment conditions [32, 77]. This comparative analysis can also be used to explore the evolutionary relationships among species, and reveal instances of horizontal gene transfer, gene duplication, and loss, offering a more comprehensive understanding of the dynamic nature of chloroplast genomes [66, 76], [12]. The evolution of next-generation sequencing (NGS) technologies is revolutionizing the field of comparative organelle genomics, allowing for a more detailed and extensive analysis of cp genomes across a wide range of the angiosperm phylogeny [79, 88, 104]. Recent studies are discovering a surprising level of diversity in cp genome structures, contradicting the earlier notion of their relative conservation. The significant variation in gene content, order, and in some cases presence of inverted repeats is suggesting that cp genomes are more dynamic and adaptable than previously understood, with evolutionary mechanisms such as rearrangements playing a crucial role in their diversity [26, 30]. Notable examples include the highly diverse chloroplast genomes found in groups such as Campanulaceae [29, 43], Geraniaceae [11, 89] and Papaveraceae [8]. Another angiosperm group with significant cp genome variability is Passiflora, known as passionflowers, which exhibits an extremely dynamic cp genome evolution, resulting in a high number of rearrangements, gene losses, and a particular case of loss of one inverted repeat [9, 10, 78]. All these events make this genus an excellent model for studying cp genome evolution.

Passiflora is the richest genus from the family Passifloraceae, comprising ca. 520 species with great diversity in size and shape of flowers. The geographic distribution of this group is particularly neotropical, being found mainly in the Americas; however, occurrences in Southeast Asia, Oceania and Australia have also been reported [87]. The first phylogenetic studies on this genus [61, 62, 97] resulted in poorly resolution of the Deidamioides clade, suggesting the paraphyletic position of this subgenus. Considering that plastid phylogenomics has been extensively used to solve the relationships in different plant groups [56, 59, 60, 94, 99], this analysis could also be useful to get new insights on the relationships wthin Passiflora. Clarfying, the current accepted taxonomic status of Passiflora is its subdivision into four subgenera, Astrophea, Decaloba, Deidamioides and Passiflora. However, with unresolved relationships within Deidamioides and the addition of Tetraphatea, it is critical to review these classifications.

In this study, we obtained 11 new plastomes for the genus Passiflora and conducted comparative genomics and phylogenetic inferences in order to investigate the evolutionary relationships between species within the genus. We also tested the robustness of the plastid phylogenetic analysis in comparison to nuclear rDNA phylogenies. Our aims were to answer: (i) What are the implications of structural variations in the cp genomes for understanding the Passiflora phylogenetic relationships? (ii) What do gene sequence comparisons in Passiflora cp genomes reveal about the impacts of gene variation and positive selection? (iii) What are the implications of positive selection in cp genes for the rapid radiation and diversification of the genus Passiflora? (iv) How the phylogeny from cp genomes compare to the one obtained from nuclear 35S rDNA in terms of resolving phylogenetic relationships within Passiflora?

Material and methods

Plant material and chloroplast (cp) DNA obtention

The present study includes 11 species of Passiflora comprising the subgenera: Astrophea (1). Decaloba (3), Passiflora (6) and Tetrapathea (1). Fresh leafs were collected from the Italian National Collection of Passiflora, Ripalta Cremasca, Italy (Supplementary Table 1), and prior to cpDNA extraction, the intact isolation of chloroplast organelles was performed based on liquid nitrogen-sucrose gradient method [85]. Briefly, in this procedure, approximately 10 g of fresh leaves were frozen in liquid nitrogen and macerated. The macerated leaves are then resuspended in 200 ml of isolation buffer, containing 50 mM Tris–HCl (pH 8.0), 0.35 M sucrose, 7 mM ethylenediaminetetraacetic acid (EDTA), 5 mM 2-mercaptoethanol, and 0.1% bovine serum albumin. This mixture was incubated in the dark for 10 min before being filtered through Miracloth, followed by centrifugation at 1,000 × g for 10 min to form a pellet. This pellet was washed and re-centrifuged under the same conditions. For further purification, the pellet was resuspended in isolation buffer and layered over a 20/45% sucrose gradient, then centrifuged at 2,000 × g for 30 min to isolate chloroplasts on the gradient interface. These chloroplasts were collected, diluted, and centrifuged at 3,000 × g to finally obtain a pellet of purified chloroplasts.

Subsequently, the resulted chloroplast pellet was lysed in 2% CTAB buffer, and the protocol was followed by two extractions using chloroform:isoamyl alcohol (24:1), an isopropanol precipitation, and finally washing the pellet in an ethanol solution (70%), that was then dried and resuspended in 40 μL of TE buffer.

CpDNA sequencing and assembly

The cpDNA samples were used for long-read sequencing on the PacBio platform, and large-insert (10 kb) libraries were constructed using the Barcode method with 150 ng of pure high molecular weight DNA. The fragment sizes, quality, and DNA concentration of the libraries were checked using a Fragment Analyzer (Agilent Technologies) and Qubit 2.0 Fluorometer (Invitrogen). Sequencing was performed on two SMRT cells using P6 polymerase with C4 chemistry on a PacBio RS II platform at the NGI Platform (Uppsala, Sweden).

The raw data were demultiplexed and filtered based on quality, removing reads with a quality score below 0.75 and length shorter than 500 bp, then, the sequences were initially assembled and corrected using CANU assembler [44]. Chloroplast contigs were extracted in Geneious by aligning the assembled contigs to the complete Passiflora cp genomes available. The complete cp genomes were subsequently constructed in Geneious by connecting the contigs using the "de novo Assembly" function, and final assemblies were verified by mapping the raw reads to the final contig using the "Map to Reference" function in Geneious.

Chloroplast genome annotation and comparative analysis

The new complete plastomes were annotated using the GeSeq (Organellar Genome Annotation) online program with default settings in order to identify protein-coding gene sequences (CDS), rRNAs and tRNAs based on cp reference sequences and BLAST homology searches [86]. The annotation was followed by manual corrections for start and stop codons, and intron positions in the GenomeView software,the OGDRAW software was used for constructing the circular and linear plastome maps [1, 28].

A multiple sequence alignment was performed to assess synteny and potential rearrangements among the cp genomes obtained. The progressive aligner from Mauve v.2.4.0 [15] was used to identify Locally Collinear Blocks (LCBs), which can reveal conserved and rearranged genomic regions across the species studied. Additionally, Mauve also put these blocks in order and orientation illustrating evolutionary and structural variations. We used CPJSdraw to illustrate the dynamic changes occurring in the inverted repeats (IRs) junction sites, in order to identify potential expansions and contractions of the IRs.

Repeat identification

For the prediction of microsatellites, also known as Simple Sequence Repeats (SSRs), we used MISA-web [6]. The parameters for SSR search were defined as follows: Motifs ranging from one to six nucleotides in length, with minimum number of repeats set at 10 for mononucleotide, 5 for dinucleotide, and 4 for trinucleotide SSRs, and three repeats for tetra-, penta-, and hexanucleotide SSRs. Additionally, we used REPuter [50] to detect direct and palindromic repeated elements in the DNA sequences, with specific parameters set for a minimum repeat size of 30 base pairs and sequence identities of at least 90% (corresponding to a Hamming distance of 3).

Plastid phylogenomic studies

Phylogenomic inferences were conducted based on 60 available plastomes from Passifloraceae (11 generated in this study, and 49 from species whose cp DNA sequences were obtained from NCBI’s Genbank). In addition, to obtain a rooted tree, the plastome of Populus trichocarpa (Salicaceae) was used as outgroup (Supplementary Table 1).

We used a set of 68 chloroplast protein-coding genes to perform the phylogenetic analysis: atpA, atpB, atpE, atpF, atpH, atpI, ccsA, cemA, clpP, matK, ndhA, ndhB, ndhC, ndhD, ndhE, ndhF, ndhG, ndhH, ndhI, ndhJ, ndhK, petA, petB, petD, petG, petL, petN, psaA psaB, psaC, psaI, psaJ, psbA, psbB, psbC, psbD, psbE, psbF, psbH, psbI, psbJ, psbK, psbL, psbM, psbN, psbT, psbZ, rbcL, rpl2, rpl14, rpl16, rpl23, rpl33, rpl36, rpoB, rpoC1, rpoC2, rps2, rps3, rps4, rps8, rps11, rps12, rps14, rps15, rps19, ycf3, ycf4.

The gene sequences were extracted from our data set (consisting in a total of 61 taxa), aligned individually at nucleotide level in MUSCLE, and then all the individual alignments were contatenated in a phylip matrix. The matrix was then analyzed in ModelFinder [37] to determine the best evolutionary model according to the Akaike Information Criterion (AIC). The Maximum Likelihood (ML) analysis was performed using RAxML version 8.2.4 [81]. The GTR + G + I model of nucleotide substitution was selected for the whole dataset, and a bootstrap analysis was performed with 1,000 replicates.

Additionally, a phylogenomic analysis using the complete plastomes of 24 species from subgenus Passiflora was performed (Supplementary Table 1). The sequences were aligned at nucleotide level in MAFFT version 7.221, using the FFT-NS-2 algorithm with default settings [39]. The substitution model for Bayesian inference (BI) was estimated in ModelFinder [37] and selected according to AIC. The Bayesian inference was conducted in MrBayes, with 10,000,000 generations, sampling one tree every 1,000 generations and discarding the first 25% of trees as burn-in. The analysis convergence was monitored using an average standard deviation of frequencies below 0.01, effective sample size (EES) above 200 for all parameters and potential scale reduction factor (PSRF) values close to 1.0. The phylogenetic trees from ML and BI inferences were visualized and edited using FigTree version 1.4.4 [70].

Selective pressure analysis

We calculated the ratio of non-synonymous (dN) to synonymous (dS) substitutions (ω = dN/dS) for 68 cp protein-coding gene sequences. The coding sequences (CDS) from each gene was separately aligned using MUSCLE [18], with manual curation in Geneious. Subsequently, in order to identify genes potentially under positive selection, we calculated the dN/dS for each CDS alignment using the CODEML program in the PAML. While ω > 1 can suggest positive selection, it is essential that this ratio be statistically significantly greater than 1 before inferring adaptive evolutionary changes. In contrast, ω = 1 would indicates neutral selection, and ω < 1 points towards evidence of purifying selection [96]. Firstly, we used several models to estimate the selection pressure on genes: M0 (one ω), which assumes a single ω ratio for all sites; M1a (neutral), which allows for two categories of sites, one with ω = 0 and one with ω = 1; M2a (selection), which adds an additional category with ω > 1 to the M1a model, allowing for positive selection; M7 (beta), which assumes a beta distribution for ω across sites (0 < ω < 1); and M8 (beta & ω), which extends M7 by adding an extra category for ω > 1. The identification of positively selected sites was combined with the Naive Empirical Bayes (NEB) and the Bayesian Empirical Bayes (BEB) methods. To compare the models, we performed likelihood ratio tests (LRT) for M1a vs. M0, M2a vs. M1a, and M8 vs. M7, calculating the test statistic as twice the difference in log-likelihoods. The p-value was obtained from the chi-squared distribution of this statistic.

Additionally, to detect selection on specific lineages within the phylogenetic tree, particularly focusing on the clade of the genus Passiflora, we used branch models. We tested the one-ratio model, which assumes the same ω ratio for all branches, and the two-ratio model, which allows two different ω ratios: one for the foreground branches (specific to the Passiflora clade) and one for the background branches. We then performed the Likelihood Ratio Test (LRT) to determine whether the likelihood of the two-ratio model is significantly different from that of the one-ratio model by comparing two times the log likelihood difference. P-values were computed using a chi-square distribution with one degree of freedom [95].

35S rDNA sequence assembly and phylogenetic analysis

The complete 35S rDNA sequence was obtained by genome skimming from the sequencing data of 31 Passiflora species: 20 of them sequenced in a previous study using the Illumina next generation sequencing strategy [9], and 11 sequenced in this study using the PacBio sequencing approach (Supplementary Table 1). The de novo assembly of the 35S for each species was performed in NovoPlasty [16] for Illumina paired-end reads, and in Geneious for filtered reads generated with Pacbio. Annotation was conducted in Geneious and sequences were extracted to create two datasets for phylogenetic inferences based on 18S and 26S complete sequences.

Sequences were processed in a local pipeline, aligned individually in MUSCLE [18], and concatenated in an interleaved matrix. The best evolutionary model was estimated in ModelFinder [37] in accordance with AIC. The phylogenetic reconstruction was performed based on maximum likelihood (ML) analysis in RA × ML v.8.2.4 [81]. The GTR + G + I substitution model was selected to infer the phylogenetic relationships using non-parametric bootstrap analysis with 1000 repetitions.

Bayesian inference was performed using MrBayes v.3.2.5 [73] with GTR + G + I evolutionary model. Markov chain algorithm (MCMC) was conducted with 10,000,000 generations, sampling one tree every 1000 generations, with the first 25% of trees discarded as burn-in. The convergence of analysis was confirmed by an average standard deviation of frequencies below 0.01, effective sample size (ESS) above 200 for all parameters and potential scale reduction factor (PSRF) values close to 1.0. All trees were visualized and edited using the FigTree v.1.4.3 [70].

Results

Chloroplast genome structural features and gene content

Our analysis revealed significant variations in cp genome size, gene content, and structural organization. All Passiflora cp genomes obtained in the present study presented the typical quadripartite structure found in most angiosperms, consisting of a large single copy (LSC) region and a small single copy (SSC) region separated by two long inverted repeats (IRs) (Table 1, Supplementary Fig. 1). The genome sizes ranged from 132,736 bp in P. adenopoda (Decaloba) to 163,292 bp in P. rusbyi (Astrophea). The large single-copy (LSC) regions exhibited substantial size differences, from 47,752 bp in P. intricata (Decaloba) to 89,229 bp in P. rusbyi (Astrophea). However, while the SSC length was approximately 13 kb for all species analyzed, a large variation in the IR regions was observed. The IR regions, which are crucial for maintaining cp genome stability, showed variation of ca. 25 kb, spanning from 19,852 bp in P. adenopoda (Decaloba) to 55,553 bp in P. intricata (Decaloba). Regarding GC content, the cp genomes presented a relatively narrow range, from 36.2% in P. tetrandra (Tetrapathea) to 37.5% in P. racemosa (Passiflora). In total, the cp genome annotation resulted in 104 to 109 unique genes representing 70 to 75 protein-coding genes, 30 tRNA and 4 rRNA genes.

Table 1 Chloroplast genomic characteristics of 11 Passiflora species from different subgenera

A comparison between the subgenera reveals notable infrasubgeneric variation in Decaloba, with an example of the IR size differences between P. adenopoda and P. intricata, suggesting significant genomic restructuring in this subgenus. Meanwhile, the subgenus Passiflora shows consistent gene content across its species, and a variability of ca. 6 kb for the IR sizes in the species analyzed.

Comparative genomics

The comparative analysis of the cp genomes across different Passiflora species revealed distinct patterns of genomic organization, with 10 synteny blocks (LCBs) or conserved regions among the cp genomes aligned. In general, rearrangements characterized by inversions, IR expansions and contractions were detected (Fig. 1).

Fig. 1

Comparative analysis of cp genomes across different subgenera in Passiflora. The alignments were obtained in Mauve, and colored bars represent syntenic blocks across the genomes, while connecting lines show the correspondence between these blocks. Below the synteny blocks, gene annotations for each chloroplast genome are displayed, with rRNA genes highlighted in red. Additionally, the inverted repeat (IR) regions are marked as light red bars in the annotation section

We found a similar structure in the cp genomes of P. rusbyi and P. tetrandra, representing the subgenera Astrophea and Tertraphatea. No major rearrangements could be observed between these two species, which only differ by an increase of ca. 2 kb in the LSC region in P. rusbyi. In contrast, subgenus Decaloba demonstrates extensive genomic rearrangements, for instance, P. intricata exhibits notable expansions and contractions in several LCBs, representing mainly an expansion of ca. 25 kb in the inverted repeat regions compared with Astrophea and Tertraphatea. Additionally, Decaloba exhibited the highest infrasubgeneric variation, with the IRs varying in ca. 35 kb among the species. Apart from that, when comparing P. intricata with the other species within the same subgenus, P. adenopoda and P. xiikzodz, we found inversions and translocations as variations in number and structure of LCBs, suggesting a dynamic organellar genome evolution in this subgenus. Within the species of the subgenus Passiflora, the cp genomes maintain a high degree of synteny, and only expansions of the IRs were observed as rearrangements. The largest IR expansion, observed in P. racemosa, is approximately 6 kb larger than the smallest IR found in P. garckei. However, in comparison with other subgenera, an inversion of approximately 30 kb in the LSC region is observed among the species of subgenus Passiflora.

When comparing the junctions of the quadripartite structure across the cp genomes in different Passiflora subgenera, distinct patterns of gene arrangement were found (Supplementary Fig. 2). In Astrophea, represented by P. rusbyi, the gene boundaries at the junctions between the inverted repeats (IRs) and the small single-copy (SSC) regions show an IR expansion, incorporating a copy of the ndhH gene. This configuration is different from most other species, where the ycf1 gene typically is in the IRs boundaries. Decaloba is the most variable subgenus in terms of gene arrangements at the junctions, a variation largely resulting from observed rearrangements, including IR expansion and contraction. P. adenopoda exhibits a contraction of the IR, leading to rps3 gene spanning from the IRb into the LSC. On the other hand, P. intricata shows a large IR expansion, incorporating several genes found in other species in the LSC region, extending up to the ndhC and ndhK genes. Additionally, P. xiikzodz displays distinct gene placement with psaI located at one of the IR borders, and ycf2 spanning between the IR and SSC junction. In contrast to the variability observed in Decaloba, species within the subgenus Passiflora generally exhibit a conserved gene order at the junctions. However, a slight difference is observed, particularly with the rps15 gene, which in some species of this subgenus spans from the IR to the SSC region.

Repeat distribution

The distribution and characteristics of repeats in Passiflora cp genomes reveal distinct patterns across different species and subgenera (Fig. 2). The identification of simple sequence repeat (SSR) types revealed that mono and dinucleotide repeats are the most prevalent types across all species (Fig. 2a), with the highest numbers of mononucleotide repeats observed in P. intricata and P. xiikzodz, both species from subgenus Decaloba. The number of trinucleotide repeats varied, with the highest number in P. intricata and the lowest in P. tetrandra. Furthermore, the occurrence of tetranucleotide repeats is significantly lower, only one repeat was identified in both P. rusbyi and P. adenopoda, and no penta- or hexanucleotide repeats were found in any of the species. Regarding the proportional representation of SSR motifs (Fig. 2b), the A/T motif is the predominant in most species showing similar proportions. Other motifs like C/G, AG/CT, AAG/CTT, and ACG/CGT show varying proportions, indicating species-specific preferences for certain SSR motifs.

Fig. 2

Distribution and characteristics of repeats in Passiflora cp genomes. a) Distribution of SSR types (mono, di, tri, tetra, penta, hexa) across different Passiflora species; b) Proportional representation of SSR motifs within each species; c) Number of repeats categorized by repeat types (F: Forward, P: Palindrome, C: Complement, R: Reverse); d) Number of repeats within specified length ranges (30–39, 40–49, 50–59, 60–69, ≥ 70) for each Passiflora species

We also identified long sequence repeats (> 30 bp) from the types: Forward (F), Palindrome (P), Complementary (C), and Reverse (R) (Fig. 2c). The most common repeat type across all species was Forward, with P. tetrandra (subgenus Tetraphatea) and P. intricata (subgenus Decaloba) showing the highest counts, 414 and 420 respectively. Palindromic repeats were also highly represented, with these same species having the highest counts, while species from subgenus Passiflora presented similar count values, ca. 30. The complementary and reverse are the least commonly found repeats, apart from the 48 reverse repeats found in P. rusbyi (subgenus Astrophea). Repeats within the range of 30–39 bp are the most prevalent across most species (Fig. 2d), with P. rusbyi and P. tetrandra showing particularly high frequencies. However, in P intricata, the 40–49 bp range is the most prevalent. Longer repeats (> 70 base pairs) are relatively uncommon in subgenus Passiflora but were found in high frequency in two species of the subgenus Decaloba and in both species from the subgenera Astrophea and Tetrataphea.

Plastid phylogenomics

The final alignment matrix of the plastid protein-coding genes was 59,350 characters long, and the ML analysis resulted in a tree with high support (BS = 100) for most of the nodes (Fig. 3). We used Populus trichocarpa (Salicaceae) to root the tree, and the base of the tree we found the other Passifloraceae species, Adenia manii, Dilkea retusa and Mitostemma brevifilis with high support (BS = 100).

Fig. 3

RAxML phylogenetic reconstruction of Passiflora evolutionary history based on 68 chloroplast protein-coding genes under the GTR + G + I substitution model. Node supports are indicated by bootstrap values (BS)

Regarding the formerly described Passiflora subgenera, our results revealed the monophyly of Astrophea, Passiflora and Decaloba. In Astrophea, a clade consisting of P. haematostigma and P. rhamnifolia species (section Pseudoastrophea) was found, grouping as sister to P. rusbyi (BS = 100) from the section Botryastrophea. On the other hand, the species from section Capreolata (P. cerradensis and P. pittieri) grouped together (BS = 100). Our results show the polyphyletic positioning of species from subgenus Deidamioides, while P. arbelaezii (section Tryphostemmatoides) formed a clade (BS = 100) which is sister to Astrophea, P. obovata, in subgenus Deidamioides, was embedded in subgenus Decaloba. In addition, another group consisted of Deidamioides species was observed: P. contracta from section Tetrastylis and P. deidamioides from section Deidamioides (BS = 100) appeared as monophyletic and sister to Decaloba.

The species from subgenus Decaloba formed a monophyletic group; however, P. microstipula grouped with P. obovata (Deidamioides). A clade with high support (BS = 100) was observed consisting of species of supersection Auriculata (P. jatunsachensis, P. rufa, P. intricata and P. auriculata). The species of supersection Cieca (P. xiikzodz, P. suberosa and P tenuiloba) grouped together and sister to supersection Bryonoides. In our study, section Decaloba appears as paraphyletic, with the group of P. filipes and P. lutea as sister to the species of section Xeragona (P. capsularis and P costaricensis). In addition, P. tetrandra, the type species of subgenus Tetrapathea, grouped as sister to subgenus Decaloba.

Our analysis included 28 species from subgenus Passiflora, confirming the monophyly of this clade. However, when comparing the distinct supersections within it, paraphyly was observed. In the tree, two major clades with high support (BS = 100) were found in subgenus Pasiflora, one consisted of species from supersection Stipulata, and the other comprising species from different supersections (Coccinea, Distephana, Laurifolia, and Passiflora). Although most species from Stipulata grouped together as a single cluster, this supersection is paraphyletic with P. actinia grouping closer to the species from supersection Laurifolia and Passiflora.

Regarding the relationships within the supersection Stipulata, the placement of species into some sections appeared paraphyletic, P. loefgrenii is from section Kermesinae, the same of P. edmundoi and P. watsoniana; however, in our findings, it was embedded in section Granadillastrum, a clade formed by P garckei, P. menispermifolia, P. oerstedii and P. retipetala. The supersection Passiflora also appeared as polyphyletic: P. cincinnata, P. edulis, P. recurva and P. serratifolia that belong to the same series in this supersection grouped each one with different sub-clades with species from supersection Laurifolia. In addition, P. vitifolia, from supersection Coccinea, grouped with P. miniata and P. cristalina from the supersection Distephana with high support (BS = 100).

To compare the results with the clustering into supersections of subgenus Passiflora, a phylogenomic analysis was performed based on whole cp genome sequences (Supplementary Fig. 3). The tree topology was congruent with the one based on the cp genes, with the supersections Passiflora and Laurifolia as paraphyletic.

Estimation of dN/dS ratios and positive selection

We estimated the ratio of non-synonymous (dN) to synonymous (dS) substitutions was performed for 68 cp protein-coding genes from 61 species (57 from genus Passiflora). Various site models (M0, M1a, M2a, M7, and M8) were employed to compare model fit using likelihood ratio tests (LRTs). Of these, only comparisons of M1a vs. M2a and M7 vs. M8 were used to specifically test for the presence of positive selection. Considering the dN/dS ratios of Model 0, most of the genes were identified under purifying selection (dN/dS ratio < 1), indicating selective pressure to maintain their functions (Fig. 4). However, when comparing the models, the results showed several genes with sites evolving under positive selection, as evidenced by the log-likelihood values and the likelihood ratio test (LRT) results (Supplementary Table 2). For instance, the atpA shows a significant result with log-likelihood values of -5370.334007 for M2a and -5416,734,687 for M1a, which resulted in a likelihood ratio test statistic of 92.801360, indicating sites under strong evidence of positive selection (Supplementary Table 2). Other genes such as atpB, atpE, atpF, and ccsA also exhibit sites under significant positive selection with p-values well below the threshold of 0.05. Functionally, the genes under positive selection belong to categories that are crucial for cellular and metabolic processes. The atp genes (atpA, atpB, atpE, atpF and atpH) are involved in ATP synthesis, which is essential for cellular energy production. The ndh genes (ndhA, ndhC, ndhF, ndhI, ndhK) are part of the NADH dehydrogenase complex, which plays a key role in the electron transport chain and cellular respiration. Additionally, genes like rpoB, rpoC1, and rpoC2, which encode components of the RNA polymerase complex, showed very high likelihood ratio test statistics and significant p-values (for example, rpoB with a test statistic of 232.414876), indicating strong positive selection in these fundamental transcriptional machinery genes.

Fig. 4

Distribution of dN, dS and dN/dS (nonsynonymous/synonymous substitution rate) for 68 cp protein-coding genes obtained from codeml under model M0. The genes were grouped per chloroplast function categories. Each box represents the interquartile range (IQR), showing the middle 50% of the data, with the horizontal line within each box corresponding to the median values for the corresponding category. Outliers are labeled with gene names

The likelihood ratio tests (LRT) for positive selection using the branch-site model with subgenus Passiflora as the foreground branch reveal strong evidence of positive selection for several genes (Table 2). Specifically, the genes atpA, atpB, ccsA, clpP, matK, petA, petD, petL, psaA, psbA, psbF, rbcL, rpl16, rpoC1, rpoC2, rps15, and rps3 show significant LRT statistics and p-values (p < 0.05), indicating that they have undergone adaptive evolution. For instance, the rbcL gene shows a highly significant p-value of 1.45e-9. Similarly, rpoC2 with a p-value of 2.57E-10 and atpB with a p-value of 6.27e-5 also indicate a robust evidence for positive selection.

Table 2 Likelihood ratio test (LRT) results for positive selection in cp genes using the branch-site model from codeml and assuming the foreground branch as subgenus Passiflora. The null model assumes no positive selection, whereas the alternative model allows for positive selection on the foreground branch. A significant LRT statistic and p-value (< 0.05) indicate evidence for positive selection

Nuclear phylogeny based on complete 18S/26S (35S rDNA cistron) gene sequences

Through genome skimming, the sequencing of Passiflora cpDNAs enabled the assembly of the complete 35S rDNA cistron, which was used to perform a phylogenetic analysis based on the complete 18S and 26S gene sequences. The phylogenetic tree using both BI and ML methods resulted in similar topologies (Fig. 5).

Fig. 5

Maximum Likelihood and Bayesian phylogenetic inference of Passiflora genus based on the nuclear 18S/26S gene sequences

The results strongly support the monophyly of subgenus Astrophea (BS = 98, PP = 1). However, despite the high values that support the grouping of subgenera Passiflora (BS = 84, PP = 1) and Decaloba (BS = 100, PP = 1), the monophyly of these subgenera was not observed, since P. adenopoda (Decaloba) was embedded within subgenus Passiflora. Additionally, the subgenus Deidamioides appears paraphyletic, with P. deidamioides emerging as sister to P. tetrandra (Tetraphathea). While the phylogeny provides high bootstrap support for higher classification levels such as subgenera, it reveals low bootstrap support for numerous nodes within subgenus Passiflora, indicating unresolved relationships. The tree also shows a polytomy and very short branches, which contributed to the lack of resolution for species such as P. chaparensis, P. alata, and P. phoenicea. Furthermore, incongruences between different trees (Figs. 3 and 5), generated by two distinct datasets, were primarily observed in the clustering within subgenus Passiflora. Unlike the nuclear phylogeny, in the cp genome tree, P. adenopoda was placed into the Decaloba clade, highlighting the differences in classification derived from the different genomes.

Discussion

Structural rearrangements in Passiflora chloroplast genomes

In the past, chloroplast genomes of angiosperms were known to be highly conservative in structure and gene order, but studies on comparative genomics have been revealing that this is not always the case [8, 89, 93]. Passiflora, a genus within the Passifloraceae family, is an example of the diversity and complexity of cp genome structure and evolution, with significant structural rearrangements, including inversions, duplications, IR expansions/contractions, and even losses of certain gene regions [9, 10, 69, 78].

Our results confirm the highly dynamic evolution of cp genomes in this genus, as we identified gene losses and large-scale inversions, such as inversions and IR expansions/contractions, in the newly assembled cp genomes. For instance, the large inverted repeat regions (IRs), which are usually quite stable in most angiosperms, show considerable variation in Passiflora species. Our analysis resulted in an IR size difference of almost 35 kb in the subgenus Decaloba, leading to large scale cp gene rearrangements in the species of this subgenus. It is hypothesized that inverted repeat regions are crucial for maintaining cp genome stability [55, 67], and a recent study by Krämer et al. [45] demonstrates that the removal of these regions can lead to a reduced number of plastid ribosomes and an increased total number of plastid genomes, highlighting the functional importance of these regions. In general, the IRs of Passiflora exhibit a huge variation in their evolutionary history, which in some cases has led to higher substitution rates for genes incorporated into IR expansions [78], and in the most extreme case has resulted in the loss of one inverted repeat copy [9]. One of the potential sources of the cp genome rearrangements is the repeat sequences, which we found to be abundant in the cp genome of some Passiflora species. Repeats can facilitate recombination, a process that could result in rearrangements, such as inversions, duplications and deletions [55, 63, 72, 93].

The study of cp genome structure could significantly help resolving the evolutionary relationships among plant species. The cp genome rearrangements have been identified as potentially valuable for plant taxonomic studies [17, 34, 71]. Particularly, the cp gene order has proven to be a good marker for phylogenetic studies in Campanulaceae [13]. Additionally, in legumes, the loss of an inverted repeat (IR) region has enabled the classification of an extensive number of papilionoid genera into what is described as the Inverted Repeat-Lacking Clade (IRLC) [51, 58, 91]. Similarly, the rearrangements we identified could be a rich source of markers for circumscribing species and tracing the evolutionary histories in Passiflora.

Passiflora plastid phylogenomics

Plastid phylogenomics has become a valuable tool in plant taxonomy, using both a set of genes [33, 94] or the whole cp genome sequences [35, 100]. The increasing number of available cpDNA sequences and the application of phylogenomics have significantly enhanced the accuracy of phylogenies, resulting for example, in higher tree resolution at low taxonomic levels [53, 68, 92, 102].

Passiflora has a complex taxonomic history, with Killip [40] first proposing the division into 22 subgenera based on morphological traits. This number was later drastically reduced to four by Ulmer & MacDougal [87], who, considering both morphological and ecological information, suggested the existence of the subgenera Astrophea, Decaloba, Deidamioides, and Passiflora. This reduction highlights the significant morphological diversity within Passiflora species, leading to considerable distinctions between the subgenera. For example, Decaloba (ca. 230 species) are herbaceous vines with small fruits and flowers, whilst Passiflora (ca. 250 species) are lianas or herbaceous vines with large, colorful flowers and long tubes. Deidamioides (14 species) are characterized by tiny stipules, petiolar glands, small bracts, and a plicate operculum [87]. Astrophea (ca. 60 species) consists of woody lianas, shrubs, or small trees lacking tendrils or presenting short spines. This morphological diversity is mirrored in their pollination strategies. For instance, species in the supersection Tacsonia are adapted for hummingbird pollination [2], while others, such as those of the series Tetrastylis (Deidamioides), are pollinated by small and large insects and bats [7, 22]. Despite their predominantly neotropical distribution, 22 species from supersection Disemma (Decaloba) are found in Southeast Asia, Australia and New Zealand [47].

To address these taxonomic complexities, Feuillet and MacDougal [21] proposed subdivisions within the four subgenera into supersections. Further supporting this effort, the first molecular phylogeny based on nuclear and cpDNA markers revealed three well-supported clades: Astrophea, Decaloba, and Passiflora, but did not resolve the relationships in Deidamioides [61]. Subsequent studies also recovered Deidamioides as polyphyletic, proposing eight subgenera [97]. Our plastid phylogenomic analysis also supports these findings, revealing high support (BS = 100) for most nodes and confirming the monophyly of Astrophea, Passiflora, and Decaloba. Further adding to the complexity, Krosnick et al. [46] proposed Tetrapathea as a new subgenus, represented in this study by P. tetrandra. In our analysis, Tetrapathea was placed as a sister to Decaloba, corroborating previous phylogenetic findings [48].

Our trees also indicate the polyphyletic nature of subgenus Deidamioides, with P. arbelaezii forming a clade sister to Astrophea, P. obovata embedded within subgenus Decaloba, and P. contracta + P. deidamioides clustering as sisters to Decaloba. These findings suggest that the current taxonomic classification of Deidamioides does not accurately reflect its evolutionary history. Given these complexities, our results highlight the need for a comprehensive revision of subgenus Deidamioides. A more detailed phylogenetic study including additional nuclear markers and broader taxon sampling could provide deeper insights into the evolutionary relationships within this group to allow a revised circumscription. Within subgenus Passiflora, our analysis including 28 complete cp genomes resulted in paraphyly observed within distinct supersections, as in supersection Stipulata and Passiflora. This pattern within subgenus Passiflora may be due to recent divergence and rapid radiations in this subgenus, which would result in incomplete lineage sorting impacting the resolution of relationships in lower taxonomic levels.

Adaptive evolution in Passiflora cp genomes

After unraveling the evolutionary dynamics of Passiflora cp genomes, our phylogenetic study supports the occurrence of a rapid radiation in this genus, an event where species diversify quickly from a single ancestor into an extensive range of species. Importantly, rapid radiation is driven by ecological opportunities, which often leads to selective pressure that shape genomic adaptations. By examining the ratio of non-synonymous (dN) to synonymous (dS) divergence (dN/dS), we found that the majority of genes exhibited overall dN/dS ratios < 1, indicating that purifying selection predominates, a pattern which does not exclude the presence of positively selected sites within these genes. Notably, certain genes, including clpP, presented overall dN/dS ratios significantly greater than 1, providing strong evidence for positive selection acting on at least a subset of their codons.

The high dN/dS ratios for clpP gene found in our results was also observed previously for other angiosperm groups, such as Oenothera, Silene [19] and Acacia [90]. In Passiflora, the evolution of clpP gene is also marked by two independent events of intron loss in subgenera Decaloba and Passiflora [9]. The chloroplasts Clp proteases are crucial for protein degradation and plant development. The disruption of clpP gene, encoding the ClpP1 subunit, caused severe developmental defects in mutant tobacco plants, including shoot system loss [49]. Shikanai et al. [80] showed that clpP gene is vital for chloroplast biogenesis, affecting early developmental stages in leaf primordia and the maturation of grana, by degrading denatured proteins and rapidly turning over regulatory factors. Additionally, evidence from a study in Arabidopsis thaliana reveals the involvement of ClpP in the latter stages of both high light and cold shifts, suggesting a role in acclimation to changing growth conditions [103].

Most surprisingly, the branch-site model test for positive selection revealed several genes under adaptive evolution in the subgenus Passiflora. The genes identified under positive selection include those involved in photosynthesis (psaA, psbA, rbcL), transcription and translation (rpoC1, rpoC2, rps15), and metabolic processes (clpP, ccsA). Particularly, subgenus Passiflora, which includes more than 230 species, has shown signs of rapid radiation in previous studies, with a diversification rate of 12 species per million year [62, 74]. Interestingly, the petL gene, which encodes a small subunit of the cytochrome b6f complex, also showed evidence of positive selection. A recent study has demonstrated that PetL is crucial for photosynthetic cold acclimation, highlighting its role in enabling plants to adapt to cold environments [27]. This is particularly relevant as some species of the subgenus Passiflora inhabit high elevations in the Andes Mountains, where cold acclimation is essential for survival. Another gene with strong positive selection in this subgenus is the rbcL, which encodes the large subunit of the ribulose-1,5-bisphosphate carboxylase/oxygenase (RuBisCO), a key enzyme in the Calvin cycle of photosynthesis [4]. Regarding positive selection of rbcL, it was demonstrated that stressful environments of the Mediterranean coast facilitated the diversification of Rubisco within the genus Limonium, resulting in different photosynthetic CO₂ assimilation rates and plant growth responses to severe water stress [25]. Additionally, the adaptive evolution of the Rubisco L-subunit, particularly in response to climatic changes during the Oligocene and Miocene, was crucial for the diversification and environmental adaptability of species in the Brassicaceae [52].

Our study highlights that although purifying selection predominates in most chloroplast genes of Passiflora, some genes harbor sites under positive selection. Such adaptive evolution may contribute to the variability observed in this genus, including habitat diversity and morphological differences. Future research could further elucidate the specific environmental or functional drivers behind the positive selection of these cp genes, enhancing our understanding of plant adaptation and evolution.

Nuclear phylogeny based on 35S rDNA sequences

We also generated a phylogeny based on the complete 18S and 26S rDNA sequences, which resulted in some well-supported clades, such as the monophyly of subgenera Astrophea. However, in general, this nuclear phylogeny showed lower bootstrap support values compared to the results from plastids. These differences can be attributed to the fact that the nuclear phylogeny was based only on two highly conserved markers, while the plastid trees were derived from a set of genes with varying levels of polymorphism. This expanded dataset in the plastids analysis likely enhanced the robustness in the phylogenetic resolution.

In plant phylogenetic studies, the ITS region of the 35S rDNA gene is a potential nuclear marker due to its high level of variation [3, 20, 36, 38]. Previously, the ITS region was used to reconstruct the first molecular Passiflora phylogenies [48, 61]. However, the very high levels of ITS variation could lead to the loss of phylogenetic signal by saturation, mainly in higher differentiated levels, as reported by Muschner et al. [61]. In this scenario, we decided to use the complete 18S and 26S rDNA genes, more conserved, which resulted in strong support for some nodes. Maia et al. [54] also suggested the potential of these gene sequences to infer phylogenies when analyzing a large sample of angiosperm species, but the low phylogenetic signal of this region resulted in low support for some clades. Aygoren Uluer et al. [5] used the 26S rDNA gene to infer interfamilial relationships in Fabales and noted a lack of support across the majority of nodes in the tree, especially for Leguminosae, as well as concerns regarding possible paralogy problems. In this way, the resulting trees of the 35S ribosomal cistron must be interpreted with caution, since this region occurs in multiple copies in the genome impacting the assembly of this region and the detection of paralogy. This is because it is challenging to determine the accurate homology for a single copy from 35S rDNA [23]. To address this issue, the use of nuclear low copy genes in plant phylogeny [31, 75, 98] would be an alternative to propose a nuclear phylogenetic hypothesis for Passiflora. For example, using target capture approaches and the Angiosperm 353 probe set, which targets single nuclear genes and have been shown to provide high-resolution phylogenetic trees [24, 57, 64].

Conclusions

Our study reveals significant structural variations in the cp genomes of Passiflora species, elucidating their complex evolutionary relationships. These variations, including inversions, IR expansions and contractions, and differences in gene content, could be used as valuable markers to study the evolution of this genus. Gene sequence comparisons indicate that while some of the cp genes are under purifying selection, some genes, such as clpP, show evidence of positive selection, suggesting adaptive evolution. Additionally, the branch-site model test for positive selection revealed several genes under adaptive evolution specifically in subgenus Passiflora, including atpA, atpB, ccsA, petL, and rbcL. These findings highlight the potential adaptive roles which these genes played in the rapid radiation and ecological success of this subgenus. The phylogenetic analysis based on cp genomes provided a highly resolved tree, confirming the monophyly of the subgenera Astrophea, Passiflora, and Decaloba. In contrast, the nuclear 35S rDNA phylogeny did not provide much resolution for the tree, showing lower bootstrap support for many nodes within these subgenera. Notably, the polyphyletic position of subgenus Deidamioides, indicated by our findings, suggests the need for further studies and possibly a taxonomic revision to accurately reflect its evolutionary history.