Evolution of phytochemical diversity in Pilocarpus (Rutaceae)

Phytochemistry 163 (2019) 132–146

Contents lists available at ScienceDirect

Phytochemistry journal homepage: www.elsevier.com/locate/phytochem

Evolution of phytochemical diversity in Pilocarpus (Rutaceae) a,∗

b

c

c,d

Daniella M. Allevato , Milton Groppo , Eduardo Kiyota , Paulo Mazzafera , Kevin C. Nixon

a

T

a

Cornell University, L.H. Bailey Hortorium, Section of Plant Biology, School of Plant Sciences, Cornell University, Ithaca, NY, USA USP Ribeirão Preto, Departamento de Biologia, Faculdade de Filosofia Ciências e Letras de Ribeirão Preto, Universidade de São Paulo, Ribeirão Preto, Brazil UNICAMP, Departamento de Biologia Vegetal, Instituto de Biologia, Universidade Estadual de Campinas, Brazil d Escola Superior de Agricultura Luiz de Queiroz, Departamento de Produção Vegetal, Universidade de São Paulo, Piracicaba, Brazil b c

ARTICLE INFO

ABSTRACT

Keywords: Pilocarpus spp. Rutaceae Jaborandi Chemosystematics Chemotaxonomy Phylogenetic signal Evolution Chemical pathways Coumarin Imidazole alkaloid

The evolution of phytochemical diversity and biosynthetic pathways in plants can be evaluated from a phylogenetic and environmental perspective. Pilocarpus Vahl (Rutaceae), an economically important medicinal plant in the family Rutaceae, has a great diversity of imidazole alkaloids and coumarins. In this study, we used phylogenetic comparative methods to determine whether there is a phylogenetic signal for chemical traits across the genus Pilocarpus; this included ancestral reconstructions of continuous and discrete chemical traits. Bioclimatic variables found to be associated with the distribution of this genus were used to perform OLS regressions between chemical traits and bioclimatic variables. Next, these regression models were evaluated to test whether bioclimatic traits could significantly predict compound concentrations. Our study found that in terms of compound concentration, variation is most significantly associated with adaptive environmental convergence rather than phylogenetic relationships. The best predictive model of chemical traits was the OLS regression that modeled the relationship between coumarin and precipitation in the coldest quarter. However, we also found one chemical trait was dependent on phylogenetic history and bioclimatic factors. These findings emphasize that consideration of both environmental and phylogenetic factors is essential to tease out the intricate processes in the evolution of chemical diversity in plants. These methods can benefit fields such as conservation management, ecology, and evolutionary biology.

1. Introduction Plants produce a vast assortment of chemical compounds ranging from primary metabolites, necessary for plant growth and survival, to plant specialized metabolites (PSM) important in growth and plant defense. Many studies have found that closely related species have highly conserved and similar metabolites (Bisby et al., 1980; Harborne, 1984; Hegnauer, 1994; Thaler et al., 2012; Rønsted et al., 2012). The similarity of PSM chemical structures in related species is what led to the development of a classification system known as chemotaxonomy or chemosystematics (C. De S. Abbott, 1886; Penfold, 1927; McNair, 1935; Hegnauer, 1994). However, various studies have shown that closely related species can also have divergent chemical traits; this could be due to disruptive selection in new environments (Latta and Linhart, 1997; Geber and Griffen, 2003; Rueffler et al., 2006), competition in the same environment (Mraja et al., 2011), overdispersion of species in the same community (Becerra 2015) or convergent evolution (Fraenkel, 1959; Pichersky and Lewinsohn, 2011).

There are many genetic mechanisms that might increase the diversity of plant specialized metabolites, from the simplest point mutation mechanisms (Dixon, 2001) to gene duplication. This duplication can revolutionize PSMs through neofunctionalization producing new PSMs (Zhang, 2003), sub-functionalization providing for specialization of a pathway (Lynch and Conery, 2000), and in the case of allopolyploidy (the combination of two different genomes) providing a diverse array of PSMs. Recent studies depict the role of small RNA's and epigenetics in chemical diversity; however, our understanding remains fragmented (Boyko and Kovalchuk, 2008, 2011; Urano et al., 2010; Parent et al., 2012). Additionally, the promiscuity of enzymes in plants (Pandya et al., 2014) has facilitated the production of a vast array of compounds without the need for highly specific enzymes (Firn and Jones, 2000, 2003; Pichersky and Lewinsohn, 2011). Variation in chemical diversity can be assessed qualitatively and quantitatively, and several studies have found that plants appear to adapt and modify the chemistry expressed in response to both time and space, encompassing many ecological roles (Moore et al., 2014). The

Corresponding author Cornell University, L.H. Bailey Hortorium, Section of Plant Biology, School of Plant Sciences, Cornell University, Ithaca, NY, USA. E-mail addresses: [email protected] (D.M. Allevato), [email protected] (M. Groppo), [email protected] (E. Kiyota), [email protected] (P. Mazzafera), [email protected] (K.C. Nixon). ∗

https://doi.org/10.1016/j.phytochem.2019.03.027 Received 23 October 2018; Received in revised form 29 March 2019; Accepted 30 March 2019 Available online 09 May 2019 0031-9422/ © 2019 Elsevier Ltd. All rights reserved.

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

complexity of alkaloid biosynthetic pathways with stereospecific roles indicates high levels of specialization of alkaloids to their environment (Waller and Nowacki, 1978; Scott, 1980). In addition, studies on Heracleum and Ruta (Rutaceae) have shown that variation in temperatures and seasonality can also affect the quantities of coumarins produced in these plants (Zobel and Brown, 1990, 1995; Zobel et al., 1991). Though many chemical ecology studies have found strong associations between chemical traits and bioclimatic factors, species are non-independent, due to phylogenetic relationships, and therefore it is important to account for the dependence of chemical traits on the shared evolutionary history between species. Chemotaxonomists have extensively studied patterns of plant specialized metabolites to infer relationships among taxa in Rutaceae, a family with 1000's of terpenes, coumarins, and alkaloids (Waterman, 1975, 1990; Gray and Waterman, 1978; Mabry and Ulubelen, 1980; Waterman and Hussain, 1983; Silva et al., 1988). In fact, Rutaceae is the plant family with the richest diversity of alkaloid classes including: acridines, furoquinolines, quinolines, quinazolines, indoloquinazolines, canthinones, benzylisoquinolines, and aromatic amides and amines (McNair, 1935; Waterman, 1975; Price et al., 1989). Of these classes, four are essentially limited to Rutaceae: imidazoles, furoquinolines, acridines, and indoloquinazolines. Major chemical markers studied in Rutaceae have also included coumarins, phenylpropanoids synthesized from the amino acid phenylalanine, due to a great diversity of furanocoumarins and pyranocoumarins (Gray and Waterman, 1978; Harborne, 1984, 1993; Murray, 2002). Pilocarpus Vahl, a genus of Rutaceae in the subtribe Pilocarpinae (tribe Galipeae), is known for its high bioactivity due to the presence of many alkaloids, terpenoids and coumarins (Santos and Moreno, 2004). Pilocarpus comprises 16 species and has a widespread distribution from southern Mexico and the Antilles, to northern Argentina (Skorupa, 1996). Of the members in the subtribe Pilocarpinae, three of the four genera have the typical furoquinolines found throughout Rutaceae, as well as a few with acridines and indolequinazolines (Silva et al., 1988). Many studies that have looked at the metabolites in the Rutaceae plant family, have found a progression (from basal species to derived species) of the prominence of certain compounds. Within Rutaceae, the prominence of alkaloids in basal lineages differs compared to the more derived lineages where there is a decline of alkaloids correlated with a rise in coumarins (oxygen derivatives), and subsequently a rise in flavonoids, lignans, limonoids, quassinoids and triterpenes (Silva et al., 1988). This could possibly be due to nitrogen economy or the evolution

of complex pathways involving many precursors. Pilocarpus has many coumarins and terpenes, none of the ubiquitous alkaloids found in the other Rutaceae, but instead has developed a unique set of imidazole alkaloids from histidine. Though there have been many studies on these imidazole alkaloids, their biosynthesis is not fully understood (Sawaya et al., 2010). It is possible that abiotic/biotic pressures led to the development of this unique set of compounds, necessary after a loss of the ancestral Rutaceae alkaloids in the genus. In this study, a phylogenetic perspective was used to gain a better understanding of the evolution of imidazole alkaloids and coumarins in species of Pilocarpus collected in Brazil (Fig. 1). We predicted that imidazole alkaloids would have a greater association with phylogenetic relationships than coumarins since imidazole alkaloids are only found in this genus. In addition, based on published work on Pilocarpus and other species in Rutaceae, we predicted that bioclimatic factors have shaped the evolutionary history of chemical traits in Pilocarpus. To test these hypotheses we 1) assessed the alkaloid and coumarin variation in the genus, evaluating the effect of phylogenetic relationships on compound presence and quantitative variation; 2) evaluated the relationship of chemical traits with bioclimatic factors to determine whether particular alkaloids and coumarins are associated with specific bioclimatic factors. Comparison of models with and without phylogenetic history helped to assess whether environmental convergence or evolutionary history could best predict chemical patterns in the genus. 2. Results and discussion 2.1. Phylogeny of Pilocarpus All three methods (Parsimony, Bayesian, Maximum Likelihood) confirmed that Pilocarpus species are monophyletic, and the genus is made up of two major clades (Fig. 2). These two major clades are geographically defined: Clade 1 (P. microphyllus Stapf ex Wardlew., P. trachyllophus Holmes, P. carajaensis Skorupa, P. jaborandi Holmes, P. racemosus Vahl, P. peruvianus (J.F. Macbr.) Kaastra) present in the tropical northern region of Brazil as well as Central America, and Clade 2 (P. grandiflorus Engl., P. riedelianus Engl., P. sulcatus Skorupa, P. spicatus A. St.-Hil., P. giganteus Engl., and P. pauciflorus A. St.-Hil.) present in mid-coastal to southern regions of Brazil. In addition, these two clades are subtended by P. pennatifolius Lem., which is determined to be the first branch of the Pilocarpus clade according to both Bayesian and Maximum Likelihood estimates. However, the Parsimony analysis

Fig. 1. Pilocarpus species collection localities. Dots on the map represent collection sites. The colors of the dots correspond to the species on the phylogeny to the right, with some species collected at multiple localities. *This phylogeny was estimated using Bayesian methods described in this paper. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.) 133

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

Fig. 2. Phylogeny of Pilocarpus estimated using three methods. (a) Parsimony in tnt, (b) Maximum Likelihood in RaxML, and (c) Bayesian in BEAST. The red star in (c) refers to the node where two fossils of Ptelea, one from the mid-Eocene and one from mid-Miocene, were used to date the Bayesian tree in BEAST. The calibrated phylogeny suggests a Miocene diversification of Pilocarpus (in yellow). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.) 134

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

