Phyllodes and bipinnate leaves of Acacia exhibit contemporary continental-scale environmental correlation and evolutionary transition-rate heterogeneity

Sclerophylly is hypothesised to be a response to nutrient-deficient soils, especially to those low in phosphorus (Loveless 1961; Beadle 1966; Eckstein et al. 1999). Sclerophyllous leaves are tougher and less prone to physical damage, both of which increase leaf longevity and so reduce nutrient loss by leaf-fall (Read and Sanson 2003). Sclerophyllous traits may be exaptations for drought tolerance (Boughton 1986; Hill and Brodribb 2001), and may have been further modified to enhance drought resistance and xeromorphy (Crisp and Cook 2013). Phyllodes are a sclerophyllous leaf type, involving the elaboration of petiolar tissue, and suppression of leaf blade expression (Gardner et al. 2005; Leroy and Heuret 2008).

Acacia Mill. (Leguminosae) is the only phyllode-bearing mimosoid legume genus in Australia, with more than 90% of the 1015 described Acacia species having phyllodes as their major adult foliage expression; the remaining species have bipinnate adult leaves (Brown et al. 2011). With three known exceptions (Acacia diaphyllodinea Maslin, A. confusa Merr. and A. glutinosissima Maiden & Blakely, Murphy et al. (2010)), all phyllodinous Acacia species are heteroblastic, that is, they experience a phase change from bipinnate juvenile to phyllodinous adult foliage as they grow (Wang et al. 2011). The ontogeny of seedling development in heteroblastic Acacia includes at least one bipinnate seedling leaf (Murphy et al. 2010; Brown et al. 2012), although the timing of the phase change can be flexible (Rose et al. 2019). In contrast, none of the Acacia species with bipinnate adult foliage experiences a heteroblastic phase change.

Acacia phyllode morphologies exhibit considerable variation (Fig. 1) and lability (Gardner et al. 2005), with the broadly circumscribed character states of uninerved versus multinerved phyllodes having multiple origins within Acacia (Ariati et al. 2006). Bipinnate leaves are thought to be the ancestral foliage type in the Acacia clade (Leroy and Heuret 2008), but this has not been tested empirically. In a broader context, bipinnate leaves are likely to be ancestral within the Caesalpinioideae, as all lineages, except Acacia, have species with bipinnate adult foliage. Within Acacia, multiple lineages express bipinnate adult foliage, and at least two independent phylogenetic reversion events are recognised (Murphy et al. 2000, 2003). Certainly, phyllodes are an apomorphy within the Caesalpinioideae, but an outstanding question is whether the most recent common ancestor of extant Acacia spp. possessed bipinnate leaves or phyllodes.

Fig. 1.  A selection of Acacia species, showing variation in shoot, flower, phyllode and leaf morphology, not to scale. A. A. baileyana. B. A. barretiorum. C. A. buxifolia. D. A. Cape le Grande. E. A. colletioides. F. A. dealbata. G. A. deltoidea. H. A. dimorpha. I. A. sp. indet. J. A. Grass Patch. K. A. gunnii. L. A. humifusa. M. A. kelleri. N. A. obtriangularis. O. A. paradoxa. P. A. platycarpa. Q. A. pulchella. R. A. tetragonophylla. S. A. translucens. T. A. Wathroo. Bipinnate adult foliage illustrated in A, D, F and Q, others phyllodinous. All photos by Russell Barrett, with permission.

Phyllodes have been hypothesised to represent a key morphological innovation facilitating the diversification of Acacia in water-limited environments (Boughton 1986; Brodribb and Hill 1993). Supporting this hypothesis, phyllodes are more drought tolerant, and have higher water-use efficiency than do leaves in conditions of water deficit (Brodribb and Hill 1993; Wright and Westoby 2002). In contrast, bipinnate leaves have higher maximum rates of carbon dioxide uptake per unit mass, and higher photosynthetic rates than do phyllodes per unit mass and per unit area. They are also likely to be more efficient at light capture, given their horizontal orientation and spread of photosynthetic tissue (Rose et al. 2019). Both leaves and phyllodes of Acacia have high water-use efficiency compared with other C3 plants (Brodribb and Hill 1993), but being capable of high water-use efficiency in times of water stress does not necessarily imply high water-use efficiency at other times. For example, in South Africa, the bipinnate Acacia mearnsii has higher water use than does woody and grassland vegetation in riparian sites where it is invasive (Dye et al. 2001). It is possible that this higher water use supports higher relative growth rates.

