Strains and plasmids
Standard methods for the growth, maintenance, and transformation of yeast and bacteria, and for the manipulation of DNA, were used throughout. The yeast
S. cerevisiae strains used in this study are BY4741 (
MATa leu2Δ
his3Δ
met15Δ
ura3Δ), BY4741-derived deletion mutants lacking
ste7,
fus3,
kss1,
ptp2,
ptp3,
msg5, or
ptp2/ptp3 (
ptp2::
URA3, ptp3::
KanMX), or BY4741 expressing
STE7,
FUS3, or
MSG5 C-terminally fused to a tandem-affinity purification (TAP) tag (Open Biosystems, Huntsville, AL). All strains are
BAR1+ and therefore do not undergo sustained arrest at the pheromone doses used. The filamentous-responsive element (FRE) transcription reporter (Ty1-lacZ) used in this study was described previously (
Maleri et al., 2004 
).
Expression plasmids encoding
STE7 (pNC752) and the feedback phosphorylation-deficient
ste7A7 mutant (pNC769) were described previously (
Maleri et al., 2004 
). Additional expression plasmids used in this study are those containing
FUS3,
KSS1,
PTP2,
PTP3, and
MSG5. Each gene was amplified using flanking PCR primers that anneal 600 base pairs upstream or 600 base pairs downstream of the open reading frame. The PCR products were then subcloned to pRS316 and/or pRS305 (for pRS305-PTP2::URA3) and/or pRS306 (for pRS306-
fus3K42R) (Invitrogen, Carlsbad, CA). Mutation of pRS316-FUS3 or pRS316-KSS1 to obtain the analogue-sensitive alleles
FUS3Q93A or
KSS1E94A was conducted with the QuikChange site-directed mutagenesis kit (Agilent, Santa Clara, CA) according to the manufacturer's directions. Expression plasmids or the corresponding empty vector control were transformed into cells and maintained in standard SCD dropout medium (Bio 101, Carlsbad, CA).
TCA acid extraction of protein for immunoblot analysis
Cells were collected to prechilled 50-ml tubes containing 10 mM NaN3 (final concentration) and centrifuged, and the cell pellets were stored at –80°C. Alternatively, in , cells were collected in prechilled 50-ml tubes containing trichloroacetic acid (TCA; 5% final concentration). For preparation of extracts, cell pellets were thawed on ice and resuspended in 250 μl of ice-cold TCA buffer (10 mM Tris, pH 8.0, 10% TCA, 25 mM NH4OAc, 1 mM EDTA). Cells were disrupted by vortexing with 100 μl of glass beads in five 1-min bursts with chilling on ice in between. Lysates were transferred to new tubes and centrifuged for 10 min at 16,000 × g at 4°C. Pellets were resuspended in 0.1 M Tris (pH 11.0) and 3% SDS, boiled for 5 min, and then centrifuged at 16,000 × g. The resulting supernatant was separated, and protein concentration was determined using the DC protein assay (Bio-Rad, Hercules, CA). Twenty micrograms of protein was used per time point.
Whole-cell protein extracts were resolved by 12% SDS–PAGE and immunoblotting with phospho-p42/44 antibody at 1:500 (9101L; Cell Signaling Technologies, Danvers, MA;
Sabbagh et al., 2001 
), G6PDH antibody at 1:10
5 (Sigma-Aldrich, St. Louis, MO), or anti-protein A antibodies (Sigma-Aldrich). Immunoreactive species were visualized by chemiluminescent detection (Chemical ECL-Plus; Pierce, Rockford, IL) of horseradish peroxidase-conjugated antibodies (sc-2006; Santa Cruz Biotechnology, Santa Cruz, CA) at 1:10,000. The signal was visualized by chemiluminescent detection with minimal exposure to x-ray film. Band intensity was quantified by scanning densitometry with ImageJ (National Institutes of Health).
Transcription reporter assay
Cell cultures bearing the Ty1-lacZ reporter were grown to A600 of 0.8, dispensed at 90 μl into a 96-well plate, and mixed with 10 μl of α-factor. After incubation at 30°C for 90 min, 20 μl of 130 mM PIPES (pH 7.2), 0.25% Triton-X100, and 0.5 mM fluorescein di-β-galactopyranoside (M0250; Marker Gene Technologies, Eugene, OR) was added to each well, mixed, and incubated at 37°C for 30 min. The reaction was stopped by the addition of 20 μl of 1 M Na2CO3. Fluorescence was quantified using a fluorescence plate reader at A750 nm (SpectraMax M5; Molecular Devices, Sunnyvale, CA).
Model equations
To investigate the mechanisms responsible for signal specificity between the yeast mating response and invasive growth pathways, we devised six differential equation models. Each model corresponds to a different mechanism of cross-inhibition. All six models assume that the total Kss1, Ptp2, and Ptp3 concentrations remain constant for the duration of the experiments. Our experimental results using TAP-tagged proteins strongly support this assumption (unpublished data). We have observed that Fus3 and Msg5 concentrations increase two- to threefold following stimulation with 3 μM of pheromone (see Figure S1C). Therefore all six models take pheromone-dependent transcriptional induction of Fus3 and Msg5 into account. We assumed that the initial values of active Kss1 and Fus3 were zero. In the fus3Δ strain, there is a significant amount of active Kss1 prior to stimulation. However, deletion of FUS3 increases the expression level of Kss1, and for this reason we did not include the fus3 deletion time series in the data set used to test the models. There is also a small increase in the amount of basal Fus3 activity in the Ste7A7 mutant. For simplicity, this slight elevation of activity was not considered in the models.
We investigated two functional forms for the temporal profile of active Ste7. Both scenarios ignore the initial activation phase of Ste7. This approximation is based on the observation that Kss1 activation is very rapid, with peak activity occurring around 10 min following pheromone stimulation. Because Ste7 activation must be at least this fast, it is reasonable to assume that Ste7 reaches maximum activation levels immediately following pheromone stimulation. In the first scenario, Ste7 activity follows a decreasing Hill function of the form [Ste7*] = [Ste7]
0/[1 + (
k5
t)
n] (see Figure S2B), where [Ste7*]
0 is the initial concentration of active Ste7 and the parameter
k5 determines the time at which active Ste7 has been reduced to half its original value. This form of Ste7 activity was motivated by our previous work, which suggests that the upstream signaling proteins in the pheromone pathway function to convert pheromone dose information into signal duration (
Behar et al., 2008 
). In the second scenario, the active Ste7 concentration [Ste7*] decreases exponentially in time. That is, the activity profile of Ste7 has the following form: [Ste7*] = [Ste7*]
0 exp(−
k5
t), where again [Ste7*]
0 is the initial concentration of active Ste7 and
k5 is the rate at which Ste7 activity decreases. To minimize the number of free parameters, we assumed a Hill coefficient of
n = 8, which produces a step-like Ste7 response (Figure S2B) consistent with the profile suggested by our recent studies. We also assume that all Ste7 molecules are active at
t = 0 and set [Ste7*]
0 = 700 in both cases. The remaining parameter
k5 in the active Ste7 time profiles is determined by fitting the models to the experimental data. In models I, III, IVa, and IVb, the unphosphorylated Fus3 concentration [Fus3] is governed by the following equation:
The first term on the right-hand side of Eq. 1 represents the phosphorylation of Fus3 by active Ste7. Our previous work demonstrated the slow phosphorylation rate of Fus3 depends on full catalytic activity (
Hao et al., 2008a 
). That is, a mutant containing a “kinase-dead” version of Fus3 displayed rapid activation kinetics similar to Kss1. Therefore we allowed the rate constant
k1′ for Fus3 activation in the inhibitor pretreated
fus3-as strain to vary from the value
k1 in strains containing wild-type Fus3. The second term in Eq. 1 models the dephosphorylation of Fus3 by the phosphatases Ptp2, Ptp3, and Msg5. Again, the symbol * indicates the active form of the phosphatases. The models make different assumptions about how the phosphatases are activated (see below). If a phosphatase is constitutively active, then the active form is equal to the total concentration of the phosphatase. The last two terms model Fus3- and Kss1-dependent induction of Fus3. We modeled transcriptional induction using Hill kinetics. The Hill coefficient was taken to be 2 and the
Km and
Vmax values were free parameters determined by fitting the models to experimental data. We did not include a term for degradation of Fus3, because its half-life was measured to greater than 2 h and increased upon stimulation with pheromone (
Wang et al., 2006 
). The equation for the phosphorylated (active) Fus3 concentration, [p-Fus3], is given by
Because we assume that the total Kss1 concentration [Kss1]Total = [Kss1] + [p-Kss1] is constant in time, we only need to consider the active concentration [p-Kss1].
Model IIa assumes the phosphatases are constitutively active, but active forms of both Kss1 and Fus3 inhibit active Ste7. In this model, we have
Model IIb assumes the phosphatases Ptp2, Ptp3, and Msg5 are constitutively active, and Fus3 and Kss1 increase their own dephosphorylation rates. In this model, the equation for unphosphorylated Fus3, phosphorylated Fus3 and phosphorylated Kss1 become
In model IVa, Fus3 inhibits the activation of Kss1, and the three phosphatases are constitutively active. This leads to the following equation for [p-Kss1]:
In model IVb, the phosphatases are again assumed to be constitutively active, and Fus3 causes an increase in the dephosphorylation rate of Kss1 through a phosphatase-independent mechanism. In this case, the equation for [p-Kss1] is
where the term α[p-Fus3][p-Kss1] models the Fus3-dependent Kss1 dephosphorylation rate. For models I and III, the equation for [p-Kss1] is
In models I, II, and IV, Msg5 is constitutively active. In this case, the equation for [Msg5*] is
The first term on the right-hand side,
k3, models the constitutive synthesis of Msg5. The second two terms model increased synthesis due to induction by Fus3 and Kss1, respectively. The Hill coefficient
n was assumed to be 2 and the
Km and
Vmax values were free parameters estimated by fitting the models to experimental data. The final term in Eq. 12 models Msg5 degradation. For model III, in which Msg5 requires activation by Fus3, [Msg5] is substituted for [Msg5*] in Eq. 12, and the term −
k8[Msg5][p-Fus3] is added to the right-hand side of this equation. The equation for the active Msg5 concentration is given by
where we have assumed that the degradation rates of the active and inactive forms of Msg5 are the same.
In model III, in which Fus3 phosphorylates and activates Ptp2 and Ptp3, the equations for the active form of these two phosphatases are
where, for simplicity, we assume that the phosphatases are not dephosphorylated during the time course of the experiment. For models I, II, and IV, in which these two phosphatases are constitutively active, we have [Ptp2*] = [Ptp2]
Total and [Ptp3*] = [Ptp3]
Total for all times.
Normalizing the data
To estimate the parameter values of the six models, we used the time-course data shown in Figure S1. Each experiment was run with wild-type cells as a control. The experimental method used to measure MAP kinase activity only allows us to measure relative changes in activity. Therefore to compare the model output with data sets from different experimental runs requires a scaling of the experimental data in a way that does not alter the shape of the time-course data or change the relative relationships between the genetically altered strain and the wild-type control. In mathematical terms, if Xij represents the experimental measurement at time point i of the jth experiment, then Xij is related to the actual concentration by an unknown scale factor kj, which will vary for different experimental runs. Our method for selecting an appropriate set of scale factors is based on the assumption that the wild-type cells should respond in roughly an equivalent manner for each experimental run. Therefore, to find the scale factor kj for each experimental time course, we constructed an optimization problem that determines k values by minimizing the variability in the wild-type results.
Assume
Xij represents the wild-type measurements and suppose we have
m sets of time series containing
n data points. To find the set of
k values, we minimize the quantity
for
j = 1 to
m and pick the set of
k values that produced the overall minimum value of
F. We then use this set of
k values to scale the time series data for both the wild-type and genetically altered cells. In this way, we preserve the shape of the time courses and the relative ratios for the wild-type and genetically altered strains. Figure S2A shows times series for wild-type Fus3 and Kss1 activity using this scaling method. The solid black line connects the mean values at each time point. For the Fus3 and Msg5 total concentration data, the models were only fit to the shape of the curves.
Data fitting
For evaluation of the ability of each model to reproduce the experimental results, the ordinary differential equations presented above were fit to the experimental data using a Monte Carlo–based approach. The method makes use of the Metropolis algorithm, in which a Boltzmann factor containing the SSD between the experimental data and model output is used to generate acceptance probabilities for suboptimal parameter sets. In this way, the algorithm samples the parameter space using a random walk that is biased toward parameter sets that produce a good fit to the data. Briefly, the algorithm works as follows. An initial set of parameter values is picked. Using these parameter values, the model equations are solved for each experimental condition. The quantity exp(−SSD1/(2σ2)) is calculated, where SSD1 is the sum of the squared differences between the model output and the experimental results. Next, a candidate parameter set is generated at random. The distribution for the candidate parameter values is taken to be a multivariate Gaussian with the means equal to the current parameter values and SDs equal to a specified percentage of the mean values. This percentage is adjusted so that roughly 20% of the candidate sets are accepted. If the randomly chosen parameter values are less than zero, the candidate set is discarded and another candidate set is generated. The model equations are solved again with the candidate parameter set, and the quantity exp[−SSD2/(2σ2)] is computed, where SSD2 is sum of the squared differences using the candidate parameter set. If exp[−SSD1/(2σ2)] × exp[SSD2/(2 σ2)] > 1, then the candidate set produces a better fit to the data, and it is saved. Otherwise, a uniform random number r is generated, and the candidate set is saved if r < exp[−SSD1/(2 σ2)] × exp[SSD2/(2 σ2)]. The variance σ2 is estimated from the variability observed in the wild-type experimental data. If accepted, the candidate parameter set becomes the current set, and the process is repeated. Because there is a finite probability of accepting a suboptimal parameter set, the method has the advantage of being able to escape from local minima and can also be used to determine how well the parameter values are constrained by the data. The accepted parameter sets were recorded and plotted for model IVa in Figure S9.