resulted in the placement of both P. pennatifolius and P. peruvianus as the first branch. This uncertainty can be confirmed in both the Bayesian and Maximum Likelihood analyses, which portray less certainty (via posterior probabilities and bootstrap support) when placing P. peruvianus at the base of Clade 1 (Clade with P. microphyllus) (Fig. 2). There has only been one other study on the phylogeny of Pilocarpus, and it was estimated using Bayesian methods and with the following genes: trnG-S, 5.8S, ITS 1 and 2 (Oliveira, 2008). In this study, P. peruvianus was placed at the base, and P. pennatifolius was nested within the monophyletic clade of Pilocarpus, as a sister clade of P. spicatus (Oliveira, 2008). For our comparative analyses we chose to use the branch lengths from our Bayesian phylogenetic tree, as we calibrated the tree using Ptelea fossils (Fig. 2c).

(Clade 2 = P. grandiflorus, P. riedelianus, P. sulcatus, P. spicatus, P. giganteus, and P. pauciflorus). Interestingly, the coumarin diversity results are reversed on the phylogeny (Fig. 3b). There are still two major clades: one with high coumarin diversity (Clade 2 = P. grandiflorus, P. riedelianus, P. sulcatus, P. spicatus, P. giganteus, and P. pauciflorus) and the other with lower coumarin diversity (Clade 1 = P. microphyllus, P. trachyllophus, P. carajaensis, P. jaborandi, and P. racemosus). However, we now see that Clade 1 (the clade with the greatest alkaloid diversity) has the lowest coumarin diversity, and Clade 2 (the clade with the lowest alkaloid diversity) has the greatest coumarin diversity. In addition, P. pennatifolius at the base contains very high numbers of coumarins and the more derived species have the opposite of escalation, a decline in numbers of coumarins, especially for Clade 1. Phylogenetic signal, a statistical measure assessing phylogenetic dependence, was calculated with Blomberg's K and we found that it was not significant for alkaloid diversity or coumarin diversity (p-value > 0.05). This implies that evolutionary history was not the main contributor for the diversity of compounds, and instead there are other factors affecting the diversity of compounds. One thing to keep in mind is that these two clades are regionally distinct; Clade 1 is mostly found in northern Brazil, the Amazon, and Central America, whereas Clade 2 is mostly found in the southern and eastern regions of Brazil. This distinction is important because as these species dispersed and developed into these different regions they were exposed to a variety of

2.2. Diversity of coumarins and alkaloids across Pilocarpus has an inverse relationship The diversity (total number) of compounds among the species in this clade for both alkaloids and coumarins was reconstructed on the phylogeny using Maximum Likelihood (Fig. 3a total alkaloids, Fig. 3b total coumarins). The color-coded branches in Fig. 3 indicate low diversity (yellow) to high diversity (blue) of total compounds. In Fig. 3a, two major clades are observed: one with a higher diversity of alkaloids (Clade 1 = P. microphyllus, P. trachyllophus, P. carajaensis, P. jaborandi, and P. racemosus) and the other with a lower diversity of alkaloids

Fig. 3. Chemical diversity (total number of compounds) reconstructed across the phylogeny of Pilocarpus is reversed when comparing (a) alkaloid diversity with (b) coumarin diversity. The trait value gradient bar depicts the number of compounds present for that taxon ranging from a low (yellow) to high (blue) number of compounds. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.) 135

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

environmental pressures. Evolutionary conservation of naturally occurring compounds comes with a metabolic cost (i.e. carbon), and this cost is a tradeoff that can lead to decreases in other processes such as plant growth (Nitao and Zangerl, 1987). Therefore, in new environments directional selection of adaptive evolution could lead to a reduction of compounds that were needed in an old environment, as they have a high cost for no reward, and instead lead to an increase in production of beneficial compounds for the new environment, as they improve fitness and survival (Olson-Manning et al., 2012). One interesting example of possible convergence is the similarity of coumarin and alkaloid diversity in P. sulcatus and P. trachyllophus, two sympatric species growing in the Caatinga of Brazil. It is also important to note that P. pennatifolius, at the base of the Pilocarpus clade, has a higher diversity for both alkaloids and coumarins.

2.4. Reconstruction of pilocarpine reveals a species with the greatest concentration The discrete (stochastic character mapping of presence/absence) and continuous (Maximum Likelihood analysis of quantitative variation) reconstructions of the important medicinal compound pilocarpine (A4) are visualized in Fig. 5. Fig. 5a depicts the presence/absence of pilocarpine, establishing its presence in all but two species: P. spicatus and P. giganteus. The pilocarpine trait is lost at the orange “/” on the phylogeny and later regained at the black “/“, this is further confirmed by the marginal probabilities at each node. Next, in our continuous trait reconstructions we determined that P. sulcatus had the greatest concentration of pilocarpine, followed by P. trachyllophus and then P. grandiflorus (Star symbol indicates P. sulcatus in Fig. 5b). It is intriguing that the majority of other species in the same clade as P. sulcatus (Clade 2) have some of the lowest values for concentrations of pilocarpine, including the absence of pilocarpine itself. As we see in the figures, pilocarpine concentration is not conserved on the phylogeny and appears to be a divergent trait. In this case, there must be other factors affecting pilocarpine expression. A comparison of the distributions of species in Clade 2 confirmed that the distribution of P. sulcatus is distinct. Pilocarpus sulcatus is the only species in Clade 2 that has a limited distribution in the Caatinga vegetative areas of southern Bahia and northern Minas Gerais, and it can be found occurring alongside P. trachyllophus from Clade 1 (Skorupa, 1996). Although P. sulcatus is one of the smallest shrubs in the genus (height = 1–2 m), its growth habit is unique with a low density of branches and simple leaf clusters present at long internodes. In addition, P. sulcatus differs from other species in Pilocarpus as it has very conspicuous veins, unique perforated pollen, reddish leaf trichomes, and few pellucid dots. This is the first instance that this species has been assessed for imidazole alkaloids, and this study found that P. sulcatus has the largest concentration of pilocarpine compared to the extractions of the 11 other species in this study.

2.3. Reconstruction of coumarin biosynthetic diversity across Pilocarpus Specialized metabolites can also be considered non-independent traits, since they can be part of the same biosynthetic pathway or network. Although the pathway of imidazole alkaloids is unknown, the relationship between coumarin compounds and pathway enzymes is mostly known or approximated (Bourgaud et al., 2006). Therefore, we were able to calculate the Sørensen distance between coumarin compounds for each species, taking into consideration the proportion of shared enzymes to account for shared biosynthetic origins, in order to determine whether there was a large correlation between coumarin compounds (Sørensen, 1948; Junker et al., 2017). Fig. 4 depicts the reconstruction of biosynthetic diversity using the mean of Sørensen distances for each pairwise compound comparison. Clade 2 mostly depicts larger Sørensen distances or few shared enzymes in the pathways for compounds present in Clade 2 species. The majority of species in Clade 1 have small Sørensen distances; therefore, they have more shared enzymes and possibly a greater correlation between compounds present. Overall the majority of species in the Pilocarpus clade have larger mean Sørensen distance values, representing fewer shared biosynthesis enzymes and a lower correlation of coumarin chemical traits (Sørensen, 1948). A larger Sørensen distance is important because the independence of traits is an important consideration for comparative methods. In addition, the mean Sørensen distance values for these species are approximately in the same range as the values found in previous studies for the correlation of chemical traits in plant species (Sørensen, 1948; Junker et al., 2017). When we calculated the phylogenetic signal of the mean Sørensen distance we found that it was not significant (λ = 7.58−0.05 and p-value > 0.05), thus there is a more random distribution of biosynthetic diversity across the genus.

2.5. Reconstruction of chemical traits across Pilocarpus In Fig. 6, the discrete character reconstructions for three specific compounds are visualized. In Fig. 6a coumarin (C25) was detected in all of Clade 2 with one subsequent loss on the P. spicatus branch. This is in contrast to Clade 1, which appears to lose this trait at the base of Clade 1, though it is regained by P. trachyllophus. As C25 is the precursor to a variety of coumarins and furanocoumarins in the genus, it is possible that its disappearance in Clade 1 could be due to its role in forming other specialized coumarins further down the pathway. If there were selection for downstream products, the precursor (C25) would be Fig. 4. Reconstruction of the biosynthetic diversity or similarity of coumarin biosynthetic enzymes in Pilocarpus. The trait value gradient represents the mean of the Sørensen distances for pairwise compounds present in each species, and this exact value is labeled to the right of each taxon name. A small Sørensen distance (yellow) signifies small differences in enzymes, therefore a greater amount of shared enzymes in biosynthesis. A greater Sørensen distance (blue) signifies few or no shared enzymes in the biosynthesis. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

136

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

Fig. 5. Pilocarpine trait reconstructions. (a) Stochastic character mapping (SCM) of presence (black)/absence (orange) for pilocarpine. Orange “/” marks loss of the trait and black “/” the trait is gained. (b) Maximum Likelihood analysis of quantitative variation in pilocarpine, with the star symbol referring to the species with the greatest pilocarpine concentration. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

utilized to form other intermediates and compounds along the path, leading to a concentration that is too low to be detected. Since these two major clades are also geographically restricted, it is not possible to determine if these differences are due to random mutation, environmental specialization, or phylogenetic relationships. Next in Fig. 6b for citropten (C27), there appears to be two losses and two gains in each clade. Finally, in Fig. 6c the reconstruction of the coumarin osthol (C29) appears to have a somewhat random distribution across the genus, with multiple gains and losses. Although some compounds were not detected in certain species, it is possible that these species could still produce the compounds. This lack of detection could be due to a variety of reasons including low compound concentrations, developmental mechanisms, and environmental plasticity. All of the continuous chemical coumarin traits were reconstructed on the phylogeny and are visualized in Supplemental Fig. 1. Interestingly it appears that C25, C27, and C29 have a similar

distribution in Clade 1, being absent in all but one species, P. trachyllophus. One possibility is that P. trachyllophus needed to adapt to the precipitation and temperature extremes, through the expression of more coumarins for defense (Zobel and Brown, 1995). This could be tested with a common garden experiment, verifying if P. trachyllophus continued to produce these compounds in a tropical humid environment. The other species in Clade 1 could also be tested to evaluate whether translocation to a drier environment could lead to expression of C25, C27, and C29. On the other hand, the absence of C25, C27, and C29 in Clade 1 species could be a trade-off. All the other species in Clade 1 can be found in tropical humid regions, facing a greater diversity of pathogens and herbivores. These biological pressures could lead to an increase of other beneficial defense compounds, which could be accompanied by a decrease of the C25, C27, and C29 coumarins.

137

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

