bri1-5/bak1-1D, bri1-5/brs1-1D and bri1-5/bri1-1D partially reconstitute bri1-5 gene expression
To better understand the molecular mechanisms of key BR signaling genes, we performed a phenotypic screening and expression analysis of bri1-5 and its three activation-tag suppressors along with their corresponding wild type. Two suppressor strains bri1-5/bak1-1D, bri1-5/brs-1D were obtained from [12] and an additional bri1-5/bri1-1D mutant was selected in the framework of the current study (see Materials and methods). Sequencing the BRI1 flanking region from the suppressor bri1-5/bri1-1D showed that the activation tag was inserted 534 bp downstream of BRI1 gene (Figure 1, A). All suppressor mutants were shown to indeed overexpress the activation tagged gene as confirmed by Real-Time qPCR (RT-qPCR) (Table S1). Phenotypically, all light-grown bri1-5 suppressors (bri1-5/bak1-1D, bri1-5/brs1-1D and bri1-5/bri1-1D) lines displayed larger seedlings and enhanced epidermal cell length as compared to the bri1-5 mutant (Figure 1-B, Figure 2). We observed that the bri1-5/bak1-1D line best approximated wild type epidermal cell length (Figure 2). To gain insight into which pathways in each of the studied lines were responsible for the restoration of the bri1-5 phenotype towards the wild type phenotype, we performed gene expression analysis. All suppressor lines, together with their corresponding WT and with the bri1-5 background were sampled atTo assess the reproducibility of the expression analysis, we measured the extent to which the expression profiles of replicate samples were similar using Principal Component Analysis (PCA): PCA indeed showed that most of the variation in gene expression between the samples could be assigned to differences in genetic background and not to differences between replicates, confirming the reproducibility. (Figure S1). In addition, microarray results were confirmed using RT-qPCR for a selected set of genes (Figure S2).
We determined for each mutant line its differential expression versus the same common reference i.e. the expression state in WS2 (Figure 3) and used these values to assess to what extent the different lines share the same/different differentially expressed genes (aberrantly expressed versus the WT control) (see vendiagram, Figure 3). In this vendiagram genes that for instance are differentially expressed in the bri1-5 mutant, but not in any of the suppressor lines represent the genes that are affected by the bri1-5 mutation and that could be restored by overexpression of any of the genes being overexpressed in the suppressor lines. Genes that are in the intersection of a suppressor line and bri1-5 represent genes of which the aberrant expression in the bri1-5 background could not be restored in the respective suppressor line with which the overlap was observed. The venndiagram represented in Figure 3 and the scatter-plots in Figure S3 (A-C) show that of all suppressor lines, bri1-5/bak1-1D could restore the largest number of genes that were affected in expression in the bri1-5 (about two-thirds of the genes that were differentially expressed in bri1-5 were no longer differentially expressed in bri1-5/bak1-1D). This is in line with its observed phenotypic behavior as indeed bri1-5/bak1-1D seems to also phenotypically best compensate for the bri15 mutation. Both the bri1-5/brs1-1D and bri1-5/bri1-1D lines could recover the same one-third of bri1-5 affected genes (see vendiagram Figure 3). The latter is also illustrated in Figure S3 panel A-C which show that bri1-5/brs1-1D and bri1-5/bri1-1D (R2 = 0.20) are the most similar from an expression impact point of view (i.e. affecting the same genes), followed by bri1-5/brs1-1D and bri1-5/bak1-1D (R2 = 0.12), whereas bri1-5/bri1-1D and bri1-5/bak1-1D display the lowest similarity (R2 =0.029) (Figure S3 D-F). This high similarity between brs1-1D and bri1-1D indicates that overexpression of BRI1 and BRS1 in a bri1-5 background restores the expression of the same genes and hence that BRI1 and BRS1 must have close roles in the BR signaling pathway.
Identification of restoring and compensatory pathways
We assumed that if the suppressor strains compensate the phenotype of the bri1-5 mutant, they could do so because they either restore the pathways disrupted in the bri1-5 mutant to wild type levels or they alter the pathways that compensate for the bri1-5 affected pathways. Both mechanisms could be potentially reflected in the expression data. To identify restored and compensatory pathways shared by all suppressor mutants, we first focused on genes with altered expression in the bri1-5 mutant, but not in any of the suppressor mutants (group A, 118 presumably restored genes) and vice versa genes that were not differentially expressed in the bri1-5 mutant, but in all suppressor strains (group B: 23 presumable compensatory genes) (Figure 3). Next to identifying genes involved in compensatory/restoring pathways shared by all mutants, we also extracted the genes that would be involved in brs1-1D specific compensatory mechanisms (groups C, 333 genes) because as compared to BAK1 and BRI1, the role of the BRS1 gene in BR signaling is yet less characterized. Another gene set we analyzed in-depth consists of the genes that are affected in the bri1-5 mutant, but that were not recovered by any of the suppressors (149 genes, group D). These genes could explain the residual phenotypic differences between the bri1-5 mutant, the suppressor strains and the wild type. To identify how the genes of the compensatory and restoring groups act together in pathways, we performed GO and network analysis (see Materials and Methods).
Restoring and compensatory pathways activated in all suppressor lines
Genes restored to wild type state by all suppressors (Group A) were enriched for ‘defense response’ and ‘hormone response’ (Additional file 1). Compensatory genes, (group B, i.e. genes that are significantly differentially expressed in all three bri1-5 suppressors, but not in bri1-5 mutant) are enriched in ‘response to abscisic acid (ABA) signaling’ including six responsive genes to ABA (DTX50, HVA22D, GRP23, COR15B, SAG113, HAB1), indicating that ABA signaling has been affected to compensate the bri1-5 deficiency in all suppressor strains (hypergeometric test, p-value =5.5e-06, Additional file 1 ). To identify the link between the restoring and compensatory pathways, we mapped the genes of both group A and B on the interaction network and used Phenetic to extract sub-networks. Phenetic identifies 6 sub-networks connecting the genes of group A and B. Note that these subnetworks contain, next to the genes of group A and B also connector genes. These are genes that are not differentially expressed themselves but that are still recovered by the network analysis because of their high connectivity with genes of group A/B. As they are needed to connect genes of group A/B, they are most likely involved in the same process as represented by genes of group A/B. Five of these subnetworks showed enrichment in known GO functions (being enriched in respectively ABA, ethylene, auxin, cytokinin and ROS signaling (Figure 4). This indicates that these pathways contribute to recovering bri1-5 signaling deficiency.
In-depth analysis shows that the subnetwork enriched in ABA signaling (Figure 4, network 1) contains two main negative regulators of ABA signaling,HAI1 and HAB1, which are up-regulated in all bri1-5 suppressor lines compared to wild type, but not differentially expressed in the bri1-5 mutant (compensatory genes). Interestingly, the list of compensatory genes (group B) contains 6 targets of the ABA signaling pathway (enrichment p-value: 5.5e-6; DTX50, HVA22D, GRP23, COR15B, SAG113, HAB1), indicating that ABA signaling has indeed been negatively affected in the suppressor strains to compensate for the bri1-5 signaling deficiency. These 6 additional genes could not be connected by Phenetic on the interaction network, implying they are quite distantly located from each other in the network and hence likely belong to different biological pathways.
The aforementioned negative regulators of ABA signaling, HAI1 and HAB1, belong to the protein phosphatase 2C (PP2C) gene family which has nine members in total (HAI2, HAB2, HAB1, HAI3, PP2CA, ABA1, AHG1, ABI2, AHI1). In our analysis, most of the PP2C genes appeared to be differentially expressed in at least two suppressors, but they failed to pass the strict FDR threshold for the third suppressor to be included in group B (Table S3). PP2C is known to repress ABI5, the main activator of ABA signaling [23]. PP2C is also known to represses BIN2 activity [6, 7]. As BIN2 activates ABI5 by phosphorylating SnRKs [6, 7], repressing ABA signaling by PP2C via blocking SnRKs phosphorylation seems to compensate for the deficiency in BRI1 mediated signaling. The subnetwork enriched in ABA signaling (network 1) also contains members of the PYR/PYL/RCAR family as connector genes (RCAR10, SRK2B, SRK2G, RCAR4, RCAR13, PYL8, RCAR14, RCAR12, RCAR11, PYL9, RCAR6, SRK2D). The PYR/PYL/RCAR family constitutes the receptor of ABA signaling and promotes the activation of SnRKs by repressing PP2C [15, 24]. The fact that the SnRKs and PYR/PYL/RCAR genes were identified as connector genes implies that they are likely involved in the pathways that connect the genes of group A/B. They are most likely not primarily regulated at expression level given their role in phosphorylation-mediated signaling [7, 25]. This explains why they were detected as connector genes and not retrieved by differential expression analysis.
We could not find any link in the literature to explain how PP2C can be up-regulated by brassinosteroid signaling in order to repress ABA signaling. Our network identified a TF (DREB1B) of which the aberrant expression in the bri1-5 mutant is restored in all suppressor lines. According to the KEGG pathway, DREB1B directly regulates PYR/PYL/RCAR [26]. Hence, this TF might the missing link between PPC2 and BR signaling. This aforementioned TF could also be a good candidate to explain the crosstalk between ABA and ethylene signaling (subnetwork 2) in response to BR signaling as this TF links the subnetwork 1 to subnetwork 2. Indeed the ethylene signaling pathway seems to be resorted by all suppressors (Figure 4, subnetwork 2).
Next to the ABA and ethylene subnetwork, also subnetworks related to other hormone signaling processes like auxin (subnetwork3), cytokinin (subnetwork5) and ROS signaling (subnetwork6) were detected. It is well known that crosstalk between these phytohormone signaling pathways exists [27, 28]. Hence interfering with one pathway e.g. ABA signaling through BR signaling might affect other phytohormone pathways as well. Overall the subnetwork enriched in ABA signaling contains mostly genes of group B (compensatory genes) whereas the subnetworks related to the other phytohormone signaling pathways belong to group A (restoring genes). Based on their differential expression behavior we can conclude that ABA signaling is mostly a compensatory pathway, whereas the other hormone signaling pathways contribute to restoring pathways that were affected in the bri1-5 strain.
Negative feedback of BR signaling on BR biosynthesis:
Group D (genes differentially expressed in the bri1-5 lines, but not fully restored by any suppressors) consisted of four cytochrome P450 genes (CYP90C1, DWF4/CYP90B1, CYP85A2 and CYP90A1). These were aberrantly up-regulated in bri1-5 and in all three suppressor mutants, although to different degrees (gradually stronger upregulated in bri1-5 than in the suppressors). This indicates that the aberrant expression of these genes in the bri1-5 mutant could only be partially compensated for in the suppressor mutants. Cytochrome P450 genes play a role in BR biosynthesis by converting the sterol ‘campesterol’ to BRs [29]. So the fact that mutations in BR signaling genes affect the expression of BR biosynthetic genes indicates that a negative feedback exists of BR signaling on BR biosynthesis.
If indeed a negative feedback exists between BR signaling and biosynthesis, this feedback should be reflected in quantitative differences in overexpression of the BR signaling and biosynthesis genes in the bri1-5 and suppressor mutants. The better the signaling can be restored in the suppressors (as reflected by the phenotype), the less we expect the expression of the BR biosynthesis to be aberrant. We indeed found that the expression of the BR-biosynthesis genes (CYP90C1, CYP90A1, CYP85A2, CYP90B1) are less affected in the strains that better mimic the wild type phenotype (see Figure 5, the best suppressor of bri1-5, bri1-5/bak1-D, shows the lowest upregulation of the biosynthesis genes). This further supports the existence of negative feedback from BR regulation on BR biosynthesis: a more sustained BR signaling results in decreased BR biosynthesis, whereas suboptimal BR signaling is compensated for by higher transcriptional activity of BR biosynthetic genes.
BRs possibly involved in providing an optimal environment for BRI1 and ligand binding
Unlike for BRI1 and BAK1, much less is known about the role of BRS1 in BR signaling. Here we had a closer look at genes of group D which are exclusively differentially expressed in bri1-5/brs1-1D mutant and hence comprise compensatory pathways specific for bri1-5/brs1-1D. GO enrichment showed that the genes uniquely involved in compensating for the bri1-5 mutant by the brs1-1D suppressor (in group C, 333 genes form Venn diagram), are overrepresented in iron ion homeostasis/ferroxidase activity (down-regulated) (p-value 3.0e-09), and glutathione transferase (up-regulated) (Figure 6). As in the ferroxidase reaction, four H+ are used to catalyzes the oxidization of Fe2+ to Fe3+, repressing this reaction results in the accumulation of H+ which can be transported to the apoplast via plasma-membrane pumping mediated by ATPase (H+-ATPase transporters) [30]. Accordingly, we found that the main inhibitor of H+-ATPase transporters, CBC1, was significantly down-regulated in bri1-5/brs1-1D (log-fold -1.7, adj p-value 8.36e-06) but not in the other suppressors. This implies that H+-ATPase transporters are more active in bri1-5/brs1-1D to export H+ from cytosol into apoplast and making the apoplast more acidic. In line with this hypothesis, the up-regulated glutathione transferase activity in the brs1-1D mutant might be essential to compensate for the more acidic environment and would be required for redox homeostasis. The observed acidification could generate a cellular environment that improves BRI1-BRs binding or BRI1-BAK1 dimerization and hence restoration of the bri1-5 mutant phenotype.
Link between stress response and BR signaling
According to literature, there is crosstalk between BR signaling and other hormone signaling in response to stress, especially via ABA and auxin signaling [27]. In the absence of BRs (or low amounts of BRs), BIN2, activates ABA and auxin signaling, resulting in the induction of stress response genes [6, 7, 15]. On the other hand, some stress-response genes are known to be targets of BZR1 and BES1 [7, 15] (Figure 7), indicating that also when BR levels are high, stress response genes can be activated. These observations show that balanced BR levels are needed for normal growth and that deviation from the optimal levels (either too high or too low) would activate stress response mechanisms. We observed that by partially recovering bri1-5 signaling deficiency by suppressors mutants, the transcript level of some stress-response genes is restored to normal but other stress-response genes become induced (Additional file 1: GO enrichment for genes exclusively differentially expressed in each suppressor from the venndiagram, “GO_only_bri1-1D”, “GO_only_bak1-1D”, “GO_only_brs1-1D”). This observation is in line with this complex effect of BR and BR signaling on stress response pathways.