Skip to main content
Download PDF
  • Research article
  • Open access
  • Published:

Exploring molecular evolution of Rubisco in C3 and CAM Orchidaceae and Bromeliaceae

BMC Evolutionary Biology volume 20, Article number: 11 (2020) Cite this article

Abstract

Background

The CO2-concentrating mechanism associated to Crassulacean acid metabolism (CAM) alters the catalytic context for Rubisco by increasing CO2 availability and provides an advantage in particular ecological conditions. We hypothesized about the existence of molecular changes linked to these particular adaptations in CAM Rubisco. We investigated molecular evolution of the Rubisco large (L-) subunit in 78 orchids and 144 bromeliads with C3 and CAM photosynthetic pathways. The sequence analyses were complemented with measurements of Rubisco kinetics in some species with contrasting photosynthetic mechanism and differing in the L-subunit sequence.

Results

We identified potential positively selected sites and residues with signatures of co-adaptation. The implementation of a decision tree model related Rubisco specific variable sites to the leaf carbon isotopic composition of the species. Differences in the Rubisco catalytic traits found among C3 orchids and between strong CAM and C3 bromeliads suggested Rubisco had evolved in response to differing CO2 concentration.

Conclusions

The results revealed that the variability in the Rubisco L-subunit sequence in orchids and bromeliads is composed of coevolving sites under potential positive adaptive signal. The sequence variability was related to δ13C in orchids and bromeliads, however it could not be linked to the variability found in the kinetic properties of the studied species.

Background

Crassulacean acid metabolism (CAM) is one of the three mechanisms found in vascular plants for the assimilation of atmospheric CO2. The CAM pathway is characterized by the temporal separation of carbon fixation: CO2 is initially fixed by phosphoenolpyruvate carboxylase at night [1,2,3]. The resulting organic acids are stored in the vacuole and, during the day, decarboxylation of these compounds provides CO2 at high concentrations for assimilation by Rubisco [1]. This mechanism makes it possible for CAM plants to close their stomata during the day when the evaporative demand is higher, so the water cost of CO2 gain is significantly reduced in CAM plants [2, 3]. This fact, along with other anatomical features that minimize water loss, increases the water use efficiency (WUE) of CAM plants several fold compared to C3 and C4 plants [4, 5]. The selective advantage of high WUE likely accounts for the extensive diversification and speciation among CAM plants, particularly in water-limited habitats [2, 6]. Indeed, CAM has been reported for 343 genera in 34 families and ca. 7% of all vascular plant species are estimated to exhibit CAM [7,8,9].

The CO2 is considered the central driving force for the earliest evolution of CAM [10, 11]. Actually, it is thought that CAM photosynthesis appeared as the result of adaptive selection related to the decline in the atmospheric CO2 concentration and progressive aridification, in a similar manner to the Miocene expansion of grasses with C4 mechanism [12,13,14]. In essence, CAM constitutes a CO2-concentrating mechanism originated through daytime malate remobilization from the vacuole and its decarboxylation with regeneration of CO2. This mechanism leads to CO2 partial pressure ranging between 0.08 and 2.50% in the leaf air spaces during CAM phase III, i.e., when Rubisco is active [15]. With up to 60-fold increase of CO2 level in the photosynthesizing organs as compared to the atmospheric CO2 partial pressure, this is the strongest known increase of internal CO2 partial pressures of CO2-concentration mechanisms [16]. This striking increase in CO2 concentration might directly impact the functioning of the enzymes in the Calvin cycle, in particular the C3 carboxylating enzyme Rubisco, by substrate-saturating its carboxylase activity.

Rubisco has evolved to optimize catalysis according to the availability of CO2 in the vicinity of its catalytic sites [17,18,19]. In principle, optimization of Rubisco to the prevailing environment has to inevitably deal with the trade-off between Rubisco affinity for CO2 (1/Kc) and the enzyme’s turnover rate [17, 18, 20, 21]. Hence, C3 species with low CO2 concentration at the site of carboxylation, such as those from dry and hot environments, tend to present Rubiscos with higher affinity for CO2, albeit with lower maximum rate of carboxylation (kcatc) [19, 22]. On the contrary, in C4 plants, with 6–10-fold increase of CO2 concentration in the bundle-sheath cells compared to the atmosphere, Rubisco has specialized towards increased kcatc [23,24,25,26,27].

In the recent years, signatures of positive selection acting on particular amino acid residues of Rubisco have been found using phylogenetic analysis of different taxonomic groups, confirming variation trends in the evolution of the Rubisco kinetics to changing intracellular concentrations of CO2 in C3 and C4 plants [22, 28,29,30,31,32,33,34]. In contrast to C3 and C4 species, the molecular evolution of Rubisco in CAM plants has been poorly investigated. While Kapralov and Filatov (2007) [28] included representatives of CAM pathway in their study, investigation of selection associated with CAM was not among their research objectives.

The exploration of the natural variation of Rubisco catalytic traits by means of molecular and biochemical approaches is far from being complete and is considered a promising way 1) to increase our understanding on how the environmental conditions shape Rubisco evolutionary fine-tuning, and 2) to find Rubisco variants with increased efficiency to use in existing engineering programs aiming at improving Rubisco performance [35, 36]. Previous results, suggesting that high selection pressure on Rubisco has particularly occurred in species from extreme environments and/or possessing innovative adaptations such as carbon concentrating mechanisms [19, 27, 30, 31, 34, 37], make CAM plants a prime subject for understanding mechanisms of evolution in Rubisco.

We undertook the comparative analysis of Rubisco evolution in closely related species from the Orchidaceae and Bromeliaceae families possessing C3 and CAM pathways. These two Neotropical plant families represent an outstanding example of adaptive radiation in plants with a striking ecological versatility, occupying habitats extremely different in the ecophysiological demands, among which epiphytic life forms predominate [13, 38,39,40]. Orchids and bromeliads contain approximately half of the total CAM plant species, and evidence of selection for weak and strong models of CAM has been reported for both families [13, 41,42,43,44]. Importantly, CAM pathway evolved several times independently within the Orchidaceae and Bromeliaceae families [3, 40, 45, 46], making them an ideal model to compare Rubisco evolution between CAM and C3 related species. Our hypothesis was that a large variability in the L-subunit exists in bromeliads and orchids, and that part of this molecular variability was positively selected to improve the catalytic performance of Rubisco according to the specific physiology of CAM and C3 species. To test this hypothesis, we characterized the chloroplast rbcL gene to explore the variability of the Rubisco L-subunit within these families and to search for specific amino acid replacements associated with CAM. Thereafter, Rubisco catalytic parameters were measured in representative species to infer the biochemical impact of amino acid replacements within the Rubisco L-subunit. Intra-molecular coevolution analysis was also conducted to further understand the importance and correlation of the Rubisco L-subunit sites under selection with the functionality of Rubisco. Finally, a decision tree (DT) model was implemented to find correlations between Rubisco L-subunit amino acid replacements and essential variables of the species, including carbon isotopic discrimination and habitat preference.

Results

Variation of leaf traits and their correlation with the photosynthetic mechanism

As for the purposes of the present study, the classification of the photosynthetic mechanism of the studied species and varieties into different CAM levels is required. Several approaches have been used to determine the presence of CAM, including some leaf morphological traits, such as the leaf thickness, the leaf mass area (LMA) or the FW/DW ratio, indicative of succulence, the leaf carbon isotopic composition (δ13C), the time-course of leaf gas-exchange and leaf titratable acidity, and the activity of specific enzymes [3, 5, 13, 42,43,44,45, 47,48,49]. The accurate categorization of species into CAM or C3 requires the combination of several of these methods. However, in the present study, due to the number of selected species and to the number and variety of analyses, we discriminated the photosynthetic mechanism on the sole basis of δ13C values, and then correlated this with leaf morphological traits in species for which the live specimens were available. The same approach has been applied before in orchids [6, 50,51,52] and bromeliads [13]. We are aware that whole-tissue δ13C alone does not provide a precise indication of the contributions of dark and light CO2 fixation to total carbon gain. This is why we adopted a conservative strategy, including the group of weak CAM as a ‘buffer’ between C3 and strong CAM.

Among the 78 orchids, the leaf δ13C values ranged from − 36.6 to − 12.4 ‰, and from − 33.3 to − 9.4‰ among the 144 bromeliads (Table 1 and Additional file 1: Table S1). The frequency diagram of species and varieties displayed a bimodal distribution with peaks around − 25 and − 15 ‰, in agreement with previous surveys [13, 41, 43, 44, 50, 54] (Fig. 1). Attending to the δ13C values, 63 orchids were classified as C3 (δ13C < − 22.9 ‰) and 15 as CAM (δ13C > − 18 ‰). Among the bromeliads, the relative presence of CAM pathway was more evident, with 74 species classified as C3, 13 as “weak CAM” and 57 as CAM (Table 1, Additional file 1: Table S1 and Fig. 1).

Table 1 List of sequenced Orchidaceae and Bromeliaceae species
Full size table
Fig. 1
figure 1

