|
nif:isString
|
-
To formulate the model, we used a rule-based approach, which simplifies the consideration of site-specific details (e.g., tracking the phosphorylation status of amino acid residues) [51]. We used the model-specification language of BioNetGen, which is called BNGL, to define the types of molecules included in our model and to write rules for interactions. The conventions of BNGL are described elsewhere in detail [52,53]. A rule provides necessary and sufficient conditions for occurrence of an interaction and a rate law that governs all reactions implied by the rule. The model accounts for seven types of molecules: (FKBP1A-bound) rapamycin, a natural product, and the proteins AMPK, MTOR, RPTOR, ULK1, EIF4EBP1, and AMBRA1. The model includes 27 unidirectional rules and 2 bidirectional rules. Below, we present and discuss the BNGL-encoded molecule type definitions and rules that comprise the model.
The focus of this section is on explaining the formal representation of molecules considered in the model. Rapamycin in complex with FKBP1A (UniProt entry P62942) binds MTOR [17,18], which prevents MTOR from associating with RPTOR [54]. For simplicity, we only consider the form of rapamycin in complex with FKBP1A, which we denote as rapamycin*. In the model, this form of rapamycin is named rapa and is taken to have a single component named mtor, which is responsible for interaction with MTOR: rapa(mtor) The level of rapamycin* is treated as an input. The level of rapamycin* is not affected by the dynamics captured in the model. We chose to study responses to rapamycin rather than physiological stimuli, such as amino acids or growth factors, because rapamycin treatment is a stimulus that can be readily manipulated in the laboratory and rapamycin directly modulates MTORC1 activity. Regulation of MTORC1 activity in response to physiological stimuli is less direct and mediated by a complex regulatory network, which integrates multiple signals.
As a simplification, we only consider the form of AMPK with activating phosphorylation in its kinase domain, i.e., the form of AMPK with phosphorylation at T172 (T183) in the PRKAA2 (PRKAA1) isoform of the α subunit of AMPK (UniProt entries P54646 and Q13131). We denote this form of AMPK as AMPK*. In the model, AMPK* interacts directly with ULK1 via an undefined region in AMPK*. Furthermore, AMPK* phosphorylates sites in RPTOR and ULK1, and AMPK* is a substrate of ULK1. The kinase domain of AMPK* is considered implicitly. We represent AMPK* as follows: AMPK(ulk1,ST∼0∼P) The ulk1 component represents the undefined region in AMPK* that mediates interaction with ULK1. The ST component represents serine/threonine residues in the α subunit of AMPK that are phosphorylated by ULK1. This component is taken to have an internal state, which is either “unphosphorylated” (0) or “phosphorylated” (P). Internal states are abstractions useful for representing local properties of sites, such as location, conformation, or post-translational modification status. In a molecule type definition, the names of all possible internal states of a site are listed after the name of that site, with each state name being prefixed by a tilde. The level of AMPK* is treated as an input. The level of AMPK* is not affected by the dynamics captured in the model.
MTOR (UniProt entry P42345) interacts directly with rapamycin* and RPTOR [54]. The interactions with rapamycin* and RPTOR, which are mutually exclusive, are mediated by the FRB domain and HEAT repeats in MTOR, respectively [17,18,54]. Furthermore, MTOR, when part of MTORC1, phosphorylates substrates in ULK1 [10,14,55] and EIF4EBP1 [56,57]. We consider the kinase domain of MTOR implicitly. Thus, MTOR is represented as follows: MTOR(HEAT,FRB) In the model, MTOR’s kinase is taken to be active when bound to RPTOR, such that only MTOR in complex with RPTOR is able to phosphorylate substrates. As a simplification, the only components of MTORC1 that we track in our model are MTOR and RPTOR. Thus, the complex of MTOR and RPTOR is taken to be equivalent to MTORC1.
RPTOR (UniProt entry Q8N122) interacts directly with MTOR, ULK1, and EIF4EBP1. The interactions with MTOR and EIF4EBP1 are mediated by WD40 repeats and the RNC domain in RPTOR, respectively. The interaction with ULK1 is mediated by an undefined region in RPTOR. Several sites in RPTOR, including S792, S855, and S859, are substrates of AMPK and/or ULK1 [11,15]. We represent RPTOR as follows: RPTOR(RNC,ulk1,WD40,S792∼0∼P,S855_S859∼0∼P∼PP) The ulk1 component represents the undefined region in RPTOR that mediates interaction with ULK1. The S792 component is allowed to have one of two internal states: 0, which is taken to represent the unphosphorylated form of S792, or P, which is taken to represent the phosphorylated form of S792. As a simplification, the S855 and S859 residues are lumped together and represented by a single component named S855_S859. This component is taken to have one of three internal states: unphosphorylated (0), singly phosphorylated (P), or doubly phosphorylated (PP).
ULK1 (UniProt entry O75385) interacts directly with RPTOR and AMPK; these interactions are mediated by undefined regions within the proline/serine-rich (PS) domain of ULK1 [10]. Several sites in ULK1, including S317, S758, and S778, are substrates of MTORC1 or AMPK, and ULK1 phosphorylates several substrates in RPTOR, AMPK and AMBRA1. We consider the kinase domain of ULK1 implicitly. Thus, ULK1 is represented as follows: ULK1(straptor,stampk,S317∼0∼P,S758∼0∼P,S778∼0∼P) The components named straptor and stampk represent the undefined regions responsible for interactions with RPTOR and AMPK. The component S758 is a substrate of MTORC1, and the components S317 and S778 are substrates of AMPK. These components are each taken to have one of two internal states: unphosphorylated (0) or phosphorylated (P). In the model, ULK1’s kinase is taken to be active when phosphorylated at both S317 and S778 and not phosphorylated at S758.
When recruited to RPTOR via its RCR domain, EIF4EBP1 (UniProt entry Q13541) can be phosphorylated by MTORC1 at several sites, including T37, T46, S65 and T70 [58]. Phosphorylation of S65 and T70 is more sensitive to rapamycin treatment than T37 and T46 [59,60]. Thus, we represent EIF4EBP1 as follows: EIF4EBP1(RCR,S65_T70∼0∼P) As can be seen, as a simplification, we represent S65 and T70 (and other less rapamycin sensitive sites) together by a single component named S65_T70. This component is taken to have an internal state: either unphosphorylated (0) or phosphorylated (P). EIF4EBP1 is included in the model as a representative of MTORC1 substrates involved in regulating translation. EIF4EBP1 represses translation when hypophosphorylated and is inert when phosphorylated by MTORC1.
AMBRA1 (UniProt entry Q9C0C7) is phosphorylated by ULK1. We represent AMBRA1 as follows: AMBRA1(ST∼0∼P) The ST component represents serine/threonine residues that are phosphorylated by ULK1. This component is taken to have an internal state: either unphosphorylated (0) or phosphorylated (P). AMBRA1 is included in the model as a reporter of ULK1 kinase activity and autophagy level. Phosphorylation of AMBRA1 by ULK1 promotes autophagy.
The focus of this section is on explaining the formal representation of interactions considered in the model. The interactions discussed here are those that seemed most relevant for understanding mutual inhibition of MTORC1 and ULK1 after a literature search aimed at identifying interactions among AMPK, MTORC1, and ULK1. Not all known interactions among this triad are included in the model; omitted interactions are considered below and also in the Discussion section. FKBP1A-bound rapamycin (rapamycin*) can interact with MTOR provided that MTOR is not bound to RPTOR. Rapamycin* binding to MTOR and RPTOR binding to MTOR are mutually exclusive. We represent reversible binding of rapamycin* to MTOR as follows: rapa(mtor)+MTOR(HEAT,FRB)<->rapa(mtor!1).MTOR(HEAT,FRB!1)a1,d1(1) where a1 and d1 are the forward and reverse rate constants for this interaction. Note that the rule of Equation (1) has both a forward and reverse direction, as indicated by the symbol “<->.” Thus, the rule can be read from left to right, or right to left; it is one of the two reversible rules included in the model. As indicated by sharing of the bond index “1” by the mtor and FRB components on the right-hand side of the rule in Equation (1), rapamycin* binds the FRB domain in MTOR. (Bond indices are prefixed by an exclamation mark.) Inclusion of the HEAT component in the left-hand side of Equation (1) without specification of a binding partner indicates that this component must be free in order for the rule to be applicable. This constraint is introduced to ensure that rapamycin* and RPTOR bind MTOR with mutual exclusivity. In the model, only RPTOR interacts with HEAT repeats in MTOR. Thus, the rule of Equation (1) indicates that MTOR must be free of RPTOR to interact with rapamycin*. We note that the reverse (i.e., right-to-left) reading of the rule of Equation (1) indicates that the HEAT repeats in MTOR must be free in order for rapamycin* to dissociate from MTOR. Writing the rule with this restriction, which is spurious (but inconsequential for the model as specified), has the benefit of allowing for concise specification of the forward and reverse rules for binding of rapamycin* to MTOR. Only a single line of BNGL code is necessary. The restriction has no consequence whatsoever because the rules of the model never allow the HEAT and FRB sites in MTOR to be bound at the same time. If the model were modified in some way to allow simultaneous binding of these sites, then the rule of Equation (1) would likely need to be rewritten as two separate, unidirectional rules.
RPTOR can interact with MTOR provided that MTOR is not bound to rapamycin*. We represent reversible binding of RPTOR to MTOR as follows: RPTOR(WD40)+MTOR(HEAT,FRB)<->RPTOR(WD40!1).MTOR(HEAT!1,FRB)a2,d2(2) where a2 and d2 are the forward and reverse rate constants for this interaction. As indicated, WD40 repeats in RPTOR interact with HEAT repeats in MTOR. Inclusion of the FRB component in the left-hand side of Equation (2) without specification of a binding partner indicates that this component must be free in order for the rule to be applicable. In the model, only rapamycin* interacts with FRB domain in MTOR. Thus, the rule of Equation (2) indicates that MTOR must be free of rapamycin* in order to interact with RPTOR.
RPTOR can interact reversibly with ULK1. This interaction, which is mediated by undefined regions in RPTOR and ULK1, is destabilized by phosphorylation of RPTOR at S792, S855, and S859. We represent reversible binding of RPTOR to ULK1 using three unidirectional rules, as follows: RPTOR(RNC,ulk1,S792~0,S855_Ser859~0)+ULK1(straptor)-> \ RPTOR(RNC,ulk1!1,S792~0,S855_Ser859~0).ULK1(straptor!1)a3(3A) RPTOR(ulk1!1).ULK1(straptor!1)->RPTOR(ulk1)+ULK1(straptor)d3(3B) RPTOR(ulk1!1,S792~P,S855_Ser859~PP).ULK1(straptor!1)-> \ RPTOR(ulk1,S792~P,S855-9~PP)+ULK1(straptor)d3max(3C) Note that unidirectionality is indicated by the symbol “->.” A unidirectional rule is only read from left to right. The first rule, Equation (3A), indicates that RPTOR, when free of EIF4EBP1 (which binds the RNC domain in RPTOR) and not phosphorylated at S792, S855, and S859, is able to associate with ULK1 with rate constant a3. (This rule requires that the cognate binding sites in RPTOR and ULK1, which are represented by the components named ulk1 and straptor, be free.) The second rule, Equation (3B), indicates that RPTOR is able to dissociate from ULK1 with rate constant d3. Dissociation can occur whenever a bond between RPTOR and ULK1 exists. In other words, according to this rule, there are no contextual constraints on dissociation. The third rule, Equation (3C), indicates that dissociation occurs with a different rate constant, d3max (>d3), when the indicated sites in RPTOR are phosphorylated. We note that a backslash (\) is used to mark a line break in BNGL.
RPTOR can interact reversibly with EIF4EBP1. This interaction occurs between the RNC domain in RPTOR and the RCR domain in EIF4EBP1. We represent reversible binding of RPTOR to EIF4EBP1 using two unidirectional rules, as follows: RPTOR(RNC,ulk1,S792~0,S855_Ser859~0)+EIF4EBP1(RCR,S65_T70~0)-> \ RPTOR(RNC!1,ulk1,S792~0,S855_Ser859~0).EIF4EBP1(RCR!1,S65_T70~0)a4(4A) RPTOR(RNC!1).EIF4EBP1(RCR!1)->RPTOR(RNC)+EIF4EBP1(RCR)d4(4B) The first rule, Equation (4A), indicates that RPTOR, when free of ULK1 (which interacts with an undefined region in RPTOR that we refer to as ulk1) and not phosphorylated at S792, S855, and S859, is able to associate with EIF4EBP1 with rate constant a4. The sites that mediate association, the RNC domain in RPTOR and the RCR domain in EIF4EBP1, must be free and available for interaction, as indicated. The second rule, Equation (4B), indicates that RPTOR is able to dissociate from EIF4EBP1 with rate constant d4.
AMPK can interact reversibly with ULK1. This interaction occurs between an undefined region in AMPK, denoted ulk1, and an undefined region in ULK1, denoted stampk. We represent reversible binding of AMPK to ULK1 using two unidirectional rules, as follows: AMPK(ulk1,T172~P)+ULK1(stampk,S758~0)-> \AMPK(ulk1!1,T172~P).ULK1(stampk!1,S758~0)a5(5A) AMPK(ulk1!1).ULK1(stampk!1)->AMPK(ulk1)+ULK1(stampk)d5(5B) The first rule, Equation (5A), indicates that AMPK, when phosphorylated at T172 and free of ULK1, is able to associate with ULK1, provided that ULK1 is not phosphorylated at S758. Recall that phosphorylation of S758 inhibits AMPK-ULK1 interaction. Association of AMPK and ULK1 occurs with rate constant a5. The sites ulk1 and stampk must be free for association to occur, as indicated. The second rule, Equation (5B), indicates that AMPK is able to dissociate from ULK1 with rate constant d5.
MTORC1 phosphorylates ULK1 at S758: We represent MTORC1-mediated inhibitory phosphorylation of ULK1 at S758 using a unidirectional rule, as follows: MTOR(HEAT!1).RPTOR(WD40!1,ulk1!2).ULK1(straptor!2,stampk,S758~0)-> \MTOR(HEAT!1).RPTOR(WD40!1,ulk1!2).ULK1(straptor!2,stampk,S758~P)p1(6) where p1 is a (pseudo first-order) rate constant. This rule requires that the stampk site in ULK1 be free (of AMPK) in order for MTORC1 to phosphorylate S758. This requirement is introduced for the following reason. It is known that phosphorylation of S758 prevents AMPK from binding ULK1. Thus, we assume that S758 is in the interface between these two proteins. If this assumption is correct, then it follows that S758 should not be accessible for phosphorylation when ULK1 is bound to AMPK (i.e., when the stampk site is bound). The rule of Equation (6) indicates that phosphorylation of S758 by MTORC1 further requires MTOR, RPTOR and ULK1 to be together in a complex. In the model, this complex forms through association of RPTOR with MTOR according to Equation (2) and through association of RPTOR with ULK1 according to Equation (3). The order of association events is of no importance. In the model, components of MTORC1 beyond MTOR and RPTOR (e.g., MLST8) are considered implicitly. We assume that MTORC1 is fully formed and competent to phosphorylate substrates whenever MTOR and RPTOR are in association. Phosphorylation of S758 in ULK1 is inhibitory because phosphorylation of this residue prevents AMPK from binding ULK1 (see Equation (5A)).
MTORC1 phosphorylates EIF4EBP1 at S65 and T70: We represent MTORC1-mediated phosphorylation of S65 and T70 in EIF4EBP1 using a unidirectional rule, as follows: MTOR(HEAT!1).RPTOR(WD40!1,RNC!2).EIF4EBP1(RCR!2,S65_T70~0)->\MTOR(HEAT!1).RPTOR(WD40!1,RNC!2).EIF4EBP1(RCR!2,S65_T70~P)p2(7) where p2 is a (pseudo first-order) rate constant. Equation (7) indicates that phosphorylation of S65 and T70 requires MTOR, RPTOR and EIF4EBP1 to be together in a complex. Recall that S65 and T70 in EIF4EBP1 are represented, for simplicity, as a single component named S65_T70. (See the molecule type definition for EIF4EBP1.) ULK1 phosphorylates RPTOR at S792, S855, and S859: We represent ULK1-mediated inhibitory phosphorylation of S792, S855, and S859 in RPTOR using three unidirectional rules, as follows: RPTOR(ulk1!1,S792~0).ULK1(straptor!1,S317~P,S778~P)-> \RPTOR(ulk1!1,S792~P).ULK1(straptor!1,S317~P,S778~P)p3(8A) RPTOR(ulk1!1,S855_S859~0).ULK1(straptor!1,S317~P,S778~P)-> \RPTOR(ulk1!1,S855_S859~P).ULK1(straptor!1,S317~P,S778~P)2*p4(8B) RPTOR(ulk1!1,S855_S859~P).ULK1(straptor!1,S317~P,S778~P)-> \RPTOR(ulk1!1,S855_S859~PP).ULK1(straptor!1,S317~P,S778~P)p4(8C) where p3 and p4 are (pseudo first-order) rate constants. As indicated, phosphorylation of S792, S855, and S859 in RTPOR by ULK1 requires that RPTOR and ULK1 be connected through the binding sites ulk1 in RPTOR and straptor in ULK1. Phosphorylation of S317 and S778 in ULK1, which makes straptor competent for interaction with ulk1 in RPTOR, is an additional prerequisite for ULK1-mediated phosphorylation of S792, S855, and S859 in RPTOR. Recall that S855 and S858 in RPTOR are represented, for simplicity, as a single component named S855_S859, which is allowed to be unphosphorylated (0), singly phosphorylated (P), or doubly phosphorylated (PP). (See the molecule type definition for RPTOR.) The rate constant in the rule of Equation (8B) is specified as twice that for the rule of Equation (8C) because the 0 to P state transition for site S855_S859 (i.e., the transition from an unphosphorylated state to a singly phosphorylated state) is expected to occur twice as fast as the P to PP transition (i.e., the transition from a singly to doubly phosphorylated state).
ULK1 phosphorylates AMBRA1 at S/T residues: We represent ULK1-mediated phosphorylation of AMBRA1 at undefined serine and threonine (S/T) residues using a unidirectional rule, as follows: ULK1(straptor,S317~P,S778~P)+AMBRA1(ST~0)-> \ULK1(straptor,S317~P,S778~P)+AMBRA1(ST~P)p5(9) where p5 is a (pseudo first-order) rate constant. Equation (9) indicates that ULK1, when phosphorylated at S317 and S778 and not bound to RPTOR (i.e., when the straptor site in ULK1 is free), is able to phosphorylate AMBRA1 at undefined S/T residues, which are represented as a single component named ST. This rule reflects the inhibitory effect of RPTOR on ULK1 kinase activity and the requirement for activating phosphorylation at S317 and S778, which is mediated by AMPK. Note that ULK1 and AMBRA1 are not required to form a complex for ULK1-mediated phosphorylation of AMBRA1 to occur. Thus, we are assuming that the enzymatic reactions defined by the rule of Equation (9) are occurring in the regime of substrate limitation, far from enzyme saturation.
ULK1 phosphorylates AMPK at S/T residues: We represent ULK1-mediated phosphorylation of AMPK at undefined serine and threonine residues using a unidirectional rule, as follows: ULK1(straptor,S317~P,S778~P)+AMPK(ST~0)-> \ULK1(straptor,S317~P,S778~P)+AMPK(ST~P)p6(10) where p6 is a (pseudo first-order) rate constant. The rule of Equation (10) is the same as that of Equation (9) except that the ULK1 substrates, which are represented by the component named ST, are in AMPK instead of AMBRA1.
AMPK phosphorylates ULK1 at S317 and S778: We represent AMPK-mediated activating phosphorylation of ULK1 at S317 and S778 using two unidirectional rules, as follows: AMPK(ulk1!1,T172~P,ST~0).ULK1(stampk!1,straptor,S317~0)-> \AMPK(ulk1!1,T172~P,ST~0).ULK1(stampk!1,straptor,S317~P)p7(11A) AMPK(ulk1!1,T172~P,ST~0).ULK1(stampk!1,S778~0)-> \AMPK(ulk1!1,T172~P,ST~0).ULK1(stampk!1,S778~P)p8(11B) where p7 and p8 are (pseudo first-order) rate constants. Each rule indicates that AMPK kinase activity requires 1) phosphorylation at T172 in the kinase domain of the PRKAA2 isoform of the α subunit (or equivalently, phosphorylation at T183 in the kinase domain of the PRKAA1 isoform of the α subunit), 2) the absence of phosphorylation at S/T residues in the α subunit outside the kinase domain, and 3) association of AMPK with ULK1, which is mediated by the ulk1 and stampk sites in AMPK and ULK1, respectively. Furthermore, AMPK-mediated phosphorylation of S317 in ULK1, but not S778, requires that ULK1 be free of RPTOR (i.e., the straptor site in ULK1 must be free). The reason for this constraint is that RPTOR appears to mask S317 when bound to ULK1, which would prevent phosphorylation of S317 by AMPK. RPTOR interacts with the PS domain in ULK1, which contains S317 [10,55].
AMPK phosphorylates RPTOR at S792: We represent AMPK-mediated inhibitory phosphorylation of RPTOR at S792 using a unidirectional rule, as follows: AMPK(T172~P)+RPTOR(S792~0)->AMPK(T172~P)+RPTOR(S792~P)p9(12) where p9 is a (pseudo first-order) rate constant. Equation (12) indicates that AMPK*, regardless of the phosphorylation status of S/T residues in the α subunit outside the kinase domain (i.e., regardless of the internal state of the ST component of AMPK), is able to phosphorylate RPTOR at S792. Note that AMPK and RPTOR are not required to form a complex for AMPK-mediated phosphorylation of RPTOR to occur. Thus, we are assuming that the enzymatic reactions defined by the rule of Equation (12) are occurring in the regime of substrate limitation, far from enzyme saturation.
Dephosphorylation of S758 in ULK1: We represent dephosphorylation of S758 in ULK1 as follows: ULK1(S758~P)->ULK1(S758~0)u0(13) where u0 is a (pseudo first-order) rate constant. In most cases, for simplicity, rules for dephosphorylation reactions are assigned the dephosphorylation rate constant u0. We assume that dephosphorylation reactions are mediated by constitutively active (i.e., unregulated) phosphatases, which are taken to be present in excess. This rule exemplifies our treatment of phosphatases, which generally are not as well characterized as kinases.
Dephosphorylation of S65 and T70 in EIF4EBP1: We represent dephosphorylation of S65 and T70 in EIF4EBP1 as follows: EIF4EBP1(S65_T70~P)->EIF4EBP1(S65_T70~0)u0(14) where u0 is a (pseudo first-order) rate constant. This rule is similar to that of Equation (12).
Dephosphorylation of S792, S855 and S859 in RPTOR: We represent dephosphorylation of S792, S855 and S859 in RPTOR as follows: RPTOR(S792~P)->RPTOR(S792~0)u1(15A) RPTOR(S855_S859~P)->RPTOR(S855_S859~0)u0(15B) RPTOR(S855_S859~PP)->RPTOR(S855_S859~P)2*u0(15C) where u0 and u1 are (pseudo first-order) rate constants. We assume that S792 in RPTOR is dephosphorylated more slowly than typical sites considered in the model. Consequently, the rule of Equation (15A) is assigned the rate constant u1 (instead of u0). The rate constant in the rule of Equation (15C) is specified as twice that for the rule of Equation (15B) because the PP to P state transition for site S855_S859 is expected to occur twice as fast as the P to 0 transition.
Dephosphorylation of S/T residues in AMBRA1: We represent dephosphorylation of S/T residues in AMBRA1 as follows: AMBRA1(ST~P)->AMBRA1(ST~0)u0(16) This rule is similar to that of Equation (13). Dephosphorylation of ST residues in AMPK: We represent dephosphorylation of S/T residues in AMPK as follows: AMPK(ST~P)->AMPK(ST~0)u2(17) where u2 is a (pseudo first-order) rate constant. We assume that S/T residues in AMPK are dephosphorylated more slowly than typical sites considered in the model, so the rule of Equation (17) is assigned the rate constant u2 (instead of u0). This rule is otherwise similar to that of Equation (13).
Dephosphorylation of S317 and S778 in ULK1: We represent dephosphorylation of S317 and S778 in ULK1 as follows: ULK1(straptor,S317~P)->ULK1(straptor,S317~0)u0(18A) ULK1(S778~P)->ULK1(S778~0)u0(18B) These rules are similar to that of Equation (13) except that the rule of Equation (18A) imposes a constraint on dephosphorylation of S317 in ULK1. Namely, dephosphorylation of this site is only allowed to occur when ULK1 is free of RPTOR (or equivalently, the straptor site in ULK1 is free). This requirement is introduced because S317 overlaps with the straptor site. We assume that S317 cannot be phosphorylated or dephosphorylated when RPTOR is in contact with the straptor site.
Dephosphorylation of S792 in RPTOR: We represent dephosphorylation of S792 in RPTOR as follows: RPTOR(S792~P)->RPTOR(S792~0)u0(19) This rule is similar to that of Equation (13). Even though the interactions considered in our model are fairly well established, information about the quantitative factors that influence these interactions (e.g., rate constants and concentrations) is scarce. To cope with this knowledge gap, we assigned values to the parameters of the model, which are summarized in Table 1, somewhat arbitrarily within ranges deemed to be reasonable, and we eschewed fine-tuning of parameter values. Thus, for each rate constant, we only specified an order of magnitude (i.e., we only used a single significant digit, 1). To offset uncertainty about parameter values, we focused on qualitative system behavior through bifurcation analyses and we also performed a systematic parameter sensitivity analysis. Our rationale for setting parameter values is further explained below.
Rate constants for dephosphorylation reactions: The model includes 10 rules for dephosphorylation reactions, which are given by Equations (13)–(19). Each of these rules is associated with a rate constant and, in accordance with the conventions of BNGL, a mass-action rate law that has this rate constant as its only parameter. We associated 8 of the 10 rules for dephosphorylation with a common rate constant, u0. We set the nominal value of u0 at 10−2 s−1 (Table 1). This setting essentially establishes the characteristic response time of the system, which is on the order of minutes [14]. Shang et al. [14] observed that inhibitory phosphorylation of ULK1 decreased over a period of minutes after activation of autophagy through starvation or rapamycin treatment. We associated the rule of Equation (15A), which characterizes dephosphorylation of S792 in RPTOR, with the rate constant u1. We set the nominal value of u1 at 10 times less than the value of u0 (Table 1). This setting is motivated by interactions of phosphorylated S792 (pS792) with 14–3–3 proteins, which are known to protect pS792 from dephosphorylation [11]. Thus, by setting u1 = 0.1 u0, we are implicitly considering the 14–3–3 binding partners of pS792 and accounting for shielding of pS792 from phosphatases by these binding partners. We associated the rule of Equation (15A), which characterizes dephosphorylation of the serine/threonine residues in AMPK that are substrates of ULK1, with the rate constant u2. We set the nominal value of at u2 at 100 times less than the value of u0 (Table 1). By setting u2 = 0.01 u0, we are assuming that the time scale of negative feedback from ULK1 to AMPK is much slower than the time required to respond to a stress stimulus. A slower time scale is necessary to allow for a period of autophagy, i.e., a period during which ULK1 kinase activity is sustained, which is presumably required for autophagosome production. Fast negative feedback would quickly shut off ULK1 activity after a stress stimulus, potentially severely limiting autophagy.
Rate constants for phosphorylation reactions: The model includes 10 rules for phosphorylation reactions, which are given by Equations (6)–(12). Each rule is associated with a rate constant. We assigned the phosphorylation rate constants p1, p2, p3, p4, p7, and p8 a common nominal value, 10 s−1 (Table 1). The value assigned to these rate constants is high (cf. the value assigned to u0, 0.01 s−1) because each of these rate constants characterizes a pseudo first-order reaction in which an enzyme and one of its substrates are colocalized within a protein complex, i.e., confined together within a small reaction volume. The remaining phosphorylation rate constants, p5, p6, and p9, characterize pseudo second-order reactions. We assume that each of these (overall) reactions corresponds to a Michaelis-Menten reaction scheme operating in the regime of substrate limitation, far from enzyme saturation. Thus, each second-order rate constant corresponds to a ratio of the form kcat/KM. We set the nominal value of p5 at 10−4 (molecule/cell)−1s−1 (Table 1). This setting is arbitrary and has little influence on predicted system behavior (Fig. 8) because p5 characterizes ULK1-mediated phosphorylation of AMBRA1. In the model, AMBRA1 phosphorylation serves as a reporter of ULK1 kinase activity. We set the nominal value of p6, which characterizes ULK1-mediated inhibitory phosphorylation of AMPK (i.e., the strength of negative feedback), at 10−6 (molecule/cell) −1s−1 (Table 1). In initial exploratory simulations, this parameter setting was found to allow for both robust activation of autophagy and nearly complete inhibition of autophagy through negative feedback from ULK1 to AMPK. It should be noted that ULK1 is responsible for phosphorylating multiple serine/threonine residues in AMPK, which in the model are lumped together. Thus, p6 is an effective rate constant. We set the nominal value of p9, which characterizes AMPK-mediated inhibitory phosphorylation of RPTOR, at 0 (Table 1). In initial exploratory simulations, we found that phosphorylation of S792 in RPTOR by AMPK limits the region of bistability (in the case without negative feedback from ULK1 to AMPK) and the region of oscillatory behavior (in the case with negative feedback from ULK1 to AMPK) in the parameter space of inputs. Thus, we set the nominal value of p9 to 0, and we gave non-zero values of p9 special attention in our evaluation of the model. In other words, we performed bifurcation analyses in which p9 was a bifurcation parameter (Fig. 7).
Rate constants for association and dissociation reactions: The model includes two reversible rules and seven unidirectional rules for association/dissociation reactions, which are given by Equations (1)–(5). The rules are associated with a total of 11 rate constants. Nominal values for the rate constants of association and dissociation reactions (Table 1) were set with the considerations noted below in mind. We selected values for a1, d1, a2, and d2 such that the affinity of rapamycin* for MTOR is 10-fold greater than the affinity of RPTOR for MTOR. In other words, we selected values for these rate constants such that a1/ d1 = 10 a2/d2. This approach is motivated by the ability of rapamycin to disrupt binding of RPTOR to MTOR [61]. In initial exploratory simulations, we found that the rate constants governing ULK1 interactions with AMPK and RPTOR affect the bias of a cell toward translation (or autophagy). An increase in the affinity of ULK1 for AMPK makes a cell more prone for autophagy. Conversely, an increase in the affinity of ULK1 for RPTOR makes a cell more prone for translation. Thus, we selected values for a3, d3, a5, and d5 that yield roughly balanced propensities for translation and autophagy. The value of d3,max was set such that d3,max is much greater than d3: we set d3,max = 100 d3. Recall that RPTOR, when phosphorylated at particular residues, dissociates from ULK1 faster than it otherwise would (cf. Equations (3B) and (3C)). The rate constant d3 characterizes slow dissociation of RPTOR from ULK1 and the rate constant d3,max characterizes fast dissociation of RPTOR from ULK1. Finally, we selected values for a4 and d4, which govern interaction of EIF4EBP1 with RPTOR, so as to avoid sequestration of MTORC1 by EIF4EBP1. Recall that EIF4EBP1 is included in the model simply to serve as a reporter of MTORC1 kinase activity.
In U2OS cells, for example, the total cellular abundances of ULK1 and AMBRA1 are on the order of 104 copies per cell, whereas the total cellular abundances of AMPK, MTOR, RPTOR, and EIF4EBP1 are far greater, roughly on the order of 106 copies per cell [62]. We assume that all proteins are present at effective concentrations of around 104 copies per cell. This assumption is introduced to account for sequestration of AMPK and MTORC1 away from ULK1 by their binding partners not considered in the model, as well as the different spatial distributions of ULK1 and MTORC1 (punctate and localized to membrane sites of autophagosome formation in the case of ULK1 vs. cytosolic and lysosomally localized in the case of MTORC1) [40,63]. Comparable effective levels of ULK1 and MTORC1 is consistent with the observation that ULK1, through phosphorylation of RPTOR, is able to relieve inhibition of autophagy by MTORC1, which seems unlikely if the entire pool of MTORC1 in a cell is available to interact with ULK1. We set concentrations or copy numbers as follows. We set the abundances of MTOR and RPTOR to be equal, at 2×104 copies per cell. This setting was guided by the 1:1 stoichiometry of MTOR and RPTOR in MTORC1. We set the abundance of ULK1 at half the abundance of MTOR, 104 copies per cell. It is important for ULK1 to be less abundant than MTORC1 so that MTORC1 is able to fully repress ULK1 while some fraction of MTORC1 remains available to interact with EIF4EBP1, which is included in the model as a reporter for translation activity. Full repression of ULK1 seems reasonable because autophagy can be repressed to undetectable levels when cells are not stressed. The levels of EIF4EBP1 and AMBRA1 were set at 104 copies per cell. These settings are largely inconsequential for predicted system behavior, because there is no feedback from these proteins to other molecules that are considered in the model. Recall that EIF4EBP1 and AMBRA1 are included in the model only for the purpose of providing readouts of MTORC1 and ULK1 activities, which presumably correlate with translation and autophagy levels. We treat the levels of AMPK* and rapamycin* as controllable inputs. Recall that AMPK* is distinct from total AMPK and active AMPK.
Our model omits several known interactions among the proteins of interest. Furthermore, AMPK, RPTOR, and ULK1 each contain more sites of phosphorylation than are included in the model. These interactions and sites, which are discussed below, were omitted to keep the model as simple as possible. Some of the omitted sites have roles similar to sites included in the model. In RPTOR, S863 and S877 are additional substrates of ULK1 that are involved in ULK1-mediated inhibition of MTORC1 kinase activity [15]. In ULK1, S479 and S556 appear to play a role analogous to that of S758; in other words, they appear to be additional substrates of MTORC1 involved in MTORC1-mediated inhibition of ULK1 kinase activity [14]. Besides S317 and S778, ULK1 contains additional sites phosphorylated by AMPK, including S467, S555, S637, and T659, which are involved in activation of ULK1 [64,65,66]. The model does not individually consider the multiple serine/threonine substrates of ULK1 that are found in the catalytic α subunit of AMPK, which are involved in negative feedback regulation of AMPK by ULK1 [24]. Instead, these sites are lumped together, i.e., they are treated as a single site. One of the omitted sites has an unclear functional role. Shang et al. [14] reported that ULK1 is phosphorylated at S638 by both AMPK and MTORC1. This finding is enigmatic because phosphorylation of other sites in ULK1 by MTORC1 has an inhibitory effect on ULK1 kinase activity, whereas phosphorylation of other sites in ULK1 by AMPK has an activating effect. Thus, the role of S638 may be distinct from that of other AMPK and MTORC1 substrates. Because the role of this site is unclear, we did not include S638 in the model. In formulating the model, we focused on interactions among AMPK, MTORC1, and ULK1, because this triad of kinases is recognized as playing a critical role in regulation of autophagy and translation [50]. However, the triad network is embedded within a larger regulatory network. Below, we call attention to complicating features of this larger network, which are beyond the scope of this study but that are, under certain conditions, likely to affect the behavior of the triad network. It is important to call attention to these complicating features to recognize the limitations of the present study. Recently, Nazio et al. [38] reported positive feedback from AMBRA1 to ULK1. In this feedback loop, ULK1-activated AMBRA1 interacts with the E3 ubiquitin ligase TRAF6, which enables TRAF6 to mediate the attachment of K63-linked chains of ubiquitin to ULK1, which stabilizes ULK1 kinase activity. Interestingly, AMBRA1-mediated ubiquitylation of ULK1 is antagonized by MTORC1-mediated phosphorylation of S52 in AMBRA1, which provides an additional (indirect) avenue for MTORC1 to inhibit ULK1 [38], which could strengthen the robustness of oscillatory behavior. In the model, we included AMBRA1 only as a reporter of ULK1 kinase activity, essentially considering a situation where, for simplicity, positive feedback from AMBRA1 to ULK1 is abrogated. Likewise, we did not consider negative feedback regulation of MTOR by amino acids generated through autophagy [28], which is similar to but distinct from negative feedback regulation of AMPK by ULK1. In addition to being a component of MTORC1, MTOR is a component of a second complex (MTORC2), which is distinguished by the subunit RICTOR. Following growth factor simulation, MTORC2 promotes full activation of AKT, which subsequently signals downstream to MTORC1 [67]. Although rapamycin is largely selective for MTORC1, on long time scales (e.g., 24 hours), rapamycin causes dissociation of MTORC2 [68]. Thus, it should be understood that in our analyses we are considering responses to rapamycin on relatively short time scales. As mentioned earlier, our model simplifies the regulation of MTORC1 by AMPK. It is well established that AMPK positively regulates TSC2 (within the TSC complex, which also consists of TSC1 and TBC1D7) via phosphorylation, which in turn suppresses the activity of the small GTPase RHEB, a direct activator of MTORC1 [45]. For simplicity, our model considers only the more direct route by which AMPK regulates MTORC1, which is achieved through phosphorylation of RPTOR.
Simulations, bifurcation analysis, and sensitivity analysis: We performed simulations using a deterministic method available within BioNetGen [53]. This method, which is an indirect method, involves first processing the rules of the model specification (S1 File) to obtain the reaction network implied by the rules, as well as the corresponding ordinary differential equations (ODEs) for mass-action kinetics. The reaction network consists of 173 chemical species and 6,581 unidirectional reactions. (The size of the reaction network reflects the number of protein phosphoforms and protein complexes that can arise from the interactions represented by the rules of the model.) BioNetGen’s built-in ODE solver, CVODE from the SUNDIALS package [69], was then used to numerically integrate the ODEs, using default settings. The steps described above were performed automatically by BioNetGen and invoked using point-and-click commands available within RuleBender [70], which provides a graphical user interface for accessing BioNetGen’s capabilities. Each one-dimensional bifurcation analysis was performed numerically as follows. We added a rule to the model specification (S1 File) that has the form 0->p k_epsilon or p->0 k_epsilon, where 0 is a source or sink, p is the bifurcation parameter, and k_epsilon is sufficiently small. The added rule has the effect of either very slowly increasing or very slowly decreasing the value of the bifurcation parameter (e.g., the level of AMPK* or rapamycin*), such that other processes are in a pseudo steady state, which slowly changes, as the bifurcation parameter is varied. Simulations were performed with the bifurcation parameter gradually increasing (from a sufficiently small value) and also gradually decreasing (from a sufficiently large value) to find stable steady states and the lower and upper bounds of stable limit cycles for values of the bifurcation parameter of interest. Two-dimensional bifurcation analyses were performed similarly. We allowed one parameter to vary while holding the other at a fixed value. We then repeated this procedure for different values of the second parameter. The numerical methods described above, in contrast to continuation methods implemented in software tools such as AUTO and MatCont, do not allow for characterization of unstable steady states or unstable limit cycles, but with these methods, it is possible to analyze the (large) system of ODEs of interest, for which continuation methods are prohibitively inefficient. In the sensitivity analysis of Fig. 8, we performed two series of one-dimensional bifurcation analyses for each of the 22 parameters considered. In the first of the two series, the bifurcation parameter was the level of AMPK*. In all cases, we considered levels of AMPK* from 0 to 106 copies per cell. In the second of the two series, the bifurcation parameter was the level of rapamycin*. In all cases, we considered levels of rapamycin* from 0 to 105 copies per cell. In each series, we considered 17 different values for one of the 22 parameters: we set Pm = Pm,0 ×10n/4 for n = −8,…,+8, where m is an index, in the range [1,…, 22], that identifies the parameter and Pm,0 is the nominal value of the parameter given in Table 1. The results from each series of bifurcation analyses were used to find the range of values for parameter Pm for which the following pattern, illustrated in Fig. C (S2 File), holds true: 1) low AMBRA1 phosphorylation and high EIF4EBP1 phosphorylation at low values of the bifurcation parameter, 2) oscillations in AMBRA1 and EIF4EBP1 phosphorylation at intermediate values of the bifurcation parameter, and 3) high AMBRA1 phosphorylation and low EIF4EBP1 phosphorylation at high values of the bifurcation parameter. The ranges so found are reported in Fig. 8.
|
![[This page as RDF]](http://data.gesis.org/somesci/static/rdf-icon.gif){kind=icon}