The physiological differences between leaves and phyllodes are significant, given Australia’s dry climate (Brodribb and Field 2010; Crisp and Cook 2012). Physiological performance studies have demonstrated that phyllodes are better capable of responding to, and recovering from, severe water stress than are bipinnate leaves (Brodribb and Hill 1993). Acacia species with bipinnate adult foliage would, therefore, be expected to occur in environments that are wetter, and not subject to severe water stress, than are species with phyllodes.

Acacia can be subdivided into two or three macroevolutionary cohorts with different diversification rates (Nge et al. 2020; Renner et al. 2020), one of which comprises a lineage occurring in mesic south-eastern Australia (Renner et al. 2020). This south-eastern mesic lineage has a higher mean diversification rate than the other extant Acacia lineages, contains both phyllodinous and bipinnate species, and may provide additional insight into the contribution that foliage type has made to the radiation of Acacia. The different diversification dynamics and high representation of bipinnate species may mean that transition rates among phyllode and bipinnate adult foliage are different in the south-eastern mesic lineage from the rest of the Acacia lineage.

The main aim of this paper is to examine why, given a progressively drying Australian climate, Acacia lineages have evidently reverted to bipinnate adult foliage. We address this question in two steps. We use geospatial and climatic data to ask whether bipinnate Acacia species occur within a continental climatic context predicted by the physiological performance of their leaves. We confirm that bipinnate adult foliage represents a phylogenetic reversion, and ask how many phylogenetic transitions to bipinnate foliage have occurred in Acacia, and at what rate. We then synthesise the results from character reconstruction and analysis of climate variables to derive a testable model that may explain these observational data, wherein winter growth provides a competitive advantage to bipinnate Acacia in mesic habitats.