Frequency diagram according to the leaf carbon isotope composition (δ13C, ‰) of the 78 orchids and 144 bromeliads studied (see Table 1 and Additional file 1: Table S1)

Full size image

FW/DW ranged between 2.7 (Puya humilis) and 15 (Pleurothalis nuda), the leaf thickness between 0.2 (Elleantus furfuraceus and Tillandsia biflora) and 3.8 mm (Puya laxa), and the LMA between 28.6 (Lycaste cruenta) and 457.3 g m− 2 (Puya stenothyrsa) (Table 1 and Additional file 1: Table S1). Among orchids and bromeliads, strong CAM plants presented significantly (p < 0.05) higher LMA than C3 plants (Additional file 2: Table S2). The leaf thickness was also significantly (p < 0.05) higher in strong CAM compared to C3 plants. Higher LMA and leaf thickness in CAM than C3 plants have been previously described in taxonomically diverse groups with CAM species [5, 42, 43, 50]. There were no significant differences in FW/DW between strong CAM and C3 plants, in agreement with previous surveys in orchids [50], but in disagreement with other surveys in phylogenetically diverse groups [11].

The δ13C values tend to be less negative with increasing leaf thickness and LMA when considering bromeliads and orchids together and separately (Fig. 2). As reported in previous findings [41, 50], most species with leaves over 1 mm thick and 100 g m− 2 had δ13C values less negative than − 18 ‰, indicative of strong CAM. A notable exception was the bromeliad Puya sanctae-crucis, with a leaf thickness of 2.6 mm and C3-type δ13C values (− 25.9 ‰). Exceptions of the correlation between leaf thickness and δ13C values have been related to the relative contribution of hydrenchyma to total leaf thickness [10, 43, 50].

Fig. 2
figure 2

Relationship between the leaf carbon isotope composition (δ13C) and the leaf thickness and the leaf mass per area (LMA) for the orchids and bromeliads listed in Table 1 and Additional file 1: Table S1. In (a) and (b) all orchids and bromeliads are plotted together, in (c) and (d) only orchids, and in (e) and (f) only bromeliads. Filled black symbols correspond to strong CAM species of orchids (▲) and bromeliads (●); symbols in grey correspond to weak CAM species of orchids () and bromeliads (); open symbols correspond to C3 species of orchids (∆) and bromeliads (○). Values are means (n = 4). Regression coefficients between parameters were performed with R [55]

Full size image

Variability in the Rubisco L-subunit sequence of orchids and bromeliads: analyses of positive selection and intra-molecular coevolution

The L-subunit sequence of bromeliads was 480 amino acids long, while orchid sequences presented diverse length: Myoxanthus exasperatus, Pleurothallis chloroleuca, P. nuda and P. cardiothallis were 483 amino acids long; Acianthera pubescens and Epidendrum paniculatum were 479 amino acids long; and the rest of orchid sequences consisted of 481 amino acids. The number of variable amino acid sites was 73 within the orchids and 38 within the bromeliads (Additional file 7: Excel S1 and Additional file 8: Excel S2).

The rbcL topologies were constructed for the 78 orchids and 130 bromeliads separately (Figs. 3 and 4). The topology was largely congruent with previously obtained phylogenies [38, 56] and accepted subfamilies divisions.

Fig. 3
figure 3

Orchidaceae topology based on rbcL sequences (species listed in Table 1 and Additional file 1: Table S1). In blue: CAM species

Full size image
Fig. 4
figure 4

Bromeliaceae topologies based on rbcL sequences (species listed in Table 1 and Additional file 1: Table S1). In blue: CAM species

Full size image

A total of 13 sites were identified under positive selection within orchids, while within bromeliads signatures of positive selection were identified in 10 sites (Table 2). Common sites under positive selection for both orchids and bromeliads were 142, 225, 251, 449, 468 and 478, while the sites 89, 265, 461, 470, 477, 479 and 481 were exclusive for orchids, and the sites 28, 91, 255 and 270 were found only in bromeliads. However, LRTs (Table 2) indicated that the models assuming positive selection on all branches were not significantly better than the models without positive selection (p = 1). For this reason, we will refer to these positively selected sites along the manuscript as candidate sites under positive selection. Moreover, no single codon was identified as evolving under positive selection in branches or clades leading to the CAM species.

Table 2 Rubisco L-subunit sites candidate to positive selection
Full size table

Within orchids, 23 amino acidic sites were identified under co-evolution in the L-subunit of Rubisco, distributed in 11 coevolution groups (Additional file 3: Table S3). The residue 475 was the one that appeared as coevolving in more groups, with a total of seven interactions. The residue 466 appeared as coevolving into six groups, and residues 26, 28, 439, 443, 449, 461, 468, 470, 477, 478 and 479 presented four interactions. Other residues with co-evolving interactions were 33, 265, 279, 328, 334, 340, 341, 353, 359 and 447.

Within bromeliads, 20 amino acidic sites were identified under co-evolution, distributed in 2 coevolution groups, being the residues 449 and 478 the ones with two interactions and the rest of co-evolving sites appeared in only one group (Additional file 3: Table S3).

Decision tree analysis applied to the observed variability in the Rubisco L-subunit

In the orchids dataset, the DT model denoted a link between the external variables (δ13C and habitat preference) and the Rubisco L-subunit variable sites 89, 224, 225, 375 (Table 3, Additional file 5: Figure S1). The xerrors calculated for each variable site were 0.96, 0.74, 0.90 and 0.51, respectively. According to the xerror, the sites that were best explained by the external variables were 375 and 224 followed by 225 and 89.

Table 3 Variable sites in the Rubisco L-subunit resolved with the DT model for bromeliads and orchids
Full size table

In the case of bromeliads, the DT pointed to a link between the variable sites 91, 142, 219, 225, 255, 407, 464 and 468, and the external variables (δ13C and habitat preference). The xerrors of 0.83, 0.96, 0.93, 0.89, 0.82, 0.66, 0.97 and 0.83, respectively indicated the site 407 as the best explained by the external variables, followed by 255, 91, 468, 225, 219, 142 and 464 (Table 3, Additional file 6: Figure S2).

For both orchids and bromeliads, the external variable δ13C was the one that better correlated with all variable sites (Table 3, Additional files 5: Figure S1 and Additional file 6: Figure S2).

Rubisco kinetics in orchids and bromeliads: trade-offs and correlation with leaf traits and L-subunit sequence

Among the 5 bromeliads examined, the Rubisco Michaelis-Menten constant affinity for CO2 (Kc) varied between 9.6 μM (T. biflora) and 27.4 μM (A. nudicaulis). Among the 6 orchids examined, Kc varied between 12.1 μM (S. macrantha) and 24.2 μM (M. cucullata) (Table 4). The range of variation of the maximum carboxylase rate (kcatc) was similar to that of Kc. Non-significant differences were found in the catalytic carboxylase efficiency (kcatc/Kc) among the selected species. Differences in the relative abundance of Rubisco over leaf total soluble protein ([Rubisco]/[TSP]) were observed among bromeliads (Table 4), with the strong CAM A. nudicaulis and T. bermejoensis presenting the lowest values.

Table 4 Rubisco kinetic parameters at 25 °C and candidate positively selected sites in the Rubisco L-subunit
Full size table

The two studied strong CAM bromeliads (T. bermejoensis and A. nudicaulis) averaged higher Kc values (25.3 ± 2.0 μM) compared to the three C3 bromeliads (N. innocentii lineatum, T. biflora and T. multicaulis, with 12.3 ± 0.8 μM) (Table 4). Non-significant differences were observed in kcatc/Kc between strong CAM and C3 bromeliads, because kcatc varied in the same proportion (6.0 ± 1.3 and 3.0 ± 0.2 s− 1 for strong CAM and C3 bromeliads, respectively).

Apparently, no single amino acid replacement in sites under positive selection (Table 2) or resolved with DT (Table 3) was correlated to the differences observed in Kc among the studied orchids and bromeliads (Table 4). However, in orchids, the species with low values for Kc and kcatc, S. macrantha and E. furfuraceus, presented the potentially positive and predicted replacements 89 V, 468 N, 470E and 478 L, while the species with the highest values for Kc and kcatc, L. amoena and M. cucullata, presented 89A, 468D, 470D and 478E.

Correlation coefficients between catalytic parameters, amount of Rubisco, leaf traits and δ13C were calculated for all the species using PIC analyses (Table 5). The trade-off between kcatc and affinity for CO2 (1/Kc) was observed in bromeliads and orchids at P < 0.001. Kc and kcatc correlated significantly with [Rubisco]/[TSP], LMA, leaf thickness and leaf δ13C. Finally, [Rubisco]/[TSP] was inversely correlated with LMA, leaf thickness and leaf δ13C.

Table 5 Phylogenetically Independent Contrasts (PIC)
Full size table

Discussion

Rubisco L-subunit amino acid replacements associated with CAM species