that they are dispersed more randomly than expected across the phylogeny. Calculations of the D-metric for all other compounds were not significantly different for the probability that D = 0 or the probability that D = 1; therefore, no other assumptions can be made for the distribution of these compounds. It is likely that the coumarin compounds would benefit from a larger scale phylogenetic study examining the phytochemical signal in the sub-tribe, tribe, or even sub-family of Rutaceae. This expansion is not possible for the imidazole alkaloids as they are restricted to the Pilocarpus clade and thus can only be assessed in this manner. 2.7. Phylogenetic signal of quantitative variation in compounds and bioclimatic variables Blomberg's K and Pagel's λ were utilized to assess the phylogenetic signal of quantitative variables (Table 2). This included the concentration of referenced compounds and the mean values of the five most important bioclimatic variables for species in the genus. There is one compound that has phylogenetic signal for compound concentration. In terms of Blomberg's K, C27 has significant phylogenetic signal (K = 1.049, p = 0.049), and considering Pagel's λ (λ = 1.151, p = 0.005) it also has significant phylogenetic signal. When running all the variables under three models (Ornstein-Uhlenbeck, Brownian motion, and Early Burst) and comparing the AIC values, the Brownian motion model performs best for C25, C27, C28, with the lowest AIC value of C27 as further confirmation of the phylogenetic signal found. Assessing Blomberg's K for the bioclimatic variables, there is a significant phylogenetic signal for Annual Mean Temperature (AMT) Table 1 Phylogenetic signal of presence/absence of chemical traits using the D-metric. Compounds A1 A2 A3 Pilocarpine (A4) A5 A6 A7 A8 A9 A10 A11 A12 A13 A14 A15 A16 A17 A18 A19 A20 A21 A22

Fig. 6. Ancestral chemical trait reconstructions. Stochastic character mapping (SCM) of discrete presence (black)/absence (orange) for each trait; C25 (one gain and one loss), C27 (two gains and two losses), C29 multiple gains and losses which appears random. (a) C25 = Coumarin (b) C27 = Citropten (c) C29 = Osthol. . (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

2.6. Phylogenetic signal of the presence/absence of individual compounds in the genus

Scopoletin (C23) Psoralen (C24) Coumarin (C25) Xanthotoxin (C26) Citropten (C27) Imperatorin B (C28) Osthol (C29) Imperatorin: C30

To determine if the presence/absence of individual chemical compounds was a shared phylogenetic trait, we calculated the phylogenetic signal, a measure of statistical dependence of species traits based on shared phylogenetic relationships. Using the D-metric for discrete variables we found that significant phylogenetic signal varied greatly among the compounds. (Table 1). Two compounds, C23 and C28 were not analyzed, as they are present in all species of our study. Compounds with significant phylogenetic signal that are more phylogenetically conserved than one would expect under Brownian motion include: A2, A8, A16, A18, A19, and A20. A few compounds, A3, A17, C7, had values for the D-metric that were significantly greater than 1, implying

D Imidazole Alkaloids 2.225038 −2.074712 1.856038 0.5762607 −0.1749972 −1.366992 −2.097808 −2.464568 1.206666 1.29005 −2.545128 −0.0485798 1.07127 0.2765816 −2.028621 −2.659028 1.997286 −4.957525 −6.42181 −2.339238 1.162835 2.94775 Coumarins * 3.481022 −0.0189366 1.264286 0.8296914 * 2.240458 0.07397704

P (D = 1)

P (D = 0)

0.704 0.04 0.854 0.384 0.18 0.151 0.172 0.025 0.495 0.514 0.169 0.202 0.48 0.167 0.148 0.027 0.8 < 0.001 < 0.001 0.04 0.487 0.677

0.108 0.833 0.017 0.375 0.585 0.609 0.632 0.899 0.247 0.265 0.641 0.516 0.196 0.532 0.639 0.891 0.054 0.802 0.818 0.892 0.257 0.224

* 0.667 0.132 0.575 0.36 * 0.957 0.254

* 0.129 0.546 0.21 0.259 * 0.007 0.549

P = probability that D-metric is equal to 0 or 1. If P is greater than 0.05 then it is not significantly different from that value of D. When D = 1 there is no phylogenetic signal, and there is a random distribution. When D = 0 there is phylogenetic signal and the trait follows Brownian motion. Bolded values represent significant p-values. * indicates compounds that cannot be analyzed since present in all. 138

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

Table 2 Phylogenetic signal of quantitative variables. Variable

K

p-value

λ

p-value

OU Model AIC

BM Model AIC

EB Model AIC

0.436 0.482 0.549 0.819 0.465 1.049 0.716 0.576 0.503

0.747 0.67 0.492 0.077 0.694 0.049 0.214 0.417 0.623

0.00008 0.28208 0.00008 0.45094 0.00008 1.15059 0.00008 0.00008 0.10594

1 0.517 1 0.463 1 0.005 1 1 0.807

130.29 118.88 128.32 117.20 135.99 109.89 146.51 111.03 146.82

134.00 121.51 129.25 114.80 139.13 106.24 145.12 111.03 148.51

137.67 125.17 132.91 118.47 142.80 109.91 148.79 114.70 152.18

1.232 0.334 0.382 1.032 1.014

0.004 0.962 0.768 0.013 0.022

0.98770 0.00008 0.00008 0.9243483 0.7885843

0.042 1 1 0.049 0.265

51.15 183.89 169.12 156.54 148.98

47.49 190.42 173.36 153.19 145.95

51.07 194.09 165.45 156.74 146.55

Compound Concentration Pilocarpine: A4 Scopoletin: C23 Psoralen: C24 Coumarin: C25 Xanthotoxin: C26 Citropten: C27 Imperatorin B: C28 Osthol: C29 Imperatorin: C30 Bioclimatic Variables AMT AP PCQ PWQ PDQ

OU = Ornstein-Uhlenbeck; BM = Brownian motion; EB = Early Burst model; AMT = Annual Mean Temperature, AP = Annual Precipitation, PCQ = Precipitation Coldest Quarter, PDQ = Precipitation Driest Quarter, PWQ = Precipitation Wettest Quarter * Bolded values represent significant p-values.

[K = 1.232, p = 0.004], Precipitation in the Wettest Quarter (PWQ) [K = 1.032, p = 0.013], and Precipitation in the Driest Quarter (PDQ) [K = 1.014, p = 0.022] (Table 2). When considering Pagel's λ there is a significant phylogenetic signal only for AMT (λ = 0.988, p = 0 .042) and PWQ (λ = 0.924, P = 0.049), but no significant signal for PDQ (λ = 0.789, p = 0.265). This variation in PDQ significance could be because Blomberg's K has been found to have greater sensitivity in detecting smaller changes in phylogenetic signal when compared to other metrics such as Pagel's λ, and Abouheif's Cmean (Münkemüller Tamara et al., 2012).

significant predictive regressions were that of the concentration of C25 with Precipitation in the Coldest Quarter, followed by C26 with Annual Mean Temperature, and C25 with Annual Mean Temperature. These three regressions all have negative slopes so we expect higher values of the bioclimatic factor to lead to lower levels of the respective compounds. Studies assessing the impact of temperature on the variation in furanocoumarins (C26) and coumarins (C25) in species of Rutaceae have found that these compounds are associated with temperature and water stress (Zobel and Brown, 1995; Zandalinas et al., 2016). In addition, many studies have found that the compounds in the phenylpropanoid pathway (coumarins, flavonoids, and lignins) are very important for the survival and response of a plant to the environment (Fraser and Chapple, 2011). We next tested for phylogenetic signal in the residual error of the OLS regressions, this included both Pagels λ and Blomberg's K. Pagel's λ estimated significant phylogenetic signal of the residuals for only one regression (C27 ∼ PWQ), while Blomberg's K found non-random and significant phylogenetic signal for the residuals of (C27 ∼ PWQ), (C28 ∼ PWQ), (C27 ∼ PWQ), (C27 ∼ PCQ), (C29 ∼ PCQ), (CD ∼ PCQ), (C27 ∼ AP), (C25 ∼ AP), and (C27 ∼ AMT). This analysis also demonstrated that non-random phylogenetic signal of bioclimatic variables themselves did not signify that phylogenetic signal for the residuals of the chemical-bioclimatic trait regressions would be significant (e.g. PCQ and AP variables had no significant phylogenetic signal but did have chemical trait relationships that were significant, whereas PDQ had significant signal but did not have any significant signal for chemical trait relationships). It is interesting to note that the chemical-bioclimatic trait regressions for all bioclimatic variables (other than PDQ which did not have any significant regressions) with C27 displayed phylogenetic signal in the residuals. For the compounds with phylogenetic signal in the residuals we next analyzed the phylogenetic generalized least squares (PGLS) regression to determine if including phylogeny would provide a better fit. Comparison of OLS and PGLS AIC scores, revealed no significant variation for most regressions; however, PGLS had a better fit for (C28 ∼ PWQ) (AIC scores for all regressions OLS and PGLS are provided in Supplementary Table 5). From these analyses, we see that the relationships between phylogenetic relationships and climate vary greatly for each compound. It is important to incorporate both bioclimatic and phylogenetic history into our analyses to better understand the evolution of chemical traits in Pilocarpus.

2.8. Association between chemical and climatic traits of species The relationship between species climatic traits (calculated from a species mean along each bioclimatic variable) and compound traits (concentration or diversity/total number of compounds) was evaluated using Ordinary Least Squares (OLS) regressions run in the R package nlme. The predictive power of the chemical-bioclimatic trait regressions is compound specific and varies greatly. Overall there are many OLS regressions that have significant R2 values ranging from relationships with moderate (R2 = 0.55) to higher variability (R2 = 0.13) (Table 3), with some regressions providing greater precision than others. Considering each bioclimatic trait we see: a) Annual Mean Temperature (AMT) has a significant association with C26 (R2 = 0.47, p ≤ 0.001), C25 (R2 = 0.36, p ≤ 0.05), C27 (R2 = 0.21, p ≤ 0.05), C30 (R2 = 0.14, p ≤ 0.05), and coumarin diversity (R2 = 0.24, p ≤ 0.05) (Table 3); b) Annual precipitation (AP) has a significant association with C26 (R2 = 0.20, p ≤ 0.05), C29 (R2 = 0.25, p ≤ 0.05), C30 (R2 = 0.13, p ≤ 0.10), and coumarin diversity (R2 = 0.16, p ≤ 0.10) (Table 3); c) Precipitation in the Coldest Quarter (PCQ) has a significant association with C25 (R2 = 0.55, p ≤ 0.001), alkaloid diversity (R2 = 0.24, p ≤ 0.01), coumarin diversity (R2 = 0.37, p ≤ 0.01), C30 (R2 = 0.19, p ≤ 0.05), C23 (R2 = 0.23, p ≤ 0.05), C29 (R2 = 0.17, p ≤ 0.10), and C26 (R2 = 0.13, p ≤ 0.10) (Table 3); d) Precipitation in the wettest quarter (PWQ) has a significant association with C28 (R2 = 0.32, p ≤ 0.01), C24 (R2 = 0.29, p ≤ 0.01), C23 (R2 = 0.21, p ≤ 0.05), C25 (R2 = 0.20, p ≤ 0.05), alkaloid diversity (R2 = 0.13, p ≤ 0.10), and coumarin diversity (R2 = 0.10, p ≤ 0.10) (Table 3); e) Precipitation in the driest quarter (PDQ) did not have any significant R2 values, and therefore these regressions cannot predict any of the compound traits (Table 3). In addition, we were unable to obtain significant regressions between A4 and the bioclimatic factors assessed in this paper. The most 139

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

