Mathematical modeling of the Candida albicans yeast to hyphal transition reveals novel control strategies

Candida albicans, an opportunistic fungal pathogen, is a significant cause of human infections, particularly in immunocompromised individuals. Phenotypic plasticity between two morphological phenotypes, yeast and hyphae, is a key mechanism by which C. albicans can thrive in many microenvironments and cause disease in the host. Understanding the decision points and key driver genes controlling this important transition and how these genes respond to different environmental signals is critical to understanding how C. albicans causes infections in the host. Here we build and analyze a Boolean dynamical model of the C. albicans yeast to hyphal transition, integrating multiple environmental factors and regulatory mechanisms. We validate the model by a systematic comparison to prior experiments, which led to agreement in 17 out of 22 cases. The discrepancies motivate alternative hypotheses that are testable by follow-up experiments. Analysis of this model revealed two time-constrained windows of opportunity that must be met for the complete transition from the yeast to hyphal phenotype, as well as control strategies that can robustly prevent this transition. We experimentally validate two of these control predictions in C. albicans strains lacking the transcription factor UME6 and the histone deacetylase HDA1, respectively. This model will serve as a strong base from which to develop a systems biology understanding of C. albicans morphogenesis.


Introduction
Candida albicans is a pleiomorphic, opportunistic fungal pathogen and an important cause of both superficial and systemic infections in humans, particularly in immunocompromised individuals. It is also responsible for 85-95% of all vulvovaginal infections resulting in doctor visits in otherwise healthy patients [1]. C. albicans forms biofilms on mucosal surfaces (e.g., oral, gastrointestinal tract, genitourinary tract, and vaginal) of the host as well as on surfaces of implanted medical devices (e.g., catheters, heart valves, and prosthetics), which are major reservoirs for infections [2,3].
Transitions between the yeast and hyphal phenotypes enable C. albicans to adapt to and persist in a wide range of environments. The yeast-form consists of single round cells that grow by forming daughter cells that bud and separate from mother cells. The hyphal form consists of long, multicellular branching tubular structures with parallel-sided walls, where the tips proliferate to elongate the hyphae [4]. C. albicans can also form an intermediate filamentous morphology called the pseudohyphal form, which consists of chains of cells with constrictions between mother-daughter cell pairs [4]. Transitioning from the yeast to hyphal phenotype is required for mucosal invasion [2,5] and biofilm formation, which are important mediators of infection [2,3,5,6]. The yeast to hyphal transition is regulated by many well-studied intracellular pathways that respond to external signals such as neutral or alkaline pH (pH > 6), farnesol levels, and temperature. These pathways converge on a handful of key transcription factors, defined as sequence-specific DNA-binding proteins, which regulate the transcription of hyphal-associated genes (HAGs). Epigenetic effects such as histone acetylation events at the promoters of HAGs also play important roles in the regulation of the expression of HAGs. The key negative regulator of the transcription of HAGs is Nrg1. The pattern of expression of Nrg1 and of its ability to bind to the promoter region of HAGs determines two phases of the yeast to hyphal transition. Hyphal initiation (the first cell division in the process that forms hyphae) requires a transient downregulation of the Nrg1 protein, whereas hyphal maintenance requires preventing Nrg1 from binding to the promoters of HAGs [7]. External signals initiate the downregulation of NRG1 transcription, while Nrg1 protein is prevented from binding to the promoters of HAGs by histone deacetylases (HDACs) such as Hda1.
Here we build a Boolean model integrating multiple extracellular signals governing the intracellular regulation of the yeast to hyphal morphological transition. We then use this model to conduct a thorough analysis of phenotype control, considering multiple possible control objectives and side effects, to rank the best and most robust control strategies.

Construction of the model
As a starting point to building a model of the intracellular network (decision-making process) that controls the yeast to hyphal transition (YHT), we focused on the most studied transcription factors that mediate hyphal initiation (Efg1 and Brg1), hyphal maintenance (Ume6), or inhibit these processes (Nrg1). The model includes four environmental cues known to regulate the YHT, namely pH, farnesol, temperature, and serum. These environmental stimuli regulate the activity of signaling pathways (e.g. cAMP/PKA and ESCRT) that inhibit the expression of NRG1, encoding the major YHT transcriptional repressor, and/or induce the expression of EFG1 and BRG1, encoding YHT activators. The transcription factors Efg1 and Brg1, aided by histone modifications due to histone acetyltransferases (HATs) and histone deacetylases (HDACs) activate the transcription of downstream HAGs important for the hyphal morphology (e.g., UME6, HGC1, HWP1, ALS3, and ECE1) [8]. These HAGs encode other transcription factors such as Ume6 that mediate hyphal maintenance, cyclins such as Hgc1 that determine polarized growth at the hyphal tips, and cell wall proteins such as Hwp1, Als3, and Rbt5 that are important for adhesion [7].
The network underlying the model contains several types of nodes, including environmental signals, mRNAs, proteins (signaling proteins, transcription factors, and epigenetic modulators), and processes. These node types are indicated with different symbol shapes and colors in Fig 1. Nrg1 and Efg1 are divided into multiple forms. For Nrg1, the model separately includes the NRG1 mRNA transcript (NRG1_T), as well as the Nrg1 protein bound to the promoter regions of hyphal-associated genes (Nrg1@HAGs). For Efg1, the model separately includes the EFG1 transcript (EFG1_T), the Efg1 protein (Efg1), and the Efg1 protein activated as a result of signal transduction (Efg1_active). The latter allows us to encode the negative feedback that active Efg1 has on the transcription of its own gene [9]. We also include three nodes that describe processes: hyphal_initiation, HAG_transcription, and hyphal_maintenance. Activation of the node hyphal_initiation indicates that external signals have impinged on the YHT core network, suppressing yeast-associated nodes, and beginning transcription of hyphal genes. Activation of HAG_transcription indicates that HAGs (e.g., HGC1 and HWP1) are transcribed. Activation of hyphal_maintenance indicates that the cell has entered a state of sustained hyphal growth and elongation as part of a multicellular hypha [7].
To describe the propagation of information in the network from external signals to the ultimate phenotypic output, we formulated a Boolean model. In a Boolean model, each node can be either ON or OFF (1 or 0), and the state of the whole system is given by the state of each node in the network. In general, ON should be interpreted as present, expressed, or active, while OFF indicates the opposite. Special cases include the signal nodes pH and Temperature. The OFF state of the node "pH" indicates an acidic environment (pH < 6), while its ON state indicates an alkaline or neutral environment. Temperature = 0 indicates an environment cooler than 37˚C, while Temperature = 1 indicates an environment at 37˚C.
The regulatory interactions between nodes are given by Boolean regulatory functions describing what the next state of the target node will become based on the current state of its regulators. The regulatory functions were determined from the literature. In cases where detailed knowledge was not available, we generally assumed inhibitory dominant regulatory functions. This means that multiple activating edges are related by the "or" operator, while inhibitory edges are related by the "and not" operator. Specific regulatory functions and evidence in the literature for these functions are summarized in Table 1, and additional notes for each function are provided in S1 Text.
Within the scope of the model, serum and temperature have identical downstream effects. For the sake of simplicity, we merge these two environmental signals into a single node, and refer to this node as "Temperature". Experimental evidence suggests that high temperature The shapes and colors of the nodes indicate their function, as indicated in the key on the right. Dashed edges represent functional relationships whose molecular mechanisms have not been determined. We translate this network into a Boolean dynamic model by characterizing each node with a regulatory function (see Table 1). https://doi.org/10.1371/journal.pcbi.1008690.g001