Because water-conserving and CO2-concentrating mechanism (CCM) in CAM plants provide an advantage in particular ecological conditions [57, 58], we hypothesized that positive selection of molecular changes promoting such physiological traits may have driven the evolution of CAM Rubisco, in a similar manner of positive adaptive signal associated with C4 Rubisco [29,30,31, 34].

Previous studies reported residues under positive selection in Rubisco L-subunit in different groups of plants [28,29,30, 33, 34, 59,60,61,62] and revealed that amino acid co-evolution is common in Rubisco of land plants [62, 63]. Kapralov and Filatov (2007) [28] reported a number of amino acid sites under positive selection in families sharing C3 and CAM species. In the present study, candidate positively selected sites 89, 225, 251 and 265 (in Orchidaceae) and 142, 225, 251 and 255 (in Bromeliaceae) coincided with those reported in their study (Table 2). Sites 142, 449, 461, 468, 470, 477, 478, 479, 481 in orchids, and 28, 91, 270, 449, 468 and 478 in bromeliads are reported in our study but not by Kapralov and Filatov (2007) [28], because of different sample design. Kapralov and Filatov (2007) [28] used different orchids and bromeliads species and fewer of them compared to the present study. Furthermore, all sites under putative positive selection found in this study have been reported in [28] if all phylogenetic groups sampled outside of bromeliads and orchids are taken into account too, confirming widespread convergent evolution within Rubisco among flowering plants [28].

The candidate positive sites 265, 449, 461, 468, 470, 477, 478 and 479 in orchids, and 28, 91, 142, 225, 251, 255, 270, 449 and 468 in bromeliads were identified as coevolving with other amino acid sites (Additional file 3: Table S3). This fact relates positive selection and coevolution within sites located in functionally important interfaces. This is the case of sites 91, 142, 225, 461, 468 and 470, involved in intra-dimer and inter-dimer interactions, interactions with the small subunits and Rubisco Activase, or near to active site (Additional file 4: Table S4).

The candidate positive sites 89 and 225 in orchids, and 91, 142, 225, 255 and 468 in bromeliads were also resolved with a DT (Additional file 4: Table S4). The DT related these sites with the isotopic discrimination, being the species leaf δ13C value the most important external variable (Table 3). The apparent discrepancy between the results of branch-site tests of positive selection (no signs of positive selection associated to CAM) and the DT model (amino acid replacements related to δ13C) may be attributed to methodological differences. While positive selection analyses were constrained by the binary classification of species into C3 or CAM (using labels # in the tree file for the CAM species), the DT model gains from less rigidity as considering numerical values of leaf δ13C (all the CAM values of δ13C > − 18 ‰ and the C3 values < − 22.9 ‰). On the basis of the huge variability in the concentration of CO2 at the sites of Rubisco among CAM plants due to the CAM mechanism it seems more appropriate the DT model approach [2, 4]. CAM plants are reported to adjust the expression of different phases in CAM pathway to boost the internal supply of CO2 to Rubisco [64, 65]. Recent evidence showed that adaptive forces may act on other regulation points of CAM metabolism, like the enzyme PEPC [66]. It is also important to remark that the δ13C values reported in the present study have been obtained from plants grown under different conditions, including greenhouse-grown plants and field data from literature. While this fact was unavoidable to warrant the feasibility of this study, we cannot discard variation in leaf δ13C values due to environmental variation.

Results suggest the existence of differences in the Rubisco kinetics among C3 orchids and between C3 and strong CAM bromeliads, but the molecular basis of these differences remains to be elucidated

Differences in Kc and kcatc at 25 °C were observed among C3 orchids but not among C3 bromeliads (Table 4). Of the three orchids with higher values of Kc and kcatc (M. exasperatus, L. amoena and M. cucullata), M. exasperatus and L. amoena exhibited large LMA (Table 1 and Additional file 1: Table S1). Higher LMA has been linked to increased mesophyll resistance to CO2 transfer and, therefore, low CO2 availability at the site of carboxylation [67]. This finding apparently contradicts previous reports suggesting that in C3 species low CO2 availability promotes Rubisco evolution towards higher affinity for CO2 (i.e., low Kc) at the expenses of low kcatc [19, 33, 37, 68]. The comparison of the particular microclimate where these species evolved may help in understanding the evolutionary causes of this variability among C3 orchids. Unfortunately, the lack of success in extracting sufficient active Rubisco in strong CAM orchids precluded the comparison of kinetics between orchids with different photosynthetic mechanism. Future attempts should consider the low concentration of Rubisco present in the leaves of these species, even those with C3 mechanism (Table 4).

In bromeliads, the two strong CAM species presented higher kcatc compared to the C3 species (Table 4), although the ratio kcatc/Kc was similar between the two groups. Higher values of kcatc at the expense of decreased affinity for CO2 (i.e., higher Kc) have been reported in C4 plants, and related to the operation of Rubisco at or close to substrate saturation [69]. This finding is in agreement with Lüttge (2011) [16], who reported that the Rubisco specificity for CO2/O2 (Sc/o) of two CAM species of Kalanchoë was at the lower end of the range given for vascular plants. Overall, our results and those by Lüttge (2011) [16] would be indicative of convergent evolution of Rubisco catalysis of C4 and CAM plants, in the sense of a retro-evolution under the influence of the internal high CO2 concentration. The lower ratio [Rubisco]/[TSP] found in the strong CAM species (Table 4) also mimics the lower content of Rubisco in C4 plants [70]. Nevertheless, other studies suggested that CAM Rubiscos retain high CO2 affinity (i.e., low Michaelis-Menten constant for CO2, Kc) similar to C3 plants and lower than C4 species [20, 24, 33, 71]. Although the data set available in this study is too small to identify any clear trend, the apparently contradictory results may be attributable to the inherent mechanism of CAM for modulating the relative proportions of Rubisco and PEPC-mediated uptake of atmospheric CO2 [64, 71]. Such mechanism determines a wide range of variation in the midday internal CO2/O2 ratio among different CAM plants [15, 71], and therefore, different degrees of suppression of the oxygenase activity of Rubisco. The apparent variability in Rubisco kinetics associated to CAM could be linked to the plasticity of CAM expression and duration of the different CAM phases and therefore to the different availability of CO2 to Rubisco, pointing to the CO2 as a driver to Rubisco kinetics evolution [71]. A wide survey on the full Rubisco kinetics including representatives of the different families with CAM species is required to shed light on the evolution of Rubisco kinetics in CAM plants.

The candidate sites under positive selection (Table 2) and resolved with DT (Table 3) in orchids and bromeliads with contrasting Rubisco kinetics did not provide any clear trend on the molecular determinants of L-subunit variability (Table 4). While the present results reveal that there are potential differences in Rubisco traits between phylogenetically related C3 and CAM species, more data are needed to confirm this trend and to link kinetic differences to amino acid replacements within the L-subunit. Although not well understood yet, the different expression of rbcS genes, encoding for the small subunit (S-subunit) of Rubisco, may allow optimizing the Rubisco performance in response to changing environmental conditions [30, 72,73,74]. In view of the phenotypic plasticity inherent of the CAM metabolism [4, 64] a role of the S-subunit in the catalysis of Rubisco may not be discarded and should be a matter of future studies. Alternatively, the fact that genes with similar kinetic properties have different amino acid sequences could mean that different lineages used different replacements to lead to the same kinetic changes.

Conclusion

This study presents an extensive analysis of Rubisco molecular and biochemical characterization in two angiosperm families with C3 and CAM photosynthetic pathways. The study includes, for the first time, analyses of closely related C3 and CAM species, in particular i) positive selection and coevolution analyses, along with a DT model for variable sites related to physiological and anatomical information, and ii) measurements of Rubisco kcatc and Kc, that permitted to explore the variability in the Rubisco L-subunit sequences and study their biochemical impact. Signal of positive selection was found in rbcL and it could be linked to CAM through the DT. The previously reported trade-off between Kc and kcatc was observed in a subset of studied species, with strong CAM bromeliads presenting high kcatc while C3 bromeliads presenting high affinity for CO2 (i.e., low Kc). In spite of the differences between C3 and CAM bromeliads, the observed variation in the kinetic properties of Rubisco from distinct photosynthetic pathways could not be related to positively selected residues in the Rubisco L-subunit. A deeper inspection of variation in the Rubisco L- and S-subunits and Rubisco biochemical traits across a larger number of families containing C3 and CAM species may help to resolve these questions.

Methods

Species selection

We selected orchids and bromeliads rbcL sequences from GenBank with δ13C data available from literature: 123 bromeliads and 58 orchids species (Additional file 1: Table S1). In addition, we included other 42 species for which rbcL was sequenced in the present study (Table 1). Therefore, the final list of species under study contained 144 bromeliads and 79 orchids.

rbcL amplification and sequencing