Table 3 Phylogenetic signal of residuals and OLS chemical-climatic trait regressions. OLS Phylogenetic Signal of Residuals

AMT

K

λ

Adj R2

Y-int

S.E.

p

Slope

S.E.

p

AIC

A4 C23 C24 C25 C26

0.43 0.42 0.62 0.56 0.39

7.61–5 7.61–5 7.61–5 7.61–5 7.61–5

−0.09 0.11 0.06 0.36* 0.47***

102.1 −120 218.5 243.2 529.1

166 95.3 143 79.4 147

0.55 0.23 0.15 0.01 0.005

−2.45 6.35 −7.97 −9.17 −20.7

7.15 4.09 6.16 3.41 6.35

0.73 0.15 0.22 0.02 0.008

114.77 103.61 111.78 99.95 112.37

C27

0.74•

7.61–5

0.15• PGLS

182.3 120.5

113 22.7

0.16 0.37

−7.16 −4.16

5 5.49

0.21 0.48

52.99 54.51

C28 C29 C30 AD CD

0.73 0.60 0.57 0.67 0.47

7.61–5 7.61–5 7.61–5 7.61–5 7.61–5

0.09 −0.08 0.14* −0.003 0.24*

544.2 43.0 590.8 −13.1 20.19

311 74.4 294 19.3 6.61

0.11 0.57 0.07 0.51 0.01

−19.2 −1.4 −21.3 0.81 −0.59

13.4 3.20 12.6 0.83 0.28

0.18 0.66 0.12 0.35 0.06

127.3 98.66 126.17 71.75 50.23

Phylogenetic Signal of Residuals K λ

AP Adj R

2

Y-int

S.E.

p

Slope

S.E.

p

AIC

A4 C23 C24

0.39 0.53 0.52

7.61–5 0.31 7.61–5

−0.02 −0.09 −0.07

87.6 16.9 10.5

48.5 31.8 46.1

0.10 0.60 0.83

−0.03 0.007 0.02

0.03 0.02 0.03

0.38 0.75 0.62

124.83 116.38 123.78

C25

0.73•

0.46

0.15 PGLS

76.6 72.8

27.5 24.5

0.01 0.01

−0.03 −0.03

0.02 0.01

0.11 0.09

113.48 114.03

C26

0.68

7.61–5

0.20*

151.2

54.4

0.01

−0.07

0.04

0.07

127.12

C27

1.05*

0.91

−0.2 PGLS

21.13 24.4

30.1 23.1

0.51 0.33

0.002 0.003

0.021 0.014

0.99 0.82

65.93 66.71

C28 C29 C30 AD

0.77 0.58 0.60 0.65

7.61–5 7.61–5 0.15 7.61–5

−0.06 0.25* 0.13• −0.03

45.4 49.3 239.1 1.31

101 18.7 89.2 5.93

0.66 0.02 0.02 0.83

0.04 −0.03 −0.09 0.003

0.06 0.01 0.06 0.004

0.60 0.05 0.13 0.44

139.62 105.74 136.99 82.77

CD

0.85•

0.51

0.16• PGLS

9.92 10.13

2.08 1.86

0.0007 0.0003

−0.002 −0.002

0.001 0.001

0.10 0.06

61.89 62.49

Phylogenetic Signal of Residuals K λ –5

PCQ Adj R2

Y-int

S.E.

p

Slope

S.E.

p

AIC

A4 C23 C24 C25 C26

0.42 0.49 0.54 0.468 0.54

7.61 7.61–5 7.61–5 7.61–5 7.61–5

0.62 0.23* −0.04 0.55*** 0.13•

54.7 6.22 45.8 59.2 82.4

22.1 11.8 20.2 8.8 25.3

0.03 0.61 0.04 0.0001 0.008

−0.03 0.07 −0.04 −0.09 −0.11

0.06 0.03 0.05 0.02 0.07

0.62 0.06 0.48 0.003 0.13

124.12 111.62 122.30 105.73 126.82

C27

1.05*

0.82

0.03 PGLS

21.6 28.15

9.3 12.7

0.04 0.05

−0.03 −0.02

0.02 0.02

0.26 0.27

106.75 106.99

C28 C29

0.63 0.79•

7.61–5 0.40

C30 AD CD

0.50 0.65 0.71

7.61–5 7.61–5 7.61–5

−0.06 0.17• PGLS 0.19* 0.24** 0.37**

120.7 23.4 30.4 159.1 1.85 8.19

45.1 8.7 10.9 38.1 2.25 0.80

0.02 0.02 0.01 0.001 0.43 1.3e-6

−0.07 −0.04 −0.06 −0.21 0.01 −0.006

0.12 0.02 0.02 0.11 0.006 0.002

0.57 0.10 0.01 0.08 0.05 0.02

138.24 105.52 105.47 134.97 78.41 57.75

Phylogenetic Signal of Residuals K λ A4 C23

0.44 0.38

–5

7.61 7.61–5

PDQ Adj R2

Y-int

S.E.

p

Slope

S.E.

p

AIC

−0.10 0.05

47.1 40.2

21.4 12.6

0.05 0.01

−0.01 −0.10

0.13 0.08

0.91 0.24

122.78 112.25

140

(continued on next page)

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

Table 3 (continued) C24 C25 C26

0.56 0.67 0.49

7.61–5 7.61–5 7.61–5

−0.04 −0.03 −0.10

21.6 21.76 51.41

19.3 12.9 27.2

0.29 0.12 0.08

0.09 0.06 −0.02

0.12 0.08 0.17

0.48 0.44 0.88

120.72 112.69 127.53

C27

0.86•

7.61–5

−0.19 PGLS

18.40 26.7

13.9 16.0

0.24 0.15

0.02 0.02

0.08 0.06

0.82 0.80

63.36 63.80

C28 C29 C30 AD CD

0.75 0.63 0.52 0.61 0.45

7.61–5 7.61–5 0.18 5.53–5 7.61–5

−0.08 0.00 −0.09 −0.09 −0.002

81.4 17.7 105.8 6.18 5.54

43.4 9.2 42.6 2.60 0.97

0.09 0.08 0.03 0.04 0.0002

0.132 −0.058 −0.069 −0.002 0.006

0.27 0.05 0.27 0.016 0.006

0.64 0.34 0.806 0.90 0.34

136.88 105.81 136.51 80.56 60.87

Phylogenetic Signal of Residuals K λ

PWQ Adj R

2

Y-int

S.E.

p

Slope

S.E.

p

AIC

A4

0.43

7.61–5

0.03

53.3

34.5

0.15

−0.02

0.08

0.80

123.58

C23 C24 C25 C26

0.30 0.49 0.53 0.42

7.61–5 7.61–5 7.61–5 7.61–5

0.21* 0.29** 0.20* −0.06

61.8 −23.7 −3.3 25.8

18.6 25.5 18.3 43.2

0.007 0.37 0.86 0.56

−0.09 0.16 0.09 0.06

0.04 0.07 0.05 0.11

0.07 0.03 0.07 0.59

111.19 117.60 110.91 128.10

C27

0.98*

1.15•

−0.16 PGLS

32.7 37.8

30.3 24.2

0.33 0.17

−0.03 −0.02

0.07 0.05

0.71 0.69

63.34 64.02

C28

0.84•

7.61–5

0.32** PGLS

−31.7 −86.6

55.4 53.3

0.58 0.13

0.36 0.44

0.14 0.00

0.03 0.001

133.09 131.77

C29 C30 AD CD

0.57 0.41 0.53 0.42

7.61–5 7.61–5 7.61–5 7.61–5

−0.05 0.03 0.13• 0.10•

19.3 26.5 11.5 4.26

15.2 65.6 3.73 1.48

0.23 0.69 0.01 0.01

−0.03 0.20 −0.015 0.005

0.04 0.16 0.01 0.003

0.53 0.27 0.13 0.16

107.20 136.16 79.14 60.65

AP = Annual Precipitation; AMT = Annual Meant Temperature; PCQ = Precipitation Coldest Quarter; PDQ = Precipitation Driest Quarter; PWQ = Precipitation Wettest Quarter; A4 = Pilocarpine concentration; AD = alkaloid diversity; CD = coumarin diversity; C23 = Scopoletin; C24 = Psoralen; C25 = Coumarin; C26 = Xanthotoxin; C27 = Citropten; C28 = Imperatorin-B; C29 = Osthol; C30 = Imperatorin; Pagel's λ and Blomberg's K: phylogenetic signal of residuals; Adj R2 adjusted R-squared; Y-int: Y intercept; S.E. standard error of intercepts and slopes, bolded p represents significant p value. A lower AIC (Akaike information criterion) value represents a better model fit. Signif. Codes for Adj R2: 0.001 ‘***’ 0.01 ‘**’ 0.05 ‘*’ 0.1 ‘•’

3. Conclusions

variation in development, silent metabolism, or variation at a larger taxonomic scale (Wink, 2003; Carrari et al., 2006; Pichersky and Lewinsohn, 2011). Further in-field and greenhouse studies would help to assess these possible hypotheses. In the case of the coumarins, additional research can be done to expand the detection of coumarins more broadly within Rutaceae or to the sub-family level within Rutoideae, as this would potentially provide greater power to detect phylogenetic signal in order to identify evolutionary patterns of coumarins. Although we found phylogenetic signal for the presence data of some imidazole alkaloids, expanding the sample of genera to increase significance is not possible in this case as imidazole alkaloids are only present at the small scale of the genus Pilocarpus. Of the quantitative traits, there was significant phylogenetic signal for annual mean temperature, precipitation in the wettest quarter, and precipitation in the driest quarter, indicating that related species could share similar niches. In terms of the quantitative chemical traits, only the concentration of C27 had significant non-random phylogenetic signal, and it was the only compound trait in the OLS regressions that had phylogenetic signal in the residual error for almost all regressions with bioclimatic traits. However, we only found one chemical-bioclimatic trait regression was improved with the inclusion of phylogeny: C28 ∼ PWQ. Overall, quantitative variation in coumarin compound traits was most significantly attributed to environmental convergence and independent of phylogenetic relationships. Compounds such as alkaloids and coumarins have been shown to differ with regard to presence and quantity based on different