PLOS COMPUTATIONAL BIOLOGY
and serum are both required to achieve sustained hyphal growth [10]. When comparing the model's results with experimental findings, we equate the ON state of the input "Temperature" with 37˚C and the presence of serum in the medium.

The model recapitulates the biological phenotypes and the trajectory of the YHT
We describe the dynamics of the YHT model using two distinct methods: general asynchronous update and stochastic propensity. With both methods, the system evolves until it reaches a stationary state or a group of states that it oscillates within. The term for such final states is "attractor". Attractors represent stable biological differentiation states or phenotypes. See Methods for more details about the update schemes and attractors.
When considering all eight combinations of the three input signals, the YHT network model has 27 attractors (13 if the value of the input signals is not considered), which we broadly categorized as one of four phenotypes: yeast, yeast-like, hyphal-like, and hyphal ( Fig  2). Phenotype identification was based on the values of the hyphal_initiation, hyphal_maintenance, and HAG_transcription nodes in the model, as well as on the expression of key transcription factors. One group of attractors corresponds clearly to the yeast state based on the activity of the YHT inhibitor Nrg1 (expressed as the ON state of NRG1_T and Nrg1@HAGs), the inactivity of Brg1 and Ume6, and the lack of hyphal initiation, hyphal maintenance and HAG transcription. The three attractors in this group, marked in blue in Fig 2, only differ in the state of the input nodes Temperature and Farnesol, while the pH must be 0 (acidic environment). We therefore named this group of attractors "yeast". Another group of attractors (marked by yellow color in Fig 2) has hyphal_initiation = hyphal_maintenance = HAG_transcription = 1. The eight attractors in this group share the activation of hyphal-associated Table 1. Boolean regulatory functions of each node in the model. Each function indicates the next state of the node as a function of the current state of its regulators. For simplicity the state of each node is represented with the name of the node. The model is provided in BooleanNet (S1 File) and SBML qual (S2 File) file formats.

PLOS COMPUTATIONAL BIOLOGY
The group of attractors marked in orange in Fig 2 fail to activate hyphal_initiation and hyphal_maintenance, yet they exhibit expression of BRG1 and UME6, as well as active HAG transcription. Due to the presence of these hyphal characteristics we named this group of attractors "hyphal-like". These attractors may describe a pseudohyphal phenotype. Unlike hyphae, which only grow at the tip, any cell within pseudohyphae can divide and branch, but, unlike hyphae, the daughter cells remain attached to the mother cells. The formation of pseudohyphae likely involves the same transcriptional core as the YHT, and involves transcription of a subset of HAGs [4]. These features are recapitulated by the hyphal-like attractors of our model.
Depending on the values of the external signals, we found up to four coexisting attractors. For example, when pH = Farnesol = Temperature = 0, there are four possible attractors, which belong to the yeast, yeast-like, hyphal-like and hyphal attractor groups, respectively. In general, the yeast-like, hyphal-like, and hyphal attractor groups contain a stable attractor regardless of the external signals while the stability of the yeast attractor requires pH = 0 and either Farnesol = 1 or Temperature = 0.
To simulate the YHT we started in a yeast attractor with pH = Temperature = Farnesol = 0, then set pH = 1. We observed trajectories that converged to any of the three other attractor groups: yeast-like, hyphal-like, and hyphal ( Fig 3A-3C). A prominent trajectory of our model reproduces the known features of complete YHT in response to alkaline pH: upregulation of BRG1, hyphal initiation, HAG transcription, hyphal maintenance ( Fig 3C). We verified that setting Temperature = 1 could also induce the YHT, in agreement with the observation that 37˚C induces the YHT [7]. We then undertook a systematic analysis of the outcomes of simulations for every environmental setting using two update schemes (general asynchronous or stochastic propensity). Table 2 indicates the probability of converging into each of the four phenotypes (attractor groups) for every environmental setting when starting from an arbitrary initial state that does not already have hyphal_initiation = 1 or hyphal_maintenance = 1 or from a yeast attractor. In both update schemes, when pH = 0 and either Farnesol = 1 or Temperature = 0 (top row of each table panel) only the yeast and hyphal-like phenotypes are reachable from an arbitrary state, and a system that starts in a yeast state stays in that state. Indeed, yeast is the dominant growth form of C. albicans wildtype strains in an acidic environment with temperature lower than 37˚C [4]; in the following, we will refer to this environmental condition as a yeast-favoring condition. When either pH = 1, or Farnesol = 0 and Temperature = 1, the yeast attractor is no longer stable, and the system converges into the hyphal-like phenotype, the hyphal phenotype, or to the yeast-like phenotype.

Stable motif analysis reveals decision points for successful and failed YHT
To understand how the system makes decisions to evolve toward a specific attractor, we performed stable motif analysis [20] on the YHT model (Fig 4). Stable motifs represent subsets of the Boolean system that, once they achieve a certain state, become locked in that state [21]. They are thus the building blocks of attractors. Expressed more precisely, the fixed state of the nodes of each stable motif determines a trap space, i.e. a region of state space that once entered cannot be exited [22]. Notably, these trap spaces are independent of the update schedule. The YHT network has a single stable motif, which does not depend on environmental conditions. This stable motif consists of hyphal_initiation = 1, which expresses the irreversible nature of hyphal initiation. In addition, there are six conditionally stable motifs. Conditionally stable motifs are only stable motifs if some external condition is met, such as a fixed state of an environmental source node or the stabilization of a parent stable motif [23]. Particularly, the conditionally stable motif outlined in blue involves the activation of the main hyphal inhibitor Nrg1 and the inactivation of the hyphal activators Brg1, Ume6 and HDACs (Fig 4). As can be seen in Fig 1, there is a mutual inhibitory relationship between NRG_T and NRG1@HAGs, on one hand, and Brg1, Ume6 and HDACs, on the other hand. This blue conditionally stable motif expresses one of the two possible states of that mutual inhibitory relationship, and is conditioned on the OFF state of both ESCRT and cAMP/PKA, which is true for the three yeastfavoring environmental conditions described by (pH = 0) AND (Farnesol = 1 OR Temperature = 0). The conditionally stable motif hyphal_initiation = 0 can also lock in under the same however, the cell fails to lock in hyphal_maintenance, and the yeast program is reestablished. This trajectory ends in the yeast-like attractor. (C) A yeast cell completes the YHT, and ends in the hyphal attractor. (D-E) Motifs from the network controlling these windows of opportunity. White nodes begin as ON at the start of the window of opportunity, while black nodes begin as OFF. (D) Once Brg1 is activated, hyphal_initiation must be activated before HATs are silenced to continue the YHT. If HDACs silence HATs first, the network cannot complete hyphal_initiation, and instead reaches a hyphal-like (pseudohyphal) phenotype. (E) NRG1_T must be silenced to begin hyphal_initiation, but can turn back ON once hyphal_initiation has started. If Nrg1@HAGs reactivates before Brg1 activates HDACs, the system reverts to a yeast-like phenotype. This window may be skipped entirely if HDACs have activated prior to hyphal initiation. https://doi.org/10.1371/journal.pcbi.1008690.g003