For the 42 sequenced species the genomic DNA was isolated from dry leaves using a DNeasy Plant Mini Kit (Quiagen Ltd., Crawley, UK) in accordance with the manufacturer’s protocol. For rbcL amplification, we designed the primers esp2F (5′-AATTCATGAGTTGTAGGGAGGGACTT-3′), B1R (5′-CAATTAGGAGAACAAAGAGGAA-3′), O2F (5′-GAGTAGACCTTGTTGTTGTG-3′) and 1925R (5′-GACACGAGATTCTACGAGA-3′), and used 1494R (5′-GATTGGGCCGAGTTTAATTTAC-3′) [75]. The BioMix Red reaction mix (Bioline Ltd., London, UK) was used to carry out the polymerase chain reaction (PCR) with the following conditions: 1 initial cycle of 95 °C, 2 min; 55 °C, 30 s; 72 °C, 4 min followed by 36 cycles of 93 °C, 30 s; 53 °C, 30 s; 72 °C, 3.5 min. PCR products were visualized on 1% agarose gels, purified using the High Pure PCR Product Purification Kit (Roche, Germany) and sequenced with an ABI 3130 Genetic analyzer and the contings were assembled using BioEdit v7.1.3 software [76]. Novel sequences have been submitted to GenBank (Table 1). Nucleotide sequences were converted into amino acidic sequences with MEGA 5 [77] and then aligned using MAFFT v5 [78].

Plant material

Among the full dataset of sequences, there were a total of 58 living specimens available at Heidelberg Botanical Garden (Heidelberg University, Germany) (Table 1 and Additional file 1: Table S1). The growth conditions in the glasshouses corresponded approximately to their natural environmental conditions. Natural daylight was supplemented by additional artificial light (photon flux density of 275 μmol m− 2 s− 1) all over the year. From May until October the glasshouses were partially shaded (approx. 65%). Day-time and night-time minimum temperatures were in the range of 18–20 °C and 14–18 °C, respectively. Relative humidity was kept within the range 70–95%. Plants were watered daily and using a conventional nutrient solution once a week. There were some special cases, e.g., Puya was cultivated under dry conditions and full sun light. All ‘grey tillandsia’ were kept outside the glasshouse (shaded as described above) from May to October with daily watering, but they were kept much dryer during the winter season, when these plants do not grow and rest.

Species and varieties were classified according to their habitat preference into epiphyte, terrestrial or lithophyte [53]. A complete documentation is accessible at [79].

Leaf traits, carbon isotopic composition and photosynthetic mechanism classification

Material for leaf traits and carbon isotopic composition determination consisted in four fully expanded leaves (replicates) in mature stage sampled from different individuals in June 2014.

After the thickness of the leaf lamina was measured between the leaf margin and midrib of the middle portions of leaves using a slide caliper (Vernier Caliper, Series 530, Mitutoyo Europe GmbH), the leaf was detached and the fresh weight (FW) immediately recorded. The leaf area of the same sample was measured after digitalizing the leaf and using the ImageJ software [80]. The dry weight (DW) was obtained after drying the leaves in a ventilated oven at 60 °C until constant weight (typically after 2 days). The leaf mass area (LMA) was calculated as the ratio between the dry weight and the area.

Values for the leaf carbon isotopic composition (δ13C) were taken from bibliography (Table 1 and Additional file 1: Table S1), except for Acineta densa, Bulbophyllum lobbii, Elleantus furfuraceus, Epidendrum rigidum, Laelia speciosa, Lycaste cruenta, Maxillaria cucullata, Pleurothallis nuda, Sobralia macrantha, Cryptanthus fosterianus and Puya humilis for which it was measured. The dried leaves used for the characterization of the leaf traits were ground into powder and subsamples of 2 mg were analyzed. Samples were combusted in an elemental analyzer (Carlo-Erba, Rodano, Italy). The CO2 was separated by chromatography and directly injected into a continuous-flow isotope ratio mass spectrometer (Thermo Finnigan Delta Plus, Bremen, Germany). Peach leaf standards (NIST 1547) were run every six samples. The δ13C was calculated as: δ13C sample (‰) = (((13C/12C) sample/(13C/12C) standard) - 1) 1000 [81] and values were referred to a Pee Dee Belemnite standard.

The photosynthetic mechanism of the species was inferred from the δ13C values, following previous surveys in orchids and bromeliads [43, 45]. Species with δ13 C > − 18 ‰ and < − 22.9 ‰ were classified as CAM and C3 respectively. In literature, species with δ13 C between − 18 ‰ and − 22.9 ‰ are considered as “weak CAM” (Fig. 1). According to Winter and Holtum (2002) [82], δ13C values below − 25 ‰ may indicate that CO2 fixation occurs exclusively in the light, while δ13C values above − 21.9 ‰ reflect that at least 50% of CO2 fixation occurs in the dark.

Detection of positive selection in rbcL

Positive selection acting on the Rubisco L-subunit was analyzed with the PAML package v4.7 [83] and PAMLX [84]. Codeml program [85] was used to calculate the non-synonymous (dN) and synonymous (dS) substitution rates across codons and the dN/dS ratio (ω). This ratio represents the selective pressures acting on the protein-coding gene with values of ω < 1, ω = 1, and ω > 1 being indicative of purifying selection, neutral evolution and positive selection, respectively.

The tree topologies based on rbcL sequences were constructed using maximum-likelihood inference conducted with RAxML version 7.2.6 [86]. It was done without species with δ13 C values between − 18.0 and − 22.9 ‰ because species with values around − 20 ‰ might be weak CAM and other may be pure C3 with no detectable CAM [41]. Therefore, the tree topologies were finally constructed with 78 orchids and 130 bromeliads (Figs. 3 and 4) and edited with Fig Tree v 1.4.0 [87].

Site models allow the ω ratio to vary among codons in the protein [88]. To identify signatures of adaptive evolution we performed two nested maximum likelihood tests: M1a vs. M2a and M8a vs. M8 [89, 90]. The null M1a model assumes purifying selection or nearly neutral evolution without positive selection and allow codons with ω < 1 and/or ω = 1, but not codons with ω > 1. The M2a model allows for codons under positive selection (ω > 1). Model M8a assumes a discrete beta distri- bution for ω, which is constrained between 0 and 1 including a class with ω = 1. Model M8 allows the same distribution as M8a with an extra class of codons under positive selection with ω > 1. Posterior probabilities for site classes were calculated with Bayes Empirical Bayes (BEB) [90].

Branch-site models allow ω to vary both among sites in the protein and across branches on the tree with the aim to detect positive selection affecting a few sites along particular branches. The branch-site A model was applied for branches leading to CAM species and for clades containing CAM species. The branch types are specified using labels in the tree file; e.g. if the dataset has CAM branch types, they are labelled using #. Species with δ13 C > − 18 ‰ were classified as CAM. The A1-A LRT compared the null model A1 with the nested model A. Both the A1 and A models allow ω ratios to vary among sites [83, 91]. The A1 model allows 0 < ω < 1 and ω = 1 for all branches and also two additional classes of codons with fixed ω = 1 along pre-specified branches, while restricted as 0 < ω < 1 and ω = 1 on background branches. The alternative A model allows 0 < ω < 1 and ω = 1 for all branches and also two alternative classes of codons under positive selection with ω > 1 along pre-specified branches, while restricted as 0 < ω < 1 and ω = 1 on background branches.

We performed three LRTs to compare the nested site models M1a-M2a, M8-M8a and branch-site models A-A1. LRTs involves the comparison of the log-likelihood values of the simple and the complex nested models and twice their difference follows a chi-square distribution with the degrees of freedom (df) being the difference in the number of free parameters between the models. For the comparison of models M1a-M2a, M8a-M8 and A1-A the df was 2, 1 and 1, respectively.

Analysis of intra-molecular coevolution in the amino acidic sequence of the Rubisco large subunit

Intra-molecular coevolution analysis was performed with the program CAPS [92, 93]. The algorithm implemented in this program identifies co-evolving amino acid site pairs by measuring the correlated evolutionary variation at these sites using time corrected Blosum values. CAPS take into account the time of sequences divergence such that correlated variation that involves radical amino acid substitutions is considered to be more likely at longer evolutionary times following a Poisson model [92, 93]. Accordingly, the transition between two amino acids at each site is corrected by the divergence time of the sequences. Synonymous substitutions per site do not affect the amino acid composition of the protein and are neutrally fixed in the gene, being the number of such substitutions proportional to the time of sequence divergence. In this respect, time since two sequences diverge is estimated as the mean number of substitutions per synonymous site between the two sequences being compared. Correlation of the mean variability is measured using the Pearson coefficient. The significance of the correlation coefficient is estimated by comparing the real correlation coefficients to the distribution of resampled correlation coefficients.

Decision tree (DT) model

DT model analysis (‘rpart’ package in R v3.1.1 [55]) was used to relate the proportion of amino acid presence in all variable sites of the L-subunit of Rubisco to species-specific traits (δ13C and habitat preference), denoted as external variables, as listed in Table 1 and Additional file 1: Table S1.

