Sampling

We primarily sampled Vincetoxicum fuscatum and V. hirundinaria subsp. adriaticum, the hypothesized parental species of V. pannonicum (Borhidi,2003,2006). To investigate the possible involvement of other parental species, we also included other European taxa (Table 1), some of which – similarly to V. hirundinaria subsp. adriaticum – are thought to be subspecies of a remarkably variable species, V. hirundinaria Medik. These taxa differ from each other only in minor traits (Markgraf,1972). Unlike V. hirundinaria subsp. adriaticum, all other subspecies are associated with Eastern European populations. We included species that are known to be diploid (2n = 22) and one species of unknown ploidy (V. maeoticum Barbar.). Samples were placed in silica-gel and were stored at room temperature until DNA extraction. Sequences of temperate Asian outgroup species (V. macrophyllum Sieb. & Zucc. and V. mongolicum Maxim.) and of one accession of V. rossicum (Kleopow) Barbar. were obtained from GenBank (Table 1).

Genetic Markers

The nuclear ribosomal internal transcribed spacer (nrITS) region is commonly used for plant phylogenetics (Baldwin et al.,1995) and has been demonstrated to preserve signs of hybridization and polyploidization (Wendel,2000); we therefore chose this nuclear marker for the current study. Additionally, to reconstruct the phylogeny of the maternally inherited plastid DNA, we screened the four commonly used markers accD-psaI, petA-psbJ, psbM-trnD, and trnH-psbA for variable positions. According to classical phylogenetic knowledge (Wendel & Doyle,1998), if diversification of V. pannonicum is truly a result of allopolyploid hybridization, incongruent placement on nuclear and plastid phylogenetic trees can help identify the parental species. In the case of autopolyploidy, we expected the two independent datasets to place the species in similar positions on both phylogenetic trees.

Laboratory Procedures

Genomic DNA of samples collected in the field was isolated as described by Laczkó et al. (2019). Polymerase chain reactions (PCR) were conducted in 25-μL reaction volumes containing 1× reaction buffer (Thermo Scientific; Carlsbad, CA, U.S.), 0.2 mM of each dNTP (Thermo Scientific), 25 mM MgCl₂, 0.25 mg bovine serum albumin (Invitrogen; Carlsbad, CA, U.S.), 0.2 μM of each primer, 0.03 U DreamTaq Green DNA Polymerase (Thermo Scientific), and 1 μL unquantified template DNA. The nrITS region was amplified using primers ITS1A (Gulyás et al.,2005) and ITS4 (White et al.,1990), and the plastid region was amplified using primers trnH^GUG (Tate & Simpson,2003) and psbA (Sang et al.,1997). Further genic markers of the plastid genome (accD-psaI, psbM-trnD, and petA-psbJ) were amplified using primers published by Small et al. (1998) and Shaw et al. (2005,2007). PCR conditions for nrITS were as follows: initial denaturation at 94 °C for 3 min, followed by 35 cycles of denaturation at 94 °C for 30 s, annealing at 51 °C for 30 s, extension at 70 °C for 40 s, and a final extension step at 72 °C for 7 min. We used the following touchdown PCR regime for all plastid markers: 94 °C for 3 min, 20 cycles of 94 °C for 30 s, 58 °C for 30 s in the first cycle and subsequently lowering the temperature by 0.5 °C in each cycle, and 72 °C for 50 s; the same cycling conditions were repeated for another 20 cycles but with constant annealing temperature of 48 °C in each cycle, followed by final extension at 72 °C for 7 min. PCR amplicons were sequenced using a commercial service (Macrogen, Korea). In three samples (“adrKam,” “pannOdv,” and “pannSzars”), we observed unambiguous double-peaks in the nrITS sequences which suggested paralogous copies (Nieto-Feliner & Roselló,2007); we therefore performed vector cloning to separate these paralogous variants. Purified PCR products were transferred and ligated into a pGEM-T Easy Vector System II (Promega; Madison, WI, U.S.), according to the manufacturer’s instructions. Transformed cells were screened with ampicillin, and recombinant plasmid DNA was isolated from white colonies by suspending them in 40 µL sterile MilliQ water followed by incubation at 98 °C for 5 min and centrifuging at 20,000 × g for 1 min. The nrITS region of eight clones per individual was amplified and sequenced as described above. For subsequent phylogenetic analyses, we chose one representative from each nrITS ribotype groups (Figure 2A).