PLOS COMPUTATIONAL BIOLOGY
set of environmental conditions. In a yeast-favoring environment, if the system starts from an initial condition in which hyphal_initiation = 0 (which is the typical case), then this value is stable, and the yeast-like and hyphal phenotypes will be unreachable (this finding is also reflected in Table 2).
A conditionally stable motif expressing the activity of Brg1 and HDACs and the inactivity of the hyphal inhibitor Nrg1@HAGs has two variants. The first variant, shown in brown outline, is a stable motif in the environmental conditions (pH = 0) and (Farnesol = 1 or Temperature = 0). The second variant, shown with a brown background, also includes the inactivity of HATs. This variant is a stable motif in the environmental conditions (pH = 1) or (Farnesol = 0 and Temperature = 1). The conditionally stable motif outlined in green is a subset of the blue conditionally stable motif and is conditioned on hyphal_initiation = 1. The conditionally stable motif outlined in pink overlaps with the brown stable motif, and is conditioned on hyphal_initiation = 0.
The sequence of which subsequent stable motifs may lock in after a given stable motif locks in is shown in the stable motif succession diagram [21] (Fig 4). The succession diagram confirms the simulation results that the yeast attractor is only reachable when pH = 0, and either Farnesol = 1 or Temperature = 0. The trajectory toward the yeast attractor involves the stabilization of the blue conditionally stable motif and the hyphal_initiation = 0 conditionally stable motif. These motifs are independent of each other and thus could activate in either order in an arbitrary trajectory; this is indicated by the bidirectional arrow in Fig 4. The other three attractors are also reachable in these environmental conditions through the successive lock-in of two conditionally stable motifs. The reason that the simulations reported in Table 2 only reach the yeast or hyphal-like attractors is that those simulations have hyphal_initiation = 0 in the initial condition, which is a stable motif under these environmental conditions, and thus is immediately locked in, restricting the allowed successions.
In the hyphal-inducing environmental conditions (pH = 1 or Farnesol = 0 and Temperature = 1) the locking in of the brown conditionally stable motif can be paired with either state of hyphal_initiation. According to the regulatory function for hyphal_initiation (Table 1), locking in the brown stable motif takes away the possibility of hyphal_initiation to turn ON if it was initially OFF. If hyphal_initiation turns on prior to the locking-in of the brown stable motif, the system converges into the hyphal phenotype. As soon as the brown stable motif locks in, hyphal initiation is prevented from turning ON, and thus it will lock in the OFF state; the system will converge to the hyphal-like attractor. The locking-in of hyphal_initiation = 1 can be followed up by the brown stable motif or the green conditionally stable motif, which expresses the state opposite of the brown stable motif. The first succession leads to the hyphal Table 2. The probability of converging to each attractor from an arbitrary initial state with hyphal_initiation = 0 and hyphal_maintenance = 0, or from a yeast state. While three attractor groups are stable in any environment, their reachability depends on the environment and on the initial state. For example, in a yeast-favoring environment (first row), if the system starts in a yeast state, it will remain in that state. The results are qualitatively the same whether general asynchronous update (top) or stochastic propensity update (bottom) is used.

PLOS COMPUTATIONAL BIOLOGY
yeast-favoring environmental conditions in that the EFG1 transcript, Efg1 protein and HATs oscillate. Efg1 oscillations in C. albicans have previously been proposed as a mechanism of maintaining cell to cell variability of Efg1 within a population [19].

Stable motif decision points correspond to YHT windows of opportunity
In the stable motif succession diagram, when there are multiple edges emerging from a single motif, the system taking one of these edges may represent an irreversible commitment to one of two mutually exclusive trajectories. Fig 4 shows two branch points that determine commitment relevant to the YHT: the edges emerging from the solid brown conditionally stable motif to either the tan or orange motifs indicate commitment to either the hyphal or hyphal-like attractor groups. Similarly, the edges emerging from the tan stable motif to either the solid brown or green motif dictate commitment to the hyphal or yeast-like attractor groups. The choice of one path versus another depends on the timing of specific events.
The first branch point depends on a sequence of events starting when Brg1 turns ON. While Brg1 remains ON, deactivation of HATs will follow via the inhibitory path Brg1 ! HDACs ─� HATs. However, activity of HATs is a requirement for hyphal_initiation to turn ON. If hyphal_initiation activates before the node HATs turns OFF, then the system will follow the path toward the hyphal attractor group. Conversely, if HATs turns OFF before hypha-l_initiation activates, the system will proceed toward the hyphal-like attractor group. For example, in the trajectory in Fig 3A, HATs turn OFF before hyphal_initiation turns ON, causing the system to proceed to the hyphal-like phenotype. This corresponds to the small incoherent feedforward loop illustrated in Fig 3D. The second branch point depends on the timing of events following the activation of hypha-l_initiation. If the brown motif (which contains Nrg1@HAGs = 0) locks in, or has locked in prior to the activation of hyphal_initiation, then the system will proceed to the hyphal attractor group. Yet, NRG1_T returns once hyphal_initiation activates, and can lead to activation of Nrg1@HAGs, which is part of the green motif. Nrg1@HAGs activity is sufficient to lock in the green motif, which expresses the deactivation of the core hyphal program, leading to the yeastlike attractor group. Thus the YHT depends on a race to exclude Nrg1 from the promoter region of HAGs following the reactivation of NRG1 transcription. This race corresponds to multiple negative regulatory pathways between hyphal_initiation and hyphal_maintenance ( Fig 3E), mediated through Nrg1@HAGs.
These decision points, and their corresponding races, reflect the documented C. albicans YHT "window of opportunity" [7,10]. This is a transient period, beginning with the downregulation of NRG1 and ending with the subsequent re-expression of NRG1, in which the hyphal program may be established. If it does not establish prior to the window closing, the cells do not complete the YHT. As described above, our model reproduces this behavior, resolving the window of opportunity into two distinct decision points (one regarding hyphal initiation and the other regarding hyphal maintenance), and describing the specific mechanisms by which the window can be missed.