For each variable site, the program builds a DT as follows. Based on the external variables (δ13C and habitat preference), the species are separated into two groups, in which the variability of that site is as low as possible. The analysis is repeated for each subgroup using again the two external variables. The process continues until the lowest xerror [94] for the entire DT is obtained. In the case of δ13C as an external variable the whole range of numerical values were considered, species with δ13 C > − 18 ‰ and < − 22.9 ‰, and in the case of habitat as an external variable, three options were possible for each species (epiphyte, lithophyte, terrestrial), so we have given a proportional value (0.34, 0.33, 0.33) for the construction of the DT.

The quality of the DT is categorized by its entropic error (xerror) as a function of the proportion of correct predictions and the complexity of the tree. The lower the xerror, the higher the correlation between the external variable and the variable site. Only DTs with xerror < 1 were selected. The relative importance of an external variable is computed as a function of the reduction of errors that the selected external variable produces on the variable site.

Rubisco kinetics measurements

For the catalytic characterization of Rubisco, a number of orchid and bromeliad species was selected as representing the different photosynthetic types and reflecting the maximum variability in positively selected residues of the Rubisco L-subunit sequence. The list of species initially selected was: Acineta densa, Bulbophyllum lobbii, Epidendrum ciliare, Epidendrum difforme, Nidularium fulgens, Nidularium innocentii var. innocentii, Nidularium regelioides, Maxillaria cucullata, Oncidium lineoligerum, Lockhartia amoena, Elleanthus furfuraceus, Myoxanthus exasperatus, Sobralia macrantha, Nidularium innocentii lineatum, Tillandsia biflora, Tillandsia multicaulis, Aechmea nudicaulis var. aureo-rosea and Tillandsia bermejoensis. Specimens of these species sent from Heidelberg were grown in the glasshouse at the University of the Balearic Islands under similar conditions described for Heidelberg.

Different protein extraction media were tested on these species. These tests determined that the most appropriate protein extraction media were buffers A and B. Extraction buffer A consisted of 100 mM Bicine-NaOH (pH 8.0), 0.1 mM EDTA, PEG4000 (6% w/v), 20 mM DTT, 50 mM 2-mercaptoethanol, 2 mM MgCl2, 10 mM NaHCO3, 1 mM benzamidine, 1 mM β-aminocaproic acid, 2 μM pepstatin, 10 μM E-64 (Sigma, USA) 10 μM chymostatin, 2 mM phenylmethylsulfonyl fluoride and 25 mg mL− 1 PVP. Extraction buffer B consisted of 350 mM HEPES-KOH (pH 8.0), 6% (w/v) PEG4000, 2 mM MgCl2, 0.1 mM EDTA, 1 mM benzamidine, 1 mM ε-aminocaproic acid, 10 mM NaHCO3. Added into the mortar: 7 μL β-mercaptoethanol, 400 μL DTT (1 M), 4 μL pepstatine, 4 μL E-64, 4 μL chymostatin, 10 μL PMSF, 75 mg PVP and 75 mg PVPP.

Extraction buffer A worked with M. cucullata, O. lineoligerum, L. amoena, E. furfuraceus, M. exasperatus, S. macrantha, N. innocentii lineatum, T. biflora and T. multicaulis. Leaf soluble protein of A. nudicaulis var. aureo-rosea and T. bermejoensis was successfully extracted using buffer B.

As for the remaining species, A. densa, B. lobbii, E. ciliare, E. difforme, N. fulgens, N. innocentii var. innocentii, N. regelioides, these two buffers yielded poor soluble protein and low Rubisco activity, and up to other four extraction buffers were tested by varying both the components and their concentration. However, none of these buffers extracted sufficient amount of active Rubisco to characterize the kinetic constants.

Leaf soluble protein was extracted on fully expanded leaves of 3–4 plants per species by grinding 0.40–0.60 g of leaf samples in a mortar with 2 mL of ice-cold extraction buffer. The proportion of leaf total soluble protein that is accounted for by Rubisco ([Rubisco]/[TSP]), the Rubisco Michaelis-Menten constant for CO2 (Kc) and the maximum rate of carboxylation (kcatc) were measured at 25 °C in semi-purified extracts following [33]. Rates of Rubisco 14CO2-fixation using the activated protein extract were measured in 7 mL septum capped scintillation vials in reaction buffer (110 mM Bicine-NaOH pH 8.0, 22 mM MgCl2, 0.4 mM RuBP and ~ 100 W-A units of carbonic anhydrase), equilibrated with nitrogen (N2). Different concentrations of H14CO3− (0, 6.7, 26.7, 53.3, 88.9, 122.2 and 155.6 μM for orchids, and 0, 6.7, 26.7, 53.3, 88.9, 122.2, 155.6 and 190 μM for bromeliads; each with a specific radioactivity of 3.7 × 1010 Bq mol− 1) were prepared in the scintillation vials as described previously [33]. Assays (1.0 mL total volume) were started by injection of activated leaf extract and stopped after 60 s with the addition of 1 M formic acid. The acidified mixtures were dried and the 14C products determined via scintillation counting. Concentrations of CO2 in solution in equilibrium with H14CO3− were calculated assuming a pKa for carbonic acid of 6.23.

Triticum aestivum cv. Cajeme was grown from seeds at the UIB under full irrigation and frequent fertilization with Hoagland’s solution [95]. Rubisco was extracted from wheat mature leaves using extraction buffer A, and the kinetic parameters measured following the same procedures as with orchids and bromeliads.

Statistical analyses

Univariate analysis of variance (ANOVA) was used to statistically examine the differences among species and photosynthetic mechanisms for the Rubisco kinetic parameters, [Rubisco]/[TSP] and leaf mass per unit area (LMA). Phylogenetic Independent Contrast (PIC) analysis was performed using R packages APE and GEIGER [96, 97]. Significant differences between means were revealed by Duncan analyses (p < 0.05) [98]. Regression coefficients between parameters were performed with R [55].

Availability of data and materials

The datasets supporting the conclusions of this article are included within the article and its additional files.

Abbreviations

ANOVA:

Analysis of variance

BEB:

Bayes Empirical Bayes

CAM:

Crassulacean acid metabolism

DT:

Decision tree

DTT:

Dithiothreitol

DW:

Dry weight

EDTA:

Ethylenediaminetetraacetic acid

FW:

Fresh weight

LMA:

Leaf mass area

LRT:

Likelihood ratio test

L-subunit:

Large subunit

PCR:

Polymerase chain reaction

PEG:

Polyethylene glycol

PEPC:

Phosphoenolpyruvate carboxylase

PIC:

Phylogenetic independent contrast

PVP:

Polyvinylpyrrolidone

PVPP:

Polyvinylpolypyrrolidone

RuBP:

Ribulose 1,5-bisphosphate

TSP:

Total soluble protein

WUE:

Water use efficiency