Our results suggest that both phylogenetic relationships and environmental factors have shaped compound variation in Pilocarpus. The genus is made up of two regionally distinct clades, Clade 1 (Northern Brazil, Amazon, Central America) and Clade 2 (Southern and Eastern Brazil), with each clade displaying a similar variation in chemistry and an opposite emphasis in compound class expression (Clade 1: high alkaloid diversity with low coumarin diversity, and Clade 2: high coumarin diversity with low alkaloid diversity). In addition, since the biosynthesis of coumarins is well known we were able to assess Sørensen distances, finding that Clade 1 has more shared conserved enzymes than Clade 2. Though most experts agree that related species share similar traits, it is also very possible that related species have divergent chemistry. This difference in chemistry could be driven by many factors including speciation due to expansion into different environments or differentiation into different niche/adaptive spaces in the same location (Mraja et al., 2011; Becerra, 2015). In this study we found that two sympatric species growing in the Caatinga biome of Brazil, P. sulcatus and P. trachyllophus, originated from different clades, but their chemical traits were more similar to each other than to their own clades. Pilocarpus sulcatus was also revealed to contain the highest concentration of pilocarpine, followed by P. trachyllophus and P. grandiflorus. A lack of correlation of specialized chemistry with phylogeny could be due to convergent evolution, induction by environmental pressures, 141

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

environmental conditions, such as changes in temperature, elevation, precipitation and soil type (Levin, 1976; Dixon and Paiva, 1995; Zobel and Brown, 1995; Martz et al., 2009). This compound production in response to the environment, is very important for the survival and evolution of plants (Fraser and Chapple, 2011). The advancement of phylogenetic comparative methods enables the discovery of relationships between the environment and phytochemical diversity. Together these complementary approaches have aided in a better understanding of the coumarins and alkaloids present in the genus Pilocarpus, as well as provided a foundation to inform studies on the biosynthesis of these compounds. Future work identifying the biosynthetic enzymes present in different species would further aid in understanding the evolution of these biosynthetic pathways. Rutaceae is a phylogenetically complex family with a rich diversity of chemistry. The combination of phylogenetic comparative methods, environmental factors, and biochemical data is beneficial for studies on the evolution of chemical diversity, ecological studies, and the conservation of bio-diverse areas.

centrifugation were pooled, filtered in 0.2 μm membranes, and analyzed. 4.3.2. Extraction protocol for coumarins Phenols were extracted from 200 mg leaf tissue in 2 ml aqueous 80% ethanol, vortexed for 1 min. The supernatant was recovered by centrifugation and filtered in 0.2 μm membrane (Durand-Hulak et al., 2015; Vialart et al., 2012). Each individual plant had a total of 3 biological replicates for each extraction. 4.4. UPLC-ESI-MS methods 4.4.1. Alkaloids Five microliters of each alkaloid extraction were injected into an Acquity UPLC coupled with a TQD triple-quadrupole mass spectrometer (Micromass-Waters, Manchester, England). Mass spectrometer conditions were: capillary 3.0 kV, cone 30 V, extractor 1 V, ion source temperature 150 °C, desolvation temperature 300 °C and column temperature 30 °C, in electrospray ionization in positive mode. The chromatographic column used was a Polaris 180 Å C18-A, 2.0 × 100 mm (Varian), 5 μm, pre-equilibrated with an ammonium acetate buffer, 10 mmol/L, pH 3.0 (solvent A) and acetonitrile (solvent B). The elution was carried out using a gradient ranging from 5 to 25% of solvent B in 8 min, then the solvent B concentration remained in 25% for 2 min and returned to 5% for 2 min, with a flow rate of 300 μL/min.

4. Experimental 4.1. General experimental procedures instruments standards purity LC-MS analyses were performed with an Acquity UPLC coupled to a TQD triple-quadrupole mass spectrometer (Micromass-Waters, Manchester, England). All standards were of analytical grade: coumarin standards (scopoletin, psoralen, coumarin, xanthotoxin, citropten, osthol, imperatorin) and the alkaloid standard (pilocarpine) were purchased from Sigma-Aldrich (St. Louis, MO). There are various statistical approaches to evaluating correlated traits in phylogenetics. Our analysis has included the quantification of phylogenetic signal using three metrics: Pagel's λ to determine the diversity in the genus, the D-metric to determine the signal of individual compounds, and Blomberg's K statistic to assess quantitative variation in compounds (Pagel, 1992; Blomberg et al., 2003). Next, ancestral state reconstruction was performed to better characterize evolutionary patterns. Finally, to assess the relationship of metabolites in the genus to extrinsic bioclimatic factors, OLS and PGLS regressions were evaluated.

4.4.2. Coumarins Coumarin sample extracts were diluted (1:4, v/v) in aqueous 80% ethanol and 4 μl of each sample were used for analysis. Mass spectrometer conditions were: capillary 3.0 kV, cone 35 V, extractor 1 V, ion source temperature 120 °C, desolvation temperature 350 °C and column temperature 30 °C, in electrospray ionization in positive mode. The chromatographic column used was a Waters Acquity C18-BEH (2.1 × 100 mm, 1.7 μm) column and the elution was carried out using Milli-Q water with 0.1% formic acid (solvent A) and acetonitrile (solvent B) in a gradient ranging from 10 to 100% of solvent B in 8 min. 4.5. UPLC-MS data scoring

4.2. Field collections of plant tissue

Chemical data on alkaloids and coumarins were scored for presence/absence, as well as concentration (Supplemental Fig. 2, Supplemental Table 2; Supplemental Table 3). Both alkaloid and coumarin samples were normalized by the exact dry mass of ground tissue extracted. In addition, each species had multiple individual replicates, as well as three biological replicates (Supplemental Table 1). For the qualitative analyses we considered a compound to be present for the entire species when the majority had the compound present. Of the 30 compounds we evaluated in this study, there were only four instances where a compound was absent from one conspecific individual (P. pennatifolius: C25 and C27; P. riedelianus: A2; P. spicatus: C27). For the alkaloid extracts, 22 alkaloids were identified for presence/absence data and pilocarpine, one of the imidazole alkaloids, was quantified (Table 4, Supplemental Table 2; Supplemental Table 3; Supplemental Figure 4a). Compound elution and m/z were used to identify the imidazole alkaloid using the method developed by Sawaya et al. (2008) and (2011). Peak areas for the alkaloid pilocarpine were determined from Total Ion Current Chromatograms (TIC) (Supplemental Figure 4a). Quantification of pilocarpine was done using a standard curve calibration with a pilocarpine standard (Sigma Aldrich: analytical grade) using a method developed by Sawaya et al. (2008) and (2011). As there are no other standards commercially available for the other imidazole alkaloids we are only presenting quantitative results for one alkaloid, pilocarpine. The standard calibration curve for pilocarpine was run with six calibrators ranging from 0.5 μg/mL – 100 μg/mL (i.e. 0.5, 1, 5, 10, 50 and 100 μg/mL). The Limit of Detection (LOD) and the Limit of Quantitation (LOQ) were calculated as LOD = 3.3σ/slope and

Field collected plants that were sampled include a total of 134 individuals (Fig. 1, Supplemental Table 1. Each individual plant underwent two types of extractions: coumarin and alkaloid. For each extraction there were three different leaf replicates, giving a total of 411 coumarin extractions and 411 alkaloid extractions. 4.3. Biochemical extraction methods To validate the biochemical extractions from field collected leaf tissue, a degradation study on two live specimens of Pilocarpus pennatifolius from the New York Botanical Garden was performed to examine degradation of alkaloids and coumarins in the leaves. Silica-dried leaf tissue on ice showed the least degradation and was not significantly more degraded relative to immediately frozen tissue. This method of preservation of leaf tissue for biochemical extractions has been used in many studies that require plant collections in remote areas, as the removal of water reduces enzymatic activity in the leaf tissue and the cold temperature helps to prevent enzymatic reactions due to handling (Kim and Verpoorte, 2010; Fine et al., 2013). 4.3.1. Extraction protocol for alkaloids Biochemical extractions were performed on approximately 10 mg of silica-dried leaf tissue following modifications on the method developed by Avancini et al. (2003). Leaf tissue was ground finely and extracted twice with 1 ml of aqueous 0.1% formic acid and placed in an ultrasonic water bath for 20 min. The supernatants recovered by 142

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

Table 4 Compounds identified in extractions. N°

Compound Name

Molecular Formula

[M+H]+ m/z

Rt (min)

Identification Source

A1 A2 A3 A4

3H-imidazole-4-ylmethyl)-3-methyl-5H-furan-2-one dehydropilocarpine Pilocarpidine/isopilocarpidine Pilocarpine/isopilocarpine

C9H10N2O2 C10H12N2O2 C10H14N2O2 C11H16N2O2

179 193 195 209

1.36 1.38 2.37 2.73

A5, A8, A11

Hydroxy-phenyl-methyl)-4-(3H-imidazole-4-ylmethyl)dihydro-furan-2-one Pilosin and isomers

C15H16N2O3

273

3.15, 4.1, 4.43

Sawaya et Sawaya et Sawaya et Sawaya et standard Sawaya et

C16H18N2O3

287

Sawaya et al. (2011)

3-Hydroxymethyl-4-(3-methyl-3H-imidazole -4-yl)-1phenylbutan-1-one 3-Benzoyl-4-(3-methyl-3H-imidazole-4-ylmethyl)dihydrofuran-2-one 3- Benzylidene-4-(3H-imidazole-4- ylmethyl)-dihydrofuran-2-one anhydropilosine 3-Benzyl-4-(3H-imidazole-4-ylmethyl) -dihydro-furan-2one Scopoletin Psoralen Coumarin Xanthotoxin Citropten Imperatorin B [271-C5H8] Osthol Imperatorin

C15H18N2O2

259

3.52, 3.52, 4.17, 4.17, 4.65 6.07, 6.12

C16H16O3N2

285

6.18

Sawaya et al. (2011)

C15H14N2O2

255

6.58, 7.49

Sawaya et al. (2011)

C16H16N2O2 C15H16N2O2

269 257

6.96, 8.95, 10.38 7.12, 7.49

Sawaya et al. (2011) Sawaya et al. (2011)

C10H8O4 C11H6O3 C9H6O2 C12H8O C11H10O4 C11H6O4 C15H16O3 C16H14O4

193 187 147 217 207 203 245 271

3.46 5.97 4.63 6.09 5.89 7.20 7.38 7.24

Analytical Analytical Analytical Analytical Analytical Analytical Analytical Analytical

A6, A7, A9, A10, A12 A13, A14 A15 A16, A20 A17, A21, A22 A18, A19 C23 C24 C25 C26 C27 C28 C29 C30

al. al. al. al.

(2011) (2011) (2011) (2011); Analytical

al. (2011)

Sawaya et al. (2011)

standard standard standard standard standard standard standard standard

-

LCMS LCMS LCMS LCMS LCMS LCMS LCMS LCMS

*N° refers to the identified compounds (A = Alkaloid, C=Coumarin).