Network control predictions
Using the Boolean YHT model, we sought to identify interventions, such as controlling the state of one or more nodes or deleting or activating an edge, that could prevent the YHT. Interventions were identified using several different control strategies and objectives, summarized in Table 3. We have applied feedback vertex set (FVS) control, stable motif control, and algebraic edge control with canalizing function analysis to our network model, as well as systematic simulations of perturbations.

Feedback vertex set control
A feedback vertex set (FVS) is a collection of nodes in a network whose removal results in a network with no cycles (no feedback loops). On a network with no feedback loops, dynamical processes described by Boolean or differential equation models have a single attractor [24,26]. FVS control thus predicts that by fixing all nodes in a given FVS, as well as all source nodes, to match a particular attractor, one can force the system from any state into that attractor [27].
Once the system achieves that target attractor, control of the FVS nodes may be relaxed, though control of the source nodes must be maintained. Unlike the other control methods, FVS only requires knowledge of the network's topology (Fig 1), that is, the collection of nodes and edges, as well as knowledge of the attractors, but it otherwise requires no specific details of the regulatory functions.
The YHT network contains a strongly connected component (feedback-rich subgraph) of 10 nodes. The FVS of the network consists of three nodes: Nrg1@HAGs, at least one node from the Efg1 feedback loop, and hyphal_initiation (which has a self-loop). FVS control predicts that a control strategy for ensuring that the system converges into the yeast attractor is to maintain Nrg1@HAGs = 1, Efg1 = 1, and hyphal_initiation = 0 to eliminate feedback sets, as well as ensure pH = 0 and either Farnesol = 1 or Temperature = 0 to control source nodes (see panel A of S1 Fig). Conversely, FVS control into the hyphal attractor group requires setting Nrg1@HAGs = 0, Efg1 = 1, and hyphal_initiation = 1 (panel B of S1 Fig). As this attractor group is reachable in any environmental condition, the source nodes do not need to be controlled.
FVS control provides a sufficient condition to maintain a given attractor. Nevertheless, it may be that a subset of nodes can still accomplish the control objective. This is especially important to identify here, as fixing nodes such as hyphal_initiation may have no obvious biological implementation.

Stable motif based control
Stable motif control seeks to determine a sequence of driver nodes that, if transiently maintained in a fixed state, will lock in a sequence of stable motifs that will force the system into a desired attractor from any initial condition [20,28]. A variation of stable motif control also identifies sequential control of driver nodes that drives the system into one of multiple target attractors [29].
Fig 5 shows driver sequences needed to drive the YHT network into an attractor with hyphal_maintenance = HAG_transcription = 0 (corresponding to one of the yeast or yeastlike attractor groups), or hyphal_maintenance = 1 (corresponding to a hyphal phenotype). Unlike FVS control, stable motif control is able to force the system to specific attractor groups Table 3. An overview of the control prediction approaches we apply to the C. albicans YHT network, their objectives, and the types of interventions they require.

Control Method
Objective Interventions by only controlling two nodes, whose identities depend on the control objective and environmental conditions. This is because the trajectory to each attractor contains a succession of two motifs (Fig 4), and each multi-node stable motif of the YHT network can be locked in by fixing the state of one node. For example, temporarily controlling hyphal_initiation = 1 (after which it locks in), followed by temporarily controlling Brg1 = 0 (until the green stable motif in

PLOS COMPUTATIONAL BIOLOGY
locks in) is sufficient to achieve hyphal_maintenance = HAG_transcription = 0 in any environment. Other alternatives for ensuring the lock-in of the green stable motif following hyphal initiation are setting HDACs = 0 or Nrg1@HAGs = 1. In all three of these cases the system will follow a trajectory toward a yeast-like attractor. Driving the system to a hyphal phenotype in any environment and from any initial condition requires temporarily holding hyphal_initiation = 1 (after which it locks in), followed by holding any node value of the brown-outlined motif (NRG1@HAGs = 0, HDACs = 1, or Brg1 = 1). The stable motif based control sets still have the shortcoming of involving direct control of hyphal initiation, which is difficult to implement experimentally.

Simulated perturbations
We systematically simulated permanent deletions (holding in the state 0) and activations (holding in the state 1) of the nine nodes that are not signals, signaling intermediaries, or phenotypic outcomes using both general asynchronous and stochastic propensity updating schemes. To identify interventions that block YHT, we began the simulations in a yeast state placed into an environment with Farnesol = Temperature = 0, pH = 1. As indicated in Table 2, the unperturbed system undergoes the YHT in about 37% of trajectories for general asynchronous update (53% for stochastic propensity update), otherwise missing one of the two windows of opportunity. We found seven node interventions (when grouping all Efg1 knockouts, and Nrg1 knockouts together) that ensure none of the trajectories converge to a hyphal phenotype (Table 4). Among these, the interventions Brg1 = 0, HDACs = 0, Ume6 = 0, or Nrg1@HAGs = 1 led to the complete elimination of the hyphal attractor. Indeed, it was found experimentally that deletion of BRG1 [30], HDA1 [10], or UME6 [31] led to impairment of the YHT. These four states are incompatible with both variants of the brown motif, whose locking-in is necessary for the hyphal attractor. Instead, the trajectories starting from the yeast attractor converge into either a hyphal-like or yeast-like attractor (in case of Ume6 = 0) or solely to a yeast-like attractor (in the other three cases). In the remaining four cases the hyphal attractor exists but it is not reachable from an initial condition corresponding to yeast. In the case of deletion of EFG1, the system stays in the yeast attractor, and for HATs = 0 it converges into a hyphal-like attractor. Indeed, it was found experimentally that deletion of EFG1 prevented hyphal formation [32] and deletion Table 4. Probabilities of completing the YHT in the condition pH = 1 and Farnesol = Temperature = 0, under single-node interventions or WT. Transition probabilities are estimated from 500 simulated trajectories starting from the yeast attractor corresponding to the given environment (see Fig 2). When a perturbed system's attractor is different from the WT attractor, it is classified into a phenotype (attractor group) by similar criteria as in Fig 2; [17].
In the case of simulated constitutive activation of Ume6, the YHT propensity decreased compared to wildtype. This happens in the model because Ume6 inhibits Nrg1@HAGs, decreasing the YHT window of opportunity. In contrast, simulated constitutive activation of the EFG1 transcript or Efg1 protein increased the YHT propensity. This result is consistent with the hyphal morphologies observed when EFG1 is overexpressed [33].
Unlike the previous control methods, which by design force the system into a pre-existing attractor, this method can introduce new attractors. For example, a simulated deletion of BRG1 prevented the YHT, consistent with experimental observations of defective hyphal elongation in mutants lacking BRG1 [30], and introduced a new yeast-like attractor in which NRG1_T and Nrg1@HAGs oscillate. All attractors for each single-node perturbation are shown in S3 File, and the new attractors observed for permanent deletions and activations of single nodes are indicated in S2 Fig. We also identified four perturbations that eliminate the yeast attractor even in yeast-favoring environmental conditions, namely Brg1 = 1, HDACs = 1, NRG1_T = 0, and Nrg1@HAGs = 0. These states are incompatible with the blue motif (Fig 4), whose locking in is necessary for the yeast attractor. Any of these node states can ensure the locking in of both versions of the brown conditionally stable motif in Fig 4. When this motif locks in, the current state of hyphal_initiation determines whether the system converges to the hyphal-like or hyphal attractor groups. Hyphal initiation following a change in environment requires that hyphal-inducing signals propagate through two parallel pathways: activating the brown motif (via NRG1_T downregulation), and activating HATs (via Efg1). Due to the proximity of the controls to the brown motif, in all simulations the brown motif locked in before HATs turned on, causing the system to converge to a hyphal-like attractor. Thus, our model predicts that these perturbations would lead to filamentation (likely pseudohyphae) even in yeast-favoring environmental conditions. Indeed, experiments indicate that deletion of NRG1 leads to pseudohyphae [34,35]. In contrast to the model prediction, C. albicans cells engineered to ectopically express BRG1 under yeast-favoring conditions stayed in a yeast phenotype [30]. One potential explanation of this discrepancy is that Brg1 translated from ectopically expressed BRG1 may not be able to recruit HDACs to the promoter region of HAGs to exclude Nrg1@HAGs; indeed the same study found that in strains ectopically expressing BRG1, the Brg1 protein could not bind to the promoter region of HAGs.