References

  1. Osmond CB. Crassulacean acid metabolism: a curiosity in context. Annu Rev Plant Physiol. 1978;29:379–414.

    Article  CAS  Google Scholar 

  2. Lüttge U. Ecophysiology of crassulacean acid metabolism (CAM). Ann Bot-Lon. 2004;93:629–52.

    Article  CAS  Google Scholar 

  3. Winter K, Garcia M, Holtum JA. On the nature of facultative and constitutive CAM: environmental and developmental control of CAM expression during early growth of Clusia, Kalanchoë, and Opuntia. J Exp Bot. 2008;59:1829–40.

    Article  CAS  PubMed  Google Scholar 

  4. Cushman JC. Crassulacean acid metabolism. A plastic photosynthetic adaptation to arid environments. Plant Physiol. 2001;127:1439–48.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Nelson EA, Sage TL, Sage RF. Functional leaf anatomy of plants with crassulacean acid metabolism. Funct Plant Biol. 2005;32:409–19.

    Article  PubMed  Google Scholar 

  6. Silvera K, Neubig KM, Whitten WM, Williams NH, Winter K, Cushman JC. Evolution along the crassulacean acid metabolism continuum. Funct Plant Biol. 2010;37:995–1010.

    Article  CAS  Google Scholar 

  7. Winter K, Smith JAC. An introduction to crassulacean acid metabolism. Biochemical principles and ecological diversity. In: Winter K, JAC S, editors. Crassulacean acid metabolism. Berlin: Springer; 1996. p. 1–13.

    Chapter  Google Scholar 

  8. Holtum JA, Winter K, Weeks MA, Sexton TR. Crassulacean acid metabolism in the ZZ plant, Zamioculcas zamiifolia (Araceae). Am J Bot. 2007;94:1670–6.

    Article  CAS  PubMed  Google Scholar 

  9. Valdez-Hernández M, González-Salvatierra C, Reyes-García C, Jackson PC, Andrade JL. Physiological Ecology of Vascular Plants. In: Islebe GA, et al., editors. Biodiversity and Conservation of the Yucatán Peninsula. International Publishing Switzerland: Springer; 2015. p. 97–129.

    Chapter  Google Scholar 

  10. Griffiths H. Carbon dioxide concentrating mechanisms and the evolution of CAM in vascular epiphytes. In: Lüttge U, editor. Vascular plants as epiphytes. Berlin: Springer; 1989. p. 42–86.

    Google Scholar 

  11. Keeley JE, Rundel PW. Evolution of CAM and C4 Carbon-Concentrating Mechanisms. Int J Plant Sci. 2003;164:S55–77.

    Article  CAS  Google Scholar 

  12. Ehleringer JR, Monson RK. Evolutionary and ecological aspects of photosynthetic pathway variation. Annu Rev Ecol Syst. 1993;24:411–39.

    Article  Google Scholar 

  13. Crayn DM, Winter K, Smith JAC. Multiple origins of crassulacean acid metabolism and the epiphytic habit in the Neotropical family Bromeliaceae. P Natl Acad Sci USA. 2004;101:3703–8.

    Article  CAS  Google Scholar 

  14. Arakaki M, Christin PA, Nyffeler R, Lendel A, Eggli U, Ogburn RM, Spriggs E, Moore MJ, Edwards EJ. Contemporaneous and recent radiations of the world's major succulent plant lineages. P Natl Acad Sci. 2011;108:8379–84.

    Article  CAS  Google Scholar 

  15. Lüttge U. CO2 concentrating: consequences in crassulacean acid metabolism. J Exp Bot. 2002;53:2131–42.

    Article  PubMed  CAS  Google Scholar 

  16. Lüttge U. Photorespiration in phase III of crassulacean acid metabolism: evolutionary and ecophysiological implications. In: Lüttge U, Matyssek R, Cánovas Ramos FM, editors. Progress in Botany 72. Berlin: Springer; 2011. p. 371–84.

    Chapter  Google Scholar 

  17. Tcherkez GG, Farquhar GD, Andrews TJ. Despite slow catalysis and confused substrate specificity, all ribulose bisphosphate carboxylases may be nearly perfectly optimized. P Natl Acad Sci. 2006;103:7246–51.

    Article  CAS  Google Scholar 

  18. Savir Y, Noor E, Milo R, Tlusty T. Cross-species analysis traces adaptation of Rubisco toward optimality in a low-dimensional landscape. P Natl Acad Sci. 2010;107:3475–80.

    Article  CAS  Google Scholar 

  19. Galmés J, Andralojc PJ, Kapralov MV, Flexas J, Keys AJ, Molins A, Conesa MÀ. Environmentally driven evolution of Rubisco and improved photosynthesis and growth within the C3 genus Limonium (Plumbaginaceae). New Phytol. 2014;203:989–99.

    Article  PubMed  CAS  Google Scholar 

  20. Badger MR, Andrews TJ, Osmond CB. Detection in C3, C4 and CAM plant leaves of a low K m(CO2) form of RuDP carboxylase, having high RuDP oxygenase activity at physiological pH. In: Avron M, editor. Proceedings of the third International Congress on Photosynthesis. Amsterdam: Elsevier; 1974. p. 1421–9.

    Google Scholar 

  21. Tcherkez G. Modelling the reaction mechanism of ribulose-1,5-bisphosphate carboxylase/oxygenase and consequences for kinetic parameters. Plant Cell Environ. 2013;36:1586–96.

    Article  CAS  PubMed  Google Scholar 

  22. Hermida-Carrera C, Kapralov MV, Galmés J. Rubisco Catalytic Properties and Temperature Response in Crops. Plant Physiol. 2016;171:2549–61.

    CAS  PubMed  PubMed Central  Google Scholar 

  23. Yeoh HH, Badger MR, Watson L. Variations in K m (CO2) of ribulose-1,5-bisphosphate carboxylase among grasses. Plant Physiol. 1980;66:1110–2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Yeoh HH, Badger MR, Watson L. Variations in kinetic properties of ribulose-1,5-bisphosphate carboxylases among plants. Plant Physiol. 1981;67:1151–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Seemann JR, Badger MR, Berry JA. Variations in the specific activity of ribulose-1,5-bisphosphate carboxylase between species utilizing differing photosynthetic pathways. Plant Physiol. 1984;74:791–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Ghannoum O, Evans JR, Chow WS, Andrews TJ, Conroy JP, Von Caemmerer S. Faster Rubisco is the key to superior nitrogen-use efficiency in NADP-malic enzyme relative to NAD-malic enzyme C4 grasses. Plant Physiol. 2005;137:638–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Kubien DS, Whitney SM, Moore PV, Jesson LK. The biochemistry of Rubisco in Flaveria. J Exp Bot. 2008;59:1767–77.

    Article  CAS  PubMed  Google Scholar 

  28. Kapralov MV, Filatov DA. Widespread positive selection in the photosynthetic Rubisco enzyme. BMC Evol Biol. 2007. https://doi.org/10.1186/1471-2148-7-73.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  29. Christin PA, Salamin N, Muasya AM, Roalson EH, Russier F, Besnard G. Evolutionary switch and genetic convergence on rbcL following the evolution of C4 photosynthesis. Mol Biol Evol. 2008;25:2361–8.

    Article  CAS  PubMed  Google Scholar 

  30. Kapralov MV, Kubien DS, Andersson I, Filatov DA. Changes in Rubisco kinetics during the evolution of C4 photosynthesis in Flaveria (Asteraceae) are associated with positive selection on genes encoding the enzyme. Mol Biol Evol. 2011;28:1491–503.

    Article  CAS  PubMed  Google Scholar 

  31. Kapralov MV, Smith JAC, Filatov DA. Rubisco evolution in C4 eudicots: an analysis of Amaranthaceae sensu lato. PLoS One. 2012;7:e52974.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Young JN, Rickaby REM, Kapralov MV, Filatov DA. Adaptive signals in algal Rubisco reveal a history of ancient atmospheric carbon dioxide. Philos T R Soc B. 2012;367:483–92.

    Article  CAS  Google Scholar 

  33. Galmés J, Kapralov MV, Andralojc P, Conesa MÀ, Keys AJ, Parry MA, Flexas J. Expanding knowledge of the Rubisco kinetics variability in plant species: environmental and evolutionary trends. Plant Cell Environ. 2014;37:1989–2001.

    Article  PubMed  CAS  Google Scholar 

  34. Rosnow JJ, Evans MA, Kapralov MV, Cousins AB, Edwards GE, Roalson EH. Kranz and single-cell forms of C4 plants in the subfamily Suaedoideae show kinetic C4 convergence for PEPC and Rubisco with divergent amino acid substitutions. J Exp Bot. 2015;66:7347–58.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Whitney SM, Houtz RL, Alonso H. Advancing our understanding and capacity to engineer nature’s CO2-sequestering enzyme, Rubisco. Plant Physiol. 2011;155:27–35.

    Article  CAS  PubMed  Google Scholar 

  36. Parry MA, Hawkesford MJ. An integrated approach to crop genetic improvement. J Integr Plant Biol. 2012;54:250–9.

    Article  PubMed  Google Scholar 

  37. Galmés J, Flexas J, Keys AJ, Cifre J, Mitchell RAC, Madgwick PJ. RP Haslam, Medrano H, Parry MAJ. Rubisco specificity factor tends to be larger in plant species from drier habitats and in species with persistent leaves. Plant Cell Environ. 2005;28:571–9.

    Article  Google Scholar 

  38. Givnish TJ, Barfuss MH, Van EB. Riina R, Schulte K, Horres R, Sytsma KJ. Phylogeny, adaptive radiation, and historical biogeography in Bromeliaceae: insights from an eight-locus plastid phylogeny. Am J Bot. 2011;98:872–95.

    Article  PubMed  Google Scholar 

  39. Stevens PF. 2013. Angiosperm Phylogeny. http://www.mobot.org/mobot/research/apweb/. Accessed 21 Sept 2013.

  40. Crayn DM, Winter K, Schulte K, Smith JAC. Photosynthetic pathways in Bromeliaceae: phylogenetic and ecological significance of CAM and C3 based on carbon isotope ratios for 1893 species. Bot J Linn Soc. 2015;178:169–221.

    Article  Google Scholar 

  41. Silvera K, Santiago LS, Winter K. Distribution of crassulacean acid metabolism in orchids of Panama: evidence of selection for weak and strong modes. Funct Plant Biol. 2005;32:397–407.

    Article  CAS  PubMed  Google Scholar 

  42. Motomura H, Ueno O, Kagawa A, Yukawa T. Carbon isotope ratios and the variation in the diurnal pattern of malate accumulation in aerial roots of CAM species of Phalaenopsis (Orchidaceae). Photosynthetica. 2008;46:531–6.

    Article  CAS  Google Scholar 

  43. Silvera K, Santiago LS, Cushman JC, Winter K. The incidence of crassulacean acid metabolism in Orchidaceae derived from carbon isotope ratios: a checklist of the flora of Panama and Costa Rica. Bot J Linn Soc. 2010;163:194–222.

    Article  Google Scholar 

  44. Borland AM, Zambrano VAB, Ceusters J, Shorrock K. The photosynthetic plasticity of crassulacean acid metabolism: an evolutionary innovation for sustainable productivity in a changing world. New Phytol. 2011;191:619–33.

    Article  CAS  PubMed  Google Scholar 

  45. Silvera K, Santiago LS, Cushman JC, Winter K. Crassulacean acid metabolism and epiphytism linked to adaptive radiations in the Orchidaceae. Plant Physiol. 2009;149:1838–47.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Givnish TJ, Barfuss MH, Van Ee B, Riina R, Schulte K, Horres R, Sytsma KJ. Adaptive radiation, correlated and contingent evolution, and net species diversification in Bromeliaceae. Mol Phylogenet Evol. 2014;71:55–78.

    Article  PubMed  Google Scholar 

  47. Winter K, Holtum JA. The effects of salinity, crassulacean acid metabolism and plant age on the carbon isotope composition of Mesembryanthemum crystallinum L., a halophytic C3-CAM species. Planta. 2005;222:201–9.

    Article  CAS  PubMed  Google Scholar 

  48. Vargas-Soto JG, Andrade JL, Winter K. Carbon isotope composition and mode of photosynthesis in Clusia species from Mexico. Photosynthetica. 2009;47:33–40.

    Article  CAS  Google Scholar 

  49. Martin CE, Mas EJ, Lu C, Ong BL. The photosynthetic pathway of the roots of twelve epiphytic orchids with CAM leaves. Photosynthetica. 2010;48:42–50.

    Article  Google Scholar 

  50. Winter K, Wallace BJ, Stocker GC, Roksandic Z. Crassulacean acid metabolism in Australian vascular epiphytes and some related species. Oecologia. 1983;57:129–41.

    Article  PubMed  Google Scholar 

  51. Earnshaw MJ, Winter K, Ziegler H, Stichler W, Cruttwell NEG, Kerenga K, Gunn TC. Altitudinal changes in the incidence of crassulacean acid metabolism in vascular epiphytes and related life forms in Papua New Guinea. Oecologia. 1987;73:566–72.

    Article  CAS  PubMed  Google Scholar 

  52. Kluge M, Vinson B, Ziegler H. Ecophysiological studies on orchids of Madagascar: incidence and plasticity of crassulacean acid metabolism in species of the genus Angraecum Bory. Plant Ecol. 1997;135:43–57.

    Article  Google Scholar 

  53. Koch MA, Schröder N, Kiefer M, Sack P. A treasure trove of plant biodiversity from the 20th century: The Werner Rauh heritage project at Heidelberg Botanical Garden and Herbarium. Plant Syst Evol. 2013;299:1793–800.

    Article  Google Scholar 

  54. Crayn DM, Smith JAC, Winter K. Carbon Isotope Ratios and Photosynthetic Pathways in the Neotropical Family Rapateaceae. Plant Biol. 2001;3:569–76.

    Article  CAS  Google Scholar 

  55. R Core Team. R: A language and environment for statistical computing: R Foundation for Statistical Computing; 2014. http://www.R-project.org/

  56. Cameron KM, Chase MW, Whitten WM, Kores PJ, Jarrell DC, Albert VA, Yukawa T, Hills HG, Goldman DH. A phylogenetic analysis of the Orquidaceae: evidence from rbcL nucleotide sequences. Am J Bot. 1999;86:208–24.

    Article  PubMed  Google Scholar 

  57. Edwards EJ, Osborne CP, Strömberg CA, Smith SA. The origins of C4 grasslands: integrating evolutionary and ecosystem science. Science. 2010;328:587–91.

    Article  CAS  PubMed  Google Scholar 

  58. Christin PA, Sage TL, Edwards EJ, Ogburn RM, Khoshravesh R, Sage RF. Complex evolutionary transitions and the significance of C3-C4 intermediate forms of photosynthesis in Molluginaceae. Evolution. 2011;65:643–60.

    Article  PubMed  Google Scholar 

  59. Iida S, Miyagi A, Aoki S, Ito M, Kadono Y, Kosuge K. Molecular adaptation of rbcL in the heterophyllous aquatic plant Potamogeton. PLoS One. 2009;4:e4633.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  60. Miwa H, Odrzykoski IJ, Matsui A, Hasegawa M, Akiyama H, Jia Y, Sabirov R, Takahashi H, Boufford DE, Murakami N. Adaptive evolution of rbcL in Conocephalum (Hepaticae, bryophytes). Gene. 2009;441:169–75.

    Article  CAS  PubMed  Google Scholar 

  61. Kato S, Misawa K, Takahashi F, Sakayama H, Sano S, Kosuge K, Kasai F, Watanabe MM, Tanaka J, Nozaki H. Aquatic plant speciation affected by diversifying selection of organelle DNA regions. J Phycol. 2011;47:999–1008.

    Article  PubMed  Google Scholar 

  62. Sen L, Fares MA, Liang B, Gao L, Wang B, Wang T, Su YJ. Molecular evolution of rbcL in three gymnosperm families: identifying adaptive and coevolutionary patterns. Biol Direct. 2011;6:29.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  63. Wang M, Kapralov MV, Anisimova M. Coevolution of amino acid residues in the key photosynthetic enzyme Rubisco. BMC Evol Biol. 2011. https://doi.org/10.1186/1471-2148-11-266.

  64. Dodd AN, Borland AM, Haslam RP, Griffiths H, Maxwell K. Crassulacean acid metabolism: plastic, fantastic. J Exp Bot. 2002;53:569–80.

    Article  CAS  PubMed  Google Scholar 

  65. Maxwell K, Griffiths H, Helliker B, Roberts A, Haslam RP, Girnus J, Robe WE, Borland AM. Regulation of Rubisco activity in crassulacean acid metabolism plants: better late than never. Funct Plant Biol. 2002;29:689–96.

    Article  PubMed  Google Scholar 

  66. Silvera K, Winter K, Rodriguez BL, Albion RL, Cushman JC. Multiple isoforms of phosphoenolpyruvate carboxylase in the Orchidaceae (subtribe Oncidiinae): implications for the evolution of crassulacean acid metabolism. J Exp Bot. 2014;65:3623–36.

    Article  PubMed  PubMed Central  Google Scholar 

  67. Flexas J, Ribas Carbo M, Diaz Espejo A, Galmés J, Medrano H. Mesophyll conductance to CO2: current knowledge and future prospects. Plant Cell Environ. 2008;31:602–21.

    Article  CAS  PubMed  Google Scholar 

  68. Delgado E, Medrano H, Keys AJ, Parry MAJ. Species variation in Rubisco specificity factor. J Exp Bot. 1995;46:1775–7.

    Article  CAS  Google Scholar 

  69. von Caemmerer S, Furbank RT. Modeling C4 photosynthesis. In: Sage RF, Monson RK, editors. C4 plant biology. San Diego: Academic Press; 1999. p. 173–211.

    Chapter  Google Scholar 

  70. Wessinger ME, Edwards GE, Ku MSB. Quantity and kinetic properties of ribulose 1,5-bisphosphate carboxylase in C3, C4, and C3-C4, intermediate species of Flaveria (Asteraceae). Plant Cell Physiol. 1989;30:665–71.

    CAS  Google Scholar 

  71. Griffiths H, Robe WE, Girnus J, Maxwell K. Leaf succulence determines the interplay between carboxylase systems and light use during Crassulacean acid metabolism in Kalanchoë species. J Exp Bot. 2008;59:1851–61.

    Article  CAS  PubMed  Google Scholar 

  72. Ishikawa C, Hatanaka T, Misoo S, Miyake C, Fukayama H. Functional incorporation of sorghum small subunit increases the catalytic turnover rate of Rubisco in transgenic rice. Plant Physiol. 2011;156:1603–11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Cavanagh AP, Kubien DS. Can phenotypic plasticity in Rubisco performance contribute to photosynthetic acclimation? Photosynth Res. 2013;119:203–14.

    Article  PubMed  CAS  Google Scholar 

  74. Morita K, Hatanaka T, Misoo S, Fukayama H. Unusual small subunit that is not expressed in photosynthetic cells alters the catalytic properties of Rubisco in rice. Plant Physiol. 2014;164:69–79.

    Article  CAS  PubMed  Google Scholar 

  75. Chen ZD, Wang XQ, Sun HY, Yin H, Zhang ZX, Zou YP, Lu AM. Systematic position of the Rhoipteleaceae: evidence from rbcL sequences. Acta Phytotaxon Sin. 1998;36:1–7.

    Google Scholar 

  76. Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucl Acid S. 1999:95–8.

  77. Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol. 2011;28:2731–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  78. Katoh K, Kuma KI, Toh H, Miyata T. MAFFT v5: improvement in accuracy of multiple sequence alignment. Nucleic Acids Res. 2005;33:511–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. WRHP database. http://scriptorium.cos.uni-heidelberg.de. Accessed 7 Apr 2018.

  80. Rasband, W.S., ImageJ, U. S. National Institutes of Health, Bethesda, Maryland, USA, http://imagej.nih.gov/ij/, 1997–2015. Accessed 9 Oct 2014.

  81. Farquhar GD, Richards RA. Isotopic composition of plant carbon correlates with water–use efficiency of wheat genotypes. Funct Plant Biol. 1984;11:539–52.

    Article  CAS  Google Scholar 

  82. Winter K, Holtum JA. How closely do the δ13C values of crassulacean acid metabolism plants reflect the proportion of CO2 fixed during day and night? Plant Physiol. 2002;129:1843–51.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  83. Yang Z. PAML 4: phylogenetic analysis by maximum likelihood. Mol Biol Evol. 2007;24:1586–91.

    Article  CAS  PubMed  Google Scholar 

  84. Xu B, Yang Z. PAMLX: A graphical user interface for PAML. Mol Biol Evol. 2013;30:2723–4.

    Article  CAS  PubMed  Google Scholar 

  85. Yang Z. Likelihood ratio tests for detecting positive selection and application to primate lysozyme evolution. Mol Biol Evol. 1998;15:568–73.

    Article  CAS  PubMed  Google Scholar 

  86. Stamatakis A. RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics. 2006;22:2688–90.

    Article  CAS  PubMed  Google Scholar 

  87. Rambaut A. Figtree 1.4.0. 2012. http://tree.bio.ed.ac.uk/software/figtree/.

    Google Scholar 

  88. Yang Z, Nielsen R. Codon-substitution models for detecting molecular adaptation at individual sites along specific lineages. Mol Biol Evol. 2002;19:908–17.

    Article  CAS  PubMed  Google Scholar 

  89. Wong WS, Yang Z, Goldman N, Nielsen R. Accuracy and power of statistical methods for detecting adaptive evolution in protein coding sequences and for identifying positively selected sites. Genetics. 2004;168:1041–51.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  90. Yang Z, Wong WS, Nielsen R. Bayes empirical Bayes inference of amino acid sites under positive selection. Mol Biol Evol. 2005;22:1107–18.

    Article  CAS  PubMed  Google Scholar 

  91. Yang Z, dos Reis M. Statistical properties of the branch-site test of positive selection. Mol Biol Evol. 2011;28:1217–28.

    Article  CAS  PubMed  Google Scholar 

  92. Fares MA, McNally D. CAPS: coevolution analysis using protein sequences. Bioinformatics. 2006;22:2821–2.

    Article  CAS  PubMed  Google Scholar 

  93. Fares MA, Travers SA. A novel method for detecting intramolecular coevolution: adding a further dimension to selective constraints analyses. Genetics. 2006;173:9–23.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  94. Breiman L, Friedman JH, Olshen RA, Stone CJ. Classification and Regression Trees. Wadsworth: Chapmann & Hall; 1984.

  95. Hoagland DR, Arnon DI. The water-culture method for growing plants without soil. In: Arnon DI, editor. Circular: California Agricultural Experiment Station; 1950.

  96. Paradis E, Claude J, Strimmer K. APE: analyses of phylogenetics and evolution in R language. Bioinformatics. 2004;20:289–90.

    Article  CAS  PubMed  Google Scholar 

  97. Harmon Luke J, Weir JT, Brock CD, Glor RE, Challenger W. GEIGER: investigating evolutionary radiations. Bioinformatics. 2008;24:129–31.

    Article  CAS  PubMed  Google Scholar 

  98. IBM Corp. IBM SPSS Statistics for Windows, Version 21.0. Armonk, NY: IBM Corp; 2012.

