Phylogeny of Aconitum Subgenus Aconitum in Europe

Phylogenetic relations within Aconitum subgen. Aconitum (Ranunculaceae) in Europe are still unclear. To infer the phylogeny of the nuclear (ITS) region and chloroplast intergenic spacer trn L (UAG) - ndh F of the chloroplast DNA (cpDNA), we analyzed 64 accessions within this taxon, 58 from Europe and six from the Caucasus Mts. Nuclear ITS sequences were identical in 51 European and two Caucasian accessions, whereas the remaining sequences were unique. cpDNA sequences could be categorized into five haplotypes, i.e., A–E , including a European-Caucasian Aconitum haplotype B . Ten cpDNA sequences were unique. A 5-bp indel distinguished the diploids from the tetraploids. None of the extant European diploids were basal to the tetraploid local group. A phylogenetic tree based on combined ITS and cpDNA sequences (bayesian inference, maximum likelihood, minimal parsimony) placed Aconitum burnatii (Maritime Alps, Massif Central) and A. nevadense (Sierra Nevada, Pyrenees) in a sister group to all other European species. A Bayesian relaxed clock model estimated the earliest split of the Caucasian species during the Late Miocene [ca. 7 million years ago (Mya)], and the divergence of A. burnatii and A. nevadense from the European genetic stock during the Miocene/Pliocene (ca. 4.4 Mya). Diploids in Europe are likely to be descendants of the Miocene European-Caucasian flora linked with the ancient Asian (arctiotertiary) genetic stock. The origins of the tetraploids remain unclear, and it is possible that some tetraploids originated from local, now extinct diploids. Both the diploids and tetraploids underwent rapid differentiation in the Late Pliocene – Quaternary period.

In a preliminary study, we found two cpDNA haplotypes in A. subgen. Aconitum from the Carpathian Mts that generally fit the cytogenetic (diploids vs. tetraploids) and taxonomic sectional division (sect. Cammarum vs. sect. Aconitum) (see Mitka et al., 2016) criteria.
Here, we aimed to resolve the complicated genetic relationships among the Aconitum taxa throughout its European range, using plastid (cpDNA) and nuclear DNA (internal transcribed spacer, ITS) sequences of species distributed across Western, Central, and Southern Europe, and in the Caucasus Mts. The primary purpose of our phylogenetic analyses was to demonstrate the relationships between the studied species and regions; thus, the pattern of clades retrieved here should not be used solely as a justification for taxonomic decisions (Hörandl, 2006). Taking the differences in the cytogenetic and ecological profiles of the European and Caucasian Aconitum into consideration, we attempted to determine if: (i)  .
Aconitum L. sectio Aconitum subsectio Aconitum (2n = 32) the European tetraploids originated in situ from the diploid genetic stock, and (ii) genetic signatures exclusive to diploid and tetraploid species exist, using phylogenetic analyses based on ITS and cpDNA sequences (trnL (UAG) -ndhF region).

Taxon Sampling
The present study included 64 accessions representing A. subgen. Aconitum in Europe, all of which were sequenced for the first in this study. These accessions were as follows: sect.

DNA Extraction, Amplification, and Sequencing
Recently collected samples (stored as silica-dried leaves) or herbarium specimens of all accessions were obtained (Table S1). Samples for DNA extraction were prepared from these materials, using ca. 2 cm 2 of the fully developed leaf blade with no symptoms of damage due to insects or fungal infections. Samples were ground in 2 mL microcentrifuge tubes with three stainless steel beads (ϕ 3 mm) by shaking in an oscillation mill Germany) for 4 min at 25 Hz. DNA was then extracted separately for each sample with Genomic Mini AX Plant DNA extraction kit (A&A Biotechnology, Poland), according to the manufacturer's protocol.
Two target fragments were used for phylogenetic reconstruction: a fragment of the maternally inherited cpDNA separating plastid trnL (UAG) and ndhF genes (positioned between sites 115,891 and 114,942 relative to the A. kusnezoffii complete plastid genome), and the biparentally inherited ITS region of the ribosomal RNA gene cluster, a tested marker in Aconitum allowing resolution of the infrageneric phylogeny within the genus (Jabbour & Renner, 2011b;Kita & Ito, 2000;Kita et al., 1995;Luo et al., 2005;Utelli et al., 2000;L. Wang et al., 2009).
Successful amplification was confirmed by agarose gel electrophoresis, and positive PCR products were purified using Clean-Up DNA purification kit (A&A Biotechnology). The purified PCR products were used as templates in the sequencing reactions.
Sequencing was performed using BigDye Terminator v.3.1 Cycle Sequencing Kit (Life Technologies, USA) in a T100 thermal cycler (Bio-Rad) and 3500 Series Genetic Analyzer (Life Technologies), using standard protocols.