LOQ = 10σ/slope. The regression equation for pilocarpine was y = 246995x + 26 with a coefficient of determination (R2) of 0.999, LOD of 11.9, and LOQ of 36.1. For the coumarin extracts, eight coumarins were identified and quantified by comparison with analytical grade standards and literature, comparing retention times and m/z fragments. (Masson et al., 2016) (Table 4). Peak areas of the coumarin compounds were determined from Selected-Ion Monitoring chromatograms (SIM) (Supplemental Figure 4b). Quantification of the coumarin concentrations was done through a regression analysis of coumarin reference standards (Sigma Aldrich: analytical grade). The standard calibration curves for the coumarin standards were run with six calibrators ranging from 0.1 μg/mL – 50 μg/mL (i.e. 0.1, 0.5, 1, 5, 10, and 50 μg/mL). A regression equation for scopoletin (C23) of y = 280644x + 712997 was acquired with a R2 of 0.995, LOD of 4.4, and LOQ of 13.4. A regression equation for psoralen (C24) of y = 333999x + 899532 was acquired with a R2 of 0.992, LOD of 5.9, and LOQ of 17.7. A regression equation for coumarin (C25) of y = 173342x + 158725 was acquired with a R2 of 0.999, LOD of 2.4, and LOQ of 7.3. A regression equation for xanthotoxin (C26) of y = 483498x + 638896 was acquired with a R2 of 0.999, LOD of 2.1, and LOQ of 6.4. A regression equation for citropten (C27) of y = 710179x + 26 was acquired with a R2 of 0.971, LOD of 10.6, and LOQ of 32.2. A regression equation for imperatorin B (C28) of y = 699690x + 36 was acquired with a R2 of 0.975, LOD of 9.7, and LOQ of 29.5. A regression equation for osthol (C29) of y = 439401x + 642153 was acquired with a R2 of 0.998, LOD of 2.5, and LOQ of 7.5. A regression equation for imperatorin (C30) of y = 39338x + 25777 was acquired with a R2 of 0.999, LOD of 4.3, and LOQ of 12.9.

Cornell University. DNA extractions were sent to University of Minnesota Genomics Center, Minneapolis, Minnesota, USA for enzyme optimization, library preparation, and RADseq/SBG Sequencing-based Genotyping. Enzymes used for digestion were PstI-HF and BtgI, and the samples were run on a Nextseq 500 High Output 150-bp SR run. Mean quality scores for all libraries were greater than or equal to Q30, indicating an accurate base call probability of 99.9%, and all barcodes were detected. 4.7. Phylogenetic methods The sample SR raw RADseq fastq files were demultiplexed and assembled de novo using ipyrad (Eaton, 2014; Eaton et al., 2017). This included quality filtering, clustering within samples, calculation of joint estimation of error and heterozygosity, consensus base calling, clustering across samples, and finally alignment. This concatenated SNP dataset was exported (tnt, phylip, nexus formats) for further analysis. Phylogenies were estimated using Parsimony in tnt, Bayesian in BEAST 2.4.8, and Maximum Likelihood in RAxML v. 8.0. (Nixon and Carpenter, 1993; Nixon, 1999; Drummond and Rambaut, 2007; Stamatakis et al., 2008; Goloboff et al., 2008). The outgroups used in the analyses are in Supplemental Table 4 (Groppo et al., 2012), and all phylogenetic trees were rooted with Ptelea as there were two Ptelea fossils available for the phylogenetic analysis (outgroups chosen according to procedures described in Nixon and Carpenter, 1993). Pilocarpus species collected in the field were used for both chemical and phylogenetic analyses, extracting as many individuals as possible. Extractions from herbarium specimens were done for species not collected in the field including P. peruvianus (J.F. Macbr.) Kaastra, P. alatus C.J.Joseph ex Skorupa, P. demerarae Sandwith, and P. manuensis Skorupa; however, RADseq library preparations for 3 of these species (P. alatus, P. demerarae, and P. manuensis) were not successful and were not used in this analysis. These species distributions include: P. alatus in the Northern Amazonian region, P. demerarae in Guyana, and finally P. manuensis in Peru and in the western Amazonian Brazilian state of Acre. The distribution of P. manuensis is similar to P. peruvianus, so it is possible that its inclusion in the phylogeny could assist with a better

4.6. DNA extraction and sequencing DNA extractions were done with a few modifications of Barry et al. (2005). Modifications include the addition of RNase A before incubation in the water bath, the use of ice-cold ethanol and ice-cold isopropanol. Genome size (2C = 5.05 pg) for optimization of sequencing was determined using Flow Cytometry in the Roeder Laboratory, 143

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

