Changes in floral shape: insights into the evolution of wild Nicotiana (Solanaceae)

Abstract

Floral shape and size play a role in plant diversification and reproductive isolation. Different floral forms can rise by selective pressures by pollinators/ecological constraints and/or genetic drift. Nicotiana (Solanaceae) has c. 82 currently recognized species grouped in 13 sections. Nicotiana forgetiana and N. alata belong to section Alatae and present different floral traits related to their primary pollinators. These species co-occur in a region of plant endemism in southern Brazil marked by a mosaic between open grasslands and Araucaria forest. Here, we conducted a population-level sampling across the range of N. forgetiana and combined geometric morphometric analyses and ecological niche modelling to shed light on the patterns underlying intraspecific floral shape variation. Corolla shape and size varied significantly across the geographical distribution of N. forgetiana and ‘rastroensis’, a putative new species. These floral shape differences were related to variations in temperature, precipitation and elevation. We also observed a range of intermediate floral traits in two populations, and our results of geometric morphometric analyses reveal morphological overlap between N. alata and N. forgetiana. Moreover, we found that habitat suitability for both species was impacted by past climatic oscillations, with severe reduction for N. forgetiana during the Last Glacial Maximum. We discuss the results to shed light on the evolution of N. forgetiana.

INTRODUCTION

Impressive diversity in flower shape and size characterizes the flowering plants (Armbruster, 2014). Corolla limb and corolla tube shape and size affect pollinator attraction (Gómez et al., 2008, 2009a), playing a role in reproductive isolation between species (Grant, 1994). However, the role of pollinators in disrupting species boundaries is also well documented (e.g. Ippolito, Fernandes & Holtsford, 2004). Floral shape can drive the tendency of pollinators to visit a specific species, implying evolutionary shifts and pollinator-mediated speciation (Castellanos, Wilson & Thonpson, 2004; Muchhala, 2007). On the other hand, such preferences can sometimes change, for example when a pollinator learns to obtain nectar rewards from non-preferential flowers (Riffell et al., 2013). Additionally, the floral traits adapted to attract different pollinators can emerge as a by-product of genetic drift, such as after dispersion or long-term population isolation (Hewit, 2011; Abbott, 2017), sometimes related to past climate oscillations, which changed habitat suitability for species. Moreover, climate fluctuation influences the emergence of new species and, consequently, interactions between them, resulting in hybridization (Hoffmann & Sgrò, 2011; Chunco, 2014). The behaviour of the pollinator on hybridization can lead to introgression and affect speciation dynamics (Tastard et al., 2012).

The family Solanaceae is important in economic and ecological terms and is species-rich in the Neotropical region (Olmstead, 2013). The ancestral area is in South America (Särkinen, Olmstead & Knapp, 2013; Dupin et al., 2016), and floral morphology varies among species in symmetry, shape, size and colour, making this family an excellent model with which to study pollinator attraction (Knapp, 2010). However, this wide variety of floral forms and pollinators cannot be explained soley by adaptive evolution. It is necessary to include ecological and other evolutionary aspects to understand the complex relationships between species pairs.

Nicotiana L. is the fifth-largest genus of Solanaceae, with > 80 recognized species grouped in 13 sections (Knapp, Chase & Clarkson, 2004; Knapp, 2020). Nicotiana has its highest species diversity and probable centre of origin in South America, in the Andean region, with long-distance dispersal events from the Andes towards North America and from the Atlantic towards Australia and Africa, including Nicotiana section Alatae Goodsp. in southern Brazil (Goodspeed, 1954; Clarkson et al., 2004). Most Nicotiana spp. occur in the Americas, with one species occurring in Africa and the remaining species in Australia, where the oldest section of polyploids, section Suaveolentes Goodsp., occurs (Clarkson, Dodsworth & Chase, 2017). The colonization of Australia by this group represents an adaptive radiation in the genus (Ladiges, Marks & Nelson, 2011). Nicotiana section Alatae (hereafter Alatae) is a lineage of related diploid species with n = 9 and n = 10 chromosomes (Clarkson et al., 2004, 2010; Lim et al., 2006) and diversity in floral traits and pollinators (Raguso et al., 2003; Kaczorowski, Gardener & Holtsford, 2005). This lineage occurs in grassland habitats in the Southern Cone (Knapp, 2020), mainly in Brazil, and comprises eight herbaceous species (Vignoli-Silva & Mentz, 2005). Alatae comprise a lineage of recent diversification (Lim et al., 2006; Clarkson et al., 2017; c. 3.0 Mya), and little is known about their mechanism of speciation. Pollinator shifts appear to be an essential driver of floral diversity for Nicotiana in general (McCarthy et al., 2015, 2019). In this regard, the hawkmoth- and hummingbird-pollinated flowers in Alatae show variation in corolla limb shape and colour (Kaczorowski et al., 2012) and volatile compounds (Raguso et al., 2006).

Species of Alatae have different distribution ranges; for example, N. alata Link & Otto (Figs 1, 2A) is widespread in southern Brazil, north-eastern Argentina and southern Paraguay, overlapping its distribution with N. forgetiana Hemsl. (Figs 1, 2B), which has a narrow distribution in the south-eastern Brazilian Atlantic Forest (Fig. 1). These two species display contrasting floral morphologies related to their pollinators: N. alata is pollinated by nocturnal hawkmoths, whereas N. forgetiana is visited primarily by hummingbirds and occasionally by small hawkmoths (Raguso et al., 2003; Ippolito et al., 2004; Kaczorowski et al., 2005). The two species are self-incompatible and have the same chromosome number (Ippolito et al., 2004).