Algebraic control with edge knockouts
Algebraic control [36] uses the polynomial form of the Boolean functions. Two types of control objectives can be formulated: node control and edge control. The identification of control

PLOS COMPUTATIONAL BIOLOGY
targets is achieved by encoding the nodes of interest as control variables within the functions; edges of interest are encoded as control variables within the inputs of the functions. Then, the control objective is expressed as a system of polynomial equations that is solved by computational algebra techniques. For node control, we considered the environmental condition pH = Rim8 = 1, Temperature = Farnesol = Cyr1 = cAMP_PKA = 0 and set the objective of finding node knockouts or constitutive activations for which there is an attractor of the system that has the hyphal_maintenance and HAG_T nodes OFF. Thus, we set our objective to find controls such that Fðx;ûÞ ¼x; x k ¼ 0; k 2 I, where the index set I corresponds to the indexes of hyphal_maintenance and HAG_T. We found the following node controls: Brg1 = 0, HDACs = 0, and NRG1@HAGs = 1. These interventions were also identified by simulations to block the YHT. For edge control, we set the control objective of destroying (or blocking) the fixed point x 0 corresponding to the hyphal state. Thus, our goal is to find controls u such that F(x 0 , u) 6 ¼ x 0 . Using the edge control approach, we identified the following edge deletions and activations, shown in Table 5, that are effective for blocking transition to the hyphal state.
For each edge control, we calculated the total effect (the total change in the state transition graph) corresponding to stochastic propensity update, as described in the Methods. Interventions with larger total effect induce greater systemic changes in the state transition graph of the system, and therefore may be associated with more side effects [25]. We likewise calculated the system's time to absorption (TTA), which corresponds to how long it takes the system to reach an attractor under stochastic propensity update, when starting from the yeast state and setting pH = 1. Controls that have a lower TTA indicate that the system will quickly reach the attractor. In the case when TTA = 0, the perturbed system has a stable yeast attractor even when pH = 1. This is the case for constitutive activation of the inhibitory edge Nrg1@HAGs ─� Brg1.
The edge interventions with smallest total effect (6.25%) are constitutive activation of the HATs ! Nrg1@HAGs edge, or constitutive deletion of the HDACs ─� Nrg1@HAGs edge. The effect of these interventions ensures that if Nrg1 is expressed, it will bind to HAGs, thus under this intervention the system misses the second window of opportunity (Fig 3E). With this intervention, the system converges preferentially into the YL attractor, or alternatively the HL attractor, with a slightly longer convergence time than the WT system would converge to a hyphal attractor.

Experimental verification of predicted control interventions
We performed two new experiments to test two model-predicted interventions that eliminate the YHT: deletion of UME6 and deletion of the HDAC HDA1. The model predicts both of these interventions make hyphal maintenance impossible, though it predicts differences in the final attractors. Setting Ume6 = 0 leads to a three-attractor repertoire: yeast, yeast-like, and a novel attractor group named hyphal-like 2 (see S2 Fig). The attractors of the Ume6 = 0 system do not achieve hyphal maintenance but do exhibit one or both of the other hyphal-associated phenotypic outcomes. Setting HDACs = 0 in the model leads to convergence to a yeast-like attractor, which is a much stronger departure from the unperturbed system's outcome. We determined the morphology of cells of C. albicans strains lacking UME6 and HDA1 at 90 minutes after inoculation. By this time, cells of the wildtype C. albicans strain have completed hyphal initiation and are in the elongation phase [7]. The experimental condition was neutral (pH = 7.0) RPMI-1640 medium at 37˚C; these mutant strains were not tested under these conditions in any prior publications we found. As shown in Fig 6A, compared to the wildtype strain at 90 minutes (0% yeast-form, 0% transitional-form, 100% true hyphal-form), the ume6 Δ/Δ strain displayed a moderate filamentation defect (4% yeast-form, 29% transitional-form, 67% true hyphal-form), while the hda1 Δ/Δ strain displayed a severe filamentation defect (60% yeast-form, 40% transitional-form, 0% true hyphal-form). Quantification of the yeast-wall protein (YWP1) and hyphal-wall protein (HWP1) transcript levels (markers of the yeast and hyphal phenotype, respectively) at 90 minutes were consistent with WT cells expressing HWP1, and the hda1 Δ/Δ strain expressing YWP1 (Fig 6B). The ume6 Δ/Δ cells expressed more YWP1 and less HWP1 than WT cells. Further, at 180 minutes, many ume6 Δ/Δ cells show constriction at the separation sites, which is a hallmark of pseudohyphal cells [4]. These results are consistent with the model predictions, especially in terms of the relative severity of the two defects, and provide evidence that the HL attractor may correspond to a pseudohyphal state.

Model validation by comparing to the published literature
As an additional verification of the model, we compared model-predicted outcomes of simulated controls with published results of corresponding experimental interventions. As shown in S1 Table, the model and experiments agree in 17 out of 22 cases. The five cases of disagreement pertain to NRG1 deletion and BRG1 deletion under filamentation-inducing conditions, Efg1 constitutive activation under conditions favoring the yeast-form, and Brg1 constitutive activation.