estimate of the Pilocarpus phylogeny. Parsimony analyses were run in TNT using the parsimony ratchet, drift, and tree-bisection-reconnection (TBR) branch swapping methods to estimate a consensus tree (Fig. 2a). RAxML v. 8.0 was used to estimate a maximum consensus tree using 10,000 bootstrap iterations (Fig. 2b) (Stamatakis et al., 2008). BEAST 2.4.8 analyses were run using the GTR + Γ +I model, implementing an uncorrelated lognormal relaxed molecular clock as well as a Yule and birth-death process of speciation (Gernhard, 2008). Two separate runs of 30 million MCMC generations (logged every 100) converged and the effective sample sizes (ESS) were assessed positively in Tracer. The independent runs were combined using LogCombiner with a 20% burn in, and then were used to estimate the maximum clade credibility tree in Tree Annotator (Drummond and Rambaut, 2007) (Fig. 2c). The BEAST tree was chosen for further phylogenetic comparative analyses since we could calibrate the tree using two Ptelea fossils. Ptelea paliuruoides, the oldest Ptelea fossil, is from the mid-Eocene and was collected in the “Middle Eocene Parachute Creek Member of the Green River Formation of Colorado and Utah” (Manchester and O'Leary, 2010). Ptelea enervosa, the second Ptelea fossil, is from the mid-Miocene and it was collected from localities on the Oregon-Idaho border (Call and Dilcher, 1995). A mid-Eocene fossil calibration was implemented as a minimum age for the Ptelea trifoliata node, (indicated by a red star on the phylogeny) using a log-normal distribution (Fig. 2c). Since P. peruvianus was not collected in the field and was an herbarium specimen extraction, we did not have any chemical data for this species. Therefore, we pruned it from the tree (leaving P. pennatifolius as the first branch), and P. peruvianus was thus not included in the rest of the analyses in this paper. The outgroups were also trimmed for subsequent phylogenetic comparative analyses using the ape package in R, and the branch lengths estimated with the full data set were saved for the subset of 12 taxa with phytochemical data.

compound present was calculated using Blomberg's K statistic, utilizing the ‘phylosig’ function in the phytools R package (Blomberg et al., 2003; Revell, 2012). K > 1 indicates that traits have more phylogenetic signal than expected under Brownian motion when K = 1. When K < 1 related species are less similar. To test whether K is significant, 1000 randomizations can be run to generate a null distribution and provide a p-value in phytools (Revell, 2012). In addition, Pagel's λ was also calculated to assess the phylogenetic signal of quantitative variation (Pagel, 1992). An assessment of evolutionary models was done through examination of AIC values for the following methods: OU = OrnsteinUhlenbeck; BM = Brownian motion, EB = Early Burst model. 4.9. Bioclimatic variables Bioclimatic variables (BIO1-BIO19) at 30-arcseconds were extracted from the WorldClim website (www.worldclim.org) using the longitude and latitude coordinates of each collection locality through the R package raster (Etten and Hijmans, 2010). The five most significant bioclimatic variables that differentiate niches in the genus were identified in an ecological niche modeling study of the genus (unpublished data Allevato). The mean value for each bioclimatic variable was calculated for each species to be used as a species bioclimatic trait. 4.10. Comparative analyses 4.10.1. Biosynthetic diversity To assess similarities of coumarin compounds due to the non-independence of pathways and shared enzymes, we first created a master matrix of enzyme presence-absence data to portray the enzymes used in the biosynthesis of each coumarin compound (Bourgaud et al., 2006) (Supplemental Figure 3). This master matrix was next used to create a new separate enzyme presence-absence matrix for each species that only included the coumarin compounds present for that species. Each species matrix was then used to calculate the Sørensen distance for each pairwise compound comparison. A small Sørensen distance represents a very similar pathway and many shared enzymes; whereas, a large Sørensen distance represents a different pathway with few or no shared enzymes. The mean of the Sørensen distances for each species was identified as the Biosynthetic Diversity of each species (Sørensen, 1948; Junker et al., 2017). We next reconstructed the Biosynthetic Diversity for each node using Maximum Likelihood through the ‘fastAnc’ function in phytools in R (Revell, 2012). And we assessed the phylogenetic signal for the Biosynthetic Diversity of compounds in each species, this was determined using Pagel's λ, calculated using the ‘phylosig’ function in the phytools R Package (Pagel, 1992; Revell, 2012). In this case high phylogenetic signal would indicate a situation where a clade of closely related species would have a similar biosynthetic diversity (i.e. a clade where all have high biosynthetic diversity, or a clade where all have low biosynthetic diversity).

4.8. Phylogenetic signal methods 4.8.1. Diversity signal Phylogenetic signal is the extent that related species resemble each other (Blomberg et al., 2003). First the phylogenetic signal for the overall diversity of compounds in each species was determined using Pagel's λ, calculated using the ‘phylosig’ function in the phytools R Package (Pagel, 1992; Revell, 2012). Pagel's λ ranges from 0 to 1, where λ is the transformation of the phylogeny that best predicts the variation in traits across the tree. λ = 0 means that the tree topology does not affect the variation in the trait i.e. no phylogenetic signal. When λ = 1 there is a perfect correlation between the trait and tree topology following a Brownian motion model, i.e. phylogenetic signal. Likelihood ratio tests are employed to test whether λ is significantly greater than 0. To examine the diversity of compounds among the species in this clade, the total number of compounds present in each species was scored, and the chemical diversity trait for both alkaloids and coumarins was reconstructed on the phylogeny using Maximum Likelihood.

4.10.2. Ancestral trait reconstruction Character reconstruction was assessed for both continuous and discrete chemical traits. Continuous traits were reconstructed for each node using Maximum Likelihood through the ‘fastAnc’ function in phytools in R (Revell, 2012). Discrete chemical traits were modeled using two methods, the continuous-time Markov chain model (Mk) (Yang et al., 1995) and Stochastic character mapping (SCM) (Nielsen, 2002). Mk model was performed using ‘equal rates of changes’ (ER) through the R function ‘rerootingMethod’ of phytools (Yang et al., 1995; Revell, 2012). SCM was performed in phytools using the ‘make.simmap’ function (Revell, 2012), summarizing of a set of 100 stochastic maps.

4.8.2. Individual compounds presence/absence signal The phylogenetic signal of individual compounds was examined by using a discrete binary trait: Compound present = 1, Compound absent = 0. The phylogenetic signal of these discrete binary traits was calculated using the D metric through the R packages caper and ape (Fritz and Purvis, 2010). In the case of the D metric, D = 1 indicates that a trait has evolved randomly across the tree (no signal), whereas D = 0 indicates that a trait has evolved in a way correlated with phylogeny by Brownian motion. In order to measure significance, 1000 permutations are run to obtain a p-value.

4.10.3. Ordinary least squares (OLS) and phylogenetic GLS (PGLS) Regressions of chemical and climatic traits were calculated in R using ape, caper, and nlme packages (Paradis et al., 2004; Pinheiro et al.,

4.8.3. Quantitative signal The phylogenetic signal of the quantitative variation of each 144

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al.

2017). As we observed regional similarities in compound variation, each compound with quantitative data (A4, C23eC30, diversity of compounds) was regressed over each bioclimatic variable to evaluate the association between species compound and bioclimatic traits. Previous work from our laboratory identified the five most significant bioclimatic variables that differentiate niches in the genus through an ecological niche modeling study (unpublished data Allevato). In addition, published work from our laboratory has found associations between the metabolomic profile of P. pennatifolius with bioclimatic factors such as temperature seasonality and the mean temperature in wettest quarter (Allevato et al., 2019). Specifically C23, C25, and C30 were associated with higher concentrations when faced with greater temperature seasonality. On the other hand C24, C26, C28, C29 were associated with lower concentrations when faced with greater temperature seasonality (Allevato et al., 2019). Species bioclimatic traits, or climatic niches, were calculated for each species by using the mean value for each bioclimatic variable. Ordinary Least Squares (OLS) regression was evaluated for each chemical and bioclimatic relationship. We next tested whether there was phylogenetic signal in the residual error to determine if it was necessary to run a Phylogenetic Generalized Least squares (PGLS) regression of each OLS relationship (Revell, 2010). In phylogenetic GLS (PGLS), the error term utilizes a variance-covariance matrix that incorporates the branch lengths in the phylogeny to account for shared evolutionary history of the species. We calculated PGLS in R using the nlme and caper packages, with branch lengths scaled using Pagel's λ (Pagel, 1992; Harmon et al., 2008; Revell, 2012; Pinheiro et al., 2017; Orme et al., 2018). OLS regressions were evaluated using R2 and p-values, and the comparison of OLS and PGLS models for each chemical trait was evaluated using the Akaike Information Criterion (AIC) to establish whether phylogeny improved the model's fit (Symonds and Blomberg, 2014). A lower AIC score denotes a more parsimonious model and improved model fit (Butler and King, 2004; Little et al., 2010; Koski and Ashman, 2016).

Appendix A. Supplementary data Supplementary data to this article can be found online at https:// doi.org/10.1016/j.phytochem.2019.03.027. Funding This material is based upon work supported by the National Science Foundation Graduate Research Fellowship under Grant No. (DGE1650441). In addition, this work was supported by Harold E. Moore Jr. and Robert T. Clausen Endowment Fund, BSA, SSB, ASPT, Garden Club of America, Einaudi Center for International Studies, Latin American Studies Program of Cornell, Cornell Graduate School, and the 48 backers from the crowdfunding online platform Experiment.com. PM thanks Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq-Brasil) for a research fellowship and Fundação de Amparo à Pesquisa do Estado de São Paulo (Fapesp - grant 2008/58035-6) for research support. PM and MG thank Cnpq-Brasil for his research fellowship and MG thanks Fapesp for research funding (grant 2016/ 06260-2). References Abbott, Helen C. De S., 1886. Certain Chemical Constituents of Plants Considered in Relation to Their Morphology and Evolution. Botanical Gazette. Allevato, D.M., Kiyota, E., Mazzafera, P., Nixon, K.C., 2019. Ecometabolomic analysis of wild populations of Pilocarpus pennatifolius (rutaceae) using unimodal analyses. Front. Plant Sci. 10. Avancini, G., Abreu, I.N., Saldaña, M.D.A., Mohamed, R.S., Mazzafera, P., 2003. Induction of pilocarpine formation in jaborandi leaves by salicylic acid and methyljasmonate. Phytochemistry 63, 171–175. Barry, C.S., McQuinn, R.P., Thompson, A.J., Seymour, G.B., Grierson, D., Giovannoni, J.J., 2005. Ethylene insensitivity conferred by the Green-ripe and Never-ripe 2 ripening mutants of tomato. Plant Physiol. 138, 267–275. Becerra, J.X., 2015. Macroevolutionary and geographical intensification of chemical defense in plants driven by insect herbivore selection pressure. Current Opinion in Insect Science 8, 15–21. Bisby, F.A., Vaughan, J.G., Wright, C.A., 1980. Chemosystematics: Principles and Practice. Academic Press Published for the Systematics Association by. Blomberg, S.P., Garland, T., Ives, A.R., 2003. Testing for phylogenetic signal in comparative data: behavioral traits are more labile. Evolution. International Journal of Organic Evolution 57, 717–745. Bourgaud, F., Hehn, A., Larbat, R., Doerper, S., Gontier, E., Kellner, S., Matern, U., 2006. Biosynthesis of coumarins in plants: a major pathway still to be unravelled for cytochrome P450 enzymes. Phytochemistry Rev. 5, 293–308. Boyko, A., Kovalchuk, I., 2008. Epigenetic control of plant stress response. Environ. Mol. Mutagen. 49, 61–72. Boyko, A., Kovalchuk, I., 2011. Genome instability and epigenetic modification—heritable responses to environmental stress? Curr. Opin. Plant Biol. 14, 260–266. Butler, M.A., King, A.A., 2004. Phylogenetic comparative analysis: a modeling approach for adaptive evolution. Am. Nat. 164, 683–695. Call, V.B., Dilcher, D.L., 1995. Fossil Ptelea samaras (rutaceae) in north America. Am. J. Bot. 82, 1069–1073. Carrari, F., Baxter, C., Usadel, B., et al., 2006. Integrated analysis of metabolite and transcript levels reveals the metabolic shifts that underlie tomato fruit development and highlight regulatory aspects of metabolic network behavior. Plant Physiol. 142, 1380–1396. Dixon, R.A., 2001. Natural products and plant disease resistance. Nature 411 (6838), 843–847. Dixon, R., Paiva, N., 1995. Stress-induced phenylpropanoid metabolism. Plant Cell 7, 1085–1097. Drummond, A.J., Rambaut, A., 2007. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evol. Biol. 7, 214. Durand-Hulak, M., Dugrand, A., Duval, T., Bidel, L.P.R., Jay-Allemand, C., Froelicher, Y., Bourgaud, F., Fanciullino, A.-L., 2015. Mapping the genetic and tissular diversity of 64 phenolic compounds in Citrus species using a UPLC-MS approach. Ann. Bot. 115, 861–877. Eaton, D.A.R., 2014. PyRAD: assembly of de novo RADseq loci for phylogenetic analyses. Bioinformatics 30, 1844–1849. Eaton, D.A.R., Spriggs, E.L., Park, B., Donoghue, M.J., 2017. Misconceptions on missing data in RAD-seq phylogenetics with a deep-scale example from flowering plants. Syst. Biol. 66, 399–412. Etten, J van, Hijmans, R.J., 2010. A geospatial modelling approach integrating archaeobotany and genetics to trace the origin and dispersal of domesticated plants. PLoS One 5, e12060. Fine, P.V.A., Metz, M.R., Lokvam, J., Mesones, I., Zuñiga, J.M.A., Lamarre, G.P.A., Pilco, M.V., Baraloto, C., 2013. Insect herbivores, chemical innovation, and the evolution of habitat specialization in Amazonian trees. Ecology 94, 1764–1775.

Declarations of interest The authors declare that the research was conducted without any commercial or financial relationships that could be potential conflicts of interest. Acknowledgments The author(s) thank various individuals and funding institutes for their assistance throughout this project. We would like to thank Professor Marcelo Caxambu (Universidade Tecnologica Federal do Parana), Professor Pedro Dias (USP-SP), Dr. Rubens Coelho (USP-SP), Alex Popovkin (Bahia, Brazil), Lukas Daneu (CEPEC Herbaria) for helping with the collection of Pilocarpus specimens; Dr. Jose Sena (Quercegen) for assistance in collection of Pilocarpus specimens in northern Brazil; Iguaçu Itaipu and Parque Estadual Lago Azul, PR for allowing the collection of Pilocarpus specimens from their property; WorldClim for bioclimatic data, New York Botanical Garden (NYBG) and Jardim Botânico do Rio de Janeiro for providing specimens, Dr. Adrienne Roeder for assistance with Flow Cytometry; Dr. Alexandra Christine Helena Frankland Sawaya for providing the updated extraction protocol for imidazole alkaloids. DA and KN thank E. Rodriguez, R. Raguso, W. Crepet, A. Gandolfo-Nixon, C. Martinez, R. Harbert, N. Jud, M. Haribal, K. Boroczky, and T. Choo for thoughtful discussions. Plants collected in public locations were collected under ICMBio SISBIO permit #52758-1, allowing for the collection of botanical, fungal and microbiological specimens. In addition, plants were collected in collaboration with and under the permits of Professor Milton Groppo, Professor Pedro Dias, and Professor Marcelo Caxambu. 145

Phytochemistry 163 (2019) 132–146

D.M. Allevato, et al. Firn, R.D., Jones, C.G., 2000. The evolution of secondary metabolism – a unifying model. Mol. Microbiol. 37, 989–994. Firn, R.D., Jones, C.G., 2003. Natural products ? a simple model to explain chemical diversity. Nat. Prod. Rep. 20, 382. Fraenkel, G.S., 1959. The raison d’Être of secondary plant substances. Science 129, 1466–1470. Fraser, C.M., Chapple, C., 2011. The phenylpropanoid pathway in Arabidopsis. Arabidopsis Book 9, e0152. Fritz, S.A., Purvis, A., 2010. Selectivity in mammalian extinction risk and threat types: a new measure of phylogenetic signal strength in binary traits. Conserv. Biol. 24, 1042–1051. Geber, M.A., Griffen, L.R., 2003. Inheritance and natural selection on functional traits. Int. J. Plant Sci. 164, S21–S42. Gernhard, T., 2008. The conditioned reconstructed process. J. Theor. Biol. 253, 769–778. Goloboff, P.A., Farris, J.S., Nixon, K.C., 2008. TNT, a free program for phylogenetic analysis. Cladistics 24, 774–786. Gray, A.I., Waterman, P.G., 1978. Coumarins in the rutaceae. Phytochemistry 17, 845–864. Groppo, M., Kallunki, J., Pirani, J.R., Antonelli, A., 2012. Chilean Pitavia more closely related to Oceania and Old World Rutaceae than to Neotropical groups: evidence from two cpDNA non-coding regions, with a new subfamilial classification of the family. PhytoKeys 19, 9–29. Harborne, J.B., 1984. Plant Chemosystematics. London. Academic Press, Orlando, Fla. Harborne, J.B., 1993. Introduction to Ecological Biochemistry. Academic Press, London; San Diego. Harmon, L.J., Weir, J.T., Brock, C.D., Glor, R.E., Challenger, W., 2008. GEIGER: investigating evolutionary radiations. Bioinformatics 24, 129–131. Hegnauer, R., 1994. Chemotaxonomie der Pflanzen. Birkhäuser, Basel. Junker, R.R., Kuppler, J., Amo, L., et al., 2017. Covariation and phenotypic integration in chemical communication displays: biosynthetic constraints and eco-evolutionary implications. New Phytol. 220 (3), 739–749. Kim, H.K., Verpoorte, R., 2010. Sample preparation for plant metabolomics. Phytochem. Anal. 21, 4–13. Koski, M.H., Ashman, T.-L., 2016. Macroevolutionary patterns of ultraviolet floral pigmentation explained by geography and associated bioclimatic factors. New Phytol. 211, 708–718. Latta, R.G., Linhart, Y.B., 1997. Path analysis of natural selection on plant chemistry: the xylem resin of ponderosa pine. Oecologia 109, 251–258. Levin, D.A., 1976. The chemical defenses of plants to pathogens and herbivores. Annu. Rev. Ecol. Systemat. 7, 121–159. Little, S.A., Kembel, S.W., Wilf, P., 2010. Paleotemperature proxies from leaf fossils reinterpreted in light of evolutionary history. PLoS One 5, e15161. Lynch, M., Conery, J.S., 2000. The evolutionary fate and consequences of duplicate genes. Science (New York, N.Y.) 290, 1151–1155. Mabry, T.J., Ulubelen, A., 1980. Chemistry and utilization of phenylpropanoids including flavonoids, coumarins, and lignans. J. Agric. Food Chem. 28, 188–196. Manchester, S.R., O'Leary, E.L., 2010. Phylogenetic distribution and identification of finwinged fruits. Bot. Rev. 76, 1–82. Martz, F., Peltola, R., Fontanay, S., Duval, R.E., Julkunen-Tiitto, R., Stark, S., 2009. Effect of latitude and altitude on the terpenoid and soluble phenolic composition of juniper (juniperus communis) needles and evaluation of their antibacterial activity in the boreal zone. J. Agric. Food Chem. 57, 9575–9584. Masson, J., Liberto, E., Beolor, J.-C., Brevard, H., Bicchi, C., Rubiolo, P., 2016. Oxygenated heterocyclic compounds to differentiate Citrus spp. essential oils through metabolomic strategies. Food Chem. 206, 223–233. McNair, J.B., 1935. Angiosperm phylogeny on a chemical basis. Bull. Torrey Bot. Club 62, 515–532. Moore, B.D., Andrew, R.L., Külheim, C., Foley, W.J., 2014. Explaining intraspecific diversity in plant secondary metabolites in an ecological context. New Phytol. 201, 733–750. Mraja, A., Unsicker, S.B., Reichelt, M., Gershenzon, J., Roscher, C., 2011. Plant community diversity influences allocation to direct chemical defence in plantago lanceolata. PLoS One 6, e28055. Murray, R.D.H., 2002. The Naturally Occurring Coumarins. Fortschritte der Chemie organischer Naturstoffe/Progress in the Chemistry of Organic Natural Products. Fortschritte der Chemie organischer Naturstoffe/Progress in the Chemistry of Organic Natural Products. Springer, Vienna, pp. 1–619. Nielsen, R., 2002. Mapping mutations on phylogenies. Syst. Biol. 51, 729–739. Nitao, J.K., Zangerl, A.R., 1987. Floral development and chemical defense allocation in wild parsnip (pastinaca sativa). Ecology 68, 521–529. Nixon, K.C., 1999. The parsimony ratchet, a new method for rapid parsimony analysis. Cladistics 15, 407–414. Nixon, K.C., Carpenter, J.M., 1993. On outgroups. Cladistics 9, 413–426. Oliveira, PD de, 2008. Filogenética de Pilocarpinae (Rutaceae). Olson-Manning, C.F., Wagner, M.R., Mitchell-Olds, T., 2012. Adaptive evolution: evaluating empirical support for theoretical predictions. Nat. Rev. Genet. 13, 867–877. Orme, D., Freckleton, R., Thomas, G., Petzoldt, T., Fritz, S., Isaac, N., Pearse, W., 2018. Caper: Comparative Analyses of Phylogenetics and Evolution in R. Pagel, M.D., 1992. A method for the analysis of comparative data. J. Theor. Biol. 156, 431–442. Pandya, C., Farelli, J.D., Dunaway-Mariano, D., Allen, K.N., 2014. Enzyme promiscuity: engine of evolutionary innovation. J. Biol. Chem. 289, 30229–30236. Paradis, E., Claude, J., Strimmer, K., 2004. APE: analyses of phylogenetics and evolution

in R language. Bioinformatics 20, 289–290. Parent, J.-S., Martínez de Alba, A.E., Vaucheret, H., 2012. The origin and effect of small RNA signaling in plants. Front. Plant Sci. 3. Penfold, A.R., 1927. The Essential Oil of Eucalyptus Micrantha (De Candolle). Pichersky, E., Lewinsohn, E., 2011. Convergent evolution in plant specialized metabolism. Annu. Rev. Plant Biol. 62, 549–566. Pinheiro, J., Bates, D., DebRoy, S., Sarkar, D., 2017. Nlme: Linear and Nonlinear Mixed Effects Models. Price, P.W., Waring, G.L., Julkunen-Tiitto, R., Tahvanainen, J., Mooney, H.A., Craig, T.P., 1989. Carbon-nutrient balance hypothesis in within-species phytochemical variation ofSalix lasiolepis. J. Chem. Ecol. 15, 1117–1131. Revell, L.J., 2010. Phylogenetic signal and linear regression on species data. Methods in Ecology and Evolution 1, 319–329. Revell, L.J., 2012. phytools: an R package for phylogenetic comparative biology (and other things). Methods in Ecology and Evolution 3, 217–223. Rønsted, N., Symonds, M.R.E., Birkholm, T., et al., 2012. Can phylogeny predict chemical diversity and potential medicinal activity of plants? A case study of amaryllidaceae. BMC Evol. Biol. 12, 182. Rueffler, C., Dooren, T.J.M.V., Leimar, O., Abrams, P.A., 2006. Disruptive selection and then what? Trends Ecol. Evol. 21, 238–245. Santos, A.P., Moreno, P.R.H., 2004. Pilocarpus spp.: a survey of its chemical constituents and biological activities. Rev. Bras. Ciencias Farm. 40, 116–137. Sawaya, A., Abreu, I.N., Andreazza, N.L., Eberlin, M.N., Mazzafera, P., 2008. HPLC-ESIMS/MS of imidazole alkaloids in Pilocarpus microphyllus. Molecules 13, 1518–1529. Sawaya, A.C.H.F., Abreu, I.N., Andreazza, N.L., Eberlin, M.N., Mazzafera, P., 2010. Pilocarpine and related alkaloids in Pilocarpus Vahl (rutaceae). In: Cassiano, N.M. (Ed.), Alkaloids: Properties, Applications and Pharmacological Effects. Nova Science Publishers, Inc., 978-1-61761-130-8, pp. 63–80. Sawaya, A.C.H.F., Vaz, B.G., Eberlin, M.N., Mazzafera, P., 2011. Screening species of Pilocarpus (Rutaceae) as sources of pilocarpine and other imidazole alkaloids. Genet. Resour. Crop Evol. 58, 471–480. Scott, T., 1980. Alkaloid Biology and Metabolism in Plants: by George R Waller and Edmund K Nowacki. Plenum Press, New York and London, pp. 294 1978. Biochemical Education 8, 30. Silva, MF.das.GF da, Gottlieb, O.R., Ehrendorfer, F., 1988. Chemosystematics of the Rutaceae: suggestions for a more natural taxonomy and evolutionary interpretation of the family. Plant Systemat. Evol. 161, 97–134. Skorupa, L.A., 1996. Revisao Taxonomica de Pilocarpus Vahl (Rutaceae) (Ph.D. Thesis). . Sørensen, T., 1948. A method of establishing groups of equal amplitude in plant sociology based on similarity of species and its application to analyses of the vegetation on Danish commons. Biol. Skr. 5, 1–34. Stamatakis, A., Hoover, P., Rougemont, J., 2008. A rapid bootstrap algorithm for the RAxML Web servers. Syst. Biol. 57, 758–771. Symonds, M.R.E., Blomberg, S.P., 2014. A primer on phylogenetic generalised least squares. In: Garamszegi, L.Z. (Ed.), Modern Phylogenetic Comparative Methods and Their Application in Evolutionary Biology: Concepts and Practice. Springer Berlin Heidelberg, Berlin, Heidelberg, pp. 105–130. Tamara, Münkemüller, Sébastien, Lavergne, Bruno, Bzeznik, Dray Stéphane, Thibaut, Jombart, Katja, Schiffers, Wilfried, Thuiller, 2012. How to measure and test phylogenetic signal. Methods in Ecology and Evolution 3, 743–756. Thaler, J.S., Humphrey, P.T., Whiteman, N.K., 2012. Evolution of jasmonate and salicylate signal crosstalk. Trends Plant Sci. 17, 260–270. Urano, K., Kurihara, Y., Seki, M., Shinozaki, K., 2010. ‘Omics’ analyses of regulatory networks in plant abiotic stress responses. Curr. Opin. Plant Biol. 13, 132–138. Vialart, G., Hehn, A., Olry, A., et al., 2012. A 2-oxoglutarate-dependent dioxygenase from Ruta graveolens L. exhibits p-coumaroyl CoA 2’-hydroxylase activity (C2’H): a missing step in the synthesis of umbelliferone in plants. Plant J.: Cell. Mol. Biol. (Oxf.) 70, 460–470. Waller, G.R., Nowacki, E.K., 1978. Alkaloid Biology and Metabolism in Plants. Plenum Press, New York. Waterman, P.G., 1975. Alkaloids of the rutaceae: their distribution and systematic significance. Biochem. Syst. Ecol. 3, 149–180. Waterman, P.G., 1990. Chemosystematics of the rutaceae: comments on the interpretation of da Silva & al. Plant Systemat. Evol. 173, 39–48. Waterman, P.G., Hussain, R.A., 1983. Chemosystematics in the Rutaceae: a chemotaxonomic survey of flavonoids and monoterpenes in leaves of Acmadenia species. Bot. J. Linn. Soc. 86, 227–235. Wink, M., 2003. Evolution of secondary metabolites from an ecological and molecular phylogenetic perspective. Phytochemistry 64, 3–19. Yang, Z., Kumar, S., Nei, M., 1995. A new method of inference of ancestral nucleotide and amino acid sequences. Genetics 141, 1641–1650. Zandalinas, S.I., Sales, C., Beltrán, J., Gómez-Cadenas, A., Arbona, V., 2016. Activation of secondary metabolism in citrus plants is associated to sensitivity to combined drought and high temperatures. Front. Plant Sci. 7, 1954. Zhang, J., 2003. Evolution by gene duplication: an update. Trends Ecol. Evol. 18, 292–298. Zobel, A.M., Brown, S.A., 1990. Seasonal changes of furanocoumarin concentrations in leaves of Heracleum lanatum. J. Chem. Ecol. 16, 1623–1634. Zobel, A.M., Brown, S.A., 1995. Coumarins in the interactions between the plant and its environment. Allelopathy J. 2 (1). Zobel, A.M., Wang, J., March, R.E., Brown, S.A., 1991. Identification of eight coumarins occurring with psoralen, xanthotoxin, and bergapten on leaf surfaces. J. Chem. Ecol. 17, 1859–1870.

146