Artificial hybrids between N. alata and N. forgetiana described as Nicotiana ×sanderae Hemsley were first obtained at the beginning of the last century and used for ornamental purposes (Watson, 1904). Crosses between N. alata and N. forgetiana in controlled experiments also produced viable and fertile hybrids (Ippolito et al., 2004; Bissell & Diggle, 2008, 2010), displaying flowers that range from white to purple–red in the F3 generation (Bissell & Diggle, 2008). A putative natural hybrid between those two species has not so far been observed in the field (Ippolito et al., 2004), although populations of two species occur near each other (Fig. 1). On the other hand, variation in floral morphology compared to N. forgetiana was reported in Serra do Rio do Rastro in Santa Catarina state in Brazil, which was recognized as a putative new species, ‘rastroensis’ (Kaczorowski et al., 2005). The name ‘rastroensis’ is in honour of the geographical area of occurrence, Serra do Rio do Rastro, a region within an area of mosaics between open field and forest characterized by high endemism of plant species (called subtropical highland grasslands; Iganci et al., 2011; Figs 1, 2C). The region is a zone of climate transition where climatic fluctuations in the Quaternary have caused floristic changes (Safford, 1999, 2007).

Here we evaluated floral shape and ecological suitability habitat in N. forgetiana, integrating geometric morphometric analyses and ecological niche modelling (ENM) to discuss the evolutionary dynamics of closely related species. We sampled N. forgetiana populations through its geographical distribution covering the rastroensis area and N. alata populations in the region of geographical overlap with N. forgetiana. We specifically investigated the following. (1) Do floral shape and size vary according to the geographical occurrence for N. forgetiana and rastroensis? (2) Is floral shape variation across the geographical range of N. forgetiana and rastroensis associated with climatic and elevational gradients? (3) Is there a putative natural hybrid between N. alata and N. forgetiana among the studied populations? (4) Did climatic oscillations during glacial and interglacial periods influence the potential distribution of these species? We discuss the results considering the evolution of N. forgetiana.

MATERIAL AND METHODS

Sampling and data collection

We analysed 177 plants from 11 N. alata populations, 95 from six N. forgetiana populations, 54 from three rastroensis populations, and 64 and 27 plants, respectively, from two populations classified a priori as atypical (AP01 and AP02), totalling 417 individuals (Supporting Information, Table S1). This sampling strategy covered the distribution range of N. forgetiana and included the occurrence area proposed for rastroensis (Fig. 1; Table S1). AP01 and AP02 were classified as atypical because plants display floral morphology, based on human eye assessment of corolla colour, shape and size, that did not fit with those for either N. alata or N. forgetiana (Fig. S1). Plants from AP01 and AP02 showed a mixture of floral traits found in N. alata and N. forgetiana. The number of individuals sampled per population ranged from two to 64 according to the number of plants flowering at the same time. We analysed adaxial limbs of the corolla and corolla tube for the same individual and one flower per individual/plant. Sampling effort for the adaxial limb of the corolla from N. alata was limited by crepuscular anthesis of this species (Table S1; Ippolito et al., 2004). Voucher specimens per population and morphotype were deposited in the herbaria ICN (Universidade Federal do Rio Grande do Sul) and BHCB (Universidade Federal de Minas Gerais) (Table S1). We collected all data over field expeditions during the flowering period of the species (September to December, Southern Hemisphere spring).

We collected photographic data in the natural habitat of the species. Flowers were removed from the plant to take pictures in a standardized way. We collected standardized photographic images using a Sony DSC-HX400V Digital Camera (Sony Co., Tokyo, Japan) featuring a 2.8–6.3/4.3–215 Carl Zeiss lens (Carl Zeiss Man. Co., Oberkochen, Germany) at a resolution of 1920 × 1080 pixels employing the macro function without flash (Teixeira et al., 2020). We obtained images of the adaxial surface of the corolla limb and profile view (corolla tube shape) per plant. We analysed all flowers just after anthesis to avoid ontogenetic effects (Gómez, Perfectti & Camacho, 2006).

Floral morphology

We used geometric morphometric analyses to characterize flower shape variation in N. forgetiana, N. alata and the putative new species rastroensis. We analysed corolla shape (corolla limb and corolla tube) and size variation by geometric morphometric analyses using the image data set from the natural habitat of the species. We analysed 264 flowers for the adaxial view of the corolla limb (AVCL = corolla limb) and 408 flowers for the profile view (PV = corolla tube shape) (Fig. 2; Supporting Information, Table S1). We organized a TPS file with all images using TPSUtil v.1.38 (Rohlf, 2015). We used the TpsDig 2 software (Rohlf, 2015) to digitize landmarks manually. For AVCL, we digitized 15 landmarks to the corolla tube opening, to the apex of each petal and at the junction of adjacent petals following the flower venation pattern (Fig. 2E) (Gómez & Perfectti, 2010; McCarthy et al., 2016). For PV, we digitized 11 landmarks positioned at the corolla tube extremities. The landmarks represent the limits between different tissues or anatomical structures (1, 3, 6, 7 and 10), and minimum or maximum curvature points (2, 4, 5, 8, 9 and 11) covering the inflated apex of the corolla tube (Bissel & Diggle, 2010; García et al., 2020) (Fig. 2F).

We performed a generalized Procrustes analysis (GPA) on the bi-dimensional coordinates of the landmarks to eliminate the effects of position, scaling and orientation (Rohlf & Slice, 1990). We used the resulting matrices to analyse floral shape independent of floral size. We ran a principal component analysis (PCA; Jolliffe, 2002) and canonical variate analysis (CVA) in MorphoJ (Klingenberg, 2011) using the matrices from GPA to summarize the shape variation in the dataset. We ran a permutation test using Procrustes distances with α = 0.05 and 10 000 permutations to assess the statistical significance of pairwise differences in mean shapes among species or populations. Corolla limb and corolla tube size were measured from photographs through the log-transformed centroid size of each image (Zelditch, Swiderski & Sheets, 2012), and the differences among populations or species were tested by analysis of variance (ANOVA). We performed multiple comparisons by a post hoc test (α = 0.05) of Tukey’s honest significant difference (HSD). For both views, we performed a multivariate analysis of covariance (MANCOVA) using shape as a response variable, size as a covariate, and species or population as predictors to account for allometric effects. We ran the ANOVAs and MANCOVAs in the geomorph R package with the advanced.procD.lm function (Adams, Collyer & Kaliontzopoulo, 2018).

Ecological niche modelling