PLOS COMPUTATIONAL BIOLOGY
The model predicts that deletion of NRG1 leads predominately to a hyphal-like attractor in every environment, while experimentally it was found that in hyphal-inducing conditions, deletion of NRG1 led to hyphal formation [34]. In the model, deletion of NRG1 often leads to early activation of HDACs, which precludes the temporary activation of HATs, which in turn is necessary for hyphal initiation (Fig 3D). This discrepancy could be mitigated by revising certain assumptions of the model. Activation of Brg1 may require the action of a specific activator in addition to Nrg1@HAGs inactivity; indeed Brg1 is documented to be part of a network of six cross-activating transcription factors [32], which also includes Efg1. We tested a modified regulatory function, Brg1 � = not Nrg1@HAGs and Efg1_active; however this change introduced a discrepancy in capturing the phenotype of NRG1 KO in a yeast-favoring environment. Further analysis of the functional requirement for the network in [32] and its relationship with Nrg1 may help to elucidate the regulatory function of Brg1. Another possibility is that the requirement for HATs activity for hyphal initiation may be less strong than assumed in the model.
Constitutive activation of Efg1 (simulated by maintaining Efg1_active = 1) in a yeast-favoring environment leads to an attractor that exhibits a low level of Efg1_T and Efg1 but otherwise is the same as the unperturbed system's yeast attractor. In contrast, a pseudohyphal phenotype was observed experimentally [33]. The reason the model does not recapitulate this result is the assumption that Nrg1_T downregulation requires active Efg1 in collaboration with cAMP/PKA or pH signaling. The discrepancy would be resolved by the alternative assumption that the active Efg1 is the mediator of the effect of the environment on Nrg1_T. Indeed, a model version that omits the direct effects of the environment on Nrg1_T (by omitting "or not (ESCRT or cAMP_PKA" from its update function) indicates a mixture of phenotypes, 18% YL, 33% HL, and 49% H, in a yeast-favoring environment. This model version does not introduce any significant deviations in results or any discrepancies with experimental results.
The simulated deletion of Brg1 leads to the system staying in a yeast state in hyphal-inducing conditions, while experimentally it was found that a BRG1 deletion strain exhibited competent germ tube formation and defective hyphal elongation. The discrepancy stems from the assumption that Brg1 activity is required for hyphal initiation. It is possible that other transcription factor(s) within the six transcription factor network mentioned above [32] could rescue hyphal initiation in the absence of Brg1. We tested adding a connection from Efg1_active into hyphal_initiation, replacing Brg1 with (Brg1 or Efg1_active). In this case, in hyphalinducing conditions, the cells went to the YL attractor 98% of the time. This is because signalinduced Efg1 activation starts the second window of opportunity (Fig 3E) earlier, before Brg1 is in position to activate HDACs and Ume6. Thus, we learned that this putative transcription factor is not Efg1_active, and it needs to activate hyphal_initiation at a similar time as Brg1.
Experimentally overexpressing BRG1 in a yeast-favoring environment caused the cells to remain in the yeast phenotype, while in a YHT-inducing environment led to sustained hyphal growth [30]. In contrast, constitutive activation of Brg1 in the model leads to a transition from the yeast attractor to a hyphal-like attractor nearly 100% of the time, regardless of environment. The cause of this environment-independent outcome in the model is that ectopically induced Brg1 activates HDACs, which has two environment-independent effects: it prevents HATs from activating, and it excludes Nrg1 from the promoter region of HAGs. The first of these effects obstructs hyphal initiation, while the second enables Brg1's activation of HAG transcription in the model. To restore agreement with experiments, ectopically induced Brg1 should not be able to activate HAG transcription in a yeast-favoring environment, and it should not prevent hyphal initiation in a YHT-inducing environment. Both discrepancies might be mitigated by a condition for the activation of HDACs that requires more than the presence of Brg1, however candidates for this option are not immediately apparent.

