Heterotrimeric G‐proteins mediate crucial and diverse signaling pathways in eukaryotes. Here, we generate and analyze microarray data from guard cells and leaves of G‐protein subunit mutants of the model plant Arabidopsis thaliana, with or without treatment with the stress hormone, abscisic acid. Although G‐protein control of the transcriptome has received little attention to date in any system, transcriptome analysis allows us to search for potentially uncommon yet significant signaling mechanisms. We describe the theoretical Boolean mechanisms of G‐protein × hormone regulation, and then apply a pattern matching approach to associate gene expression profiles with Boolean models. We find that (1) classical mechanisms of G‐protein signaling are well represented. Conversely, some theoretical regulatory modes of the G‐protein are not supported; (2) a new mechanism of G‐protein signaling is revealed, in which Gβ regulates gene expression identically in the presence or absence of Gα; (3) guard cells and leaves favor different G‐protein modes in transcriptome regulation, supporting system specificity of G‐protein signaling. Our method holds significant promise for analyzing analogous ‘switch‐like’ signal transduction events in any organism.
Heterotrimeric G‐proteins, composed of α, β, and γ subunits, participate in a wide range of signaling pathways in eukaryotes (Morris and Malbon, 1999). According to the typical, mammalian paradigm, in its inactive state, the G‐protein exists as an associated heterotrimer. G‐protein signaling begins with ligand binding that results in a conformational change in a G‐protein‐coupled receptor (GPCR). Once activated by the GPCR, the Gα separates from the associated Gβγ dimer and the freed Gα and Gβγ proteins can then interact with downstream effector molecules, alone or in combination, to transduce the signal. Subsequent to signal propagation, Gα re‐associates with the Gβγ dimer to reform the G‐protein complex.
There are several classical routes for signal propagation through heterotrimeric G‐proteins that have been categorized in mammalian systems (Marrari et al, 2007; Dupre et al, 2009). One route, which we designate classical I, requires the presence of both subunits, and can invoke one of two distinct mechanisms. In one mechanism, on GPCR activation, freed Gα and Gβγ each interact with downstream effectors to elicit the downstream response. In a related mechanism, Gα but not Gβγ interacts with downstream effectors, but the Gβγ dimer is nevertheless required to facilitate coupling of Gα with the relevant GPCR (Marrari et al, 2007). In a second route, which we designate classical II, it is solely the Gβγ dimer that interacts with downstream effectors; in this case, sequestration of Gβγ within the heterotrimer prevents signal propagation. In addition, a few non‐classical G‐protein regulatory modes have also been implicated in some systems, for example signaling by the intact heterotrimer in yeast (Klein et al, 2000; Frank et al, 2005). Observations such as these lead to a fundamental question, namely, which of all the theoretical regulatory modes of G‐protein signaling are realized biologically. Our study answers this question in the context of the model plant Arabidopsis thaliana, and in addition analyzes the manner in which G‐protein signaling couples with signaling by the plant hormone abscisic acid. The Arabidopsis genome encodes only one canonical Gα subunit, GPA1, and one canonical Gβ subunit, AGB1, and knockout mutants are available for each of these, allowing clear dissection of Gα‐ and Gβ‐related phenotypes.
Abscisic acid (ABA) is a major plant hormone, which inhibits growth and promotes tolerance of abiotic stresses such as drought, salinity, and cold. ABA signaling is known to interact with heterotrimeric G‐protein signaling in both developmental and stress responses in a complex manner, causing, for example, ABA hyposensitivity of guard cell stomatal opening in gpa1 and agb1 single mutants as well as agb1 gpa1 double mutants (Fan et al, 2008), but ABA hypersensitivity of the inhibition of seed germination and post‐germination seedling development in the same mutants (Pandey et al, 2006). These experimental observations implicate G‐proteins as one of the components of ABA signaling, but to date no systematic study has been conducted in either plant or metazoan systems to define the co‐regulatory modes of a G‐protein and a hormone.
In this study, we conduct genome‐wide gene expression profiling in G‐protein subunit mutants of A. thaliana guard cells and leaves, with or without treatment with ABA. By introducing one or more mediators acting downstream of the G‐protein and ABA to control transcript levels, we propose nine G‐protein/ABA signaling pathways including ABA‐independent G‐protein signaling pathways, G‐protein‐independent ABA signaling pathways, and seven distinct ABA–G‐protein‐coupled signaling pathways (Figure 1). We develop a Boolean modeling framework to systematically enumerate 14 possible theoretical regulatory modes of the G‐protein and 142 co‐regulatory modes of the G‐protein and ABA, and then use a pattern matching approach to associate target genes with theoretical regulatory modes.
Our analysis shows that the G‐protein regulatory mode that requires the presence of both Gα and Gβγ subunits (consistent with classical I mechanisms), is well represented in both guard cells and leaves. The G‐protein regulatory mode that requires a freed Gβγ subunit (classical II G‐protein regulatory mechanism) is well supported in guard cells and somewhat less so in leaves. In addition, a G‐protein regulatory mode representing a non‐classical regulatory mechanism is prevalent in guard cells but less so in leaves (Figure 5). In this regulatory mode, signaling by Gβ(γ) occurs, and this signaling is not regulated in any way by Gα.
By relating the target genes with the nine proposed G‐protein/ABA signaling pathways, we are able to gauge the plausibility of regulatory modes of the G‐protein and ABA at the pathway level. We find that G‐protein‐independent ABA signaling pathways are prevalent in both guard cells and leaves. The existence of an ABA‐independent regulatory activity of the G‐protein is well supported in guard cells, but not supported in leaves. Additive regulation by G‐protein signaling plus G‐protein‐independent ABA signaling is rare in both guard cells and leaves. In addition, combinatorial cross‐talk between G‐protein signaling and ABA signaling and additive cross‐talk between ABA–G‐protein signaling and G‐protein‐independent ABA signaling are observed in both guard cells and leaves. Our transcriptome analysis indicates that in some cases, ABA definitely does not influence G‐protein signaling, though it may do so in some other cases.
To investigate whether previously observed hypersensitivity or hyposensitivity of developmental and dynamic transient responses to ABA in G‐protein mutants is recapitulated at the level of transcriptional regulation, we compare gene regulation by ABA in guard cells and leaves of the G‐protein mutants versus wild type. We find that in guard cells, equal ABA hyposensitivity of all mutants combined is significant, although hyposensitivity in individual mutants is not. There is also a separate group of genes in guard cells that show ABA hypersensitivity in the gpa1 mutant, suggesting complex interactions between ABA and G‐protein signaling in gene regulation in this cell type. In leaves, ABA hyposensitivity of gene expression in the three individual mutants and equal hyposensitivity in all mutants are strongly supported. In addition, several of the functional categories identified by our analysis of G‐protein regulatory modes have been implicated in previous physiological analyses of G‐protein mutants, providing validation to the biological interpretation of our results.
In summary, by conducting a genome‐wide gene expression profiling study in G‐protein subunit mutants of A. thaliana guard cells and leaves and developing a Boolean modeling framework, we systematically evaluate the biological utilization of mechanisms of G‐protein regulatory action and reveal novel regulatory modes of the G‐protein. The results generate empirical evidence and insights regarding molecular events of G‐protein signaling and response at the physiological level in both plants and mammals.
Classical mechanisms of heterotrimeric G‐protein signaling are observed to function in regulation of the transcriptome. Conversely, many theoretical regulatory modes of the G‐protein are not manifested in the transcriptomes we investigate.
A new mechanism of G‐protein signaling is revealed, in which the β subunit regulates gene expression identically in the presence or absence of the α subunit.
We find evidence of cross‐talk between G‐protein‐mediated and hormone‐mediated transcriptional regulation.
We find evidence of system specificity in G‐protein signaling.
Heterotrimeric G‐proteins, composed of α, β, and γ subunits, participate in a wide range of signaling pathways in eukaryotes (Morris and Malbon, 1999; Jones and Assmann, 2004; Perfus‐Barbeoch et al, 2004). According to the typical, mammalian paradigm, in its inactive state, the G‐protein exists as an associated heterotrimer. G‐protein signaling begins with ligand binding that results in a conformational change in a G‐protein‐coupled receptor (GPCR). Once activated by the GPCR, the Gα protein, which possesses a GDP/GTP‐nucleotide‐binding site and GTP‐hydrolase activity, changes its form to a structure that allows exchange of GDP for GTP. The GTP‐bound Gα separates from the associated Gβγ dimer and the freed Gα and Gβγ proteins can then interact with downstream effector molecules, alone or in combination, to transduce the signal. Subsequent to signal propagation, the intrinsic GTPase activity of Gα eventually results in hydrolysis of bound GTP to GDP, which inactivates Gα and allows its re‐association with the Gβγ dimer to reform the inactive G‐protein complex.
There are several classical routes for signal propagation through heterotrimeric G‐proteins that have been categorized in mammalian systems (Morris and Malbon, 1999; Marrari et al, 2007; Oldham and Hamm, 2008; Dupre et al, 2009). On the basis of the requirement for both Gα and Gβγ subunits or for only the Gβγ subunit to propagate the signal, we group these routes into classical I and classical II categories, respectively. In one mechanism of the classical I type (designated here as classical Ia), on GPCR activation, freed Gα and Gβγ both interact with downstream effectors, either additively or synergistically, to elicit the downstream response. In a second mechanism (designated here as classical Ib), Gα but not Gβγ interacts with downstream effectors; in this case, there still is a requirement for the Gβγ dimer to facilitate correct targeting of Gα to the plasma membrane and coupling of Gα with the relevant GPCR (Marrari et al, 2007). In the classical II route, well established in mammalian and yeast cells, it is solely the Gβγ dimer that interacts with downstream effectors; in this case, sequestration of Gβγ within the heterotrimer prevents signal propagation.
Each one of these classical mechanisms leads to explicit predictions of how genetic knockout of either Gα or Gβ would affect production of the downstream phenotype. In classical Ia, knockout of either Gα or Gβ would lead to an impaired response. Either there could be a stronger impairment in the case of a double knockout of Gα and Gβ, or a single knockout could suffice to eliminate the response, if the single knockout reduced the signal to below a requisite activation threshold. In classical Ib, knockout of Gβ would prevent correct Gα localization and Gα‐GPCR coupling, whereas knockout of Gα would eliminate Gα‐effector coupling. Thus, each single knockout as well as a double knockout of Gα and Gβ would eliminate the response. In the classical II mechanism, knockout of Gα would free Gβ(γ) from sequestration, thereby eliciting the phenotype even in the absence of agonist activation of the GPCR; conversely, knockout of Gβ would defacto eliminate the possibility of eliciting the phenotype, even in the presence of the appropriate ligand.
Direct genetic tests of the above predictions are complicated in mammalian systems, due to the abundance of G‐protein subunits (e.g. 23 Gα subunits (including splice variants), 5 Gβ subunits, and 12 Gγ subunits encoded in the human genome (McCudden et al, 2005)) and their potential for combinatorial associations and partially redundant function. However, the presence of only one canonical Gα protein (GPA1), one Gβ protein (AGB1), and two identified Gγ proteins (AGG1 and AGG2) in the model plant Arabidopsis thaliana (Jones and Assmann, 2004; Assmann, 2005), and the availability of viable genetic knockout mutants lacking Gα (gpa1), Gβ (agb1), or both (agb1 gpa1 double) (Lease et al, 2001; Jones et al, 2003; Ullah et al, 2003; Chen et al, 2004), offers an excellent opportunity to experimentally test these predictions. Although biochemical analyses of plant G‐protein receptor and subunit‐effector coupling are just beginning, plant G‐proteins, like those of animals, have been shown to participate in multiple signaling and developmental processes, and phenotypic analysis of G‐protein mutants suggests that the above classical mechanisms also exist in plants. For example, some phenotypes, such as rounded rosette leaves, are exhibited similarly by both gpa1 and agb1 knockout mutants (Assmann, 2005), consistent with the classical Ia and Ib mechanisms. Other phenotypes are opposite in gpa1 and agb1 mutants, supporting a classical II mechanism. For example, agb1 mutants exhibit increased numbers of lateral roots compared with wild type, whereas gpa1 mutants show decreased lateral root production (Chen et al, 2006), and agb1 mutants exhibit impaired resistance to some plant pathogens, whereas gpa1 mutants exhibit enhanced resistance (Trusov et al, 2006).
In addition to the two classical mechanisms discussed above, a few non‐classical G‐protein regulatory modes have also been implicated in some systems, for example signaling by the intact heterotrimer in yeast, the possibility of varying extents of heterotrimer dissociation in mammalian cells (Klein et al, 2000; Frank et al, 2005; Digby et al, 2008), and a suggestion that Gα in Arabidopsis exists primarily in a GTP‐bound, dissociated state (Johnston et al, 2007; Temple and Jones, 2007). Observations such as these lead to a fundamental question, namely, which of all the theoretical regulatory modes in G‐protein signaling are biologically possible, exemplifying a more general question of how we can best model the effects of switch‐like signaling mechanisms that have multiple active states. It is these two questions that are addressed here.
To facilitate the discovery of non‐classical mechanisms, which arguably occur more rarely than well‐established classical mechanisms, here we generate microarray data from wild‐type, gpa1, agb1, and agb1 gpa1 mutant plants and use transcriptome analysis, in which thousands of outputs (i.e. levels of individual transcripts) can be monitored simultaneously. To assess cell/tissue specificity of G‐protein signaling mechanisms, we perform transcriptome analysis in two types of samples, stomatal guard cells and rosette leaves.
We also assay these transcriptomes in the presence or absence of the phytohormone abscisic acid (ABA), a major plant hormone that both inhibits growth and promotes tolerance of abiotic stresses such as drought, salinity, and cold (Leung and Giraudat, 1998; Finkelstein et al, 2002; Acharya and Assmann, 2009). Although a few dozen candidate plant GPCRs with predicted 7TM structure have been computationally identified (Moriyama et al, 2006; Gookin et al, 2008), and several of these have been shown experimentally to interact with GPA1 (Gookin et al, 2008), to date none of these proteins has an identified ligand. We chose ABA as a variable because ABA signaling is known to interact with heterotrimeric G‐protein signaling in both developmental and stress responses in a complex manner (Wang et al, 2001; 2006, 2009; Fan et al, 2008). For example, ABA inhibition of stomatal opening, which promotes water conservation under stress conditions by reducing water vapor efflux through microscopic stomatal pores at the leaf surface, is impaired in gpa1 and agb1 single mutants as well as agb1 gpa1 double mutants, exemplifying ABA hyposensitivity of guard cell processes (Wang et al, 2001; Coursol et al, 2003; Fan et al, 2008). By contrast, seed germination and post‐germination seedling development are hypersensitive to inhibition by ABA in G‐protein complex mutants (Pandey et al, 2006). These experimental observations suggest G‐proteins as one of the components of ABA signaling, but to date no systematic study has been conducted to define the regulatory modes of a G‐protein or the co‐regulatory modes of a G‐protein and a hormone. Further, regulation of gene expression in G‐protein complex mutants has rarely received genome‐wide analysis in either plant or mammalian systems (Ullah et al, 2003; Okamoto et al, 2009), and there is only one replicated study of the guard cell transcriptome (Leonhardt et al, 2004).
As there is no evidence for ABA or heterotrimeric G‐proteins directly interacting with chromatin or DNA to influence gene expression, in this study we posit one or more mediators acting downstream of these inputs to control transcript levels. For simplicity, we use the minimum number of mediators necessary to explain the regulatory modes. In reality, the number of components involved may be more. In the two simplest cases, transcript levels would be (1) controlled by ABA independently of the G‐protein, or, conversely (2) controlled by the G‐protein independently of ABA, that is on activation by a signal other than ABA (Figure 1A and B). The simplest case of a combined G‐protein–ABA effect on a gene's expression would be that in addition to the G‐protein regulation, ABA has an additive G‐protein‐independent regulatory effect on the gene through another signaling mediator, as shown in Figure 1C. It is also possible that after the G‐protein is activated by a signal other than ABA, the G‐protein subunit(s) and ABA combinatorially co‐regulate a mediator, which then directly or indirectly regulates the gene, as shown in Figure 1D. There is also the possibility that ABA serves as a signal that activates the G‐protein. Combined with the regulatory modes illustrated in Figure 1, a total of nine possible G‐protein and/or ABA signaling pathways can be distinguished (Supplementary Figure S1). On the basis of nine G‐protein/ABA signaling pathways and the combinatorial relationships between G‐protein subunits, we use a Boolean modeling approach to systematically enumerate all possible theoretical regulatory modes of the G‐protein and co‐regulatory modes of the G‐protein and ABA (a total of 142 regulatory modes). In this study, we determine which signaling pathways and regulatory modes are experimentally supported by transcriptome data.
The Boolean modeling approach is used to formulate the theoretical regulatory modes of the G‐protein and ABA. We then apply a pattern matching approach to associate gene expression profiles with Boolean models. Our method allows us to address the following questions: (1) Are all of the classical modes of G‐protein signaling observed and are any non‐classical modes of G‐protein signaling revealed? (2) Does the heterotrimeric G‐protein have ABA‐independent regulatory activity? (3) Is it possible that ABA and the G‐protein co‐regulate some processes, but ABA is not the signal activating the G‐protein? (4) Is there system specificity of G‐protein regulation of the transcriptome? (5) Does regulation at the level of the transcriptome recapitulate previously observed hypersensitivity or hyposensitivity of developmental and dynamic transient responses to ABA in G‐protein mutants?
According to the characteristics of the conditions under which we sample the transcriptome (G‐protein subunit mutants with or without ABA), we use Boolean variables to denote the states of GPA1, AGB1, and ABA. Specifically, GPA1=0 and AGB1=0 represents the agb1 gpa1 double mutant, GPA1=0 and AGB1=1 represents the gpa1 mutant, GPA1=1 and AGB1=0 represents the agb1 mutant, and GPA1=1 and AGB1=1 means wild type. ABA=1 indicates the presence of ABA (i.e. ABA treatment) and ABA=0 represents the absence of ABA (i.e. solvent control). The expression profile of a target gene in different genotypes or treatments is the result of the (indirect) regulatory activity of the G‐protein and/or ABA, which can be seen as a truth vector of a Boolean function determined by the states of GPA1, AGB1, and/or ABA. Statistical analysis of gene expression data verifies the validity of this Boolean framework (Supplementary information 3). The regulatory modes of the G‐protein are modeled by Boolean functions A(GPA1, AGB1), reflecting the effect of the G‐protein on the expression level of a target gene (see Materials and methods). For example, the classical I mechanism of the G‐protein can be represented by A2(GPA1, AGB1)=GPA1 and AGB1. A target gene regulated by this G‐protein mode is strongly expressed when both GPA1 and AGB1 are present (i.e. in the wild type) and weakly or not expressed in the mutant genotypes. The classical II mode of the G‐protein corresponds to A5(GPA1, AGB1)=not GPA1 and AGB1. A target gene regulated by this mode is strongly expressed in the gpa1 mutant and shows weak or absent expression in other genotypes.
The co‐regulation of a gene by the G‐protein and ABA can be described as F(ABA, GPA1, AGB1)=CABA+B(ABA, A(GPA1, AGB1)). CABA is a constant that is non‐zero if there is a G‐protein‐independent regulatory effect of ABA on a gene. B(ABA, A(GPA1, AGB1)), shown in Table I, is a nested Boolean function representing the possible combinatorial regulation of ABA and the G‐protein. The complement of B(ABA, A(GPA1, AGB1)) functions determines a total of 72 non‐trivial regulatory modes of the G‐protein and/or ABA (see Materials and methods). Together with CABA, F(ABA, GPA1, AGB1) is able to identify 142 regulatory modes of the G‐protein and ABA, which can be classified into five main regulatory categories (see Materials and methods): ABA‐only regulation of transcript levels, G‐protein‐only regulation, G‐protein–ABA additive regulation, ABA–G‐protein combinatorial regulation, and ABA–G‐protein mixed regulation. For example, B2(ABA, A2(GPA1, AGB1))=ABA and (GPA1 and AGB1) is one example of combinatorial regulation by the G‐protein and ABA in which both G‐protein subunits and ABA are required for high gene expression. B4(ABA, A(GPA1, AGB1))=ABA exemplifies regulation by ABA independent of the status of the G‐protein. CABA+B6(ABA, A6(GPA1, AGB1))=CABA+AGB1 is one example of G‐protein–ABA additive regulation, namely regulation by the G‐protein β subunit with a separate additive input from ABA.
Each regulatory mode represented by A(GPA1, AGB1) or B(ABA, A(GPA1, AGB1)) corresponds to an idealized Boolean expression pattern. To avoid potential biases in binarizing expression levels, instead of correlating absolute expression profiles, our method is based on the differential expression patterns of genes. Figure 2 illustrates four examples of idealized expression patterns and idealized differential expression patterns of genes governed by Boolean functions. To associate genes with regulatory modes of the G‐protein and/or ABA, we construct the real differential expression pattern of each gene and adopt a correlation‐based pattern matching approach to examine whether the real differential expression pattern of the gene is consistent with one of the idealized differential expression patterns given by the Boolean regulatory modes of the G‐protein and/or ABA (see Materials and methods). By looking at the number of genes belonging to each regulatory mode, we are able to gauge the plausibility of each mode.
Biologically plausible modes of G‐protein regulation
Our first aim is to determine which of the theoretically possible G‐protein regulatory modes, as reflected by the 14 Boolean functions A2(GPA1, AGB1) to A15(GPA1, AGB1), are experimentally supported by the transcriptome data. We separately consider four combinations of tissues and treatment, that is guard cells or leaves with or without ABA treatment. Using stringent thresholds for differential expression and correlation with idealized patterns (see Materials and methods), a total of 107 G‐protein‐regulated genes are identified in guard cells without ABA treatment, and 51 G‐protein‐regulated genes are identified in guard cells with ABA treatment. There are 103 G‐protein‐regulated genes in leaves without ABA treatment, and 143 G‐protein‐regulated genes in leaves with ABA treatment. These genes are listed in the Supplementary file Sup1.xls.
We observe a number of G‐protein‐regulated genes in the absence of ABA in both guard cells and leaves, indicating that the G‐protein has regulatory activity independent of ABA treatment. Only the two knocked‐out genes, GPA1 and AGB1, are common to the G‐protein‐regulated gene sets in guard cells and leaves in either of the two ABA conditions, suggesting system specificity of G‐protein action. The distribution of genes in each G‐protein regulatory mode in guard cells and leaves, with and without ABA treatment, is given in Figure 3. Note that it is possible that the genes in Ai and A17‐i are regulated by the G‐protein through the same signaling pathway, for example for A2 and A15, GPA1 and AGB1 could synergistically regulate a mediator (e.g. a transcription factor or other regulator of gene expression such as a microRNA), which activates the genes in A2 and inhibits the genes in A15. To gain insight into the overall abundance of possible G‐protein regulatory mechanisms, we plot the cumulative number of genes in each G‐protein regulatory group by merging the genes present in A17‐i with those in Ai, as shown in the inset of Figure 3 (see also Supplementary information 9 for the interpretation of the relationship between an A(GPA1, AGB1) function and its negation).
From Figure 3, we can see that A2=GPA1 and AGB1 (and/or A15=not A2=not GPA1 or not AGB1) is a popular G‐protein regulatory mode in both guard cells and leaves. This regulatory mode requires the presence of both Gα and Gβγ subunits, and thus represents the classical Ia and Ib mechanisms, in which both Gα and Gβ(γ) are required. We note that, formally, this mode is also consistent with regulation by a non‐dissociated heterotrimer (Klein et al, 2000; Frank et al, 2005; Digby et al, 2008). The regulatory mechanism A5=not GPA1 and AGB1 (and/or A12=not A5=GPA1 or not AGB1) is also strongly supported in all conditions except in leaves without ABA treatment. This mode is consistent with the classical II G‐protein regulatory mechanism, based on signaling by the freed Gβγ subunit, as it incorporates the fact that in the gpa1 mutant, Gβγ will be available to regulate downstream genes. The G‐protein regulatory mode A6=AGB1 (and/or A11=not AGB1) is prevalent in guard cells but less so in leaves. This Boolean function represents a non‐classical regulatory mechanism in which signaling by Gβ(γ) occurs, and this signaling is not regulated in any way by Gα; in other words, at no point during this signaling mechanism does Gβ(γ) assemble with Gα into a heterotrimer. The functions A3=GPA1 and not AGB1, A14=not A3=not GPA1 or AGB1 are not supported, suggesting that Gα does not regulate downstream effectors in the absence of Gβγ. The functions A4=GPA1, A13=not GPA1 are also very weakly supported, suggesting that Gα signaling unaffected by the presence or absence of Gβγ is rare. The function A8=GPA1 or AGB1 (and/or A9=not A8=not GPA1 and not AGB1) is weakly supported in the absence of ABA, but its low abundance suggests that regulation of a target equally well by either free Gα or free Gβγ is relatively rare. Interestingly, the function A7=(GPA1 or AGB1) and not (GPA1 and AGB1) (and/or A10=not A7=(GPA1 and AGB1) or not (GPA1 or AGB1)) is well represented in leaves with ABA treatment. The complex combinatorial relationship between Gα and Gβγ underlying this function suggests the involvement of more than one G‐protein effector.
Co‐regulatory modes of the G‐protein and ABA
We next combine the data with and without ABA treatment and determine which of the theoretically possible (co)regulatory modes of G‐protein and ABA are experimentally supported by the transcriptome data. By using stringent thresholds for association of genes with regulatory modes (see Materials and methods), 71 ABA–G‐protein combinatorially regulated genes, 29 ABA–G‐protein mixed regulated genes, 28 G‐protein‐only regulated genes, and 3 G‐protein–ABA additively regulated genes are identified in guard cells. In contrast, in leaves 471 ABA–G‐protein combinatorially regulated genes, 78 ABA–G‐protein mixed regulated genes, and only 3 G‐protein‐only regulated genes are identified. We also obtain 467 ABA‐only induced genes and 246 ABA‐only repressed genes in guard cells, and 309 ABA‐only induced genes and 94 ABA‐only repressed genes in leaves. All of these genes are listed in the Supplementary files Sup2.xls for guard cells and Sup3.xls for leaves.
To focus on the combinatorial G‐protein–ABA regulatory modes encapsulated in the B(ABA, A(GPA1, AGB1)) function, we cumulate ABA–G‐protein combinatorially regulated genes and ABA–G‐protein mixed regulated genes that are governed by the same Boolean function. Similarly, the genes governed by the same Boolean function in G‐protein‐only regulation and G‐protein–ABA additive regulation are also merged. We find that 33 out of the 72 theoretically possible regulatory modes are unsupported (have 0 genes) in guard cells and 34 regulatory modes are unsupported in leaves. The distribution of genes in each supported regulatory mode of the G‐protein and ABA in guard cells and leaves is illustrated in Supplementary Figures S10, S11 and S12. In Table II, we summarize the regulatory modes that govern >10% of the genes in the relevant category (ABA–G‐protein combinatorial/mixed regulation, or G‐protein‐only/G‐protein–ABA additive regulation).
We observe that in guard cells, A5 (and/or A12)‐related regulatory modes are the most representative: B12(ABA, A5)=ABA or not (not GPA1 and AGB1)=ABA or not A5, B15(ABA, A5)=not ABA or not A5, and B5(ABA, A5)=not ABA and A5. These co‐regulatory modes of the G‐protein and ABA are consistent with Figure 1D where the mediator M1 is regulated by the classical II G‐protein regulatory mechanism dependent on the Gβγ subunit on its release from Gα, and the mediator M2 is combinatorially regulated by ABA and M1. Figure 4A uses the commonly used ‘heat map’ portrayal to illustrate one of these well‐supported classical II modes, B12(ABA, A5), in which the absence of Gα frees Gβγ to downregulate downstream genes, in this instance with the effect of Gα knockout on gene expression dependent on the absence or presence of ABA. In guard cells, G‐protein‐only signaling and G‐protein–ABA additive regulation, where the relative effect of subunit knockout is the same regardless of the presence or absence of ABA, are also possible. The classical II mode B11(ABA, A5)=not (not GPA1 and AGB1) (Figure 4B) and the non‐classical mode B6(ABA, A6)=AGB1 (Figure 5) are the two most supported G‐protein‐only/G‐protein–ABA additive regulatory modes.
In leaves, A2 (and/or A15)‐related regulatory modes are the most prevalent: B2(ABA, A2)=ABA and (GPA1 and AGB1)=ABA and A2, B15(ABA, A2)=not ABA or not A2, B5(ABA, A2)=not ABA and A2, and B12(ABA, A2)=ABA or not A2. These regulatory modes are consistent with Figure 1D where the mediator M1 is synergistically regulated by GPA1 and AGB1 (classical I G‐protein regulatory mechanisms) and the second mediator M2 is combinatorially regulated by M1 and ABA. Figure 6 illustrates one of these well‐supported modes, B2(ABA, A2)=ABA and (GPA1 and AGB1), in which gene expression is highest in wild type plus ABA, due to a requirement for both a wild‐type G‐protein complement and ABA treatment to elevate transcript levels. There are only three genes in the G‐protein‐only category in leaves, and two of them are the G‐protein subunits GPA1 and AGB1 included by default; the third is At3G12730, a MYB family transcription factor. This suggests that in leaves there are very few genes that are regulated by the G‐protein independently of ABA signaling and instead the co‐regulatory mode of the G‐protein and ABA is prevalent.
To test the robustness of the above conclusions, we repeated the analysis with a looser genotype threshold and found no new representative gene groups, only more genes for the originally representative regulatory modes (see Supplementary Figures S10, 11 and S12). These results indicate that these regulatory modes are biologically plausible and there indeed exists a coupling between transcription regulation by the G‐protein and ABA in both guard cells and leaves.
Relating co‐regulated genes to the G‐protein and ABA signaling pathways
The regulatory modes and ABA additive effects identified in our analysis also enable determination of the most plausible G‐protein and ABA signaling pathways (see Figure 1; Supplementary Figure S1). ABA‐only regulated genes (i.e. B4(ABA, Ai)=ABA, B13(ABA, Ai)=not ABA, i=2,…,15), correspond to the pathway in Figure 1A, indicating that there is a G‐protein‐independent effect of ABA, and the status of the G‐protein is irrelevant. This regulatory mode is strongly supported in both guard cells and leaves: it describes 713 genes in guard cells (467 upregulated and 246 downregulated, see Sheet 5 of the Supplementary file Sup2.xls) and 403 genes in leaves (309 upregulated and 94 downregulated, see Sheet 5 of the Supplementary file Sup3.xls). This result indicates that there are G‐protein‐independent ABA signaling pathways in both guard cells and leaves. The signaling pathway of Figure 1B means that the G‐protein is turned on by an ABA‐independent signal, and there is no G‐protein‐independent effect of ABA on gene expression. This regulatory mode corresponds to G‐protein‐only regulated genes (i.e. B6(ABA, Ai)=Ai, B11(ABA, Ai)=not Ai, i=2,…,15 with non‐significant CABA), and is exhibited by 28 genes in guard cells (see Sheet 3 of Sup2.xls) but only three genes (of which two are GPA1 and AGB1) in leaves (see Sheet 3 of Sup3.xls). These results corroborate the existence of an ABA‐independent regulatory activity of the G‐protein in guard cells, through both classical (B11(ABA, A5)) and non‐classical (B6(ABA, A6)) regulatory modes. Interestingly, these ABA‐independent G‐protein regulatory modes are not supported in leaves.
The signaling pathway of Figure 1C means that in addition to ABA‐independent G‐protein regulation, there is also a G‐protein‐independent effect of ABA. This regulatory mode corresponds to G‐protein–ABA additive regulated genes (i.e. B6(ABA, Ai)+CABA=Ai+CABA, B11(ABA, Ai)+CABA=not Ai+CABA) and represents an additive cross‐talk, converging on the target gene, between ABA and the signal that activates the G‐protein. This regulatory mode has three genes in guard cells (see Sheet 4 of Sup2.xls) and no genes in leaves, indicating that additive regulation by G‐protein signaling plus G‐protein‐independent ABA signaling is rare.
The signaling pathway shown in Figure 1D represents the case where the G‐protein is turned on by an ABA‐independent signal, and there is a combinatorial G‐protein–ABA effect. This pathway corresponds to ABA–G‐protein combinatorially regulated genes with non‐significant CABA (i.e. B2(ABA, Ai)=Ai and ABA, B8(ABA, Ai)=Ai or ABA, B5(ABA, Ai)=Ai and not ABA, and B14(ABA, Ai)=Ai or not ABA, i=2,…,15) and represents a combinatorial cross‐talk that converges on the mediator M2 (see Figure 1D). This regulatory mode is observed for 71 genes in guard cells (see Sheet 1 of Sup2.xls) and 471 genes in leaves (see Sheet 1 of Sup3.xls). A similar pathway (Supplementary Figure S1E) is that in addition to the combinatorial G‐protein–ABA regulatory effect, there is also a G‐protein‐independent effect of ABA (i.e. (Ai and ABA)+CABA, (Ai or ABA)+CABA, (Ai and not ABA)+CABA, and (Ai or not ABA)+CABA, i=2,…,15). This regulatory mode corresponds to ABA–G‐protein mixed regulated genes and is supported by 29 genes in guard cells (see Sheet 2 of Sup2.xls) and 78 genes in leaves (see Sheet 1 of Sup3.xls), indicating that cross‐talk between ABA–G‐protein signaling and G‐protein‐independent ABA signaling occurs.
In mammalian systems, it is well established that a given Gα can interact with more than one specific GPCR and thus can be regulated by more than one ligand (Leaney and Tinker, 2000). In Arabidopsis, if the signal that catalyzes the dissociation of Gα and Gβγ and activates the G‐protein were ABA (see Supplementary Figure S1F–I), then without ABA treatment, the G‐protein would be off and would have no regulatory activity (described by B2(ABA, Ai)=Ai and ABA, B14(ABA, Ai)=Ai or not ABA, i=2,…,15). In addition to such regulation of the G‐protein by ABA, there could be a G‐protein‐independent effect of ABA, or/and a G‐protein–ABA combinatorial regulatory effect. The genes regulated by these modes would form an unknown subset of ABA–G‐protein combinatorially or mixed regulated genes, that is all or a fraction of the genes that have no differential expression with respect to genotypes in the absence of ABA (see Supplementary information 11 for a detailed description). The maximum possible number of genes regulated by a mode in which ABA is the signal activating the G‐protein is 38 in guard cells and 287 in leaves. In other words, these genes are consistent with the possibility of ABA serving as a ligand for a GPCR. On the other hand, for the specific case where the G‐protein functions according to the classical mechanism II, we can conclude that our analysis does not validate the scenario in which ABA is the signal that activates the G‐protein (see Supplementary information 11 for the reasoning behind this statement). Finally, the 62 genes in guard cells and 262 genes in leaves consistent with ABA–G‐protein combinatorial or mixed regulation but not consistent with ABA activation of the G‐protein are a clear indication that the G‐protein can be activated by a signal other than ABA.
ABA hypo‐/hypersensitivity of gene regulation in G‐protein mutants
For guard cell responses published in the literature to date, namely for ABA inhibition of stomatal opening and inward K+ channel regulation, gpa1 and agb1 mutants are equally hyposensitive to ABA (i.e. less sensitive than wild type) (Wang et al, 2001; Coursol et al, 2003; Fan et al, 2008). By contrast, in seed germination and post‐germination seedling growth, G‐protein subunit mutants are un‐equally hypersensitive to ABA (i.e. more sensitive than wild type): the agb1 single and the agb1 gpa1 double mutants are more hypersensitive than the gpa1 single mutant (Pandey et al, 2006). There are no extant studies regarding mature leaf sensitivity to ABA in the G‐protein mutants versus wild type, nor any transcriptome‐wide studies on the sensitivities of gene regulation to ABA in the G‐protein mutants in any cell or tissue type. Here, we use our transcriptome data to assess and compare gene regulation by ABA in guard cells and leaves of the G‐protein mutants.
ABA‐only regulated genes by definition have equal responses to ABA in all mutant genotypes. G‐protein‐only regulated genes and G‐protein–ABA additively regulated genes by definition have no or equal response to ABA in all mutant genotypes; thus, these genes are not useful for revealing ABA hyposensitivity or hypersensitivity. Therefore, we investigate those genes that exhibit ABA–G‐protein co‐regulation, comparing each ABA–G‐protein‐regulated gene's expression change in response to ABA in the wild‐type genotype versus each mutant genotype. We define hyposensitivity to ABA as the case when the gene has lesser up/downregulation in response to ABA in a mutant as compared with wild type. Similarly, we define hypersensitivity to ABA as the case when the gene has greater up/downregulation in response to ABA in the mutant as compared with wild type. We also determine the numbers of genes that have similar responses in combinations of two mutant genotypes, or all three of them.
The fraction of regulated genes exhibiting hyposensitive or hypersensitive response to ABA in guard cells and leaves of the mutants are listed in Table III. To determine the significance of these results, we compare these fractions to those generated using a random background (given in Supplementary Table S3), in which the genes are randomly scattered over all Boolean rules that reflect ABA–G‐protein co‐regulation. Significance values are computed by a binomial cumulative distribution. We find that in guard cells, equal hyposensitivity in all mutants (the last row) and hypersensitivity in the gpa1 mutant are significant (P=5 × 10−5 and P=3 × 10−8, respectively). In leaves, ABA hyposensitivity of gene expression in the three individual mutants (the first three rows) and equal hyposensitivity in all mutants (the last row) are highly significant (P<10−50 for all these cases). Using a less stringent differential expression threshold (Pgeno<0.05) leads to very similar results (data not shown).
Functional preferences in G‐protein regulatory modes
To examine whether specific G‐protein regulatory modes tend to be involved in specific aspects of plant physiology, we perform functional analysis on the co‐regulated genes in each G‐protein regulatory mode using the FunCat functional catalogue of the MIPS database (Mewes et al, 2006). The P‐values are calculated by the hypergeometric cumulative distribution and Bonferroni correction was applied for multiple testing. We applied a highly stringent cutoff of P<10−4 and identified 12 functional categories significantly enriched in G‐protein regulatory modes (acknowledging the caveat that absolute gene numbers are small in some categories), shown in Table IV. The genes in these significant functional categories are listed in Supplementary Tables S4 and S5 and in Supplementary information 13.
As we have described, there are multiple theoretical signaling pathways and regulatory modes involving the G‐protein and/or ABA. Through our identification of the Boolean rules governing genes (co)‐regulated by the G‐protein and ABA, we can relate these genes to 72 theoretical regulatory modes (Table I) and the proposed nine signaling pathways (Supplementary Figure S1), providing insights into the biological plausibility of G‐protein and ABA regulation and allowing us to address the questions posed in the Introduction.
Are all of the classical modes of G‐protein signaling supported and are any non‐classical modes of G‐protein signaling revealed?
Phenotypic analysis of Arabidopsis G‐protein mutants has indicated that classic mammalian G‐protein regulatory paradigms are also used during plant G‐protein signaling. Our analysis strongly supports the conclusion that such classical paradigms also function in the regulation of the plant transcriptome. Our analysis indicates that most of the 70 theoretical regulatory mechanisms involving the G‐protein (i.e. excluding two ABA‐only regulatory modes) are unsupported and thus may be biologically untenable, at least in the biological systems investigated here. For example, 33 out of these 70 theoretical regulatory modes are unsupported (have 0 genes) and only five are strongly supported in guard cells; the corresponding numbers for leaves are 34 and 4 (Supplementary Figure S10, S11 and S12; Table II). These results lend credence to our analysis method, as it would be unexpected if G‐protein regulation of the transcriptome were to proceed through entirely different mechanisms than the well‐established mechanisms revealed by physiological analyses.
Notably, our analysis also reveals a novel G‐protein regulatory mechanism: positive regulation by the Gβ(γ) subunit independent of the presence or absence of the Gα subunit is well supported in the guard cell transcriptome. According to this rule, AGB1, agb1 plants exhibit a mutant phenotype, whereas gpa1 plants are wild type for the same phenotype. Interestingly, among the numerous physiological and developmental responses that are G‐protein regulated in Arabidopsis as judged by mutant analysis (Perfus‐Barbeoch et al, 2004), there are a few phenotypes consistent with this pattern: agb1 mutants are resistant to the protein glycosylation inhibitor, tunicamycin, whereas gpa1 plants exhibit wild‐type sensitivities (Wang et al, 2007). This response has been associated with AGB1 resident in the endoplasmic reticulum (Wang et al, 2007), suggesting that subcellular localization may be one mechanism through which Gβ(γ) subunits are partitioned into heterotrimer‐dependent versus heterotrimer‐independent functions. Another plausible AGB1 phenotype is altered root waving and skewing, which is observed in agb1 mutants and not in gpa1 mutants (Pandey et al, 2008).
To our knowledge, Gβ(γ) function completely independent of any input from Gα (i.e. Gβ(γ) never assembles into a heterotrimer, even in the absence of agonist) has not been previously implicated in mammalian systems, except possibly for non‐Gγ‐related signaling modes of the unusual Gβ5 subunit (Dupre et al, 2009). This may be because the mechanism is unique to plants or it may be because the mechanism only occurs in specialized mammalian cell types or systems and thus has escaped detection. Alternatively, it may be that the complexity and redundancy of the G‐protein complement in mammals has simply precluded discovery of this mechanism to date, whereas our use of a simpler model system coupled with the ability to monitor thousands of phenotypes simultaneously enabled detection. Thus, future investigation of this mechanism in both plants and metazoans is warranted.
Conversely, our analysis indicates that the regulation by Gα independent of the presence or absence of the Gβγ subunit is not supported by the microarray data. It has been proposed that GPA1, owing to fast kinetics of GDP/GTP exchange and a slow rate of GTP hydrolysis, would persist constitutively in the GTP‐bound form, and thus could function independently from Gβ(γ) (Johnston et al, 2007). In our own biochemical analyses (Pandey et al, 2009), we also observed slow rates of GTP hydrolysis by GPA1 but rates of GTPγS binding were comparable to those for bovine Gα assayed by the same method. Our analysis here does not indicate independent Gα function, at least with regard to transcriptome regulation in guard cells and leaves.
Does the heterotrimeric G‐protein have ABA‐independent regulatory activity? Is it possible that ABA and the G‐protein co‐regulate some processes, but ABA is not the signal activating the G‐protein?
Previous physiological analyses of guard cell function have focused almost exclusively on the roles of heterotrimeric G‐proteins in ABA signaling (Wang et al, 2001; Fan et al, 2008). Our transcriptome analysis now reveals that in this cell type, the heterotrimeric G‐protein also has clear ABA‐independent activity. Thus, when ABA–G‐protein co‐regulatory modes are evaluated, two G‐protein‐only functions are seen in guard cells: GPA1 or not AGB1 (classical II) and AGB1. Guard cells actually sense a wide variety of signals, including light, CO2, humidity, several hormones, and pathogens (Roelfsema and Hedrich, 2005; Israelsson et al, 2006; Shimazaki et al, 2007; Melotto et al, 2008; Acharya and Assmann, 2009). Our results now suggest that potential roles of G‐proteins in these other signaling processes should also be evaluated. Indeed, there is already some evidence for heterotrimeric G‐protein involvement in guard cell response to pathogens (Zhang et al, 2008b), as well as in guard cell development (Zhang et al, 2008a).
In guard cells, there is also evidence that ABA and the G‐protein co‐regulate sets of genes, consistent with known co‐regulation by ABA and the G‐protein at the physiological level. For the 62 genes corresponding to the functions B8(ABA, Ai)=Ai or ABA, B5(ABA, Ai)=Ai and not ABA with or without a significant CABA, ABA can be excluded as the activating signal for the G‐protein. An additional 38 guard cell genes exhibit regulation consistent with ABA being the signal that activates the G‐protein; however, in the absence of additional biological information, these patterns are also consistent with ABA regulating signal transduction downstream of the G‐protein: our current data do not allow us to distinguish between these two possibilities.
Our analysis does not support G‐protein functions independent of ABA (G‐protein‐only regulatory modes) in mature leaves. There are no genes in the leaf transcriptome that exhibit an identical pattern of expression relative to genotype in both the presence and the absence of ABA. There is strong support, however, for both mixed and combinatorial regulation by ABA and the G‐protein. To date, ABA responses of rosette leaves have not been studied in the G‐protein mutants, and our results strongly suggest that such studies would be useful for deciphering novel G‐protein/ABA signaling pathways.
Is there system specificity of G‐protein regulatory mechanisms?
To account for the diversity of physiological and developmental phenotypes observed in knockout mutants of plant G‐proteins, despite the small numbers of G‐protein subunits encoded in plant genomes, it has been proposed that there is specificity in how G‐protein signaling components are used in different cells and tissues. Our analyses corroborate this hypothesis, through several lines of reasoning. First, the regulatory modes that are strongly supported by our analysis differ in guard cells versus leaves. This is true both when we evaluate G‐protein regulatory modes without considering ABA effects (Figure 3) and when we combine all the data together to analyze G‐protein/ABA co‐regulation (Table II). For example, in both the former and the latter analyses, AGB1 is an observed regulatory mode in guard cells but not in leaves and, in fact, for leaves there is no G‐protein‐only regulatory mode that is supported when G‐protein/ABA co‐regulation is considered. These results suggest that heterotrimeric G‐proteins are actually more involved in ABA signaling in leaves than in guard cells, an interpretation that is consistent with the stronger ABA hyposensitivity exhibited by the leaf transcriptome than the guard cell transcriptome of the G‐protein mutants (see next section). This new finding was unexpected, given that ABA has been integrally linked to guard cell physiology for decades, but has received little attention as a regulator of mesophyll cell biology. It is also important to note that this conclusion could not have arisen from extant physiological data, since to date no leaf responses to ABA have been evaluated at the physiological level in the G‐protein mutant lines.
Although both guard cells and leaves do exhibit ABA–G‐protein co‐regulation of gene expression, the two systems favor different co‐regulatory modes. For example, it is striking that at the transcriptome level, the classical II mechanism (not GPA1 and AGB1) is significantly more common in guard cells than in leaves, whereas the classical Ia and Ib mechanisms (GPA1 and AGB1) are possible in both tissue types, but more strongly so in leaves (Table V; see Supplementary information 14 for statistical analysis). In addition, the specific genes that are co‐regulated by ABA and the G‐protein differ between guard cells and leaves. Thus, among the 100 genes in guard cells and 549 genes in leaves that are ABA–G‐protein co‐regulated in some manner, only five genes (At1g30290, At4g29750, At1g31150, At1g74720, and At1g11210) appear in both guard cell and leaf data sets, in contrast to 106 ABA‐only regulated genes common to both systems.
Does regulation at the level of the transcriptome recapitulate previously observed hypersensitivity or hyposensitivity of developmental and dynamic transient responses to ABA in G‐protein mutants?
Transcriptome level responses to ABA show some interesting differences compared with the existing physiological literature on ABA sensitivity of G‐protein mutants (Wang et al, 2001; Pandey et al, 2006). For example, whereas ABA inhibition of stomatal opening and ion channel regulation in guard cells are ABA hyposensitive in gpa1, agb1, and agb1 gpa1 double mutants (Wang et al, 2001; Fan et al, 2008), our analysis shows that ABA hyposensitivity of the mutant guard cell transcriptome, at least at the time point sampled here, is not as strong as the hyposensitivity evidenced by these rapid responses to ABA, possibly indicating differences between protein level and transcriptionally mediated responses. Further, a set of guard cells transcripts that shows ABA hypersensitivity in the G‐protein mutants is observed. There is precedence for the same guard cell signaling element functioning in both ABA hyposensitivity and ABA hypersensitivity; for example, overexpression in guard cells of an inositol polyphosphate 5‐phosphatase, which terminates inositol phosphate signaling, results in hyposensitivity to ABA inhibition of stomatal opening and hypersensitivity to ABA promotion of stomatal closure (Perera et al, 2008). Nevertheless, transcriptome results do not necessarily mean that ABA hypersensitivity of a portion of the transcriptome results in physiological outcomes that are ABA hypersensitive as, for example, the ABA‐hypersensitive genes could encode elements that function in repression of ABA response. In either case, these results do indicate that a diversity of signaling pathways exists between the G‐protein and gene regulation in guard cells. Interestingly, the majority of genes that show ABA hypersensitivity in the gpa1 mutant (50/56) are governed by the classical II mechanism and the majority of genes that show equal hyposensitivity in all mutants (19/22) are governed by the classical I mechanisms.
Another novel observation is the strong ABA hyposensitivity of the G‐protein‐regulated mature leaf transcriptome, with all mutants showing equal ABA hyposensitivity. Previously, seed germination, seedling growth and gene expression were shown to be hypersensitive to ABA in these G‐protein mutants (Pandey et al, 2006), leading to the hypothesis that ABA hyposensitivity might be a distinctive characteristic found exclusively in the highly specialized guard cells, where ABA signaling has an integral function in cellular function. This hypothesis is not supported by our microarray data, which instead suggest the alternative hypothesis that G‐protein involvement in ABA signaling is defined at the cell or organ level.
Functional categories in G‐protein regulatory modes
If the biological interpretations indicated by our model are valid, then we would expect to find at least some level of correlation between these interpretations of our transcriptome data and results from previous physiological analyses of G‐protein mutants. Indeed, several of the functional categories identified by our analysis have been implicated in previous physiological analyses of G‐protein mutants, providing validation to the biological interpretation of our results. Three functional categories related to plant pathogen response appear in guard cell regulatory modes, consistent with known roles of the G‐protein in plant pathogen response (Trusov et al, 2006). For example, the functional category ‘disease, virulence, and defense’ is enriched in the regulatory mode B11(ABA, A5)=not (not GPA1 and AGB1)=GPA1 or not AGB1in guard cells. Specifically, At1G63750, At1G72520, LCR69, NHL10, NHL1, and AtPCB are all genes implicated in defense against pathogens. This regulatory mode is consistent with a regulatory mechanism in which the Gβγ subunit, freed from its association with Gα, downregulates transcript levels of these defense genes, and is also consistent with our earlier physiological observations that gpa1 guard cells are hyposensitive to the bacterial elicitor molecule, flg22 (Zhang et al, 2008b). Given this information, we might speculate that, in guard cells, ligands coupling to GPCRs upstream of the heterotrimer may be, or include, pathogen‐derived or ‐dependent molecules, that is elicitors. Associated with the rule A6=AGB1 in guard cells are the significant functional categories ‘de/phosphorylation’ and ‘phosphate metabolism’. The latter category contains nine kinases and two genes encoding PP2C‐type protein phosphatases. A number of kinases and PP2C phosphatases already have been shown to be integral to guard cell signaling pathways (Gosti et al, 1999; Li et al, 2000; Merlot et al, 2001; Mustilli et al, 2002; Yoshida et al, 2002; Ma et al, 2009; Park et al, 2009; Rubio et al, 2009); our results implicate new candidate proteins of this type. In guard cells, we also find that genes associated with metabolism of aromatic compounds are components of the ABA–G‐protein co‐regulated transcriptome. Although not yet studied in guard cells, there are several reports at the whole seedling level that implicate G‐proteins in the regulation of aromatic amino‐acid synthesis (2008, 2006).
In leaves, we find that the functional category ‘chloroplast’ is highly significant, consistent with earlier observations linking GPA1 to chloroplast development (Zhang et al, 2009). Interestingly, the chloroplast protein THF1 (Zhang et al, 2009), which has been shown to physically interact with GPA1 (Huang et al, 2006), appears as one of the regulated transcripts in this category. Representation of functional categories related to signaling (‘calcium binding,’ ‘cellular signaling’) is also consistent with the known signaling roles of G‐proteins. The category ‘nucleotide/nucleoside binding’ is typified by kinases, ATPases, and ATP‐binding G‐proteins, linking G‐proteins to well‐established secondary messengers and transport functions. Thus, our transcriptome analyses are consistent with previously identified biological roles of G‐proteins in Arabidopsis while also providing specific new gene targets for investigation.
One point to note is that, considering the data as a whole, the numbers of genes regulated by the ABA‐only mode (713 in guard cells and 403 in leaves) is much greater than the numbers of genes regulated by G‐protein‐only modes (28 genes in guard cells and only At3G12730 in leaves, excluding GPA1 and AGB1 themselves). These results, along with the observation that transcription factors are not found to be overrepresented among G‐protein‐only or G‐protein–ABA co‐regulated genes (data not shown), suggests that the phenotypes of G‐protein mutants do not arise due to a ‘resetting’ of the entire cellular profile but rather are likely to stem from disruption of specific G‐protein‐coupled signaling events.
In this study, we generate transcriptome data from stomatal guard cells and rosette leaves and investigate G‐protein signaling in these two systems. There may be a cellular complexity difference between guard cell and leaf samples, as guard cells are a single cell type, whereas leaves are composed of several types of cells. However, the majority of the cell types contributing to the leaf microarray are mesophyll cells; guard cells and other non‐mesophyll cell types comprise a very small percentage of the leaf. Our transcriptome measurements confirm that minor cell types such as guard cells do not significantly influence the leaf microarray. For example, in our wild‐type transcriptome data sets, 361 transcripts were detectable in all three guard cell chip hybridizations but were not present in any of the three leaf chip hybridizations. In particular, some genes that are very highly expressed in Col guard cells, for example ATHSP23.6‐MITO (AT4G25200), ATATH4 (AT3G47760), DDF1 (AT1G12610), EPF1 (AT2G20875), SP1L2 (AT1G69230), and AtSIP1 (AT1G55740) are not present in any of the Col leaf microarrays.
In addition, in this study, we obtain leaves and guard cells at one developmental time point and measure gene expression at one time point after ABA treatment. Although we acknowledge these caveats, the former is not an issue as germination and leaf development throughout the vegetative stage is quite synchronous among all four genotypes under our conditions (Nilson and Assmann, 2010). For the latter, by our 3‐h time point, the first (rapid) phase of guard cell response will have reached steady state, based on earlier studies of stomatal aperture response to ABA (Israelsson et al, 2006). Our study identifies key genes for which kinetics of gene expression can be assessed in more detail in subsequent studies.
Relevance of our Boolean framework
Although Boolean networks have been used previously to extract gene relationships from time course expression profiles (Akutsu et al, 1999; Martin et al, 2007), this is, to our knowledge, the first study using systematic Boolean rules to deduce regulatory mechanisms based on analysis of mutant transcriptomes. The microarray studies in this work offer an unprecedented complexity by studying four genotypes and two signaling conditions, in two different tissue types. With such complex data sets, pairwise differential expression analysis usually leads to results that are difficult to interpret because of the combinatorial complexity of regulation. The novelty of our method is that we provide an intuitive and reasonable way to analyze mutant expression data, which not only determines gene groups/modules, but also indicates the regulatory mechanisms connecting the genes in each group to the G‐protein and/or to ABA. This makes our method very helpful for analyzing signal transduction pathways given individual signaling components with known properties, for example a transcription factor known to vary in phosphorylation status.
When associating genes with regulatory modes, we avoid biases and binarization of the expression values by correlating differential expression patterns instead of correlating absolute expression profiles. We use a systematic pairwise comparison procedure (Figure 2) to construct the differential expression vectors of genes. Not all of these pairwise comparisons are independent or orthogonal, but these comparisons determine a vector uniquely corresponding to a Boolean function from a comprehensive perspective and give a reliable result. We also tested the use of orthogonal contrasts to construct the differential expression patterns of genes, and obtained similar results (Supplementary information 15).
Clustering methods comprise another approach commonly used to predict co‐regulated genes in microarray data sets (Thalamuthu et al, 2006; Liu et al, 2008). These methods are most useful when only rough knowledge such as upregulation or downregulation is expected as a result. The clusters do not offer information on the biological regulatory modes, whereas our method provides a putative mechanism for each group of co‐regulated genes. In addition, clustering methods cannot distinguish the extent of a regulatory effect on a gene. In our method, each gene has correlation scores indicating the strength of the gene's association with each regulatory mode. We have tested the popular k‐means clustering method to analyze our mutant expression data (Supplementary information 16) and find that most of the clusters do not have coherent expression patterns like ours (Supplementary Table S8). Moreover, these clusters lack biological interpretations such as those provided by our models.
In conclusion, our application of a Boolean framework for the analysis of microarray data sets has discovered new mechanisms of G‐protein and hormonal control. We posit that application of this approach to transcriptomic data sets from other systems will similarly provide new perspectives regarding other switch‐like signal transduction mechanisms, such as regulation by reversible post‐translational events or by combinatorial association of transcription factors.
Materials and methods
Plant materials, ABA treatment, and microarray experiments
All mutants are in the Arabidopsis Col background and have been described earlier (Jones et al, 2003; Ullah et al, 2003; Chen et al, 2004). These are T‐DNA insertional mutants that fail to produce full‐length GPA1 or AGB1 transcripts; gpa1‐4 is also confirmed as a null mutation at the protein level (Fan et al, 2008). Microarray data were generated from four genotypes (wild type, gpa1‐4 mutant, agb1‐2 mutant, and agb1‐2 gpa1‐4 double mutant) with or without ABA treatment. Arabidopsis plants were grown in growth chambers with an 8 h light/16 h dark, 20°C/18°C cycle, with light intensity of 120 μmol m−2 s−1. Three hundred Arabidopsis leaves excised from 60 to 70 five‐week‐old plants were used as the starting material for each guard cell microarray (see Supplementary information 2 for plant growth details and guard cell isolation protocol). Ten mature leaves taken from three to four plants grown side‐by‐side with the plants for guard cell isolation were used for each leaf sample. Excised leaf and isolated guard cell samples were treated with ABA (50 μM) or EtOH (solvent control) for 3 h. RNA was isolated from each sample using Trizol reagent (Invitrogen). Trizol‐isolated RNA was further purified using the Plant RNeasy kit (Qiagen, Valencia, CA). RNA samples exhibiting a 25S/18S ratio >1.4 and no significant degradation in Bioanalyzer profiles (Agilent Technology, Palo Alto, CA) were used for cDNA and then cRNA synthesis, followed by hybridization to Affymetrix ATH1 ‘whole genome’ chips, all according to standard Affymetrix protocols. For each type of sample (guard cells or leaves), three independent biological replicates were performed, resulting in a total of 48 microarray hybridizations (2 sample types × 4 genotypes × two treatments × 3 replicates).
Microarray data processing and validation
The transcriptome data reported in this paper have been deposited in the Gene Expression Omnibus database (http://www.ncbi.nlm.nih.gov/geo) with accession no. GSE19520. We applied several different approaches to assess quality control and the overall quality of the hybridizations is very high (Supplementary information 17). The average pairwise correlation among all biological replicates is 0.982, indicating that the variation arising from biological replicates is very small. The widely used RMA method implemented in the affy package in the bioconductor project was used to normalize the data for all probe sets using default parameters (Gautier et al, 2004; Gentleman et al, 2004). The RMA does not provide a way to calculate present calls and absent calls, so GeneChip suite 5.0 algorithm implemented in the affy package was used to obtain presence/absence calls for each probe set in each independent sample with default parameters. In total, there are 22 810 probe sets on the Affymetrix ATH1 Arabidopsis chips. We eliminated only those genes that were absent in all 24 samples of a given tissue type (a very loose filtering); this left 17 581 probes from guard cell samples and 17 827 probes from leaf samples for further analysis. The subsequent analysis of variance (ANOVA) step discards those genes that show no differential expression in any of our microarray comparisons.
Data obtained from the microarray hybridizations were confirmed by real‐time quantitative reverse transcriptase PCR (QPCR) of selected genes (Supplementary information 18). For each gene, the differential expression pattern of all genotype × treatment combinations was compared between microarray data and QPCR data (8 comparisons per gene, for a total of 64 comparisons). Microarray and QPCR patterns of gene differential expression exhibit an 86% (=55/64) match (see Supplementary Figure S16).
Boolean model of G‐protein regulatory modes
We represent GPA1 and AGB1 by Boolean variables that can have two states: 1 denoting ‘on’ (not knocked out) and 0 denoting ‘off’ (knocked out). Each of the four genotypes considered in our study corresponds to one of the 22=4 possible combinations of the states of GPA1 and AGB1, from GPA1=0 and AGB1=0 for the agb1 gpa1 double mutant to GPA1=1 and AGB1=1 for wild type. The regulatory modes of the G‐protein are modeled by a Boolean function A(GPA1, AGB1). There is a total of 24=16 possible Boolean functions fitting A(GPA1, AGB1), which we denote as Ai(GPA1, AGB1), i=1,2,…,16. Except A1 and A16, each of these corresponds to one regulatory mode of the G‐protein and reflects the activity of a G‐protein‐regulated mediator that further regulates the expression level of a target gene. All possible regulatory modes of the G‐protein and the associated expression patterns of their potential target genes can be found in Supplementary Table S1.
Boolean model of G‐protein–ABA co‐regulatory modes
We also represent ABA by a Boolean variable with 1 indicating the presence of ABA (i.e. ABA treatment) and 0 representing the absence of ABA (i.e. solvent control). We assume that co‐regulation by the G‐protein and ABA can be described as F(ABA, GPA1, AGB1)=CABA+B(ABA, GPA1, AGB1), where CABA is a constant representing a G‐protein‐independent regulatory effect of ABA on a given gene, and B(ABA, GPA1, AGB1) is a Boolean function representing the possible combinatorial regulation of ABA and the G‐protein. A value for CABA that is not vanishingly small indicates that there is a G‐protein‐independent signaling pathway by which ABA regulates the gene. Although a general B(ABA, GPA1, AGB1) function allows a variety of combinatorial regulatory schemes by individual G‐protein subunits together with ABA, we have not found evidence of ABA participating in contrasting combinatorial regulation with individual G‐protein subunits (see Supplementary information 5 for a detailed analysis). Thus, in the following analysis, we use the simplified form B(ABA, A(GPA1, AGB1)) that excludes these cases.
In total, there are 24=16 Boolean functions fitting the B(ABA, A(GPA1, AGB1)) form, denoted as B1, B2,…, B16. Each of these represents a set of Boolean rules due to the fact that A(GPA1, AGB1) is already a Boolean function. Because of logical equivalence, in total, there are 72 non‐trivial unique Boolean functions that we can use to distinguish between (co)regulatory modes of the G‐protein and ABA (see Supplementary information 6 for the derivation of these Boolean functions). Together with CABA, F(ABA, GPA1, AGB1) is able to identify 142 (i.e. (72−2) × 2+2) regulatory modes, which can categorize the regulation by the G‐protein and ABA into five types: ABA‐only regulation (B4, B13), G‐protein‐only regulation (B6 (ABA, Ai)=B11 (ABA, A17‐i), i=2,…,15, combined with a non‐significant CABA), G‐protein–ABA additive regulation (B6 (ABA, Ai)=B11 (ABA, A17‐i), i=2,…,15, combined with a significant CABA), ABA–G‐protein combinatorial regulation (B2(ABA, Ai)=B3(ABA, A17‐i), B5(ABA, Ai)=B9(ABA, A17‐i), B8(ABA, Ai)=B12(ABA, A17‐i), B14(ABA, Ai)=B15(ABA, A17‐i), i=2,…,15, combined with a non‐significant CABA) and ABA–G‐protein mixed regulation (B2(ABA, Ai)=B3(ABA, A17‐i), B5(ABA, Ai)=B9(ABA, A17‐i), B8(ABA, Ai)=B12(ABA, A17‐i), B14(ABA, Ai)=B15(ABA, A17‐i), i=2,…,15, combined with a significant CABA), which further allow us to examine the representation of the nine proposed G‐protein and ABA signaling pathways (Supplementary Figure S1).
Association of genes with regulatory modes
As a first step, we use ANOVA (Sahai and Agell, 2000) to identify the genes that are differentially expressed with respect to genotypes or ABA. Within the control or ABA treatment categories, the expression of each gene is affected merely by genotypes, so we perform one‐way ANOVA and denote the corresponding P‐value as Pgeno. The expression of each gene in all samples is affected by both genotypes and ABA. To extract the ABA effect, we use two‐way ANOVA and denote the ABA‐associated P‐value as PABA.
For association of genes with G‐protein regulatory modes, we use Pgeno<0.01 for each individual treatment × tissue combination. For association of genes with G‐protein–ABA co‐regulatory modes, we designate P>0.05 as no differential expression and P<0.02 as differential expression. Specifically, there are 6975 genes in guard cells and 6299 genes in leaves that have PABA>0.05 and Pgeno>0.05 (no differential expression with respect to genotypes and also no differential expression with respect to ABA) and these genes are not considered further. Those genes that have no differential expression with respect to genotypes but are differentially expressed with respect to ABA (PABA<0.02 and Pgeno>0.05 for both ABA=0 and ABA=1) are candidates for ABA‐only regulation (Figure 1A). For those genes with differential expression in at least one condition of ABA=0 and ABA=1, we determine and subtract CABA from the expression values of the genes in the condition (presence or absence of ABA) corresponding to the higher level. We choose CABA=0.6 as the lower threshold of significance (Supplementary information 7). Those genes that have genotype‐dependent differential expression in both ABA=0 and ABA=1 (Pgeno<0.02 for both ABA=0 and ABA=1) are candidates for G‐protein‐only regulation (if their CABA is not significant) (Figure 1B) or G‐protein–ABA additive regulation (if their CABA is significant) (Figure 1C). Those genes that exhibit differential expression with respect to genotypes for ABA=0 or ABA=1, but not both (Pgeno<0.02 in one case and Pgeno>0.05 in the other case) are candidates for ABA–G‐protein combinatorial regulation (if their CABA is not significant) (Figure 1D) or ABA–G‐protein mixed regulation (if their CABA is significant) (Supplementary Figure S1E).
Each regulatory mode of the G‐protein and ABA corresponds to a unique three‐valued vector, which we call an idealized differential expression pattern, through comparison of the gene expression values between pairs of genotypes or within the same genotype in the absence or presence of ABA (Figure 2; Supplementary information 8). The real differential expression pattern of a gene is constructed using a signal‐to‐noise ratio metric (Golub et al, 1999). We introduce a measure to quantify the correlation between the real differential expression pattern of a gene with a Boolean function's three‐valued vector (see Supplementary information 8). Each gene is associated with that Boolean function for which it has the maximum correlation score. We classify a gene as regulated by a particular mode if its correlation score is higher than a threshold.
False discovery rates and statistical significance evaluation
Permutation tests are used to determine the statistical significance of each gene's association with a given function, and to estimate the false discovery rate (FDR) in putative gene groups of particular type (see Supplementary information 8). The correlation threshold is chosen in such a way to control the FDR at an acceptable level (see Supplementary information 9 and 10). Specifically, for association of genes with G‐protein regulatory modes, in all four data sets, differentially expressed genes are associated with Boolean functions A(GPA1, AGB1) if the correlation between their differential expression patterns and the idealized pattern exceeds the threshold R0=1.1, which controls the FDR within 0.02 (see Supplementary Figure S6). For association of genes with G‐protein–ABA co‐regulatory modes the correlation threshold R0=0.9 is used, which controls the FDR of ABA–G‐protein combinatorially or mixed regulated genes within 0.05 and the FDR of G‐protein‐only or G‐protein–ABA additively regulated genes within 0.005 (Supplementary Figures S7 and S8). As the number of ABA‐only regulated genes is large, we set 2.0 as the threshold of ABA‐only regulated genes, which controls the FDR within 0.0005 (Supplementary Figure S9). These thresholds and FDRs apply to both guard cell data and leaf data.
This research was supported by NSF grants MCB‐0209694, MCB‐0618402, CCF‐0643529, and USDA grant 2006‐35100‐17254. We thank Dr Craig Praul for outstanding technical advice and assistance on sample preparation and microarray hybridization.
Conflict of Interest
The authors declare that they have no conflict of interest.
Supplementary information, Supplementary figs S1–16, Supplementary tables S1–9 [msb201028-sup-0001.pdf]
This file contains G‐protein regulated genes in guard cells and leaves with or without ABA treatment. [msb201028-sup-0002.xls]
This file contains genes regulated by the G‐protein and/or ABA in leaves. [msb201028-sup-0003.xls]
This file contains genes regulated by the G‐protein and/or ABA in guard cells. [msb201028-sup-0004.xls]
This is an open access article under the terms of the Creative Commons Attribution‐NonCommercial‐NoDerivs License, which permits use and distribution in any medium, provided the original work is properly cited, the use is non‐commercial and no modifications or adaptations are made.
- Copyright © 2010 EMBO and Macmillan Publishers Limited