Download references

Acknowledgements

We thank the Botanical Gardens of Heidelberg University for providing the plant material. We thank M Truyols at the UIB Experimental Field and Greenhouses (UIBGrant 15/2015) and T García at the radioisotope service at UIB. We thank CN Schröder for help screening exhaustive relevant material to determine life history traits. Funding from Klaus-Tschira-Foundation gGmbH (Heidelberg, Germany) awarded to MAK allowed for taxon sampling set-up and respective historical data collection. We thank C Douthe and J Cifre for improving different parts of the manuscript. Prof. U Lüttge is acknowledged for his stimulating discussions on CAM during the Methods in Plant Ecophysiology course held in Mallorca in June 2014. We would like to dedicate this article to our colleague Mario Fares, who participated enthusiastically in this project. Rest in peace, Mario.

Funding

Projects AGL2009–07999 and AGL2013–42364-R (FEDER/Spanish Ministry of Economy, Industry and Competitiveness – Research State Agency) awarded to J. Galmés contributed to the design of the study, collection, analysis and interpretation of data, and in writing the manuscript. Project PGC2018-096956-B-C43 (FEDER/Spanish Ministry of Economy, Industry and Competitiveness – Research State Agency) awarded to A. Mir and J. Rocha contributed to analysis and interpretation of data. Funding from Klaus-Tschira-Foundation gGmbH (Heidelberg, Germany) awarded to M.A. Koch allowed for taxon sampling set-up and respective historical data collection. C. Hermida-Carrera was supported by a FPI Fellowship BES-2010-030796.