Sequence Alignment
Individual sequencing reads were examined carefully and compiled into full contigs with ChromasPro 1.7.6 software (Technelysium, Australia). As a relatively high level of nucleotide ambiguity was detected in the ITS sequences, the two independent reads for each contig had to be compared. A given nucleotide position was deemed ambiguous when two peaks were detected at the same position in the sequencing chromatogram, and the weaker peak was at least one third as high as the stronger peak in both independent reads (Fuertes-Aguilar & Nieto-Feliner, 2003). Sequences with ambiguous positions, encoded according to the IUPAC nucleotide code, were used for all downstream analyses. Both ITS and trnL (UAG) -ndhF contigs were aligned using the Clustal W algorithm (Thompson et al., 1994) of the MEGA 6 software package (K. Tamura et al., 2013), followed by manual adjustment.

Phylogenetic Analyses
Separate analyses of the trnL (UAG) -ndhF and ITS data sets produced no significant topological discordance for incongruent nodes with Bayesian inference (BI) and maximum likelihood (ML) bootstrap proportions >70%, and the datasets were therefore concatenated, yielding a matrix of 1,531 characters and 16 accession combinations (haplogroups), plus two accessions of the outgroup (Table 2). Substitution model parameters were estimated separately for each partition, using the GTR+G model (with four rate categories) for both the trnL (UAG) -ndhF and ITS regions. The model was selected by FindModel (https://www.hiv.lanl.gov/content/sequence/findmodel/findmodel.html), which uses the ModelTest script (Posada & Crandall, 1998).
Tree searches were based on a BI method (Rannala & Young, 1996) implemented in MrBayes v.1.10.4 (Huelsenbeck & Ronquist, 2001;Ronquist et al., 2005). The analysis was carried out by sampling every hundredth generation for 5 million generations, starting with a random tree. The first 1,250 million generations were excluded as burn-in after convergence of the chains, which was evaluated by the average standard deviation of the splitting frequencies below 0.01.
The ML analysis was performed for the concatenated data set in PhyML 3.0 (Guindon et al., 2010). The GTR model of nucleotide substitution was used.
Parametric bootstrap values for ML were based on 400 replicates.
DNA sequences were also analyzed using the maximum parsimony (MP) optimality criterion (Felsenstein, 2004;Fitch, 1971) in PAUP* version 4.0.b10 (Swofford, 2002). A heuristic search was conducted with random addition, tree bisection-reconnection (TBR) branch swapping, and the MULTREES option on. The consistency index (CI) and retention index (RI) were calculated with PAUP*, excluding uninformative characters. The strict consensus tree and support for its branches were evaluated by bootstrapping (BS) (Felsenstein, 2004), with 174 bootstrap replicates, each with 10 random stepwise additions performed using the same settings as above, and no more than 100 trees were retained per replicate.

Molecular Clock Analyses
We used TempEst v.1.5.1 to test the clock-like behavior of the concatenated dataset . Divergence dating was performed in Beast v.1.10.4 (Drummond et al., 2006;Drummond & Rambaut, 2007), which employs a Bayesian Table 2 Sixteen haplogroups (based on cpDNA haplotypes and ITS ribotypes) of Aconitum subgen. Aconitum in Europe and the Caucasus (including the 10 unique haplotypes) used in the phylogenetic analyses and the accessions within each haplogroup (species codes are given in Figure 1 and Table S1).

ITS ribotype
No. of accessions Markov chain Monte Carlo (MCMC) approach to coestimate topology, substitution rates, and node ages. All dating runs relied on the GTR+G model, a Yule prior, with uncorrelated and log-normally distributed rate variation across branches.
Several estimations of the divergence time between the subgenera Aconitum (ingroup) and Lycoctonum (outgroup) are available, considering the lack of any reliable Neogene Aconitum fossils. All these estimates were based on the generally accepted substitution rates, and served as secondary calibration points in Beast MCMC analyses, verified by cross-validated calibration approaches (Jabbour & Renner, 2011a

Phylogenetic Networks
To visualize genealogical relations among the cpDNA haplotypes, we used the TCS algorithm of Clement et al. (2000), implemented in POPART software (Leigh & Bryant, 2105). It is based on the concept of statistical parsimony and aims at producing an unrooted haplotype phylogenetic network, in which two haplotypes are joined by an edge only if the "probability parsimony" exceeds 0.95 for that edge (Huson et al., 2010).

Characterization of Nucleotide Data
The aligned ITS matrix included 18 unique sequences (16 ingroup + two outgroup) and a total of 632 positions, of which 557 were constant, 60 (9%) were parsimonyinformative, and 15 were parsimony-uninformative.
The combined (cpDNA + ITS) matrix consisted of 18 unique sequences and 1,531 positions, including 1,408 constant, 35 (2%) potentially parsimony-informative, and 88 parsimony-uninformative positions. Further information on the datasets and tree statistics from MP analyses of the nuclear and chloroplast regions and concatenated data is summarized in Table 3.

Chloroplast DNA Variation and Geographic Distribution
The 64 accessions of Aconitum subgen. Aconitum (excluding the outgroup accessions) could be categorized into five cpDNA haplotypes, i.e., haplotype A (24 specimens), B (19), C (four), D (two), and E (two), whereas the remaining 10 sequences were unique ( Figure 1A-D, Table 2). The trnL (UAG) -ndhF region could not be amplified in three accessions, namely accessions 73, 79, and 108 (Table 2). Haplotype A was distributed across Europe, hapl. B in Europe and the Caucasus, hapl. C was absent in the Carpathians but present in the Alps, Sudetes, and West Balkans, and hapl. D and hapl. E occurred exclusively in the South Carpathians or the Caucasus, respectively (Figure 1).

Nuclear ITS Genotype Variation
We observed extremely low nucleotide variation in Aconitum ITS sequences. The 64 accessions (excluding two of the outgroup) were arranged in two ITS ribotypes: R1 (51 specimens) and R2 (two specimens); the remaining seven sequences represented specific ribotypes (ribotypes R3-R9; Table 2). The ITS region of four accessions, namely 02, 46, 50, and 101, could not be amplified. R1 was distributed across Europe and Caucasus, and R2, R6, and R7 occurred only in the Caucasian Mts (Table 2).

Phylogenetic Analysis
The BI tree, based on the combined DNA plastid and ITS dataset arranged into 16 haplogroups (Table 2), is presented in Figure 2. It suggested the basal, statistically supported position of the Caucasian species ( B5), was also included in this clade (Figure 1).  Table 2 and Figure 1).

Haplotype Network
The TCS haplotype network of cpDNA haplotypes ( Figure 3)

Molecular Clock Estimations
Divergence time estimates for Aconitum in Europe and the Caucasus Mts are shown in Figure 4. The Bayesian analysis showed that the earliest split of the Caucasian genetic stock occurred around 7.3 Mya (Late Miocene). The earliest divergence in Europe was between Aconitum burnatii and A. nevadense, at the Miocene/Pliocene break approximately 4.4 Mya, and the remaining European diploids and tetraploids started to differentiate ca. 2.6 Mya. Diversification within the diploid and tetraploid sections appeared at the beginning of the Quaternary 1.8 Mya and continued till 0.5 Mya (Figure 4).

Geographic-Historical Background
The occurrence of Aconitum in Central Europe can be traced back to as early as the Late Miocene, as suggested by the Aconitum pollen deposits found in the Central Paratethys realm (Central Europe) (Stuchlik & Shatilova, 1987). The Caucasian and European lines diverged in the Late Miocene, and internally diversified mainly in  Table 2 and Table S1. the Quaternary, similarly to Ranunculus s. s. (Paun et al., 2005), Syringa (Kim & Jensen, 1998), and Wulfenia (Surina et al., 2014), highlighting the significance of this period for the evolution of the European mountain and high-mountain flora.
During the Late Miocene, the temperate forests along the southern coasts of Central and Eastern Parathetys (spread to Western Asia) were continuously replaced by open woodlands. The aridization trend corresponded to forest fragmentation and appearance of open landscapes, the development of grasslands and xerophytic plant communities, and disappearance of subtropical species from the fossil flora (see Dénes et al., 2015;Ivanov et al., 2011). The process was accompanied by a remarkable shift in the composition of fossil mammal assemblages from the Early/Middle Pannonian to the Late Pannonian, reflecting an increase in the seasonality and aridity in the Pannonian Basin area (Harzhauser et al., 2004). The ongoing fragmentation of forests could have disrupted the continuous Aconitum distribution along the southern coast of the Parathetys, contributing to its geographic isolation and evolutionary divergence. This process is illustrated by the Aconitum haplogroup B1 and A. nasutum in the Caucasus and Europe.

The Role of the Caucasus and Central Asia in the European Alpine System
The Caucasus represents a spatial and evolutionary link for many European genera of Asian origin (Ozenda, 2009), for example the genera Trollius L. (Després et al., 2003), Acer L. (Grimm & Denk, 2014), and Prunus (Volkova et al., 2020). The Asian genetic stock underwent further evolutionary migration to Europe (via the Caucasus and Balkans), i.e., phylogenetic divergence leading to the origin of sister taxa (Bräuchler et al., 2004;Dumolin-Lapègue et al., 1997;Song et al., 2016). The relationships between these regions appear to be older than the Quaternary (Hantemirova et al., 2016). This scenario may have applied to only the diploid line of Aconitum, and the current links between Europe and the Caucasus have been preserved in the diploid cpDNA of haplotype B. In a study on Aconitum in Bela Krajina (Slovenia), Starmühler (1996) discovered a Caucasian species, A. × tuscheticum (N. Busch) N. Busch (see Luferov, 2000), another putative relict of the South European-Caucasian floristic links.

Independent Evolution of Diploid and Tetraploid Lines
The origin and monophyly of the core European Aconitum subgen. Aconitum remains elusive. Molecular clock analysis dated the split of the tetraploids from the diploid stock at the beginning of the Quaternary (ca. 2.6 Mya). The sister position of the diploid and tetraploid lineages could be misleading, as they could not have originated in situ from a common ancestor and might represent independent genetic lineages in Europe. In this context, A. subgen. Aconitum in Europe could be a nonmonophyletic group.
The simplest "monophyletic" scenario is that the group originated in situ from an ancient, local diploid stock. Molecular analyses did not retrieve any extant diploid species as basal to the tetraploid group in Europe, as was observed for the Japanese tetraploids, where a diploid species, A. volubile Koelle, formed a monophyletic group with all East Asian tetraploid taxa, strongly suggesting it as their ancestral species (Kita & Ito, 2000).
However, some extant European tetraploids could have originated in situ from the local, possibly extinct, diploid genetic stock (Mitka et al., 2007), e.g., A. firmum and A. superbum, presently placed in the diploid clade. Their current position among the diploids is probably a relic of their initial diploid state and subsequent tetraploidization or horizontal gene transfer via intersectional hybridization (see below). Whole-genome duplication followed by diploidization in the ancient lineages support the hypothesis of Aconitum palaeoploidy (Park et al., 2020).
Excessive accumulation of 5S rDNA clusters in Aconitum chromosomes (FISH) in the tetraploid species (A. firmum and A. plicatum), followed by a reduction of the basal genome size (Joachimiak et al., 2018), likely occurred during diploidization, which is one of the stages of the cyclical process described as the "wondrous cycle of polyploidy" in plants. It could be a nonrandom process, as suggested by the retention of the original diploid ancestral progenitor genomes (Wendel, 2015), at least partially responsible for the paraphyly of the tetraploid and polyphyly of the diploid clades.
Weak support of the European diploid clade and links with the Caucasian genetic stock (haplogroup B1) might indicate its origin from multiple ancestor species that disappeared between 4.4-2.6 Mya. If this is the case, their roots could trace back to arctiotertiary temperate forest elements of Asian origin (Baskin & Baskin, 2016;Deng et al., 2015;Popov, 1983;Zhang et al., 2014). Some of them disappeared completely, whereas some underwent evolutionary divergence, including genome doubling (tetraploidization), when the global temperatures dropped markedly towards the end of the Pliocene (Abbott, 2008;Hultén, 1937).

Palaeonedemic Status of Aconitum burnatii/nevadense
Aconitum burnatii and A. nevadense represent the oldest genetic line in Europe, dating back to ca. 4.4 Mya. Their present position at the base of the entire European genetic stock could be a result of their initial diploid status and further palaeoploidisation or speciation by ploidy (see Brochmann et al., 1998;Favarger, 1960;Verlaque et al., 1997). According to this hypothesis, both A. burnatii and A. nevadense are autotetraploids, arising from conspecific Tertiary diploid parents, which are now extinct. It may have occurred at the time of the Neogene cooling phases, which culminated in the onset of major glaciation in the Northern Hemisphere (Pearson & Palmer, 2000). It is widely accepted that environmental stress resulting from the climatic cooling episodes was the driving force behind the widespread formation of polyploids. These species often occupy habitats different from those of their diploid parents and have been proposed as superior colonizers (Baduel et al., 2018;Soltis & Soltis, 2000;Stebbins, 1984). This relict group exhibited an independent evolutionary trajectory from the European Aconitum since the Miocene/Pliocene break. Another hypothesis states that the oldest diploid genetic lineage in Europe originated from the extant Central/East Asian diploid species, and this warrants further investigation.
A mechanism underlying the origin of such a pattern could be explained using the example of Silene ciliata Pourret, whose ancestral populations in the Mediterranean Basin might have been forced to migrate northward at the onset of climatic oscillations during the Late Tertiary and the Quaternary periods, resulting in the gradual taxonomic and phylogenetic splitting of the once monophyletic group (Kyrkou et al., 2015).

Reticulation Among European Aconitum
We believe that intersectional hybridization and subsequent genetic introgression are the most relevant factors responsible for the paraphyly of the tetraploid clade ( Figure 2) (Mitka et al., 2007Sutkowska, Boroń, et al., 2017;. Hybridization is frequent in Aconitum (Kita & Ito, 2000). Present horizontal transfer of the cpDNA gene between diploid and tetraploid Aconitum species (via reverse "triploid bridge") has been reported in the Tatra Mts (Sutkowska, Boroń, et al., 2017, Zieliński, 1982a, 1982b. Such horizontal gene transfer could be responsible for the observed interchange of cpDNA between the different sections of A. subgen. Aconitum in Europe.
The TCS algorithm showed ancient reticulation among the hypothetical ancestors of the A. burnatii/nevadense group and diploid/tetraploid haplotypes. All these observations indicate the ancient and complicated evolutionary history of the subgenus in Europe, including palaeoploidization, and recent and historical reticulations.

Conclusion
The diploid and tetraploid lines of Aconitum in Europe form independent phylogenies. The links of the European and Caucasian diploid species represented by haplotype B indicate its ancient history in the region and arctiotertiary Asian origins. Paraphyly in the tetraploid clade could have been caused by ancient and present horizontal gene transfer at the section level. High-mountain European tetraploids likely originated from unknown ancestors in the Miocene age, presumably of Asian origin, as early as ca. 2.6 Mya, which is the estimated divergence time for the diploid and tetraploid lines in Europe. Similarly, presumed ancestral diploids, presently extinct, could be ancestral to the extant tetraploid A. burnatii/nevadense line, independent at least since the Late Miocene/Pliocene (4.4 Mya), which may have undergone tetraploidization and evolutionary divergence at least ca. 2.3 Mya.

Supporting Material
The following supporting material is available for this article: • Table S1. List of the Aconitum accessions from Europe and the Caucasus used in this study.