To estimate changes in suitable habitats for N. forgetiana and N. alata across time, we ENM and calculated distribution models using present and past climatic conditions. We also measured niche overlap for the different periods. For the analyses, we obtained the occurrence data from our collection points and the SpeciesLink database (http://www.splink.org.br/), manually checked for taxonomic inconsistencies using specimen information and randomly filtered points closer than 1 km² from each other (Supporting Information, Table S2). We used the maximum entropy algorithm as implemented in MAXENT v.3.4.1 (Phillips, Anderson & Schapire, 2006) to calculate a species ecological niche. We implemented auto features with seven and ten cross-validation replicates for N. forgetiana and N. alata, respectively, 5000 iterations, 10 000 background points in each run and logistic output (Phillips & Dudik, 2008). To model suitable habitat for species, we used the 19 available WorldClim bioclimatic layers (v.1.4) at smaller arc-second resolution (Hijmans et al., 2005). We used Present variables (1960–1990), and we used the Mid-Holocene (c. 6 kya), Last Glacial Maximum (c. 21 kya; LGM), and Last Interglacial (c. 120–140 kya; LIG) conditions to model past distributions.

To crop bioclimatic layers, we used DIVA-GIS v.7.5 (http://www.diva-gis.org) from 16 to 35°S and 39 to 64°C. We implemented the pairwise Pearson correlation coefficients in ENMTools v.1.4 (Warren, Glor & Turelli, 2010) to identify highly correlated variables (R > 0.75), removing the correlated variables from the analyses (Peterson, 2007). We used the following variables: BIO6 (min temperature of coldest month), BIO7 (temperature annual range), BIO9 (mean temperature of driest quarter), BIO10 (mean temperature of warmest quarter), BIO12 (annual precipitation), BIO13 (precipitation of wettest month), BIO14 (precipitation of driest month), BIO15 (precipitation seasonality – coefficient of variation), BIO16 (precipitation of wettest quarter) and BIO19 (precipitation of coldest quarter). We checked the quality of the models through the area under the curve (AUC; Pearce & Ferrier, 2000).

Additionally, we estimated the environmental niche overlap between the species for each period from the mean results of Maxent models using ENMtools to calculate Schoener’s D (Schoener, 1968). Values of Schoener’s D range from 0 (species have wholly discordant ENMs) to 1 (species have identical ENMs). The identity test was performed in ENMtools with 100 replicates to assess the statistical significance of niche overlap (Warren et al., 2008, 2010).

Association between floral shape/size and climatic variables

To evaluate the covariance between climatic variables and elevation with the shape and size of the corolla between N. forgetiana and rastroensis, we used a two-block minimum square analysis (PLS) of two blocks (Rohlf & Corti, 2000) in geomorph (Adams et al., 2018). The PLS approach treats the two blocks of variables as symmetric, maximizing the covariance between data sets without assuming one block as dependent on the other (Monteiro, Duarte & Reis, 2003). In this case, we used the climatic variables from WorldClim bioclimatic layers (v.1.4) as block 1 and shape or size from PC scores as block 2.

RESULTS

Phenotypes and field observations across the distribution of N. forgetiana

Although it is challenging to distinguish flower colours by human eye, we observed a distinct variation in colour of the corolla limb and corolla tube across the distribution of N. forgetiana compared to the description of this species (pink to red flowers, Fig. 2D). For N. alata and rastroensis populations, we cannot verify heterogeneity in corolla colour (Fig. 2A, C). However, corolla pigmentation analyses are needed to give a more accurate pattern to this variation. We found two populations with visually mixed phenotypes (AP01 and AP02), suggesting hybridization (Supporting Information, Fig. S1). The flowers of all rastroensis individuals displayed a white ring at the mouth of the corolla tube. This trait was present in 88% of the N. forgetiana flowers, 77% of AP01 flowers and 22% of AP02 flowers. Although we did not make systematic observations, we recorded small bees collecting pollen from rastroensis flowers and bees collecting nectar from N. forgetiana flowers.

Floral shape variation in N. forgetiana and rastroensis

We performed geometric morphometric analyses to examine the dimensions of corolla limb and tube shape variation between N. forgetiana and rastroensis. We did not consider two groups a priori; instead, we grouped flowers based on where they were collected for the analyses. Geometric morphometric analyses generated 26 and 18 PCs for AVCL and PV, respectively. For AVCL, the first three PCs explained c. 55% of the variation in floral shape (Supporting Information, Table S3). PC1 (22.6% of the variation) and PC2 (20.1% of variation) showed an overlap between corolla limb shape found in N. forgetiana and rastroensis (Fig. S2A). For PV, the first three PCs accounted for c. 79% of the floral variation (Table S4). PC1 (22.6% of the variation) and PC2 (20.1% of the variation) showed that N. forgetiana and rastroensis occupy different morphospace (Fig. S2B).

In general, N. forgetiana and rastroensis showed more pronounced differences in corolla tube shape (PV) than in corolla limb shape (AVCL) (Fig. 3A, B). The first two CVs explained 58.7% and 81.7% of the differences between groups concerning corolla shape for AVCL and PV, respectively. The morphospace for AVCL (Fig. 3A) seems to reflect the degree of dissection of the floral limb outline (from more to less stellate) and the relative size of the corolla tube aperture. CV1 varies both in how stellate the floral outline is and in the relative size of the corolla tube opening. Positive values indicate fewer stellate flowers with relatively smaller tube openings, and negative values indicate more stellate flowers with relatively larger tube openings. CV2 seems to vary mainly in the rotation of the flowers and some variation in relative size of the corolla tube opening. For positive values, the left hand lateral and lower petals are closer together and the corolla tube opening is smaller, whereas for negative values, the right hand lateral and lower petals are closer together and the corolla tube opening is larger. Figure 3A shows a gradient of floral shape from N. forgetiana, with more stellate petals and larger corolla tube aperture, to rastroensis, with less stellate petals and smaller corolla tube aperture. For PV, CV1 was mainly associated with throat width, with N. forgetiana showing a change in throat shape in the region connected with the two lower petals (a downward twist), and rastroensis showing a smaller inclination of the throat and larger tube (Fig. 3B). Those results revealed two main groups corresponding to N. forgetiana and rastroensis, respectively. Pairwise comparisons of corolla shape between N. forgetiana and rastroensis populations were significantly different for AVCL and PV (Table 1). Nicotiana forgetiana showed high levels of corolla variation; Forg01 was the most different population (Fig. 3A, B).

Floral shape variation between N. alata, N. forgetiana and atypical populations

Geometric morphometric analyses from all populations sampled (N. fogetiana, rastroensis, N. alata and AP populations) generated 26 and 18 PCs for AVCL and PV, respectively. For AVCL, the first three PCs explained c. 64% of the variation in floral shape (Supporting Information, Table S5). PC1 and PC2 accounted for 40.6% and 14.7% of the total variation, respectively, showing a clear difference in floral limb shape between N. alata and the remain species or populations (Fig. S2C). For PV, the first three PCs accounted for c. 90% of the floral variation (Table S6). PC1 and PC2 accounted for 78.5% and 6.6% of the total variation, respectively (Fig. S2D). Nicotiana alata, N. forgetiana and rastroensis occupied different morphospace. The population AP01 overlapped with N. forgetiana, whereas AP02 overlapped with N. forgetiana and N. alata.

Corolla shape differed significantly among species and AP populations (Table 2). The first two canonical variates (CV1 and CV2) from AVCL explained 54.5% and 14.8% of the among-group differences, respectively (Fig. 3C). For PV, the first two CVs explained 57.9% and 10.6% of the differences, respectively (Fig. 3D). For AVCL, CV1 represents changes in the angles between petals of the floral limb. Negative values for CV1 corresponded to greater angles between petals, resulting in a pentagonal floral outline and relatively larger corolla tube aperture. Positive values from CV1 corresponded to smaller angles between petals, resulting in a stellate floral outline and a relatively smaller corolla tube aperture. CV2 seems to represent both the floral outline and relative size of the corolla aperture, similar to CV1, but to a lesser extent. Moreover, more stellate flowers are associated with relatively larger corolla tube apertures (Fig. 3C). For PV, CV1 represents changes in the corolla tube width and throat width and shape (Fig. 3D). Negative values of CV1 conformed to larger corolla tubes and a straight apex. Positive values of CV1 corresponded to narrow corolla tubes and a curved corolla tube apex. CV2 was mainly associated with differences in the orientation of the corolla mouth compared to the tube. Negative values of CV2 corresponded to greater curvature of the apex of corolla tubes, and positive values corresponded to smaller curvature of the apex of the corolla tubes (Fig 3D).

Flowers of N. alata and N. forgetiana formed two different groups in morphospace, and no overlap of individuals was observed between them with regard to corolla limb shape (Fig. 3C). Moreover, population Ala05 was the closest to N. forgetiana in the morphospace. However, for corolla tube shape (Fig. 3D), an overlap of individuals of these two species was observed. This overlap is mainly due to populations Ala04 and Ala05, which showed considerable intrapopulation variation (Supporting Information, Fig. S3). Nicotiana forgetiana populations had higher intrapopulation variation than N. alata for corolla limb and tube shape. Considering all populations collected throughout the distribution of N. forgetiana (Fig.3C), rastroensis overlaps in corolla limb shape mostly with population Forg04. The separation between N. forgetiana and rastroensis in the morphospace was more evident for corolla tube shape, and overlap was primarily because of population Forg01 (Fig. S3).

Atypical populations showed different patterns of overlap in corolla limb and tube shape with the two species. For AVCL, the AP01 population showed an intermediate corolla shape and coincided with both N. forgetiana and N. alata in the morphospace. The overlap with N. alata was with population Ala05 (Fig. 3C; Supporting Information, Fig. S3). However, for PV, AP01 showed more forms found in N. alata and N. forgetiana (Fig. 3D) but did not overlap with rastroensis, unlike what was observed for AVCL (Fig. 3C, D; Fig. S3). Population AP02 showed lower variability in floral shape than AP01. The greatest degree of overlap was in AVCL morphospace with populations Forg04 and AP01 (Fig. 3C), whereas, for PV (Fig. 3C), AP02 showed the most forms found in N. forgetiana and AP01 (Fig. 3D; Fig. S3).

Floral size variation

Size variation estimated from centroid size (log) by ANOVAs differed considering all data for corolla limb size (F = 464.2; P < 0.01) and for corolla tube size (F = 322.9; P < 0.01). Moreover, all pairwise comparisons showed significant differences in corolla limb size and corolla tube size (P < 0.01). In general, N. alata presented the largest corolla limb and corolla tube size (Fig. 4A, B). The two atypical populations presented intermediate size between N. alata and N. forgetiana for corolla limb and tube size (Fig. 4A). In rastroensis, corolla size was smaller than in N. forgetiana, mainly due size the corolla tube (Fig. 4B).

Allometric effects

The MANCOVAs indicated that all groups had significant ontogenetic allometry (P < 0.01 in all cases). Size predicted 30% (R² = 0.299; F = 122.16; P < 0.01) and 8% (R² = 0.0776; F = 37.86; P < 0.01) of shape variation within groups for corolla limb shape and corolla tube shape, respectively. After considering size as a covariate, the shape variation between groups was 6% (R² = 0.0599; F = 8.14; P < 0.01) and 7% (R²² = 0.0746; F = 12.12; P < 0.01) for AVCL and PV, respectively. The groups presented differences, despite weak, in allometric trajectories for AVCL (1%; R² = 0.0127; F = 1.73; P = 0.014) and PV (3%; R² = 0.0275; F = 4.48; P < 0.01). For AVCL, the allometric trajectory of N. alata was different from the other groups (Fig. 5A). For PV, there was geographical variation in floral allometry related to N. forgetiana and rastroensis (Fig. 5B). The regression scores, considering all samples, showed that larger sizes presented stellate petals and longer and narrower corolla tubes; smaller sizes presented pentagonal corolla limb outlines and shorter and wider corolla tubes. There were two distinct groups in rastroensis with plants from all three populations (Supporting Information, Fig. S4).

Ecological suitability across time

The variable BIO6 (min temperature of coldest month) contributed most to the models of both species followed by BIO12 (annual precipitation) for N. forgetiana and BIO7 (temperature annual range) for N. alata (Supporting Information, Fig. S5). AUC values were > 0.9, and the models can be considered good. In general, the predictive distribution maps recovered the current records for both species (Fig. 6), also considering our observations during field expeditions. Projections for the past suitable habitat revealed that the current distribution is remarkably similar for Holocene suitability projections for the two species. The LIG period showed a broader suitable range of habitat for the two species, whereas the LGM was the period that most impacted the suitable habitat for both species, presenting a smaller distribution. According to the projections, N. forgetiana could may have had a fragmented area during the LGM (Fig. 6). The two species showed similar levels of niche overlap through modelled times, and Schoener’s D values ranged from 0.42 during the LGM to 0.55 during the LIG. The values of Schoener’s D obtained are outside the 99.9% intervals from an identity test and thus are statistically significant. The LGM period exhibited the lowest overlap between species (Fig. 6).

Association between floral shape/size and climatic variables

For N. forgetiana and rastroensis, analysis of the linear association between climatic variables and the AVCL form revealed a strong covariance between the two variables in the first PLS axis (R = 0.59, P < 0.01), as also seen in the PV view (R = 0.26, P < 0.01). The PLS graph (Fig. 7A, B) revealed that the shapes vary according to the climatic variables, with elevation being the most important (Supporting Information, Fig. S6). The relationship between size and climatic variables was also significant for AVCL view (R = 0.49, P < 0.01) and PV view (R = 0.74, P < 0.01), in which it is noted that the size varies mainly with elevation (Fig. 7C, D; Fig. S6).

DISCUSSION

This study aimed to characterize the floral shape and size variation at N. forgetiana population levels throughout its natural habitat to shed light on species evolution. Our results showed high variability in floral shape and size across the geographical distribution of N. forgetiana (Figs 1, 3, 4). The variation in floral shape and size was related to the putative new species, rastroensis, not yet formally described. Flowers of rastroensis have smaller corolla limbs with a pentagonal shape with relatively wider and shorter corolla tubes than those of canonical individuals of N. forgetiana. The floral form of N. forgetiana and rastroensis was associated with climatic variables and elevation. We also found increased variability in floral shape and size related to two mixed populations based on corolla limb colour. We suggest that they are probably a result of hybridization between N. forgetiana and N. alata (Figs 3, 4) in their contact area. We also showed that the species share a climatic envelope, at least in part, and probably had their suitable occurrence area impacted by the Pleistocene climatic oscillation (Fig. 6; Supporting Information, Fig. S5). These results shed light on a more complex evolutionary history for N. forgetiana in the southernmost Brazilian Atlantic forest. We discuss the findings of this study in the section below.

Population-level variability in floral shape and size between N.forgetiana and rastroensis

We found significant differences in corolla limb shape and corolla tube shape and size among the geographical area where N. forgetiana and rastroensis grow (Figs 3A, B, 4). Individual-based morphometric analyses are essential to disentangle morphological variability of plant species and contribute to recognize taxonomic entities (Gardere et al., 2019; Aristizábal, Figueiredo & André, 2020). These studies can shed light on interactions between intraspecific floral variation and differences in environmental conditions and pollinators (García et al., 2020), which drive plant diversification.

Although we observed overlapping floral limb shape among N. forgetiana and rastroensis, in general, flowers of rastroensis had a more pentagonal shape (Fig. 3A). Corolla limb shape in section Alatae varies with main pollinator groups (hawkmoth- vs. hummingbird-pollinated species) and is less variable within pollinator groups (Kaczorowski et al., 2012). Studies investigating the evolution of floral shape in allopolyploids of Nicotiana (McCarthy et al., 2016, 2019) revealed that corolla limb shape is a more labile trait than floral tube width. Floral limb shape also can be variable in other taxa of Solanaceae (e.g. N. glutinosa L.: Goodspeed, 1954; Petunia Juss.: Turchetto et al., 2014; Teixeira et al., 2020; Giudicelli et al., 2019) and may be impacted by different pollinators visiting the flowers (Venail, Dell’Olivo & Kuhlemeier, 2010). Floral limb shape is an essential trait for attracting pollinators (Gómez et al., 2008), and the rapid change in this shape can facilitate a rapid response to changes in selective pressures.

Floral tube shape was divergent between N. forgetiana and rastroensis (Fig. 3B), with rastroensis having a shorter corolla tube with a relatively wider corolla tube opening than N. forgetiana. The width of the corolla tube aperture may be related to changes in throat shape associated with the two lower petals. This throat shape reflects changes in the angle of inclination of the flower in N. forgetiana and rastroensis. The floral morphology of N. forgetiana predicts its primary pollinator (hummingbird; Kaczorowski et al., 2005). Among other traits, flowers adapted to hummingbirds are pendulous (Grant & Grant, 1968), generally due to the flexibility of floral pedicel. The shape of the corolla tube in N. forgetiana probably causes the corolla limb to face down owing to the shape of corolla tube opening in a position where it is connected to two inferior petals. On the other hand, rastroensis has a corolla tube shape with a different corolla aperture angle, changing the presentation of the two lower petals, which could serve as an airstrip (Fig. 2B) for some kinds of pollinators, such as bees. An elegant study by García et al. (2020) provided evidence for pollinator-driven selection for floral tube shape in N. glauca Graham occurring in different environments. The authors found an association between the curvature of the corolla tube and the form of the bills of sunbirds (García et al., 2020). Selection for floral tube length and shape by pollinators has long been documented, since Darwin’s experiments (Darwin, 1862) and in many taxa with different pollination syndromes (e.g. Nilsson, 1988; Whittall & Hodges, 2007; Pauw, Stofberg & Waterman, 2009; Knapp, 2010). Such studies indicated coevolution between floral tube and pollinator guilds and pollinator behaviours, a significant evolutionary force in angiosperm radiations (Hu et al., 2008; Soltis, Folk & Soltis, 2019).

McCarthy et al. (2016) suggested that allopolyploids in Nicotiana tend to have shorter and broader corolla tubes, indicating a more generalist pollination system. The wide and short corolla tube of rastroensis highlight ideas already presented in previous studies on the importance of other pollinators (e.g. bees) in the evolutionary history of section Alatae (halictid bees collecting pollen and bumblebees collecting nectar from these Nicotiana spp.; Kaczorowski et al., 2005). The putative new species rastroensis has been described in the literature as hummingbird-pollinated because of its similar floral traits to N. forgetiana (Kaczorowski et al., 2005); it had a higher nectar concentration than expected for hummingbird-pollinated flowers (as high as 57% solids, Kaczorowski et al., 2005). A high nectar concentration is a relevant floral trait attracting bees (Bolten & Feinsinger, 1978). In addition, rastroensis emits a floral volatile not found in N. forgetiana but present in species from allopolyploid Nicotiana section Suaveolentes (methyl nicotinate; Raguso et al., 2006). Methyl nicotinate has been recorded as a scent component in moth-pollinated species, suggesting their potential role as a moth attractant in daytime emissions (Raguso et al., 2003; Van der Niet, Jürgens & Johnson, 2015). These findings shed light on the role of other floral traits, beyond flower colour, in attracting pollinators in rastroensis. Field experiments are needed to characterize the pollination biology of rastroensis. Moreover, studies of volatile compounds in pollen and petals at the population level will provide more accurate information on this subject.

Studies with species belonging to section Alatae have demonstrated that corolla tube, stamen and carpel length are correlated (Bissel & Diggle, 2010). Lee et al. (2008), through experimental interspecific crosses with species of section Alatae, observed that all members could intercross; however, pollination success is highest in fertilizing ovules from species with pistil lengths like their own. They found a positive correlation between tube growth rate and pollen donor pistil length, suggesting a hybridization barrier in species with different pistil sizes. Here we demonstrated, by analysing many individuals from different populations, that rastroensis has a shorter tube length than N. forgetiana (Fig. 4). This trait, at least in part, could serve as a prezygotic barrier between them. Moreover, the wider corolla tube found in rastroensis can be related to a more generalized pollination mechanism, as suggested by McCarthy et al. (2016, 2019) for allopolyploid Nicotiana spp.

We also found that size influenced the corolla limb shape more than the corolla tube shape, and rastroensis showed evolutionary changes in allometric trajectories of the floral tube compared to other species (Fig. 5). The distinction of rastroensis involved the colonization of a new region in allometric space, especially with the corolla tube occupying two different areas (Supporting Information, Fig. S4B). We found plants from different populations in these two groups. Changes in allometric spaces are associated with differences in pollination system (Hazle & Canne-Hilliker, 2005; Summers, Hartwick & Raguso, 2015; Strelin et al., 2017). Additionally, abiotic environmental factors also play a role in the evolution of allometric trajectories. For example, it could protect from UV-B radiation at high elevations that affects the fitness of plants (Ackermann & Weigend, 2006; Zhang, Yang & Duan, 2014) and animals (Wilson & Sánchez-Villagra, 2010; Wilson, 2013). Here we found a latitudinal and elevational gradient of floral shape variation associated with N. forgetiana and the putative new species rastroensis. These patterns may be related to changes in floral form and species diversification.

Variation in floral shape and size in N. alata and N. forgetiana – hybrids?

Hybridization has been a critical driver of diversification in Nicotiana, with many allopolyploid species of different ages and homoploid hybrid species (Chase et al., 2003; Kelly et al., 2010, 2013; Clarkson et al., 2010, 2017; McCarthy et al., 2015). Species of section Alatae readily hybridize in experimental crosses, but no hybrids have been found in natural populations (Goodspeed, 1954; Ippolito et al., 2004; Lee et al., 2008). Species in this section have been studied regarding many traits related to pollinator divergence and nectar composition (Ippolito et al., 2004; Kaczorowski et al., 2005), volatile compounds (Raguso et al., 2003; 2006) and floral shape (Kaczarowiski et al., 2012). Nicotiana alata and N. forgetiana display floral morphologies mainly associated with their pollination syndromes (hawkmoth and hummingbird, respectively; Kaczorowski et al., 2005). As demonstrated by Kaczorowski et al. (2012) and in our individuals-based analyses, the floral shape of these species is significantly different. We found that the two atypical populations displayed corolla colour ranges from white to purple–red (Supporting Information, Fig. S1), intermediate floral size (Fig. 4) and a range of floral shape (Fig. 3C, D). The floral shape variation found in atypical populations, mainly in AP01, resembled the phenotype mosaic reported for artificial F2 and F3 hybrid generations between N. alata and N. forgetiana (Bissel & Diggle, 2008, 2010). Geometric morphometrics is a helpful technique to detect such mosaics resulting from hybridization between species with porous genomes (Lexer et al., 2009; Teixeira et al., 2020).

Through experimental conditions, Ippolito et al. (2004) showed that, although the primary pollinators of N. alata and N. forgetiana forage preferentially in line with the floral syndrome, both pollinators (hawkmoths and hummingbirds) could visit flowers of both species. However, hybrid seed set mainly came from N. forgetiana, indicating that corolla tube length acts as a reproductive barrier (Ippolito et al., 2004). The authors also argued that floral size differences, not just floral shape, were significant for flower discrimination by hawkmoths, and the presence of hybrid individuals could facilitate introgression (Ippolito et al., 2004). Manduca sexta, the main pollinator of N. alata (Kaczorowski et al., 2005), is not a particularly specialized pollinator, as demonstrated by some diversity of pollen present on the moth (Alarcón, Davidowitz & Bronstein, 2008). Additionally, N. forgetiana is primarily pollinated by hummingbirds (Kaczorowski et al., 2005), and it has a purple–red corolla that seems to play a significant role in attraction of birds and deterrence of bees (Cronk & Ojeda, 2008). One hypothesis that arose in this scenario is the possibility of bees serving as a vector for introgression in N. forgetiana, as variation in colour, floral shape and volatile compounds has been reported in this species (Raguso et al., 2006; and results presented here). However, other analyses, including floral biology studies and genome-wide genetic variability, will be fundamental to investigate this question.

Manduca sexta could exploit the same plant species at both life-history stages (Alarcón et al., 2008), being both a herbivore and a pollinator. A trade-off between pollination and protection against the attack of larvae of M. sexta has been reported for N. attenuata Torr. ex S.Watson, a nocturnal hawkmoth-pollinated species (Kessler, Diezel & Baldwin, 2010). The authors showed that plants attacked by hawkmoth larvae change the amount of benzyl acetone at night and then open in the morning, being pollinated by hummingbirds (Kessler et al., 2010). Gómez et al. (2009b) also showed associations between floral shape variation and regional differences in herbivores in Erysimum ediohispanicum Polatschek, revealing a balance between mutualistic and antagonistic selection. During our fieldwork, we also observed N. alata flowers open during the day, and Raguso et al. (2006) found intraspecific differences in volatile compound emissions in N. alata from different geographical origins (Rio Pelotas collection point). Those flowers opening during the day are more likely to be visited by hummingbirds, which may facilitate hybridization between N. alata and N. forgetiana. This discussion sheds light on questions still to be investigated. (1) What is the role of generalist pollinators in the evolution of Nicotiana? (2) Could strategies to help escape from herbivores increase the rate of gene flow between pairs of species with different floral syndromes?

Climatic suitability through time

Through geographical distribution and molecular phylogenetic analyses (Goodspeed, 1954; Clarkson et al., 2004, 2017), the Andean region was suggested as the centre of origin for Nicotiana. Clarkson et al. (2004) hypothesized that, at least, one dispersal event to Brazil gave rise to section Alatae, a recently diversified lineage (c. 2.5–3.0 Mya; Lim et al., 2006; Clarkson et al., 2017). The phylogenetic relationships among the species from section Alatae remain uncertain (Clarkson et al., 2004, 2010, 2017; McCarthy et al., 2019). A sister group relationship of N. forgetiana and N. bonariensis Lehm. has been suggested based on nuclear and plastid markers (Clarkson et al., 2004, 2010). At the same time, N. alata appears more closely related to N. langsdorffii Weinmann based on nuclear markers (Clarkson et al., 2010) and plastid markers (Clarkson et al., 2017), or to N. mutabilis Stehmann & Semir based on concatenated datasets from plastid and nuclear markers (McCarthy et al., 2019), which is the most recently described species in this section (Stehmann, Semir & Ippolito, 2002). Pollinator shifts have been documented in section Alatae and could be a crucial driver for cladogenesis and maintaining the species boundaries (Kaczorowski et al., 2005; Lee et al., 2008). Additionally, Lim et al. (2006) suggested that chromosome morphology could serve as a post-zygotic barrier between N. forgetiana and N. bonariensis (a white-flowered, moth-pollinated species) because the latter has a substantially different karyotype structure, whereas N. alata and N. forgetiana have a similar karyotype structure. However, the impact of climatic changes during the Pleistocene has not yet been investigated in this lineage of Nicotiana (the n = 9 clade). Here we showed that the most critical variables for N. alata and N. forgetiana were annual precipitation (BIO12) and temperature (BIO07) and the temperature in the coldest month (BIO06) and that the two species had some niche overlap mainly during the interglacial periods. Suitable habitat for N. forgetiana during the LGM was drastically reduced (Fig. 6) and potentially suffered from habitat fragmentation. Nicotiana forgetiana and rastroensis are found only in open areas from the southernmost Brazilian Atlantic Forest biome. Barros et al. (2015) investigated the main bioclimatic variables associated with diversity patterns in genera that diversified in the same region. They found a high correlation between the number of species per area and annual temperature, and some groups (e.g. Petunia) showed high haplotype diversity related to annual mean temperature. Additionally, Barros et al. (2020) demonstrated that temperature, precipitation and elevation changed with haplotype diversity in two closely related plant genera of Solanaceae occupying grasslands in South America (Calibrachoa La Llave & Lex. and Petunia). The changes in climate during the LGM, making the climate colder and drier (Behling et al., 2002), could therefore have limited N. forgetiana to small, fragmented areas.

During the Pleistocene, oscillations in climatic conditions may have influenced the population dynamic in N. alata and N. forgetiana and their pollinators. Clarkson et al. (2017) estimated the age of allopolyploid Nicotiana spp., and most originated during the Quaternary (c. 0.4, 0.6, 0.7, 1.4, 4.3 and 6.0 Mya) marked by major climatic shifts. Climatic oscillation promotes the appearance of many hybrid lineages (Hewitt, 2011) due to changes in the species distribution with the formation of refuges and colonization routes. According to our projections of suitable habitats for N. alata and N. forgetiana, these species could have experienced changes in their geographical distributions over time. They could have been restricted to potential refuge areas, e.g. during the LGM, and colonization routes in the most recent time, in the Holocene (Fig. 6). Differences in environmental conditions through the geographical area could influence flower shape. These differences could be associated with locale differences in the selection pressures affected by pollinators and ecological conditions (García et al., 2020). Here, we found differential floral shape between N. forgetiana (most north-western distribution) and rastroensis (most north-eastern distribution), forming a gradient of variation across geography and elevation. Moreover, N. forgetiana could have had a fragmented area during the LGM. These findings support the idea that the selection pressure by different pollination systems and climatic conditions influenced the formation of these two morphotypes.

The geographical area of N. forgetiana is one of high plant endemism (Iganci et al., 2011; Carnaval et al., 2014) and is identified as southern Atlantic Forest climatically stable area over the last 30 kyr (Costa et al., 2017), presenting high levels of genetic diversity (Raposo do Amara et al., 2013). This region is a mosaic of open grassland fields and Araucaria Juss. forest (Overbeck et al., 2007), sharing floristic similarities with the Andean region (Safford, 2007). Evidence from pollen records (Behling, 2002) have shown a predominance of grasslands where Araucaria forests occur today. There is evidence that grasslands extended 750 km from southern to south-eastern Brazil during glacial periods (Behling, 2002). Whereas the Holocene period (starting between 11.0 and 9.7 kyr BP) was marked by a warmer and moist climate, allowing migration of the flora surviving in glacial refugia (Leonhardt & Lorscheitter, 2009). Thus, the flora adapted to colder and drier conditions experienced an expansion during the LGM (Behling, 2002). Based on this, we expect that species from open fields, such as Nicotiana spp., experienced an expansion of their area of suitability during the LGM. However, we found contrasting patterns for Nicotiana spp., with a reduction of suitability in the LGM, especially for N. forgetiana. This pattern agrees with what was observed for Araucaria angustifolia (Bertol.) Kuntze surviving in glacial refuge areas and post-glacial expansion towards the south (Lauterjung et al., 2018; Stefenon et al., 2019). These results showed that species growing in open fields in southern South American environments may have had different responses to climate change, according to their adaptability to climatic conditions. Studies with plants growing in open fields in South America have shown evidence that even lineages within a species responded differently to climatic fluctuations during the Pleistocene, shwoing different population dynamics reflecting the pattern shaping genetic variation (Turchetto et al., 2014; Giudicelli et al., 2019).

The dynamics of contraction and expansion have been associated with speciation events in species growing in open fields in this region (Lorenz-Lemke et al., 2010), from the fragmentation of ancestral geographical areas and isolation of populations. Here we suggested that the dynamics of the interactions with pollinators and the impact of climate oscillations during the Pleistocene (in plant and pollinator responses) could lead to lineages and hybridization between N. forgetiana and N. alata. Molecular analyses at population levels with these species will be necessary to support these findings.

CONCLUSIONS AND PERSPECTIVES

Population-level studies are essential to understand the microevolutionary processes underlying the origin and diversification of plant species complexes from the Neotropical region (Pinheiro et al., 2018). This is especially important for diploid Nicotiana spp. (Knapp, 2020). Our morphometry analysis showed that N. forgetiana and rastroensis could be readily distinguished by the shape and size of the corolla tube. These findings, associated with their different geographical distribution and elevational gradient, corroborates the acceptance of the populations of Serra do Rio do Rastro and surroundings as representing a distinct species to be formally described. Ecological factors and selective pressures for pollinators could be related to putative new species of Nicotiana from highland grasslands in southern Brazil. This study was able to identify mixed phenotypes suggesting natural hybridization between N. alata and N. forgetiana. Population genomic analyses will be beneficial for investigating the levels and outcomes of hybridization in section Alatae.

SUPPORTING INFORMATION

Additional Supporting Information may be found in the online version of this article at the publisher’s web-site:

Table S1. Sampled flowers of Nicotiana species used in this study for geometric morphometric analyses. N = sample size for adaxial view of the corolla limb (AVCL), n = sample size for profile view (PV).

Table S2. Occurrence data of the species used for ecological niche modelling (ENM) obtained from collection points by the authors and the SpeciesLink database (http://www.splink.org.br/).

Table S3. Principal component analysis (PCA) for Nicotiana forgetiana and rastroensis for AVCL view.

Table S4. Principal component analysis (PCA) for Nicotiana forgetiana and rastroensis for PV view.

Table S5. Principal component analysis (PCA) of AVCL view for all populations collected in the field, N. forgetiana, N. alata, rastroensis and atypical populations (AP).

Table S6. Principal component analysis (PCA) of PV view for all populations collected in the field, N. forgetiana, N. alata, rastroensis and atypical populations (AP).

Figure S1. Phenotype mosaic found in the two atypical populations AP01 and AP02.

Figure S2. Principal component analyses (PCA). Floral shape variation in N. forgetiana and rastroensis species. A, adaxial view of the corolla limb AVCL. B, profile view PV. Floral shape variation among all populations collected in the field. C, adaxial view of the corolla limb – AVCL. D, profile view – PV.

Figure S3. Canonical variance analysis (CVA) showing the position of each population in the morphospace, for AVCL (A) and PV (B). Floral shape variation in N. forgetiana and rastroensis species. A, adaxial view of the corolla limb AVCL. B, profile view PV. Floral shape variation among all populations collected in the field. C, adaxial view of the corolla limb – AVCL. D, profile view – PV.

Figure S4. Multivariate analysis of covariance (MANCOVA) using shape as a response variable, size as a covariate and species as predictors. A, adaxial view of the corolla limb – AVCL. B, profile view – PV. Initial shape and resulting shape can be visualized in the grey and black schematic shapes on the graphic.

Figure S5. Percentage contribution of environmental variables to the ecological niche modelling (ENM) for each species. BIO6 = min temperature of coldest month, BIO7 = temperature annual range, BIO9 = mean temperature of driest quarter, BIO10 = mean temperature of warmest quarter, BIO12 = annual precipitation, BIO13 = precipitation of wettest month, BIO14 = precipitation of driest month, BIO15 = precipitation seasonality (coefficient of variation), BIO16 = precipitation of wettest quarter and BIO19 = precipitation of coldest quarter.

Figure S6. Percentages of the contribution of environmental variables to the analysis of partial least-squares (PLS). The x-axis gives climate variables BIO1–BIO19 and elevation. The y-axis gives shape differences.

ACKNOWLEDGEMNTS

We thank Carolina K. Schnitzler and Analu C. Souza for assistance with the field collection.

FINANCIAL SUPPORT

This work was supported by the Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq - 427575/2018-4 and 306086/2017-4), and the Universidade Federal do Rio Grande do Sul (UFRGS). M.A. was supported by PhD fellowship grants of Fundação de Amparo à Pesquisa do Estado de Minas Gerais and Coordenação de Aperfeiçoamento de Pessoal de Nível Superior - Brasil (CAPES) - Finance code 001.

CONFLICT OF INTEREST

The authors declare that they have no conflicts of interest.

REFERENCES

Author notes