Author information

Authors and Affiliations

  1. Research Group on Plant Biology under Mediterranean Conditions, Universitat de les Illes Balears-INAGEA, Ctra. Valldemossa km. 7.5, 07122, Palma, Illes Balears, Spain

    Carmen Hermida-Carrera, Marcel Font-Carrascosa, Arántzazu Molins, Miquel Ribas-Carbó & Jeroni Galmés

  2. Integrative and Systems Biology Group, Department of Abiotic Stress, Instituto de Biología Molecular y Celular de Plantas (CSIC–UPV), 46022, Valencia, Spain

    Mario A. Fares

  3. Department of Genetics, Trinity College Dublin, University of Dublin, Dublin 2, Ireland

    Mario A. Fares

  4. School of Natural and Environmental Sciences, Newcastle University, Newcastle upon Tyne, NE1 7RU, United Kingdom

    Maxim V. Kapralov

  5. Department of Biodiversity and Plant Systematics, Centre for Organismal Studies (COS) Heidelberg, Heidelberg University, Im Neuenheimer Feld 345, 9120, Heidelberg, Germany

    Marcus A. Koch

  6. Computational Biology and Bioinformatics Research Group, Department of Mathematics and Computer Science, Universitat de les Illes Balears, 07122, Palma, Illes Balears, Spain

    Arnau Mir & Jairo Rocha

Authors
  1. Carmen Hermida-Carrera
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Mario A. Fares
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Marcel Font-Carrascosa
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Maxim V. Kapralov
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Marcus A. Koch
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Arnau Mir
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Arántzazu Molins
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Miquel Ribas-Carbó
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. Jairo Rocha
    View author publications

    You can also search for this author in PubMed Google Scholar

  10. Jeroni Galmés
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

JG conceived the project; CH-C and JG designed the experiment; CH-C, MF-C, AMir, AMolins, MR-C. and JR performed the experiments; CH-C, MAF., MVK, MA.K, AMir, JR and JG analyzed the data and wrote the article. All authors have read and approved the manuscript.

Corresponding author

Correspondence to Jeroni Galmés.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Additional file 1: Table S1.

List of Orchidaceae and Bromeliaceae species downloaded from GenBank, accession number, leaf carbon isotope composition (δ13C, ‰), the ratio of leaf fresh mass to dry mass (FW/DW), the leaf thickness (mm), the leaf mass per area (LMA, g m− 2) and the habitat preference according to [53].

Additional file 2: Table S2.

Mean ± S.E. (n = 4) for the leaf mass per area (LMA), the leaf thickness and the leaf fresh to dry weight ratio (FW/DW) of C3, weak CAM and strong CAM for orchids and bromeliads. Values for the individual species are shown in Table 1 and Additional file 1: Table S1. Different letters denote statistically significant differences among metabolic types through Duncan test (p < 0.05).

Additional file 3: Table S3.

Coevolving groups of residues detected within the L-subunit of Rubisco within orchids and bromeliads.

Additional file 4: Table S4.

Integrative view of the Rubisco L-subunit variable sites under positive selection, coevolving and resolved with DT model as a function of external variables (δ13C and habitat preference).

Additional file 5: Figure S1

. Decision trees (DT) resolved for each Rubisco L-subunit variable site as a function of the external variables leaf δ13C (‰) and habitat preference for the orchids dataset.

Additional file 6: Figure S2

. Decision trees (DT) structure resolved for each variable site as a function of the external variables leaf δ13C (‰) and habitat preference based on the bromeliads dataset.

Additional file 7: Excel S1

. Rubisco L-subunit variable sites among orchids. Dots: the same amino acid as the first species in the list. Dash: not available.

Additional file 8: Excel S2

. Rubisco L-subunit variable sites among bromeliads. Dots: the same amino acid as the first species in the list. Dash: not available.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Hermida-Carrera, C., Fares, M.A., Font-Carrascosa, M. et al. Exploring molecular evolution of Rubisco in C3 and CAM Orchidaceae and Bromeliaceae. BMC Evol Biol 20, 11 (2020). https://doi.org/10.1186/s12862-019-1551-8

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12862-019-1551-8

Keywords

BMC Ecology and Evolution

ISSN: 2730-7182

Contact us