Discussion
Here we have developed a new Boolean model of the C. albicans yeast to hyphal transition (YHT), and demonstrated that it recapitulates several known behaviors. Our model has attractors corresponding to the yeast and hyphal phenotypes, as well as two types of attractors that correspond to different ways the system can fail to complete the YHT, by failing to pass socalled "windows of opportunity" (Fig 2). The hyphal-like attractors exhibit expression of BRG1, UME6 and hyphal-associated genes, but fail to activate hyphal initiation. Our model predicts that the YHT can arrest if histone acetyltransferases (HATs) inactivate prior to the transcription of hyphal-associated genes. This state is a point attractor in our model, which may correspond to a pseudohyphal phenotype, or it may not be a stable biological state. Further experimental investigation of this system may support or contradict the stability of this state. It may also reveal the conditions under which the transition can resume or alternatively the system can reset to a yeast state. The yeast-like attractors exhibit active hyphal initiation, but have active Nrg1 and inactive Brg1 and Ume6, which are characteristic of yeast cells. We interpret these states as the YHT arrested after hyphal initiation, followed by a resetting into the yeast-like state of the transcriptional regulators.
Through stable motif analysis, we showed that the previously described YHT window of opportunity corresponds to two branch points in stable motif succession (Fig 4). The choice of which branch the system takes depends on the specific timing of events within two small subnetworks we identified as the window of opportunity motifs (Fig 3). Both of these subnetworks contain an incoherent feed-forward loop, which is the coexistence of two short paths of opposite signs between a pair of nodes. The timing of the events in each path determines the outcome. In the first window of opportunity (Fig 4D), hyphal_initiation must be activated before HATs are silenced by HDACs to continue the YHT. In the second window of opportunity ( Fig  4E), HDACs must be activated by the time the Nrg1_T expression starts to increase again, to avoid the reactivation of Nrg1@HAGs, otherwise the YHT cannot be sustained. Further experimental investigation of these epigenetic regulatory processes will be able to determine their timing and regulation. For example, one may elucidate the timing requirement of the HAT activity by engineering a deletion of YNG2 (encoding the active subunit of the HAT NuA4) or by including an inducible promoter for a constitutively acetylated (and thus non-responsive to HDACs) mutant (K175Q) YNG2 in a WT YNG2 strain.
We probed the dynamics of the YHT using two stochastic update methods with parsimonious assumptions. General asynchronous update assumes that only one node can change state at any given moment, and each node has the same probability of being updated. Stochastic propensity update with propensity of 0.5 assumes that each node has an equal probability of following its regulatory function or maintaining its current state. Notably, in this update method multiple nodes can change state simultaneously. The emerging theoretical understanding of Boolean dynamical systems indicates that both perfect synchrony and complete lack of synchrony can limit possible state transitions in a way that leads to artificial oscillations [29,37]. These limitations are avoided in the recently proposed Most Permissive Boolean Networks [38]. We found very similar results for the two update methods: the phenotypes are the same and the probabilities to converge into each phenotype are similar (see Tables 2 and 4). More specifically, the point attractors are the same and the complex attractor possible in five input combinations resides in the same trap space (i.e., in the same branch of the stable motif succession diagram). The agreement of the two update methods on the complex attractor is consistent with the fact that this complex attractor originates in a negative feedback loop. The agreement of the two update methods on predicted control strategies increases the confidence on the applicability of these control strategies for the biological system. We applied four methods to predict control strategies that might block the YHT in C. albicans. Each method searches for different types of control strategies, with slightly different control objectives (Table 3). Based on a previous analysis of FVS and SM control [39], we expected these methods to give similar control strategies, since they both have the same control objective (forcing the system into one of its natural attractors) and are based on the control of feedback loops in the network. Based on [40], we expected the node control simulations to give interventions not found by the other methods, since the simulations use as an initial state-the yeast state-while the other control methods do not have a way to use a given initial state to make their predictions. It was not clear how much commonality algebraic edge control would have with the other node control methods, since algebraic edge control distinguishes itself by its ability to (1) disrupt individual edges and (2) destroy the target attractors. Nevertheless, some common key driver nodes and interactions emerged across the different control strategies. Nrg1@HAGs, for example, was identified as a key node in all 4 control methods. In FVS control, all but two feedback loops were broken by controlling Nrg1@HAGs (S1 Fig). Nrg1@HAGs was a participant in almost all predicted stable motif control sequences (Fig 5). Further, simulations found that constitutive activity of Nrg1@HAGs could prevent the YHT, while constitutive inactivity of Nrg1@HAGs eliminated the yeast attractor. Lastly, more than half of the algebraic edge control predictions involved constitutive activity or silencing of edges involving Nrg1@HAGs or NRG1_T, including the interventions that led to the lowest total effect (Table 4). These control results agree with the well-known central gatekeeping role that Nrg1 plays in regulating the YHT [7].
Our edge control results predict effective parsimonious interventions that target interactions as opposed to nodes. For example, we predict that intervening in the acetylation properties of the promoter regions of HAGs in such a way to decrease the ability of Nrg1 to bind there would decrease the YHT. Such intervention offers a potentially more practical alternative than genetic deletion of NRG1. We predict that disabling the capacity of Brg1 to recruit HDACs, or disabling the capacity of HDACs to block HATs, would also disable the YHT.
While one may view these control predictions as possible targets for externally applied perturbations, they also reveal much about how the system may control its own repertoire of behaviors in different environmental conditions. Our analysis reveals the importance of positive feedback loops containing an even number of inhibitory edges (i.e., inhibiting an inhibitor); these feedback loops form conditionally stable motifs that can lock-in and restrict the system's repertoire. The other regulatory motif important in this system is the incoherent feed-forward loop, which underlies both windows of opportunity discussed earlier. The most unexpected result concerns the incoherent feedback-mediated role of HATs in the YHT. First, HATs are required for hyphal_initiation [17], whereas they must be degraded for sustained hyphal maintenance [10]. The exact mechanisms mediating HAT activity, and their timing, requires further study.
There remain some experimental observations our model does not recapitulate, revealing gaps in knowledge regarding regulation of the YHT. For example, our current model does not recapitulate certain experimental observations pertaining to Brg1's ability to recruit HDACs and to regulate HAGs, potentially questioning Brg1's primacy among multiple interacting transcription factors that co-regulate HAGs [32]. Future genetic epistasis experiments, for example, that combine deletion of BRG1 with deletion of genes encoding these interacting regulatory partners, may shed new light on the hierarchical roles of these transcription factors in regulating HAGs. The discrepancy regarding the phenotype of Efg1 constitutive activation suggests that Efg1 may be the mediator of the environment's effect on NRG1 transcription. We hope that our model will lead to follow-up experiments that eliminate these gaps of knowledge. Ultimately, given the importance of the yeast to hyphal transition in modulating C. albicans virulence, understanding how cellular decisions are made to undergo this morphological transition will ultimately allow us to better understand how C. albicans causes disease in humans and how this process can be altered to prevent disease.

