Deep evolutionary comparison of gene expression identifies parallel recruitment of trans-factors in two independent origins of C4 photosynthesis

Abstract

With at least 60 independent origins spanning monocotyledons and dicotyledons, the C4 photosynthetic pathway represents one of the most remarkable examples of convergent evolution. The recurrent evolution of this highly complex trait involving alterations to leaf anatomy, cell biology and biochemistry allows an increase in productivity by ∼ 50% in tropical and subtropical areas. The extent to which separate lineages of C4 plants use the same genetic networks to maintain C4 photosynthesis is unknown. We developed a new informatics framework to enable deep evolutionary comparison of gene expression in species lacking reference genomes. We exploited this to compare gene expression in species representing two independent C4 lineages (Cleome gynandra and Zea mays) whose last common ancestor diverged ∼ 140 million years ago. We define a cohort of 3,335 genes that represent conserved components of leaf and photosynthetic development in these species. Furthermore, we show that genes encoding proteins of the C4 cycle are recruited into networks defined by photosynthesis-related genes. Despite the wide evolutionary separation and independent origins of the C4 phenotype, we report that these species use homologous transcription factors to both induce C4 photosynthesis and to maintain the cell specific gene expression required for the pathway to operate. We define a core molecular signature associated with leaf and photosynthetic maturation that is likely shared by angiosperm species derived from the last common ancestor of the monocotyledons and dicotyledons. We show that deep evolutionary comparisons of gene expression can reveal novel insight into the molecular convergence of highly complex phenotypes and that parallel evolution of trans-factors underpins the repeated appearance of C4 photosynthesis. Thus, exploitation of extant natural variation associated with complex traits can be used to identify regulators. Moreover, the transcription factors that are shared by independent C4 lineages are key targets for engineering the C4 pathway into C3 crops such as rice.

Figure 1. The C₄ maturation gradient in leaves of Cleome gynandra.

Venation, bundle sheath cell (BS) size, mesophyll (M) cell size and abundance of C₄ transcripts and proteins in the base, middle and tip of 3 mm leaves as well as fully mature leaves of C. gynandra. (A) Leaves of 3 mm length possess a gradient in venation density from base to tip, whereas in mature leaves (B) this gradient is no longer visible, insets show representative images of samples used for RNA isolation. (C) Quantification of venation density and complexity. (D) Transverse sections and quantification of BS and M cell size. (E) Quantitative RT-PCR for the CA4, PPC2, NAD-ME2 and PPDK of genes important in the C₄ cycle. (F) Abundance of carbonic anhydrase, phosphoenolpyruvate carboxylase, NAD-dependent malic enzyme and pyruvate,orthophosphate dikinase proteins from the base (B), middle (M), tip (T) and mature (Mat) leaves. Scale bars in A and B represent 0.3 mm and 3 mm respectively, while 1 mm gradations are shown within the insets.

Figure 2. Overview of the workflow and results of the conditional orthology assignment method.

Identification of homologues and quantification of gene expression after de novo assembly, for full details see Text S1. (A) Correlation in quantification derived from reciprocal best BLAST (RBB) hits in the de novo assembly and reference summed over all transcript isoforms per reference gene locus. (B) The Spearman correlation in transcript abundances between the reference guided estimation and estimates generated using different transcript orthology assignment methods on the same de novo assembled transcriptome. “RBB only” means that only the reciprocal best BLAST transcripts were selected. E-value cut-offs (e.g. 1e-5) indicate the fixed value at which sequences were determined to be homologues. OrthoMCL indicates that OrthoMCL was used to cluster and identify orthologous transcript groups. Finally, the black bar indicates the effect of varying the percentile cut-off on the abundance estimate accuracy of the conditional orthology assignment method. (C) Conditional orthology assignment method begins by performing all versus all BLAST searches of the assembled transcripts against a reference proteome. (D) The reciprocating hits (indicated by blue lines) are selected for self-training. (E) The reciprocating hits are binned according to assembled transcript length and a quadratic model is fit to the e-value and length data. (F) Non-reciprocating hits which fall above the curve are accepted as putative homologues, non-reciprocating hits which fall below the curve are rejected. (G) Correlation in quantification derived from conditional assigned transcripts using species own reference genome. (H) Correlation in quantification derived from conditional assigned transcripts using intermediary reference genome. For full details, validation and explanation please see the supplementary methods (Text S1).

Figure 3. Convergence in patterns of gene expression in leaf gradients of C. gynandra and maize.

(A) Venn diagram indicating numbers of shared and unique transcripts to each type of C. gynandra leaf tissue. (B) Major bin categories identified using Wilcoxon test implemented in Pageman tool that alter between the base, middle, tip of 3 mm and mature C. gynandra leaves. (C) Number of genes with ascending (red) and descending (grey) behaviours as leaves of C. gynandra (Cg) and maize (Zm) mature. (D) Venn diagrams depicting the total number of transcript homologues that increase or decrease in abundance as leaves of both C. gynandra and maize mature. The number of genes common to the two gradients is shown in blue, with the number of transcription factors shown in parentheses. Red circles and numbers correspond to genes that increase in abundance, while grey circles represent genes that show reduced abundance.

Figure 4. Classification of gene expression in the two C₄ species C. gynandra and maize.

As leaves of C. gynandra (A) and maize (B) mature, transcripts were classified into twenty-six behaviours, thirteen ascending (A&B) and thirteen descending. Statistically significant differences between neighbouring tissue types are delineated by red circles in ascending filters. The total number of genes within each behaviour is presented in parentheses and behaviours containing photosynthesis-related genes are annotated by red boxes around each plot (eg A3, A4 and A6). Genes of the core C₄ cycle occupy six and five of the thirteen ascending filters in C. gynandra and maize respectively (transcripts in green). (C) Venn diagram representing transcription factors showing the same behaviours as C₄-related genes in the maize and in C. gynandra leaf gradients. (D) Behaviour of homologous genes in C₃ A. thaliana.

Figure 5. Convergence of mesophyll and bundle sheath transcriptomes in C. gynandra and maize.

(A) Schematic showing M or BS accumulation of transcripts involved in the C₄ cycle. Shared parts of the pathway are annotated in red, while differences between the species are shown in grey. CA, carbonic anhydrase; PPC, phosphoenolpyruvate carboxylase; PEPC Kin, phosphoenolpyruvate carboxylase kinase, ASPAT, aspartate aminotransferase; ALAAT, alanine aminotransferase; PPDK, pyruvate-orthophosphate dikinase; TPI, triose phosphate isomerase; PGK, phosphoglycerate kinase; FBA, fructose-bisphosphate aldolase; SBP, sedoheptulose-bisphosphatase; TKL, transketolase; PRK, phosphoribulokinase; RbcS, RubisCO small subunit; RCA, RubisCO activase; FBP, fructose 1,6-bisphosphate phosphatase; RPE, D-ribulose-5-phosphate-3-epimerase; NAD-ME, NAD-dependent malic enzyme, MDH malate dehydrogenase. (B) Venn diagrams representing transcripts expressed in M (left panel) and BS (right panel) of C. gynandra and maize. Cell-specific maize data represents the overlap between two independent experiments , . (C) Venn diagrams of transcription factors expressed in M or BS in maize and C. gynandra. (D–G) Expression in M and BS cells of the 18 homologous transcription factors showing co-ordinated induction with C₄ photosynthesis genes during leaf maturation of both maize and C. gynandra. Abbreviations: Cg data from C. gynandra (this study), while Zm1 data are from Li et al (2010) and Chang et al (2012) respectively.