Phylogenetic Analyses

We tested sequences for sequencing errors and aligned them using MUSCLE software v.3.8.31 (Edgar,2004). For the trnH-psbA region, gaps were reduced to singletons and were coded as characters by applying the rules of simple gap coding proposed by Simmons and Ochoterena (2000) as implemented in SeqState software v.1.4.1 (Müller,2005). GenBank accession numbers of sequences generated in the present study are shown in Table 1.

We analyzed the nuclear and plastid datasets separately using maximum parsimony (MP), maximum likelihood (ML), and Bayesian inference (BI) approaches. Heuristic searches using MP were conducted in Paup software v.4.0b*10 (Swofford,2003) and with the tree bisection-reconnection swapping algorithm. Ten trees were held in each iteration with 1,000 random replicates, and statistical support of the tree topology was estimated by 1,000 nonparametric bootstrap (bs) pseudo-replications. Branching events were considered to be supported if the bootstrap value was >70.

For the ML reconstruction, we used the PhyML online server (Guindon et al.,2010) with the smart model selection (SMS) algorithm (Lefort et al.,2017) and default parameters. Branch support was tested using approximate likelihood ratio tests (SH-aLRT). Nodes were accepted if they produced SH-aLRT values >0.8.

For phylogenetic reconstruction based on BI, we used a Markov chain Monte Carlo algorithm as implemented in MrBayes v.3.2.6 (Ronquist et al.,2012). We set the substitution model to GTR+I (nst = 6; rates = invgamma) for both datasets, the parameters of which were adopted from the SMS results. We analyzed nrITS and plastid regions as one partition in separate analyses. Four chains were run for 25 million generations using three independent runs, saving every ten-thousandth generation. We assessed convergence and effective mixing of runs using Tracer software v.1.6 (available from http://beast.bio.ed.ac.uk/Tracer). Only branches with a posterior probability >0.95 were evaluated.

DNA sequences within the ingroup revealed a rather shallow differentiation (see “Results” section), which may cause conventional phylogenetic methods to perform poorly (Crandall,1996). Thus, to investigate genetic distances between samples, we reconstructed the NeighborNet of ribotypes using Kimura-2-parameter genetic distances in SplitsTree software v.4.14.4 (Huson & Bryant,2006) which can also reveal possible past reticulation events. Additionally, we applied the phi statistics to the nrITS dataset to test for evidence of recombination. We carried out the haplotype network reconstruction of the plastid sequences using the R package PEGAS (Paradis,2010).

Flow Cytometry

We aimed to test ploidy levels of our samples and used flow cytometric analyses to assess DNA content of cells. We measured DNA content of both V. pannonicum and V. hirundinaria subsp. adriaticum using 30 measurements per species. For these analyses, we used fresh leaves and a Becton Dickinson FACScan flow cytometer (BD Biosciences; San Jose, CA, U.S.) fitted with an Ar laser lamp (488 nm). Leaf tissues of the samples and of an external standard (Bellis perennis L., 2C = 3.38 pg) (Schönswetter et al.,2007) were chopped using a razor blade in a plastic Petri-dish containing 2 mL ice-cold Galbraith’s buffer (Doležel et al.,2007). The suspension was filtered through a 40-μm nylon mesh to remove tissue debris and was then incubated on ice for at least 10 min. Isolated nuclei in the filtered suspension were stained using propidium iodide, and relative fluorescence intensity of 10,000 particles was recorded. Sample/standard ratios were calculated using the means of fluorescence histograms visualized using FCS Express 4 software (De Novo Software; Pasadena, CA, U.S.).

Sequence Variation

The aligned nrITS data matrix consisted of 700 characters, 26 of which were variable and 20 were parsimony-informative. We found evidence of multiple nrITS copies underlying the direct sequences in three cases: two samples of V. pannonicum (samples “pannOdv” and “pannSzars”) and one sample of V. hirundinaria subsp. adriaticum (sample “adrKam”) displayed additive polymorphic sites as defined by Fuertes-Aguilar & Nieto-Feliner (2003). Eight vector-cloned products of each of these samples were sequenced, and these sequences were included in the analysis based on genetic distance; identical sequences were collapsed to a single clone ribotype for phylogenetic analyses. Interestingly, the nrITS sequence of the third sample of V. pannonicum (“pannJen”) did not show signs of intraindividual sequence polymorphism.

Out of the four tested plastid DNA regions, three appeared to be monomorphic within the ingroup. The only variable locus was trnH-psbA, thus, we present a plastid phylogeny based exclusively on this region. The trnH-psbA matrix consisted of 266 nucleotides and 10 gap characters. We found eight single-nucleotide polymorphism characters to be variable in the matrix, four of which were parsimony-informative. Nine of the 10 gap characters appeared to be polymorphic within the ingroup, and five were apparently parsimony-informative.

Phylogenetic Reconstruction

The parsimony analysis of the nrITS dataset produced two equally parsimonious trees (MPT) at 40 steps with a consistency index (CI) of 0.7750, a homoplasy index (HI) of 0.2250, and a retention index (RI) of 0.9011. These trees differed from each other only in unsupported branches based on bs. The first MPT with bs support is shown in Figure 1A. The resulting topology of ML and BI analyses were in line with MP after collapsing unsupported branches into polytomies; therefore, only support values these analyses are shown in Figure 1. The two outgroup species (V. macrophyllum and V. mongolicum) were separated with high support. Within the ingroup, we identified four clades that appeared to be metaphyletic, i.e., phylogenetic relationships between the otherwise supported clades could not be resolved unequivocally (Figure 1A). The first clade comprised only V. hirundinaria subsp. cretaceum, and the second clade comprised V. fuscatum, V. maeoticum, and one clone of V. pannonicum (“pannOdv” clone A01). The widespread Eurasian species V. hirundinaria s. str. (i.e., in a strict taxonomic sense) shared one clade with the Eastern European species (V. ucrainicum Ostapko, V. rossicum, and Vincetoxicum flavum Ostapko). All other accessions of V. hirundinaria subsp. adriaticum and V. pannonicum formed the fourth clade. Within this fourth clade, cloned sequences of each sample were shared between the two groups.

Heuristic searches using the plastid dataset resulted in a single MPT at 21 steps (CI = 0.9524, HI = 0.0476, RI = 0.9815). ML and BI analyses produced the same topology after collapsing unsupported branches using the thresholds mentioned above (Figure 1B). Based on trnH-psbA, the outgroup taxa formed a separate clade. We distinguished two groups within the European taxa, one of which comprised species from the Adriatic coast of Croatia (V. hirundinaria subsp. adriaticum), from the Balkans (V. fuscatum), and from Central Europe (V. pannonicum). Vincetoxicum hirundinaria s. str. was again assigned to a separate group with V. ucrainicum, V. rossicum, and V. flavum. These analyses also placed V. maeoticum in this clade. Within the two groups above, we exclusively identified metaphyletic lineages (Figure 1B).

Figure 1

Phylogenetic relationships reconstructed using the nrITS sequences (A) and trnH-psbA sequences (B) presented on a phylogram. Support values on the branches refer to bootstrap support by MP/SH-aLRT probability calculated using phyML/Bayesian posterior probability as estimated by MrBayes. Identification codes used here refer to Table 1.

NeighborNet analyses of nrITS sequences (Figure 2A) identified the same groups, in general agreement with the phylogenetic tree reconstruction. Vincetoxicum macrophyllum was placed on a relatively long edge with respect to the ingroup, but close to V. mongolicum. Within the ingroup, three groups of samples were identified which were separated from each other by few splits and considerable genetic distances. The first group comprised V. fuscatum which clustered with V. maeoticum and two cloned sequences of V. pannonicum (both from sample “pannOdv”; clones A01 and F01) but was connected through dense reticulate edges. Vincetoxicum hirundinaria subsp. cretaceum was placed approximately equidistantly from this group and the second group of Eastern European taxa (V. ucrainicum, V. rossicum, and V. flavum) including a sequence of the widely distributed V. hirundinaria. The third group of sequences, which contained accessions of V. pannonicum and V. hirundinaria subsp. adriaticum, was separated from the other sequences by several splits but included a long edge indicating substantial genetic distance within this sequence group. Both direct nrITS sequences and vector-cloned sequences of these species clustered in this group. Within this group, two subgroups were identified which were connected through a long edge representing a relatively deep genetic split between samples; however, they were associated with two taxa (V. hirundinaria subsp. adriaticum and V. pannonicum) originating from different sampling sites (Table 1). Although numerous reticulations occurred within the network, the phi statistics produced no evidence of recombination (p = 0.70).

Figure 2

NeighborNet of nrITS ribotypes reconstructed in SplitStree based on K2p distances (A) and haplotype network of trnH-psbA haplotypes reconstructed using PEGAS (B). Small unlabeled dots in (B) represent hypothetical haplotypes that were not sampled. Identification codes used here refer to Table 1.

The haplotype network reconstruction of trnH-psbA sequences placed V. macrophyllum equidistantly between V. mongolicum and all accessions of V. hirundinaria subsp. adriaticum and V. pannonicum, all of which shared an identical haplotype, apart from “adrJez” (Figure 2B). The latter was located one hypothetical haplotype away from other V. hirundinaria subsp. adriaticum samples and was connected to V. fuscatum through two hypothetical haplotypes. The probability of V. fuscatum connecting to the “central” haplotype of V. hirundinaria subsp. adriaticum and V. pannonicum through three steps was relatively low, however, this configuration cannot be excluded as possible placement. Vincetoxicum hirundinaria s. str. was located six steps away from this haplogroup and was connected through one step to all accessions of Eastern European samples (V. ucrainicum, V. rossicum, and V. flavum) that shared the same haplotype. Vincetoxicum maeoticum was placed one step away from the Eastern European samples but could also be connected to V. hirundinaria s. str. through an alternative link with two hypothetical haplotypes.

Flow Cytometry

The genome size of the diploid V. hirundinaria subsp. adriaticum appeared to be 2C = 1.76 ± 0.06–2.39 ± 0.02 pg (mean ± SD; samples “adrJez” and “adrKam,” respectively; mean of two samples 2C = 2.15 ± 0.34 pg), which is the approximate diploid genome size in the genus that corresponds to 2n = 22 in this taxon (Serra et al.,2001). The genome size of V. pannonicum (sample “pannOdv”) was 2C = 5.59 ± 0.05 pg, which is almost threefold larger than that of the maternal parent.

Hybridization and Allopolyploidy

Our results show the hybrid origin of V. pannonicum. Flow cytometric measurements question the tetraploidy of the species (as discussed below). Cloned ribotypes of nrITS sequences of V. pannonicum were placed into three groups (Figure 1A, Figure 2A). Two of the ribotpyes were genetically closest to the ribotype of V. fuscatum, a putative parent of the Pannonian species. The remaining ribotypes of V. pannonicum appeared to be closely related to V. hirundinaria subsp. adriaticum and clustered into two groups (Figure 1A, Figure 2A). This variation within species is most likely affected and may be explained by concerted evolution of nrITS paralogs (Álvarez & Wendel,2003; Bailey et al.,2003) within hybrid species, i.e., the nrITS array favors a given ribotype and converts different copies to this particular ribotype over time. However, not all nrITS copies were converted to the preferred ribotype, and several thus retained the parental ribotype (“pannOdv” clones A01 and F01). Therefore, paralogs within the nrITS region of V. pannonicum show direct evidence of the hybrid origin of this species as ribotypes from both parental species were isolated from the hybrid by vector-cloning. Alternatively, V. maeoticum may be a parental species, based exclusively on the topology of our nrITS tree (Figure 1A). However, this seems very unlikely as the distribution area of V. maeoticum is substantially more restricted and is distant from that of V. pannonicum; therefore, biogeography renders V. fuscatum a more likely parental candidate.

The trnH-psbA haplotypes of V. pannonicum were identical to those of V. hirundinaria subsp. adriaticum regarding a 276-bp long, fast evolving DNA sequence. The confirmed polyploidy of V. pannonicum in combination with the presented genetic structure corroborates the hypothesis of Borhidi (2003,2006) suggesting that the species emerged as a result of allopolyploid hybridization between the taxa V. fuscatum and V. hirundinaria subsp. adriaticum, where – based on the genetic relation of V. pannonicum to the latter taxon based on plastid data – V. fuscatum should have been the pollen donor. The shallow genetic differentiation between V. pannonicum and the parental species together with incomplete concerted evolution suggest relatively recent (postglacial) emergence of this species, as postulated by Borhidi (2003,2006); however, further research is needed to assess the divergence time of the species. It should be noted that all V. pannonicum samples used in the present study share a plastid haplotype with V. hirundinaria subsp. adriaticum, which would make the latter taxon an exclusive maternal parent (i.e., the ovule donor) of our focal species. This may suggest unidirectional hybridization where V. hirundinaria subsp. adriaticum was consistently an ovule donor to V. pannonicum.

We identified a remarkable discordance between the two datasets with respect to the position of V. maeoticum. The nrITS sequences placed these samples close to V. fuscatum (Figure 1A, Figure 2A), whereas trnH-psbA sequences suggested that the species were closely related to the V. hirundinaria group. This incongruence suggests another hybridization event most probably between V. fuscatum and a species from the V. hirundinaria group, which may be a further example of allopolyploidization in this genus. Unfortunately, insufficient resolution of the plastid dataset and unknown ploidy of V. maeoticum did not allow identification of the parent species.

Flow Cytometry

The genome size of V. pannonicum was almost threefold larger than that of the diploid maternal parent, V. hirundinaria subsp. adriaticum. The general pattern of genome downsizing (Leitch & Bennett,2004) would thus suggest hybridization between a diploid and a tetraploid parent. Indeed, Strid & Franzén (1981) reported V. fuscatum to be tetraploid. Thus, we cannot directly corroborate the ploidy level of our V. pannonicum sample from the same locality as originally reported based on chromosome counts (Borhidi,1968). One possible explanation for this paradox may be genomic reorganization that led to cytological diploidization which is a prerequisite of successful allopolyploid formation (Ma & Gustafson,2005). Such diploidization may have caused our focal species to appear tetraploid even though it was originally hexaploid. Our results thus warrant further research in this regard and call for revising the chromosome count of V. pannonicum using modern cytological methods.

Taxonomic Consequences

Based on our results, we suggest some taxonomical consequences. The Adriatic endemism V. hirundinaria subsp. adriaticum is commonly accepted to be a subspecies of V. hirundinaria, as listed in the Flora Europaea (Markgraf,1972). However, we observed substantial genetic distance between this taxon and V. hirudinaria s. str., under which it is classified, in the plastid and in the nuclear dataset. The genetic distance between these taxa is comparable to that between one of the outgroup species (V. macrophyllum) and the ingroup, in both datasets. Although we used only a single sample to represent V. hirundinaria in our analyses, the congruent phylogenetic placement of V. hirundinaria as closely related to Eastern European species in the study of Liede-Schumann et al. (2012) makes us confident that our sample is a reliable representation of what we believe to be typical of this species. Even though there was substantial phylogeographic structuring within V. hirundinaria, the large phylogenetic distance observed in both phylogenetic markers (Figure 1, Figure 2) suggests unequivocal separation of V. hirundinaria subsp. adriaticum from V. hirundinaria. We therefore propose that the Adriatic species be accepted as a separate species, which was originally referred to as V. adriaticum Beck.