Update methods
The state of the system at time t is given by an n-dimensional vector x t , where n is the number of nodes (15 in this case). For each node i in the network, there is a corresponding Boolean function f i (x t ) that specifies the regulation of node i. To analyze the network, we used two different types of stochastic update for the Boolean functions: general asynchronous, and stochastic propensity.
In general asynchronous updating, at each time step one node (i) is selected at random and its value is updated as We also used the stochastic propensity framework described by Murrugarra et al [41]. Briefly, for each node i, we consider two parameters p " i ; p # i 2 ½0; 1� called activation and degradation propensities. p " i is the probability that in a situation when f i calls for an update from 0 to 1 the state change happens (the state will remain 0 otherwise, with probability 1 À p " i ) and p # i is the probability of an update from 1 to 0 according to f i . If applying f i has no effect on x i , then x i keeps its current value with probability 1. In summary, In the stochastic propensity framework, each node may have different propensities to change from 0 to 1 (p " i ) versus changing from 1 to 0 (p # i ), and these propensities may be different for different nodes. Unlike general asynchronous update, it is possible for multiple nodes to simultaneously change state in the stochastic propensity framework. Because of these differences, the timescales for simulations using the two frameworks are different. Due to the lack of information on the kinetic rates or timescales of the nodes in this network, we choose all the propensities to have value 0.5.

State transition graph and attractors
Both Boolean update methods described above define a state transition graph (STG). Each possible state of the system corresponds to a node of the STG, and each directed edge indicates a possible next state after a single time step. In the general asynchronous update framework, only a single variable is updated each time step, so only states that differ by the value of a single variable can be connected in the STG. Conversely, in the stochastic propensity framework, multiple variables may simultaneously update within a single time step, and so distant nodes may be connected in the STG. This difference is partly responsible for the different timescales discussed above.
The state transition graph with all nodes and edges completely describes all possible dynamics of the system. Of special interest are the attractors of the system, which are individual states, or collections of multiple states, that have transition edges into them, but no transition edges out of them. In graph theoretical terms attractors correspond to terminal strongly connected components of the STG. Once the system enters an attractor, it cannot leave using the dynamics of the network alone. However, perturbations to the regulatory network, such as fixing node values or deleting edges, can change the underlying STG, possibly changing the set of attractors.
The distance between any two states x,y 2 [0,1] N can be described using the Hamming distance, defined as and which has the property that H(x,y) 2 {0,1,2,. . .,N}.
Given any state x, the state y = (f 1 (x), f 2 (x),. . .,f N (x)) describes the state achieved by applying each regulatory function f i (x) to each node i simultaneously. Let d x = H(x,y) count the number of nodes whose values can change by applying one of the regulatory functions f i (x). Each state x then has d x possible transitions out under general asynchronous update. Specifically, the probability of transitioning between a pair of states x and y = (f 1 (x), f 2 (x),. . .,f N (x)) under general asynchronous update is Probðx i ! y i Þ ¼ 0 whenever y i 6 ¼ f i ðxÞ: All transition edges present in the general asynchronous update STG are also present in the stochastic propensity STG, though they may have different probabilities. Conversely, the stochastic propensity STG may have many edges not present in the general asynchronous STG. Nevertheless, the two frameworks share many common dynamical behaviors. If a state x is a point attractor (i.e., an attractor comprising a single state) of the system under general asynchronous update, then by definition of an attractor, d x = 0. Thus there are 2 0 -1 = 0 transitions from x under stochastic propensity updates, indicating x is an attractor of the stochastic propensity system. Conversely, if x is a point attractor of the system under stochastic propensity updates, then 2 d x À 1 ¼ 0, so d x = 0. Thus, point attractors are preserved between both update orders. Indeed, point attractors are independent of the implementation of time [42].
In contrast, complex attractors are update-dependent [43], thus they are not guaranteed to be preserved between the two update methods. The C. albicans YHT network has five related complex attractors under general asynchronous update (one for each hyphal-inducing environmental condition) (Fig 2). Simulations of these attractors using stochastic propensity reveals that all nodes oscillating in the general asynchronous attractor also oscillate in the stochastic propensity attractor. Thus, for the YHT network, all attractors of each method are also attractors of the other.

Stable motif analysis
Stable motifs correspond to positive feedback loops in the network that, once they achieve a certain state, become locked in. A stable motif succession diagram describes paths that may be taken once a given stable motif has locked in. Edges indicate subsequent stable motifs or conditionally stable motifs that may become locked in after the edge's source stable motif is established. We created a stable motif succession network using the StableMotifs python package [29]. Several branches of the stable-motif network contained redundant information that enabled simplifications. For example, the network has three inputs: pH, temperature, and Farnesol, each of which forms both ON and OFF stable motifs (e.g., pH = 1 is a stable motif, as is pH = 0). However, many of these lead to identical successions. Such redundancies were removed to derive a parsimonious version containing all the information (Fig 4).

The total effect of deleting an edge
We implement the activation or deletion of a directed edge from node s to node t by replacing x s in the regulatory function of node t by 1 or 0, respectively. To define the total effect of deleting an edge for stochastic propensity, we compute the number of changes in the state space before and after an edge deletion. If an update function f t (x) is written in canalizing layers format [44], that is, then the percentage of change from the initial state space upon the deletion of the edge from node s to node t is at most:

Computation of the time to absorption
Consider the transition matrix of the Boolean network (as a Markov chain) in canonical form , where 0 is the zero block matrix and I is the identity submatrix. The fundamental matrix N is defined as the inverse of (I-Q). That is, N = (I − Q) −1 . The time to absorption for a transient state j is defined as the expected number of steps before absorption and can be calculated as the sum of the jth row of N (see Theorem 11.5 of [45]).

Filamentation assay
C. albicans cells were grown at 30˚C on yeast extract peptone dextrose (YPD) agar plates (2% Bacto peptone, 2% dextrose, 1% yeast extract, 2% agar) for two days. Single colonies were picked and inoculated into YPD liquid medium (2% Bacto peptone, 2% dextrose, 1% yeast extract) and grown at 30˚C overnight. Strains were inoculated from the overnight YPD culture into RPMI-1640 medium at pH = 7.0 (with glutamine and phenol red and without bicarbonate, buffered with MOPS) (MP Biomedicals, catalog# 0910601) at an OD600 = 0.2. RPMI-1640 cell cultures were incubated at 37˚C for 90 minutes and imaged by light microscopy. A minimum of 70 cells were counted to quantify the percentage of cells categorized as yeast-form cells, transitional cells (including pseudohyphal cells), and true hyphal cells for each strain.

Gene expression analysis by quantitative real-time PCR (qPCR)
RNA was extracted using the RiboPure Yeast RNA Purification Kit (Thermo Fisher Scientific, catalog# AM1926). cDNA was synthesized using the High-Capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific, catalog# 4368814). Expression levels of YWP1 and HWP1 were measured in the WT, ume6 Δ/Δ and hda1 Δ/Δ strain backgrounds by real-time quantitative PCR (qPCR) with the following primer pairs: CJNO3815 (5' ACAAATGTCAAGAAAC CACCGT 3') and CJNO3816 (5' ATCGCAAGCAACAACAGTGATA 3') for YWP1; CJN O3813 (5' GGTCAAGGTGAAACAGAGGAAG 3') and CJNO3814 (5' AATCACAAGGT TCTTCCTGCTG 3') for HWP1. Expression levels were assessed under the filamentation assay conditions described above. Normalized relative gene expression values were calculated by the ΔΔC T method using ADE2 as a reference gene and normalized with respect to the WT strain. Results are the means of two determinations.
Supporting information S1 Fig. Feedback vertex set control. FVS control strategies to drive the system into a target attractor (group). Bold edges participate in feedback loops that are broken by controlling the values of the FVS. Nodes and edges that are irrelevant to feedback vertex control are shown in light grey, while nodes of the FVS are shown with a blue outline, and the nodes are colored based on the values they require for FVS control. (A) FVS control predicts fixing environmental conditions as pH = 0 and either Farnesol = 1 or Temperature = 0, and then fixing Nrg1@HAGs = Efg1 = 1 and hyphal_initiation = 0 will force the system into the yeast attractor. Instead of Efg1, EFG1_T or Efg1_active would also be suitable targets. (B) FVS control predicts fixing Nrg1@HAGs = 0, and Efg1 = hyphal_initiation = 1 will force the system into a hyphal attractor. (TIF) S2 Fig. Novel attractors following single-node control. Single-node deletions or activations lead to new attractors. All attractors reached under YHT inducing conditions (pH = 1, Temperature = Farnesol = 0) are shown. Phenotype classification was performed as before, except the values of environment source nodes (pH, Farnesol, and Temperature), signaling intermediaries (Rim8, ESCRT, Cyr1, cAMP/PKA), and the individual perturbed node are ignored. One attractor emerged which did not fit the phenotype classifications defined previously, which we here call hyphal-like 2 (HL2). It corresponds to deletion of UME6, and has hyphal_initiation = 1, HAG_transcription = 1, and hyphal_maintenance = 0. (TIF) S1 Table. Model validation by comparison with literature. Compilation of published experimental intervention results and comparison with the relevant model results. The first column describes the intervention. The second column indicates the environmental condition used in the experiments. The composition of the various media is the following: YPD (also denoted YEPD) medium: yeast extract-peptone dextrose, pH = 7 at the beginning of culture (it decreases as the yeast breaks down dextrose); B-medium: 0.67% yeast nitrogen base, 2% Nasuccinate, pH = 6.5; RPMI-1640 supplemented with L-glutamine and buffered with morpholinepropanesulfonic acid (MOPS), pH = 7.0; Spider medium: nutrient broth, mannitol, K2PO4, agar, pH = 7.2. Because all of these experiments start with a dilution of cells into fresh medium, the farnesol level is expected to be very low (equivalent with Farnesol = 0 in the model).