Our methods followed those of Renner et al. (2020). We used the dataset published by Mishler et al. (2014), comprising two nuclear (nrITS and ETS) and four chloroplast markers (psbA–trnH, trnL–trnF, rpl32–trnL, matK) from 510 species, for phylogenetic reconstruction (see supplementary table 1 in Mishler et al. (2014) at https://doi.org/10.1038/ncomms5473, for voucher specimen details and GenBank accession numbers). Owing to the geographic location of the research groups working on this phylogenetic dataset for Acacia, the sampling of species for this dataset was skewed towards the eastern coast of Australia. This means that some groups of Acacia are over-represented in the data, for example, section Botrycephalae for which 82% of its constituent species were sampled, in contrast to an overall sampling rate of 47% for all species of Acacia (Renner et al. 2020). However, this is not as straightforward to assess as these numbers indicate, because section Botrycephalae is not monophyletic, similar to other sections of Acacia, so it is difficult to translate a geographic bias into an estimation of phylogenetic bias, especially as the relationships of unsampled species can be difficult to estimate on the basis of morphological data. We excluded three Acacia species with three or more missing markers that were non-overlapping within the data matrix, leaving 503 Acacia species and two outgroup taxa. The alignment was checked, and some sequences were manually re-aligned.

Incongruence between molecular markers is an important consideration in phylogenetic tree reconstruction. We concatenated chloroplast markers into a single alignment since the chloroplast is a predominantly non-recombining organelle, and concatenated the two nuclear markers into a single alignment because of their flanking positions in eukaryotic ribosomal cistrons. Incongruence among markers was assessed by comparing scores of Akaike information criterion for small sample sizes (AICc), similar to the approach of Walker et al. (2018). The optimal substitution model and the optimal partitioning scheme for each dataset were selected with PartitionFinder (ver. 2.1.1, wee http://www.robertlanfear.com/partitionfinder/; Lanfear et al. 2014, 2017), with GTR+G being selected for all partitions (see supplementary table 2 in Renner et al. 2020). Maximum-likelihood trees were estimated for separate nuclear and chloroplast alignments, and for an overall concatenated alignment (chloroplast + nuclear) by using IQtree (ver. 1.6.12, see http://www.iqtree.org; Nguyen et al. 2015), with each partition having a separate GTR+G model. The log-likelihood scores for trees from IQtree were used in AICc calculations. AICc scores suggested congruence between chloroplast and nuclear markers, and supported our decision to concatenate genetic markers (see supplementary table 1 in Renner et al. 2020). We then estimated a maximum-likelihood tree using RAxML (ver. 8, see https://cme.h-its.org/exelixis/web/software/raxml/; Stamatakis 2014), for use as a starting tree in divergence-time estimation.

To estimate a time-calibrated phylogeny, we employed BEAST (ver. 2.3.2, see https://www.beast2.org; Bouckaert et al. 2014). The substitution models for each marker were unlinked, and among site-rate variation modelled using a gamma distribution with four rate categories. Rate variation among lineages was modelled using a lognormal relaxed clock. The clock model was linked across molecular markers. The third codon position within matK was allowed to scale to a rate different from that of the first and second codon positions. A birth–death speciation model was used to generate a prior on the distribution of branch lengths and node depths, and the monophyly of the outgroup (Paraserianthes) was enforced. Markov-chain Monte Carlo (MCMC) chains were run for 100 million generations, logging trees sampled every 10 000 generations, and parameters sampled every 5000. A single time calibration point was used, based on fossil pollen with distinctive pseudocolpi on their surface (Macphail and Hill 2001), recorded in the late Oligocene (23 Ma), following Miller et al. (2013) and Renner et al. (2020). The use of a lognormal prior can be justified against the fossil pollen record itself; although Acaciapollenites pollen is present from the late Eocene onward (37.2–33.9 Ma; Macphail and Hill 2001), only younger records of this fossil genus can be confidently assigned on the basis of the possession of distinctive pseudocolpi on the pollen surface, characteristic of modern Acacia (Miller et al. 2013). Bounding the prior probability distribution of the Acacia origin against the younger ages of the Acaciapollenites record is appropriate because it reflects the origin of pseudocolporate pollen between these bounds, although it is likely closer to the younger limit. Therefore, to realistically model this prior, we applied a lognormal distribution with an offset of 23 Ma, median of 27.2 and 95% of the prior density between 23 and 34 Ma, following Miller et al. (2013), to the Acacia crown node. Maximum clade credibility trees were calculated from both runs after excluding the first 25% of each run as burnin, with TreeAnnotator (ver. 1.8.4, see https://beast.community/treeannotator).

Adult foliage type was scored from literature or herbarium specimens, with species having either bipinnate leaves or phyllodes. Species whose adult foliage comprises a mix of phyllodes and bipinnate leaves (heterophylly), such as A. rubida, were scored as phyllodinous (Forster and Bonser 2009). We compared four models of trait evolution using BayesTraits (ver. 3.0, see http://www.evolution.rdg.ac.uk/BayesTraitsV3.0.5/BayesTraitsV3.0.5.html; Pagel et al. 2004), which described whether transition rate parameters were linked or free, and whether transition rate parameters across the phylogeny were homogeneous or heterogeneous. Several phylogenetic transitions to phyllodes are localised within the bipinnate lineage that is part of the south-eastern mesic lineage, which was inferred to have different diversification dynamics from the rest of Acacia (see ‘event subtree 2’ in Renner et al. 2020). Therefore, we estimated transition rate parameters for that lineage separately from the rest of Acacia to test for heterogeneity in morphological transformation rates across Acacia. Under the heterogeneous models, the south-eastern mesic lineage was allowed its own transition rate matrix, separate from the rest of Acacia. Transition rate matrices were estimated using Bayesian inference with a multistate continuous-time Markov model of character evolution for discrete data implemented by BayesTraits. For each Bayesian analysis, we seeded runs with an exponential prior with a mean of 10, because we expected rate parameters to be relatively small, and scaled the tree to have a mean branch length of 0.1, so the estimates would not be too small (for an explanation of what is too small, see the BayesTraits manual; Meade and Pagel 2016). MCMC chains were run for 10 million generations, sampling every 1000 generation. To accommodate phylogenetic uncertainty, we integrated estimates across a sample of 1000 trees from the posterior probability distribution sampled by BEAST. A burnin of 10 000 generations was discarded from each run. We estimated log Bayes factors (Kass and Raftery 1995) from the marginal likelihood estimated by stepping-stone analysis in BEAST. For all models, analyses used 100 stepping stones and 10 000 replicates. Each analysis was run independently three times to confirm convergence and consistency, and parameter traces and diagnostics were estimated and inspected with Tracer (ver. 1.6, see https://beast.community/tracer; Drummond et al., 2012). We also used the variable-rates model implemented by BayesTraits to query the tree for regions having higher or lower rates of transformation. The variable-rates model uses reversible jump MCMC to identify parts of the tree in which the rate of character transformation differs. Because this model takes a single tree as input, we used the MCC tree from the BEAST analysis. The variable-rates model was run for 10 million generations, with a burnin of 1 million generations and scaled branch lengths, and the analysis was replicated three times to confirm run convergence.

We reconstructed ancestral states using maximum likelihood and Bayesian inference methods. Likelihood-based ancestral-state reconstructions were completed using the ace function in the R package Ape (ver. 5.3, see https://cran.r-project.org/web/packages/ape/index.html; Paradis et al. 2004; Popescu et al. 2012; Paradis and Schliep 2019). These analyses also presented a visual summary of state evolution on the most credible clade tree. To address outgroup sampling, we manually created trees with 13 additional outgroup taxa following the topology (branching order and branch lengths) of the phylogeny presented by the Azani et al. (2017), which, although based only on matK, does give an indication that early diversification of the Caesalpiniodeae was rapid, with lineages being separated by short branches. The ancestral-state reconstruction was repeated on this manually expanded outgroup tree. We also fixed the state at the Acacia crown node and compared marginal likelihoods and associated 2 × log Bayes factors (2ln BF) for models with phyllodes and leaves as the ancestral state under a homogeneous free model of trait evolution, estimated by stepping-stone sampling of the prior, using 100 stepping stones and 10 000 replicates, with BayesTraits. Morphological data, the MCC tree, the sample of 1000 trees from the posterior probability distribution (PPD) and control file for the heterogeneous analysis are available in the Supplementary material (Tables S1–S4 of the Supplementary material) associated with this article.

We used a set of 132 295 geo-referenced herbarium specimens of Acacia originally downloaded from the Australian Virtual Herbarium (see https://www.ala.org.au) and curated by González-Orozco et al. (2011), and later analysed by Mishler et al. (2014) and Renner et al. (2020). Thirty-five bioclimatic (BioClim) variables from the ANUCLIM (ver.6, β, see https://fennerschool.anu.edu.au/research/products/anuclim) layers, available through the Atlas of Living Australia’s Spatial Portal (spatial.ala.gov.au), were downloaded for these ~132 000 voucher specimens. For each environmental variable, the species mean was calculated; the resulting means represent the realised climatic niche centroid for each species, and were used in the regression analyses. These species means are available in the supplementary material of Renner et al. (2020).

Correlation and covariation among environmental variables were assessed using the spectral decomposition method of principal-component analysis, accomplished using the princomp function of the stats package (ver. 4.0, R Foundation for Statistical Computing, see https://stat.ethz.ch/R-manual/R-devel/library/stats/html/stats-package.html) in R (ver. 3.5.6, R Foundation for Statistical Computing, see http://www.R-project.org). The variance–covariance matrix for Z-scaled data was calculated with the cov function. The kappa value suggested a high degree of multicollinearity within the environmental dataset, and variance inflation factors estimated by the VIF package (ver. 1.0, see https://CRAN.R-project.org/package = VIF; Lin et al. 2011) indicated the same result, and retained 28 variables. Because species with bipinnate adult foliage tend to be closely related to one another, we fitted a single multivariate phylogenetic binomial generalised linear model (Ives and Garland 2010) to examine the correlation between possession of bipinnate leaves or phyllodes and the 28 climate variables retained by VIF within a phylogenetic context, with the function phyloglm in package phylolm (ver. 2.6, see https://www.rdocumentation.org/packages/phylolm/versions/2.6, accessed 21 January 2020).

Reversions to bipinnate foliage have occurred at least three times, and in three different lineages within Acacia according to maximum-likelihood ancestral-state reconstruction (Fig. 2). Within the south-eastern mesothermic radiation (event Subtree 2), there have been four inferred transformations to and from bipinnate adult foliage. The 2 × log Bayes factors estimated by BayesTraits returned very strong support for the homogeneous free rates model as the best descriptor of state transformation in Acacia foliage (Table 1). Under the homogeneous free model of trait evolution, the transition rates among phyllodinous and bipinnate adult foliage are asymmetrical (Table 2); estimated median transition rate of bipinnate leaves to phyllodes was two orders of magnitude higher than the transition rate for phyllodes to bipinnate leaves (Fig. 3). The variable-rate estimation using reversible jump MCMC recovered no evidence for significant rate heterogeneity across the phylogeny (results not shown), in accord with the inference that the homogeneous model provides a better fit with the databased on Bayes factor comparison.

Fig. 2.  A. Likelihood ancestral-state reconstruction for phyllodes and bipinnate adult foliage on the most credible clade time-tree from the time-calibration analysis; states show adult foliage. Red, bipinnate; black, phyllode. The time tree root-node median age is 26.82 Ma, the scale bar is 5 Ma. B. Distribution of bipinnate Acacia on the Australian continent. C. Distribution of phyllodinous Acacia in Australia. Records sourced from Australian Virtual Herbarium’s Atlas of Living Australia.

Table 1.  Bayes factor comparison of models of trait evolution, using 2ln Bayes factors
Positive values indicate support for the model in the column, negative values for the model in the row of the comparison. The model supported in all comparisons was the homogeneous free-parameter model, shown in bold

Table 2.  Parameters describing the rate of transformation among phyllodinous and bipinnate adult foliage in Acacia
ESS, effective sample size; s.e., standard error of the mean; s.d., standard deviation

Fig. 3.  Posterior probability distributions for transition rates among leaves and phyllodes, under homogeneous and heterogeneous models of trait evolution across the Acacia phylogeny. The homogeneous model was strongly favoured in a Bayes factor comparison.

Likelihood-based ancestral character reconstruction returned a high probability of the phyllodinous character state at the Acacia crown node, at all nodes along the Acacia backbone, and also at the MRCA of Paraserianthes plus Pararchidendron and Acacia (Fig. 2). Results from the phylogeny with an additional 13 outgroup taxa (and 13 nodes below the Acacia MRCA) had bipinnate adult foliage reconstructed along all nodes along the phylogeny backbone, and at the MRCA of Paraserianthes plus Parachidendron and Acacia. The MRCA of Acacia was again reconstructed as phyllodinous, as in the original limited-outgroup dataset. The median age for the oldest node inferred to have been leafy is ~15 million years old, but this is an outlier; the median ages of most leafy nodes are less than 9 million years old and more than half are less than 3 million years old (Fig. 4).

Fig. 4.  Boxplot of age distribution for Acacia ancestral nodes inferred to have had bipinnate adult foliage.

Bayesian inference of ancestral states under a homogeneous free model returned positive evidence for bipinnate over phyllodinous adult leaves at the Acacia crown node, whereas the heterogeneous free model was indecisive (Table 3).

Table 3.  Bayes factor comparison of different ancestral-state reconstructions at the Acacia crown node, under homogeneous and heterogeneous free transition-rate models

In all, 2 of the 28 bioclimatic variables were significantly correlated (after correction for phylogeny) with bipinnate or phyllodinous adult foliage, ‘moisture seasonality’ (Bio31) and ‘radiation during the coldest quarter’ (Bio27; Table 4), with bipinnate species means being clustered around lower values of both variables (Fig. 5). There was a trend towards ‘moisture during the coldest quarter’ (Bio35), ‘moisture index during the highest period’ (Bio29) and ‘precipitation during the warmest period’ (Bio18), correlating with the occurrence of phyllodinous and bipinnate foliage, albeit without passing our statistical significance threshold. In plots of the distribution of bipinnate and phyllodinous adult foliage along moisture variables, bipinnate species exhibited compression towards the lowest moisture seasonality values, representing those soil environments with the least difference between their wettest and driest moisture contents. This contrasts with the relatively wide distribution of bipinnate species on the axis of soil moisture values during the lowest quarter. Along the axis of radiation during the coldest quarter, bipinnate species clustered towards lower values (Fig. 5).

Table 4.  Regression statistics returned by phylogenetically corrected binary logistic regression assessing the correlation between phyllodinous and bipinnate Acacia species, and species mean environmental parameters
Significantly correlated variables are in bold. s.e., standard error; CI, confidence interval

Fig. 5.  Violin plots showing the distribution of phyllodinous (1) and bipinnate (0) species along a selection of climate variables. Significant correlations were returned for radiation in the coldest quarter and moisture seasonality. Note the contrasting distribution of bipinnate species along the moisture coldest-quarter, and moisture lowest- and warmest-quarter axes.

Why phyllodes are so widespread across the Australian continent is a research question that has been addressed by many studies, including physiological studies demonstrating higher water use efficiency for phyllodes. Reversions to bipinnate adult foliage have occurred in three different Acacia lineages, which is surprising in the Australian context because Australia has become a dry continent (Martin 2006) subject to episodic severe water stress (Murphy and Timbal 2008), and phyllodes are more resistant to drought than are bipinnate leaves. Why have some Acacia reverted to a more risky foliage type, given this continental climatic context? The answer may be found in the ecological setting wherein bipinnate species grow. With lower soil-moisture seasonality, and higher soil-moisture values generally, bipinnate Acacia grows in settings accessible to a range of other broadleaf and mesophyllous species, including rainforest pioneers. A ‘grow fast, die young’ strategy involving a winter growth period facilitated by high cold-quarter soil moisture may explain the relative restriction of bipinnate Acacia along that climatic axis. A seasonally dependent competitive growth advantage has not been posited as an explanation for the relative success of bipinnate Acacia in south-eastern Australia previously, but our results are consistent with this explanatory model, which we elaborate with reference to heteroblasty, character-state transition rates, physiological efficiencies, geographic occurrence, and environmental space occupancy in the following discussion.

Heteroblasty is a key innovation enabling sessile organisms to maximise performance in response to environmental change (Chitwood and Sinha 2016). Heteroblastic phase change can be developmentally programmed but environmentally induced (Ostria-Gallardo et al. 2016), and the timing of phase change can be plastic. Low light levels delay the phase change to phyllodinous leaves in Acacia koa and A. implexa (Forster et al. 2016; Rose et al. 2019). In A. koa, and other Acacia, bipinnately compound leaves contribute to rapid early seedling and sapling growth in forest gaps where light is, or will become, potentially limiting at low carbon cost (Brodribb and Hill 1993; Craven et al. 2010; Pasquet-Kok et al. 2010). The existence of plasticity in the timing of phase changes means that different environments could act to hasten, delay or increase plasticity in the timing of phase transition (Pigliucci 1997). Differences in the timing of phase change are heritable, and both A. koa and A. melanoxylon exhibit population-level differences in the timing of phase transition that are genetic (Farrell and Ashton 1978; Rose et al. 2019). In populations of A. koa exposed to regular drought, the phase change from bipinnate leaves to phyllodes occurs early in seedling growth, regardless of weather conditions experienced by seedlings, and so contributes to maximising survivorship in harsh, drought-prone sites (Rose et al. 2019).

MicroRNAs are evolutionarily conserved regulators of vegetative phase change, and play important roles in plant development (Lauter et al. 2005; Wu and Poethig 2006; Chuck et al. 2007; Wu et al. 2009; Cho et al. 2012). Wang et al. (2011) studied the transition from bipinnate to phyllodinous foliage during seedling ontogeny, and identified a likely causative signalling mechanism involving microRNA, wherein a reduction in microRNA expression correlates with a phase change to phyllode expression. The phase transition from bipinnate leaves to phyllodes in seedlings of A. confusa and A. colei is governed by opposing expression patterns of the microRNAs miR156 and miR172, which regulate a conserved framework of phase changes in many, if not all, angiosperms (Wang et al. 2011). miR156 represses the expression of squamosa promoter binding-like (SPL) transcription factors, shows high correlation with juvenile-like vegetative leaf traits, and its expression decreases during the juvenile to adult transformation (Wang et al. 2008, 2011). It follows, then, that overexpression of miR156 in transgenic Populus was shown to prolong a juvenile phase (Wang et al. 2011).

A related observation in flowering plants, that of a prolonging of a juvenile phase in shade, is also consistent with a role for sugar-signalling in miR156 expression, with high sugar levels signalling a down-regulation in the expression of miR156 (Yang et al. 2013; Yu et al. 2013). Nutritional status, therefore, plays a demonstrable role in phase change in plants, as it does in animals (Chitwood and Sinha 2016). Wang et al. (2011) did not study microRNA expression in bipinnate Acacia; however, if miR156 expression does govern the transition to phyllodes, then simply fixing the overexpression in Acacia seedlings could prolong indefinitely the bipinnate phase of ontogeny. This possibility is also consistent with fixed overexpression of miR156 in transgenic Populus × canadensis, which results in a prolonged juvenile phase (Wang et al. 2011).

Fixed overexpression of miR156 may alone explain leafy bipinnate adult foliage in bipinnate Acacia. However, whether patterns of miR156 expression in other bipinnate species of Caesalpinioideae, such as, for example, Paraserianthes or Albizia, are the same as in Acacia is not yet known. Given the flexibility in evolutionary transitions among leafy and phyllodinous adults in the Acacia phylogeny, the regulation of miR156 expression may have fundamentally changed in bipinnate species. Measuring miR156 expression in bipinnate Acacia species would be a first step in testing the role of this micro-RNA in governing adult leaf morphology. If it is overexpressed throughout the life cycle of bipinnate species, then silencing it may allow the expression of ‘hidden’ phyllodinous phenotypes in normally bipinnate species.

Presumably, the underlying capacity for phyllode expression is maintained by species with compound leaves. Acacia species with bipinnate adult foliage may represent a special case where selection has acted to delay phase transitions entirely; however, the genetic and regulatory mechanism for this is, as yet, incompletely understood. The transition to bipinnate adult foliage is predominantly a relatively recent phenomenon, with the majority of leafy ancestral nodes belonging to the most diverse and predominantly south-eastern Australian leafy lineage (comprising the bulk of section Botrycephalae), having median node ages less than 5 million years, and half having median ages between 3 and 1.5 Ma. This radiation was possibly promoted by wetter climatic conditions during a humid paleoclimate interval between 5.5 and 3.3 Ma that was identified to have occurred over Australia’s north-west (Karatsolis et al. 2020), and by the continental margin, more generally (Miller et al. 2012; Christensen et al. 2017).

Given that heteroblasty confers multiple benefits across different developmental stages in many plants (Jaya et al. 2010; Zotz et al. 2011), including Acacia (Pasquet-Kok et al. 2010; Wang et al. 2011), a key question is under what circumstances might heteroblasty be advantageously dispensed with? Occurrence of bipinnate species is significantly correlated with low moisture seasonality and low sunlight levels during winter (Table 4, Fig. 5), implying that most bipinnate species grow on soils that not only exhibit low moisture seasonality, but are also wetter on average than many Australian soils in absolute terms. This is as expected from physiological studies of bipinnate leaves and phyllodes (Brodribb and Hill 1993), with bipinnate leaves performing poorly when subject to severe water limitation. These physiological performance parameters correlate with continental-scale occurrence of bipinnate species within the totality of Australia’s climate space, across all Acacia and after correction for phylogeny. The significant correlation with winter light levels is also interesting, and what we might expect given the photosynthetic efficiency of leaves. Bipinnate leaves are more efficient photosynthetic organs than are phyllodes in terms of their maximum photosynthetic rate per unit mass (Brodribb and Hill 1993), which may confer a significant advantage to seedling and sapling growth rates (Morris et al. 2011). The higher mass of phyllodes is due in part to their investment in water-delivery tissues; primary nerve densities correlate with hydraulic conductance, and mass per area and bundle sheath extensions. Higher nerve density decreases the distance water travels within the phyllode, and phyllodes with high mass per area are more common in arid environments (Sommerville et al. 2012). Acacia species with bipinnate leaves have faster growth rates than do phyllodinous species, at least in part owing to the lower resource investment into water-delivery and support tissues in leaves (Atkin et al. 1998).

Bipinnate leaves should also be more efficient at capturing photosynthetically active radiation in low-light environments, given their horizontal orientation and spread of photosynthetic tissue (Brodribb and Hill 1993). This is not to say that phyllodinous species do not also occur in areas of low coldest-quarter radiation; they do, and, indeed, have the lowest mean values for this incident variable, below the bipinnate species clustered in the lower half of this variable. However, several observations from experimental studies of growth, development and ecology attest to the light-capturing efficiency of bipinnate leaves. Seedlings of Acacia implexa produce more leafy nodes under lower light levels (Forster et al. 2011), as do seedlings of Acacia koa (Rose et al. 2019), responses that reflect the nutritional status of the seedlings. Studies of invasive populations of the bipinnate species Acacia dealbata have shown that this species has high photosynthetic performance under low light, and considerable photosynthetic plasticity, both of which contribute to the shade tolerance of this species and, in turn, its invasiveness under established canopies outside its native range (Aguilera et al. 2015). Acacia dealbata invades forests in Europe (Rodríguez et al. 2017), Chile and Argentina (Spalazzi et al. 2019).

Phyllodinous Acacia species also occur at sites of low winter radiation, but they are not also restricted to sites of high coldest-quarter soil moisture. The occurrence of both phyllodes and bipinnate adult foliage in similar environmental space may represent a solution to simultaneous trade-offs along many environmental variables. For example, light and water availability impose challenges that interact within the Australian climate context.

Within the climatic envelope in which bipinnate species are found, they are compressed towards the high coldest-quarter soil-moisture values (Fig. 5, Bio35). This contrasts with their scattered distribution along the driest quarter soil-moisture axis, for which the correlation between bipinnate species mean values and lowest-quarter moisture index was not significant (Fig. 5). Perhaps, surprisingly, this suggests that many bipinnate species have an innate capacity to tolerate drought, possibly by mechanisms of avoidance, as suggested by resource investment in root biomass under water stress (Forster et al. 2016; Rose et al. 2019). An innate drought tolerance ought not preclude the occupancy of lower coldest-quarter moisture environments by bipinnate Acacia species. Yet, on average, they do not occur in these lower coldest-quarter moisture environments. This restriction may have an explanation in biotic factors that affect how competitive bipinnate Acacia species are with other woody species. For example, a period of active winter growth during times of relatively low evapotranspiration stress may provide a competitive advantage for bipinnate Acacia over other woody species. Conversely, there are three bipinnate species among those sampled that do occur in environments with high seasonal variation in soil moisture. Acacia epacantha, A. fagonioides and A. guinetti are all spinescent shrubs from south-western Western Australia, with bipinnate leaves that are smaller in size than most other species of Acacia, comprising a single pair of pinnae, each with 2–5 pairs of pinnules. These three species are nested within a group of bipinnate species that are sister to most of the Acacia phylogeny, and were traditionally recognised in section Pulchellae (Benth.) Taub. Each of these three species deploys bipinnate foliage in quite different ecological contexts from others in the current study, and a comparative investigation of their performance response to water deficit and stress could be rewarding.

Phyllodes and bipinnate adult foliage are ecologically important traits that are optimised to different climatic contexts (Renner et al. 2020). Physiological and life-history studies have previously demonstrated their optimisation within different life-history stages at the individual level, and these studies, along with our current results, suggest that the different performances of leaves and phyllodes have ecological consequences from the individual (as demonstrated by previous physiological studies) through to the macroecological and continental scale (as found here).

Within Acacia, bipinnate adult leaves have at least three separate origins within the genus. Multiple origins for bipinnate foliage were previously identified when the leafy section Botrycephalae (Benth.) Taub. was resolved as polyphyletic (Brown et al. 2006). Multiple reversals to bipinnate foliage have been inferred (Murphy et al. 2003), with neoteny being the likely developmental pathway for bipinnate adult leaves (Murphy et al. 2010). However, in the present study, it was newly discovered that transformations from phyllode to leafy adult foliage may have occurred in as few as three lineages within Acacia, and that each transformation to leaves is linked with many subsequent transformations back to phyllodes. This is consistent with Bayesian estimation of transition rates under the optimal model of trait evolution, in which the rate of transitions from phyllodinous to bipinnate adult foliage across the Acacia phylogeny was low.

Given the regulation of phyllode expression and transition from juvenile to adult foliage by miRNA156, it is likely that paedomorphosis describes each evolutionary transition to bipinnate adult foliage in Acacia; however, whether post-displacement or neoteny (sensu Alberch et al. 1979) is involved will not be resolved until the underlying gene-regulation pathway is better characterised. The existence of plasticity in the timing of phase-change initiation, pre- and post-displacement at the level of individual (see Pryer and Hearn 2009 for an explanation of terminology), suggests that evolutionary lability in phase change is possible as a result of fixation of plasticity in timing that exists within individuals. However, whether the same gene regulatory mechanism(s) are involved in each case of paedomorphosis in Acacia is currently unknown. We may hypothesise that it is not necessarily the same pathway for all bipinnate Acacia. The diversity of adult bipinnate leaf morphologies, phylogenetic distance and evolutionary time separating bipinnate lineages, and, in particular, the reduced leaves of species in section Pulchellae compared to those species in south-eastern Australia, suggest that more than one paedomorphic mode may well be involved.

Bipinnate adult foliage represents a derived character in extant Acacia, even though it is a phylogenetic transition to a symplesiomorphic state within the broader Caesalpinioideae. Bayesian estimation of ancestral states suggests that the most recent common ancestor of Acacia had phyllodinous adult foliage, and the transition to phyllodes preceded the radiation of Acacia. However, this reconstruction is based on the assumption that all phyllodes are homologous, and there is some evidence that suggests that they are not (Gardner et al. 2005). Studies of ontogeny, such as those on A. mangium, inform homology relations between phyllodes and bipinnate leaves (Leroy and Heuret 2008); however, as yet, no decisive arbitration among the three competing hypotheses of the Acacia phyllode origin has been achieved (Gardner et al. 2008). Character coding that better reflects the morphological diversity, and possible homology relations among phyllodes, may result in different ancestral-state reconstruction at the Acacia most recent common ancestor, as hinted at by Bayesian comparisons of the alternative states fixed at the Acacia crown node. This could be worth pursuing on the basis of a detailed appraisal of phyllode diversity.

No new data were generated for this study; however, the molecular sequence data used in this study are available from GenBank https://www.ncbi.nlm.nih.gov/genbank/ and see the voucher information provided at https://doi.org/10.1038/ncomms5473; the bioclimatic data are available at the Atlas of Living Australia http://www.ala.org.au, the Supplementary material associated with this publication, or from the authors by request, the morphological data are as reported in regional flora treatments and worldwide wattle, see http://worldwidewattle.com, the character matrix for which is also available in the Supplementary material associated with this publication, of from the authors by request.

This study received no specific funding.

Dr Daniel Murphy is an editor for Australian Systematic Botany but did not at any stage have editor-level access to this manuscript while in peer review, as is the standard practice when handling manuscripts submitted by an editor to this journal. Australian Systematic Botany encourages its editors to publish in the journal and they are kept totally separate from the decision-making processes for their manuscripts. The authors have no further conflicts of interest to declare.