Genetically Determined Variation in Lysis Time Variance in the Bacteriophage φX174.

Researchers in evolutionary genetics recently have recognized an exciting opportunity in decomposing beneficial mutations into their proximal, mechanistic determinants. The application of methods and concepts from molecular biology and life history theory to studies of lytic bacteriophages (phages) has allowed them to understand how natural selection sees mutations influencing life history. This work motivated the research presented here, in which we explored whether, under consistent experimental conditions, small differences in the genome of bacteriophage φX174 could lead to altered life history phenotypes among a panel of eight genetically distinct clones. We assessed the clones' phenotypes by applying a novel statistical framework to the results of a serially sampled parallel infection assay, in which we simultaneously inoculated each of a large number of replicate host volumes with ∼1 phage particle. We sequentially plated the volumes over the course of infection and counted the plaques that formed after incubation. These counts served as a proxy for the number of phage particles in a single volume as a function of time. From repeated assays, we inferred significant, genetically determined heterogeneity in lysis time and burst size, including lysis time variance. These findings are interesting in light of the genetic and phenotypic constraints on the single-protein lysis mechanism of φX174. We speculate briefly on the mechanisms underlying our results, and we discuss the potential importance of lysis time variance in viral evolution.

A key opportunity in evolutionary genetics is the decomposition of beneficial mutations into their proximal, mechanistic determinants. To date, much of this work has applied methods from molecular biology to dissect the basis of protein evolution (Harms and Thornton 2013; Dean and Thornton 2007), but life history theory (Stearns 1992; Roff 2002) similarly decomposes population growth rate into constituent life history traits. It can thus predict which mutations influencing life history will be favored by natural selection (i.e., increase population growth rate) (Stearns ). The simple life history of the lytic bacteriophage (phage) makes it particularly well suited to this sort of work. Each generation begins with a free phage adsorbing to a host cell and injecting its genome. Transcription and translation facilitate the assembly of progeny phages within the host, and, finally, host cell lysis releases the next generation of free phages. Assuming the density of host cells is not limiting, Wang provided an expression for phage growth rate in terms of mean burst size and mean time to lysis. Later, Bull extended this treatment by incorporating adsorption and phage death rates. Other studies have further developed this formalism (Rabinovitch , 2002), with many applying it to explore the mechanistic basis of beneficial mutations in lytic phages (Pepin ; Shao and Wang 2008; Bull ). Several have demonstrated that a mutation’s fixation probability depends critically on which component of life history it affects (Wahl and DeHaan 2004; Hubbarde ; Patwa and Wahl 2008, 2009).

At steady state (i.e., when host cell infections are at their stable age-of-infection distribution), population growth rate is entirely determined by mean life history traits (Bull ). But the moment a novel mutation appears, it establishes a lineage that, by definition, is far from this equilibrium. In populations deviating from a stable age distribution, variance in both lysis time and adsorption rate may play a significant role in determining growth rate (Storms and Sauvageau 2014), since it may affect how long the deviation from stable age distribution persists (Bull ).

These findings motivate interest in the possibility of genetic control of higher moments in phage life history traits. Previously, Dennehy and Wang 2011 used a microscopy-based approach to demonstrate heritable variation in lysis time mean and variance among a panel of phage λ holin protein mutants.

To explore heritable variation in phage life history traits, including lysis time variance, we studied a panel of eight genetically distinct clones of the lytic bacteriophage φX174. φX174 is an ubiquitous (Afshinnekoo ) and industrially important phage (Labrie ). It was chosen for this research in large part due to its small, well-characterized genome (Hayashi ), the first of its kind ever sequenced (Sanger ). In addition, it infects a well-studied host (Escherichia coli C), has a short generation time, and is robust to a wide range of temperatures, making it amenable to evolution experiments and other lab manipulations (Wichman and Brown 2010).

The gold standard for measuring phage lysis time and burst size is the one-step phage growth assay (Hyman and Abedon 2009). However, this method affords only limited access to information about lysis time variance and the relationship between burst size and time for a given phage. Therefore, the panel of eight clones was assessed using an updated version of a venerable experimental method first proposed by Burnet 1929, and later extended by Delbrück 1945. We simultaneously inoculated each of a large number of replicate host cell volumes with, on average, less than one phage particle of identical haplotype. These volumes were then plated at 15-sec intervals on lawns of host cells, and the resulting plaques were counted after incubation.

From this assay, we obtained data that could help determine mean and variance in both lysis time and burst size (Figure 1). However, the information we wanted was obscured by a suite of hidden variables. In wells showing lysis events, we did not know when lysis actually occurred. We also did not know the number of phage particles originally in these wells and, if that number was greater than one, how many had lysed their hosts. To infer lysis time and burst size from the two observed metrics (time of sampling and plaque count), we somehow had to circumvent these hidden variables. We overcame this challenge by developing a novel inferential apparatus that complemented the experimental framework with a maximum likelihood-based approach.

Figure 1

Typical serially sampled parallel infection assay data (pos5, = 0.71). Here, sampling was performed at 15-sec intervals between 5 min and 28.75 min. At each time point, three outcomes were possible: zero plaques, a “few” plaques (here, fewer than eight; see section Statistical methods), and more than a few. Assuming 100% plating efficiency, these outcomes implied respectively that no phage particles were added to the corresponding host cell volume; that the visible number of phage particles was added, but none lysed its host by the time of sampling; and that at least one infected cell had been lysed by the time of sampling. In these data, the first sample with more than a few plaques occurred at 15.5 min, giving an upper bound for the time to lysis (). Inset illustrates plaque counts prior to the first observed lysis event. The absence of corresponding stratification among plaque counts after 15.5 min implied high variance in burst size. Similarly, the two plaques observed at 23.5 min (arrow) corresponded to two phage particles that had not yet lysed their hosts by that time, implying high variance in lysis time. Thus, these data contain information about higher moments in both burst size and lysis time.

In this article, we describe our finding that, under the tested conditions, genetic heterogeneity could give rise to significant heterogeneity in φX174 lysis time and burst size, including lysis time variance. These findings are interesting in light of the genetic and phenotypic constraints on the single-protein lysis mechanism of φX174, and they point to promising avenues for future theoretical and experimental work.

Materials and Methods

Media, host cells, and bacteriophages

Difco LB Lenox Broth (Detroit, MI), pH 7.50, and supplemented with 2 mM CaCl2 was used throughout. We studied the performance of φX174 on E. coli C [DSMZ 13127, kindly provided by Olivier Tenaillion, Institut National de la Santé et de la Recherche Médicale (INSERM), University Denis Diderot Paris 7], the phage’s natural host. A stock of E. coli C was stored in 40% glycerol at –80°. We used this primary stock to seed bacterial colonies on agar plates stored at 4°. Colony plates were created anew every 2–4 wk. To generate experimental bacterial cultures, 10 mL LB was inoculated with a single colony and shaken in a 50-mL flask at 200 rpm and 37° for about 16 hr. The resulting stationary-phase cells were then diluted somewhere between 33- and 80-fold into 40 mL LB in a 250-mL flask. The diluted culture was shaken at 37° for about 3 hr until it reached an optical density of 0.6. Exponential-phase cells prepared in this way were used for each assay.

To synchronize phage infections, we treated cells using a starvation protocol adapted from Denhardt and Sinsheimer 1965 and kindly provided by Darin Rokyta, Florida State University. The treatment starved cells so that they “shut down,” allowing for phage adsorption but not entry. (It is likely that the treatment also had the fortunate side effect of reducing variance in time to cell division compared to an exponential-phase culture, since it was bound to cause a lag in the cell cycle. See Discussion for the cell cycle’s relevance.) First, we transferred 1 mL cells to a 1.5-mL centrifuge tube and pelleted them by spinning in an Eppendorf centrifuge (5417R) for 2 min at 20,800 rcf. After aspirating the supernatant, we resuspended the cells in 1 mL HFB-1 solution [0.06 M NH4Cl/0.09 M NaCl/0.1 M KCl/0.1 M Tris-HCl (pH 7.4) /1.0 mM MgSO4/1.0 mM CaCl2] by vortexing for 8 sec. Then, we centrifuged the solution, removed the supernatant, and resuspended in HFB-1 twice more. After removing the supernatant from the third and final HFB-1 rinse, we resuspended cells in 1 ml HFB-2 solution (HFB-1/10 mM MgCl2/5 mM CaCl2). This served as our culture of starved bacterial cells, ready for phage adsorption.

We studied a panel of eight genetically distinct clones, chosen based on reports of lysis profiles different from the wild type. These included three Epos mutants carrying point mutations in the E gene, whose product is the chief player in φX174 lysis (pos4b, pos5, and pos6, first described in Bernhardt ; kindly provided by Ry Young and Rohit Kongari, Texas A&M University); four carrying point mutations in the D-promoter region, which exerts regulatory control over transcription of the E gene (mut319, mut321, mut323, and mut324, previously described in Brown and Brown ; kindly provided by Celeste Brown and Amber Stancik, University of Idaho); and a wild type (also provided by Ry Young and Rohit Kongari) (Table 1 and Figure 2). Clones are available from D.M.W. on request. Presence of the mutations of interest was verified by Sanger sequencing of the 65–730 base region of each clone’s genome (nucleotides numbered as in Sanger ). The clones’ sequences in this region (GenBank accession KU646482–KU646588) were identical to the wild-type sequence (GenBank accession J02482.1) apart from the noted mutations. The wild type employed was the ancestor of the Epos mutants, which differed from the wild-type ancestor of the D-promoter mutants (GenBank accession AF176034) at five sites: G833A (nonsynonymous), G1650A (nonsynonymous), T2811C (synonymous), A4518G (nonsynonymous), and T4784C (synonymous) (Pearson ). The G833A difference entailed an amino acid change in the E (but not the D) gene, though it fell well outside the region coding for the E protein’s essential transmembrane region (Figure 2) (Tanaka and Clemons Jr 2012).

Table 1

Characteristics of φX174 D-promoter and E mutants

Type	Clone	Base Change(s)	Amino Acid Change(s)	Lysis Kinetics Compared to Ancestral WT
				Lysis Time on E. coli C (37°)	Via Transfected Plasmid
D-Promoter	mut319	G319Ta	V63F	+2 min	N/A
	mut321	T321C	None	+2 min	N/A
	mut323	A323G	N64G	+0–2 min	N/A
	mut324	C324T	None	+0–2 min	N/A
Epos	pos5	G624T	L19F	+10 min	Slightly faster
	pos6	G575A	R3H	Same	Much faster
	pos4B	G575A, G624T	R3H, L19F	Same	Much faster

All nonsynonymous D-promoter changes are within the C gene (Brown , 2013). L19F is an amino acid change in the transmembrane region of the E protein, while R3H is a change in E’s periplasmic N-terminus. The same mutation that changes leucine to phenylalanine in E confers a A79S change in the D protein. The five background nucleotide differences between the D-promoter and Epos mutants are not included in this table (see Materials and methods). The rightmost column refers to the bulk lysis kinetics of the mutated E genes after they had been cloned onto plasmids and inserted into E. coli K-12 cells (Bernhardt ). Relevant accession numbers are provided in Materials and methods.

A change adjacent to the D promoter’s sigma factor binding site. The three other D promoter changes are within the binding site.

Figure 2

Positions of mutations in φX174 genome. The figure shows a 832-bp region of interest, including genes C, D, E, and J. Gene C is involved in DNA replication, gene D is an external scaffolding protein required for procapsid morphogenesis, gene E is responsible for lysis of the host, and gene J is required for DNA packaging (Fane ). Gene E is shown in three sections, corresponding to three domains of the protein. The transmembrane domain binds the E protein’s substrate, MraY (Mendel ; Zheng ). Mutations are indicated with pins at the position, and the site and nucleotide change are indicated to the right of the pin head. The clone(s) carrying each mutation is indicated in parentheses. Note that one of the five background nucleotide differences between the Epos and D-promoter mutants (G833A) is included in this figure.

Phage titers were determined by combining a known volume of phage stock with 200 μL exponentially growing host cells (described above), and plating in top agar held at 42°, consisting of 7.0 g/L Difco BactoAgar (Detroit, MI) added to LB medium (as described above). After several minutes on the bench at room temperature, we incubated plates agar-side-up at 37°. Preliminary results (not shown) demonstrated that accurate plaque counts required at least 4 hr of incubation. Occasionally, for logistical reasons, we interrupted incubation with overnight refrigeration at 4°.

Serially sampled parallel infection assay

We began by adding 90 μL starved cells to each well of a TempPlate nonskirted 96-well PCR plate (USA Scientific item no. 1402-9598) at 15°. The PCR plate was mounted on a thermocycler to control its temperature precisely. Using a 12-channel pipettor, we then added 10 μl 30-PFUs/mL phage dilution to each well in each row of 12 wells, allowing for a 30-sec interval between rows. PCR plates were held at 15° for at least 20 min, during which time phages could adsorb to the starved cells. We estimate that 10–20% of phage particles adsorbed to hosts during this incubation period (Celeste Brown, personal communication; Newbold and Sinsheimer 1970).

Then, infections were initiated by raising the PCR plate temperature to 37° and adding 100 μL warm LB (incubated at 37° for at least 45 min) to each well in a row with a 12-channel pipettor, again allowing 30-sec intervals between rows. Consequently, the infective process in each row was staggered by half a minute, allowing ample time for the sampling procedure that followed. Note that the temperature increase alone did not trigger infection. The addition of rich medium, such as LB, was required for infection to advance (Newbold and Sinsheimer 1970). Based on past reports of φX174 attachment kinetics under conditions similar to our assay (Brown ), unattached phage particles likely adsorbed to hosts very quickly at this time, with an average adsorption time of less than 30 sec [assuming adsorption followed an exponential decay process (Abedon )].

It is important to note that, based on our knowledge of the mutated genes’ functions (Fane ) and past results (Brown ), the studied mutations likely had little or no effect on phage adsorption rate. This includes the nonsynonymous difference (G1650A) between the D-promoter mutants and other clones in coat protein F. While a change in the F protein could theoretically have affected adsorption (Bernal ), it is unlikely to have done so here, because the difference was merely between one basic polar residue (His) and another (Arg).

Five minutes after infection was initiated, we began sampling wells at 15-sec intervals. Sampling was performed by transferring the entire 200 μL contents of a well to a sterile 96-well PCR plate, held at 4° to halt the infective process. We assume that this sampling procedure did not artificially induce lysis in sampled infections. We staggered sampling between rows such that, factoring in the 30-sec gaps in infection between each row, we sampled one well every 15 sec between t = 5 min and t = 28.75 min. The endpoint of 28.75 min was chosen on the basis of preliminary experiments that showed most infections ending in lysis by that time. Even though not all phages would lyse their hosts by 28.75 min, we found no evidence that φX174 undergoes T4-like lysis inhibition (Abedon 1992), so we believed that samples taken between 5 min and 28.75 min would give us a reasonably complete picture of the clones’ lysis phenotypes.

At the conclusion of this process, the contents of each sampled well—now at 4°—were plated and incubated, as described, and then counted. Typical results are shown in Figure 1. Four to six replicate assays were performed for each clone.

Statistical methods

The model:

We began with an overview of the model that we assumed generated the data. We then developed methods of estimation from the model and applied them to real experimental data.

Each time an assay was performed, we assumed a small, Poisson-distributed number of phage particles had been added to each well. The mean of the Poisson was allowed to vary between assays, which is to say, we assumed that phage concentrations in the stocks differed by clone and could change over time.

Within the wells, we assumed phage particles adsorbed instantaneously to starved cells upon the initiation of infection. (We discuss the rationale underlying this assumption in the Results section). We assumed that each phage particle infected and then lysed a cell T minutes after the addition of warm media, where T could not be less than a biologically imposed minimum (the latent period), and lysis was equally probable every moment thereafter (though we also considered a model in which lysis probability changed with time, as described below). We assumed that burst size grew linearly with lysis time, a relationship generally assumed and supported by past work (Hutchison and Sinsheimer 1966; Abedon 1989; Wang ; Abedon ; Chao ; Bull ; Wang 2006; but see Hutchison and Sinsheimer 1963 and Storms ). The fate of each infecting phage particle within a well was assumed to be independent of other particles in that well. Wells were sampled and plated at 15-sec intervals from t = 5 min through t = 28.75 min, and we assumed 100% plating efficiency, implying that the plaque count on a plate revealed the true number of phage particles in a well at the time of sampling. A plaque count of zero meant no particles were in the well, a count of one to a few implied that as many particles as visible plaques were added but that no host cells had been lysed yet, and a count of more than a few represented the sum of the burst sizes of lysis events prior to the sample time, plus any infected yet unlysed cells. Note that, in the third case, we did not know how many phage particles were in the well to begin with, we did not know how many lysis events had occurred by the time of sampling, and, of the lysis event(s) that had occurred, we did not know when they occurred, except that they must have occurred prior to the moment of sampling. From a statistical perspective, this suite of hidden variables made the problem both interesting and challenging.

Notation:

Developing estimators under this model required substantial notation. Table 2 summarizes this notation for the reader’s reference.

Table 2

Summary of notation used in statistical modeling

Symbol	Meaning
β_a	The number of phages in a well from assay a is Poisson with mean β_a
p	Probability a random phage in a randomly selected well in the experiment has not lysed its host
L_a	# of wells examined in assay a
R_a	# of wells examined in assay a that show lysis events
M_a	# of wells examined in assay a that do not show lysis events but have counts ≥1
Z_a	# of wells examined in assay a that have zero phages in them
N_ar	# of phages in well r from assay a. We reorder wells such that r=1,2,⋯,Z_a have 0 phages, wells r=Z_a+1,Z_a+2,⋯,Z_a+M_a do not show lysis events (but are nonzero), and r=Z_a+M_a+1,Z_a+M_a+2,⋯,L_a are wells showing lysis events
N_a¯	Mean count in wells not showing lysis events; ∑r=Z_a+1M_aN_ar/M_a
T_iar	Time to lysis of phage i in well r from assay a
D_T	Burst size of a phage that lysed its host at time T
C_tar	Observed count in well r sampled in assay a at time t
B_tar	The number of phages in well r, assay a, that have lysed their hosts by time t when well is sampled; B_tar=∑i=1N_arI{T_iar≤t}
Y_tar	An indicator function that is 1 when B_tar>0 and 0 when B_tar=0
U_tar	The number of phages in well r, assay a, yet to lyse their hosts when sampling occurs at time t; U_tar=N_ar−B_tar
X_tdr	Sum of lysis times for all phages in well r, assay a, that have lysed their hosts prior to sampling at time t; X_tar=∑i=1N_arI{T_iar≤t}
α, μ, σ	Parameters describing the assumed linear relationship between burst size and time C_T=μ+αT+ε where ε∼N(0,σ^2)
t_0	Latent period. For t<t_0, the lysis rate = 0; for t>t_0 it is constant per unit time
w_t	The probability a phage will have lysed its host by time t
λ	Parameter relating lysis probability to time; w_t=1−e^−(t−t_0)/λ

Number of phage particles per well:

Our first task was to estimate , the Poisson parameter that governs the probability of having phage particles in each well for assay a. This was somewhat complicated by the fact that, before the latent period ended and lysis began, we observed the actual Poisson counts, but after this period, we only observed the true Poisson counts for “zero” wells and wells that do not show lysis events. In wells showing lysis events, we knew only that the Poisson count must have been greater than zero. To derive an estimator tailored to the situation, we ignored the time element of the data and defined p as the probability that a randomly selected phage particle from a randomly selected well would not have lysed its host. We let n be the Poisson number of phage particles initially in a well (i.e., before lysis occurs). Then, we noted that each sampled well fell into one of three categories: samples with zero phage particles, samples with phage particles but no lysis events, and samples showing lysis events. The probability of a well belonging to each category was, respectively, , , and .

We defined wells not showing lysis events as those having plaque counts greater than zero, and less than or equal to seven. When the Poisson parameter β is less than or equal to one, counts greater than seven are highly improbable in wells where lysis has not occurred (). Of course, it is possible that wells showing fewer than eight plaques could in fact be the results of infections. But, ultimately, we had to set the threshold for distinguishing lysis events from nonlysis events somewhere. This threshold-setting exercise was not straightforward. Complicating factors included the following:Therefore, given these intermediate plaque counts’ rarity, their potential origins in experimental error, and their potential to distort our estimates of clones’ lysis phenotypes, we chose to set the threshold for lysis events at greater than or equal to eight. It is worth noting, however, that when our data instead were analyzed using a threshold of greater than or equal to five, our results did not change qualitatively (analysis not shown).

Poisson average numbers of particles per well varied for each assay, so setting the threshold based purely on the Poisson probabilities was not trivial.

There was a nonzero chance of slight inaccuracies in some plaque counts. Rare errors could have occurred in the counting of plaques that would have increased counts on average (e.g., plaque-like air bubbles). In addition, certain experimental manipulations (e.g., mechanical stress from pipetting) could feasibly have induced very premature lysis of some cells, resulting in plates with extremely low plaque counts that otherwise would have had just one plaque.

The few plaque counts of six and seven, a couple of which were sampled at unusually early time points, could have exerted disproportionate influence on estimates of clones’ lysis phenotypes if included as lysis events. This is because, as discussed below, an artificially early lysis event would force the value of the parameter to shift to some time before that event.

We let the number of wells in assay a be and the number of wells that showed lysis events be . The likelihood of the data, then, wasWe then took the log of the likelihood, the partial derivatives with respect to and p, set each equal to zero, and did a considerable amount of algebra (see Appendix 1 for details). This resulted in the following estimator:We validated this estimator by simulating data under and 2, and we compared it to two reasonable alternatives: the “zero-class” estimator [, where is the proportion of sampled wells that had plaque counts of zero], and the “mean-count” estimator (, where X is the mean among observed plaque counts, applied only to wells sampled before we witnessed any lysis events). The results (which we will not cover in detail) revealed that Equation 2 had a very slightly positive bias ( when , falling to at ) but a lower mean squared error than the other estimators across conditions. Valuing the nearness of the estimate to the truth above all else, we adopted the new estimator of throughout.

Time to lysis:

We considered two nested models describing time to lysis. Both assumed that phages do not lyse cells until the latent period of ends. The exponential model assumed that lysis occurred with equal probability each moment after , while the Weibull model allowed the probability of lysis to change with time. Under the Weibull, the time to lysis, T, had the following distribution:where α was a hazard parameter that dictated how lysis probability changed with time. When , the Weibull model simplified to the exponential model. Because they are nested, we employed a likelihood ratio test to determine if the improved fit under the Weibull justified the additional parameter. To do so, we analyzed each genetically distinct clone separately. We fitted the data to each model using maximum likelihood (as described below), estimated parameters, obtained the log-likelihood, and took the difference: . We then simulated 100 datasets under the null (exponential) parameter estimates. For each simulated dataset, we fitted the data to both models and calculated the difference in their log-likelihoods: . We then approximated the P-value by the proportion of greater than or equal to . Large P-values (i.e., ) would indicate the simpler exponential model described the data sufficiently well. This test returned P-values greater than 0.3 for all clones, and we therefore assumed the exponential model for the duration of the study. One important implication of this assumption was that, under the exponential model, the lysis time variance was equivalent to .

Defining as the number of phages in well r, assay a that had lysed their hosts by time t, followed the Poisson distribution with rate . If we let be an indicator function that is one when is greater than or equal to one, and zero when is zero, and we let be the number of phages that had not lysed their hosts in well r, assay a at time t, we show in Appendix 2 that the log-likelihood of the data were proportional toNote that directed us to sum over all assays of a particular clone. We analyzed each clone separately. We have omitted a clone subscript throughout to avoid further complicating notation. In practice, we obtained joint maximum likelihood estimates of λ and by establishing a grid of joint values at narrow intervals, calculating the log-likelihood of each using Equation 3, and selecting those values that together give the best likelihood. Confidence intervals on λ and were obtained by parametric bootstrap. To do this, we used , , and to simulate one suite of datasets to match (in number and well count) the real datasets for the clone under analysis. We analyzed the simulated suite of datasets using Equation 4 and a grid of potential values to obtain and . The one condition we imposed in this estimation step was that must be less than or equal to the minimum observed lysis time from the real data. If we did not do this, the bootstrap would have returned values and hence confidence bounds for greater than this minimum lysis time, and we knew such values would not be valid confidence bounds when the probability of the real data under them is zero. We repeated 500 bootstrap replicates. We then sorted the values, and defined the 95% confidence interval as the central 95% of the values. We independently sorted the values and used the central 95% of these. Thus, the confidence intervals on and were not joint. We did not validate that these intervals capture the truth greater than or equal to 95% of the time, but we did confirm that the means of and were close to and , respectively, suggesting that the parametric bootstrap was a valid technique.

We also were interested in inferring mean lysis times from our data. Given the unique properties of exponential distributions, we could derive each clone’s mean lysis time simply by summing its estimated and λ. A confidence interval for this estimate could be obtained by repeating this procedure using the upper and lower bounds on the individual parameter estimates.

We now turn our attention to how to determine whether lysis time parameters λ and are different for different clones. We conducted several formal tests for differences between clones using a likelihood ratio framework. To test whether both λ and clones differ, we note that, under the null hypothesis, a single shifted exponential function would describe the data from both clones being compared; the alternative hypothesis asserted the data for each clone came from a distinct function. We calculated the log-likelihood of the data under the null () by pooling both clones’ data, fitting as just described to obtain and , and then using Equation 4. We calculated the log-likelihood under the alternative () by doing the same thing, except we fitted each clone’s datasets separately, calculated their log-likelihoods separately, and then summed. We define . To determine whether is significantly large, we needed the distribution of under the null. We obtained this by simulation. We simulated data analogous in size and sample times to the pooled data using and , repeated the fitting process, and calculated log-likelihoods and their difference to yield . We repeated this process 500 times, and took the approximate P-value of as the proportion of bootstrap simulations where was greater than .

To test whether just λ differed between pairs of clones, we used a similar likelihood ratio test. Under the null, both clones shared a common value of λ but could have different values of ; under the alternative, both λ and could differ between clones. Note that the null did not fix λ at a single specific value; rather, it merely stipulated that λ, whatever value that may be, be shared by the clones being compared. We fitted the data to the null by finding the maximum likelihood estimates for each clone under the constraint of a shared value of λ and unconstrained values of and then summing their log likelihoods (). We fitted the alternative by imposing no constraints, and summing the two log-likelihoods (). We defined . To obtain the distribution of under the null, we simulated 500 datasets using the null parameter estimates, repeated the model-fitting exercise for each, and took the difference in the log-likelihood: . The proportion of the bootstrap replicates where was greater than was the approximate P-value for this comparison. To test for significant differences in , we repeated this procedure, except was constrained and λ was allowed to vary.

To test whether the distribution of P-values calculated for each set of pairwise parameter comparisons diverged from the null expectation (i.e., a uniform distribution), we binned the P-values into four categories [0.00,0.25), [0.25,0.50), [0.50,0.75), and [0.75,1.00] and conducted an exact multinomial test against expected frequencies in each bin under the null. The null bin frequencies were equal across the four categories (i.e., 25% in each bin).

Burst size as a function of time:

We modeled the potential for burst size to change with time as a linear function, , where was the expected burst size when T minutes had elapsed since , μ was the burst size at , α was the rate of increase in burst size (i.e., the slope parameter), and ε captured the normally distributed variation from expectation (mean 0, variance ). Note that none of the terms in this model was observable. If a well sampled at time t showed a lysis event, we knew that lysis occurred before t, but we did not know exactly when. In addition, we did not know how many phage particles initially were in the well. Thus, both lysis time and burst size were hidden variables. In Appendix 3, we show that the observed count (rather than the burst size) could, like burst size, be expressed as a linear function:where was the count observed in well r sampled at time t from assay a, was the number of phage particles in well r (assay a) that had lysed their hosts by time t, was the sum of their lysis times, and δ was normally distributed with mean zero and variance . Our goal was to estimate the parameters μ, α, and σ, but because Equation 5 depended on the unobservable and , we could not use standard regression. Instead, we employed an expectation maximization (EM) algorithm, through which the problem was circumvented by replacing and with their expected values before doing regression, and then iterating (Meng and Rubin 1993). More specifically, our EM algorithm involved the following steps: (1) guess values for and ; (2) estimate the parameters μ, α, and σ by least-squares regression (equivalent to maximum likelihood) assuming the guesses of and are correct; (3) assuming the parameter estimates are correct, calculate the expected values of and , and update the guesses to these values; and (4) repeat steps 2 and 3 until convergence occurs. It turns out that both parameter estimation (step 2) and calculating expectations (step 3) involved a fair amount of math; we present this in Appendix 3.

To obtain confidence intervals and test for differences between clones in the regression parameters, we used a parametric bootstrap approach. For each genetically distinct clone, we simulated 1000 datasets using the estimated parameters and, each time, estimated α, μ and σ. We then sorted the bootstrap set of parameter values individually, using the 25th and the 975th values (i.e., the central 95%) to define the confidence intervals. To keep our analysis simpler, we did not conduct formal likelihood ratio tests for significant differences between clones for these parameters. Instead, we compared confidence intervals and applied the conservative criteria that when confidence intervals were nonoverlapping, the two clones were significantly different.

Data and code required to repeat this analysis can be found in the Supplemental Material File S1.

Data availability

Strains are available from D.M.W. upon request. File S1 contains a description of all supplemental material, including raw experimental data and the R code that was used to analyze them. Sequences of the WT (Epos mutant ancestor) and D-promoter ancestor are given by GenBank accession numbers J02482.1 and AF176034, respectively. Regions of the mutants' genomes were sequenced, and these sequences also were archived in GenBank (accession numbers KU646482 to KU646588).

Results

Thirty-seven assays were conducted for the eight clones (four to six each). (Raw data from the assays are available in Supplemental Material, File S1.) The estimated Poisson parameter ranged from 0.20 to 0.91, with a mean of 0.46. It follows that approximately 63% of all wells had no phage particles, 29% had one particle, 6.5% had two particles, 1% had three, and 0.1% had four. Observed plaque counts ranged from zero to 1270, with eight being the smallest plaque count classified as representing a lysis event (see section Statistical methods).

Across the clones, the earliest times at which lysis was observed ranged from 11.0 min (pos6) to 14.5 min (mut321, mut323). Estimates of are shown in Table 3. This range of latent periods was consistent with past results, as reported in Bernhardt (one-step growth assays in figure 2A in Bernhardt , about 12 min for wild type, pos6, and pos4B infections, though about 22min for pos5) and Brown (about 12–14 min for the D-promoter mutants and wild type).

Table 3

Estimates of λ and for each clone

			Clone
Clone	λ (min)	t_0 (min)	WT	319	321	323	324	5	6
WT	11.2 (7.5, 15.2)	12.2 (11.0, 13.2)
319	5.8 (4.0, 7.8)	12.7 (12.0, 13.0)	0.028*
321	4.8 (3.5, 6.2)	14.5 (14.0, 14.8)	0.002*	0.006*
323	6.2 (4.5, 8.5)	13.9 (13.2, 14.2)	0.062**	0.132**	0.448***
324	6.0 (4.2, 8.2)	13.4 (12.5, 13.5)	0.028*	0.548***	0.134**	0.698***
5	5.7 (4.0, 7.5)	13.2 (12.5, 13.5)	0.010*	0.768***	0.082**	0.438***	0.778***
6	8.9 (6.2, 12.2)	10.0 (9.0, 10.8)	0.012*	0.008*	< 0.002*	< 0.002*	< 0.002*	0.008*
4B	9.5 (7.0, 12.5)	12.7 (11.8, 13.0)	0.800***	0.058**	0.010*	0.192***	0.140**	0.076**	< 0.002*

95% confidence intervals are in parentheses. The diagonal matrix to the right gives the P-values (based on 500 bootstrap replicates) of tests of the null hypothesis, which holds that the clones being compared have the same joint λ and values. * P < 0.05, ** P = 0.051–0.150, *** P > 0.151.

The smallest observed value of λ among the clones was 4.8 min (mut321), meaning that the cumulative probability of lysis after rose most rapidly for this clone under the experimental conditions. The wild type was found to have the largest λ (11.2 min). Values of λ for other clones, along with confidence intervals, can be found in Table 3.

Figure 3 illustrates how the observed increase over time in the proportion of wells showing lysis events was used to estimate the cumulative lysis probability function for each clone. The figure also includes pairwise comparisons of each clone with the wild type, demonstrating clear differences between them. Note that Figure 3 does not visualize the β parameter, which, in addition to the lysis events shown, was important in the estimation of lysis time phenotypes. Taken together, the curves form three rough clusters: (1) pos6, (2) wild type and pos4B, and (3) the remaining clones. This clustering can be observed clearly when the parameters composing those functions, λ and , are plotted against each other (Figure 4).

Figure 3

Estimated cumulative lysis probability by time functions, pairwise comparisons of each clone with the wild type. The wild type panel (top left) provides an example of fitting a cumulative lysis probability by time function to wild-type data. The gray line is the best-fit curve. Dotted lines represent 100 of the 500 bootstrap replicates. Each bubble’s position shows the proportion of wells sampled at a given time in which we observe a lysis event; the size of each bubble shows the number of samples on which the proportion is based. In all other panels, the wild type is represented with a gray line and gray bubbles.

Figure 4

Estimates of parameters λ and that define the lysis probability function for wildtype and each mutant. Bars represent the approximate 95% confidence intervals for each parameter independent of the other (based on 500 bootstrap replicates).

These visual comparisons were formalized into pairwise statistical tests. Along with parameter estimates, Table 3 shows the pairwise P-values associated with the null hypothesis that the clones being compared came from the same function (i.e., that jointly, both λ and are the same for the clones). If either parameter differed, this test would reject the null. By contrast, Table 4 shows the results of two different sets of tests: the set of pairwise tests of the null in which λ was held to be the same between clones, but was allowed to differ, and a similar test in which was held to be the same, but λ was allowed to differ. In other words, the first of these two sets of tests gauged differences between clones’ values of λ, whereas the second set of tests gauged differences between clones’ values of .

Table 4

Approximate P-values of pairwise comparisons of individual parameters in lysis probability by time functions

	Clone
Clone	WT	319	321	323	324	5	6	4B
WT		0.010*	< 0.002*	0.028*	0.018*	0.006*	0.322***	0.506***
319	0.478***		0.404***	0.888***	0.870***	0.844***	0.082**	0.028*
321	< 0.002*	0.006*		0.260***	0.314***	0.464***	0.010*	< 0.002*
323	0.050*	0.132**	0.246***		0.836***	0.702***	0.136**	0.092**
324	0.084**	0.334***	0.044*	0.486***		0.752***	0.078**	0.042*
5	0.210***	0.478***	0.014*	0.364***	0.698***		0.068**	0.022*
6	0.038*	0.002*	< 0.002*	0.002*	< 0.002*	0.004*		0.780***
4B	0.552***	0.910***	0.004*	0.142**	0.276***	0.478***	0.006*

Above the diagonal are results of tests of the null in which λ is held to be the same between clones but is allowed to differ. Below the diagonal are results of similar tests in which is held to be the same but λ is allowed to differ. In other words, the first of these two sets of tests gauges differences between clones’ values of λ (above the diagonal), whereas the second set of tests gauges differences between clones’ values of (below the diagonal). Based on 500 bootstrap replicates each. * P < 0.05, ** P = 0.051–0.150, *** P > 0.151.

Notice that we have not corrected for multiple tests. Our P-values are based on 500–1000 bootstrap replicates, whereas the amount of replication needed to achieve significance under a Bonferroni correction would be in the 5000–10,000 range and computationally infeasible. We note, however, that if the joint null hypothesis were universally true (i.e., all clones followed the same function), we would expect the P-values of the 28 pairwise tests to follow a uniform distribution, and we would expect just one or two false positives (28  0.05 = 1.4). Indeed, we would be unlikely to have more than three false positives, since the binomial probability of greater than three false rejections is 0.0491. Instead, the distribution of P-values for the joint λ- comparisons skewed strongly from a uniform distribution [exact multinomial test, P < 0.0001 (McDonald 2014)], and we observed 13 pairwise results less than or equal to 0.05 (Table 3). [The binomial probability of 13 or more false rejections is .] We suggest, therefore, that most, though perhaps not all, of the differences that appeared as significant were truly significant. Similarly, the tests for pairwise differences in only λ and only yielded 10 and 13 values less than or equal to 0.05, respectively, which is far more than expected under the null of a uniform distribution (exact multinomial test, P = 0.0068 and P < 0.0001, respectively). (The probability of 10 or more false rejections is , and the probability of 13 or more is ; Table 4.)

We also note that our omission of adsorption time from our statistical model (see Materials and Methods) had little effect on these findings. By ignoring adsorption time, we treated the time to lysis as a single exponential process, when in reality it involves the sum of two exponential waiting times (Abedon ). However, additional analysis (not shown) demonstrated that, because adsorption happens at least an order of magnitude faster than lysis (Brown ), the distribution of lysis time under the one-exponential model is very similar to that generated by a more complex two-exponential model. Furthermore, the variance of our one-exponential model is slightly larger than that of the two-exponential model, which implies that tests underlying the pairwise phenotypic comparisons were conservative. To see why, recall that the variance of an exponential with mean λ is . The variance of the sum of two independent exponentials with means and is the sum of their individual variances: . But the variance of a single exponential with mean is . Therefore, although we might have gained marginally more power in our phenotypic estimates by adopting a more complex two-exponential statistical model (supplemented with data on adsorption) or following the preattachment protocols described in Brown and Cherwa Jr et al. 2009, the methods we employed were adequate.

The regression analysis of burst size by time revealed that burst size increased with time for all clones (Figure 5). The parameter estimates for each clone are provided in Table 5 along with confidence intervals. Figure 5 illustrates these regression parameter estimates graphically, showing that clones fell into four clusters of decreasing slope: (1) pos6, (2) pos4B and pos5, (3) wild type and mut319, and (4) mut321, mut323, and mut324. Formal testing of parameter differences was not conducted. However, where confidence intervals were nonoverlapping, there was good evidence for a significant difference. For slope, the confidence intervals within the four clusters were generally overlapping while intervals between clusters were not.

Figure 5

Estimated burst size by time relationships for wild type and each mutant.

Table 5

Parameter estimates and confidence intervals for burst size by time regression analysis

Clone	Slope (α^) (PFU/Min)	Intercept (μ^) (PFU)	σ^ (PFU)
WT	16.7 (15.1, 18.4)	28.4 (19.9, 36.5)	19.9 (15.6, 23.5)
319	19.3 (19.0, 19.6)	12.6 (11.4, 13.8)	3.2 (2.6, 3.7)
321	11.3 (10.7, 12.0)	13.9 (11.6, 16.1)	7.4 (6.2, 8.3)
323	12.4 (11.7, 13.1)	14.8 (12.1, 17.5)	7.2 (5.8, 8.3)
324	8.7 (7.6, 9.8)	27.5 (23.0, 32.3)	13.0 (10.7, 15.2)
5	26.7 (25.0, 28.3)	34.6 (28.1, 41.0)	17.5 (14.4, 20.2)
6	31.6 (31.1, 32.2)	8.0 (5.3, 10.7)	7.2 (5.8, 8.4)
4B	25.1 (24.5, 25.6)	16.7 (13.9, 19.2)	6.6 (5.3, 7.7)

Each confidence interval based on 1000 bootstrap replicates.

Table 6 summarizes the calculated mean and variance in lysis time for each clone.

Table 6

Estimates of lysis time mean and variance of φX174 clones

Clone	Lysis Time
	Mean (min)	Variance (min2)
WT	23.4 (18.5, 28.4)	125.4 (56.3, 231.0)
319	18.5 (16.0, 20.8)	33.6 (16.0, 60.8)
321	19.3 (17.5, 21.0)	23.0 (12.3, 38.4)
323	20.1 (17.7, 22.7)	38.4 (20.3, 72.3)
324	19.4 (16.7, 21.7)	36.0 (17.6, 67.2)
5	18.9 (16.5, 21.0)	32.5 (16.0, 56.3)
6	18.9 (15.2, 23.0)	79.2 (38.4, 148.8)
4B	22.2 (18.8, 25.5)	90.3 (49.0, 156.3)

95% confidence intervals are in parentheses. Mean lysis time = . Variance in lysis time = .

Discussion

In recent years, it has been found that the probability that a beneficial mutation will fix in a lytic phage population depends sensitively on the life history parameter it affects (Wahl and DeHaan 2004; Hubbarde ; Patwa and Wahl 2008, 2009). Furthermore, analytical work has demonstrated that variance in life history traits, not just the trait means, can influence growth rate predictions, particularly in iteroparous populations far from a stable age distribution (Bull ). Inspired by this body of work, we wished to learn more about whether bacteriophage life history traits, including higher moments (e.g., lysis time variance), could be subject to the phage’s genetic control.

To this end, we used a serially sampled parallel infection assay to interrogate the lysis time and burst size phenotypes of a panel of eight genetically distinct φX174 clones. Among these clones were four with mutations in the D promoter that have been reported to downregulate the transcription of several genes, including the lysis gene (E); three Epos clones with mutations in E that have been found to upregulate its expression; and the wild-type ancestor of the Epos mutants. (As noted in the Materials and Methods, these mutations were introduced on slightly different genetic backgrounds.)

Using a novel statistical approach that allowed us to overcome the challenge of hidden variables, we found significant variation in both lysis time and burst size among these clones, including evidence that φX174 exerts genetic control over lysis time variance under our experimental conditions.

The inferred differences among the clones’ burst size phenotypes are generally consistent with our limited understanding of the mutations’ biological effects (Table 5 and Figure 5). In accordance with their depressed rates of transcription at the D promoter (Brown ), burst size increased more slowly than the wild type in three of the four D-promoter mutants. (Burst size increased with lysis time at about the same rate as the wild type in the fourth D-promoter mutant, mut319.) By contrast, burst size increased more quickly than the wild type in the three Epos mutants. This may be related to the findings of Bernhardt that the R3H and L19F mutations lead to markedly faster E protein synthesis than the wild type, perhaps by increasing mRNA translatability. Given the compact genome of φX174 and the location of E within the D gene (Figure 2), it is possible that these mutations somehow upregulate gene expression more globally, thereby precipitating more rapid accumulation of progeny virions and larger burst sizes.

As with burst size, we inferred significant variation in lysis time phenotypes among the eight clones under the tested experimental conditions. The clones seem to differ in their overall cumulative lysis probability functions (Table 3) and in both and λ individually (Table 4).

Within these differences, we observe an inverse relationship between and λ (Figure 4): under the tested conditions, the mutations appear to delay the onset of lysis while causing the cumulative lysis probability to rise more quickly compared to the wild type (Figure 3). (The most noteworthy exception is pos6, the only clone with a value of less than that of the wild type.) A consequence of this trend is that, despite substantial differences in and λ, the clones’ mean lysis times (computed as the sum of and λ) are rather similar (Table 6). A mechanistic investigation of the apparent maintenance of mean lysis time is beyond the scope of this study, but we speculate that it may reflect constraints on lysis timing imposed by the host’s cell cycle (for φX174: Witte ; Young ; and other phages as well: Hadas ; Storms ). While mutations in E and related genes could potentially bias lysis timing earlier or later, or make the process more or less noisy, the host’s physiology and cell cycle may prevent extreme divergence from a common mean.

In fact, given the genetic and phenotypic constraints on φX174 life history, it is notable that burst size and lysis time phenotypes appear to diverge as much as they do. A primary genetic constraint is the phage’s compact genomic architecture. The E gene is entirely embedded within the structurally essential D gene (Barrell ; Figure 2), and the D promoter simultaneously regulates the transcription of E, D, and four other structural genes downstream (i.e., J, F, G, and H) (Fane ). The work presented here stands with and expands upon previous studies (e.g., Brown , 2013; Bernhardt ) in demonstrating that, despite the relatively limited genotype space accessible to the organism, the genomic architecture of φX174 permits a significant degree of flexibility in lysis time and burst size phenotypes. Still, the genetic constraints on these phenotypes are clearly visible in the form of pleiotropic effects on burst size and lysis time in our data.

The inferred disparity among our clones’ lysis time parameters is also interesting given the aforementioned phenotypic constraints imposed by the host’s cell cycle and physiology. Whereas λ and many other larger DNA phages effect lysis via a tightly regulated, two-protein lysis strategy, φX174 and other Microviridae lyse their hosts by way of a single protein (Young ; Young and Wang 2006). The lysis timing of φX174 is thought to be dictated largely by the timing of host cell division and, therefore, compared with phages like λ, even less subject to the phage’s genetic control (Witte ; Young ; Zheng ; Dennehy and Wang 2011; Singh and Dennehy 2014; Storms ). But our demonstration of genetic control of both lysis time mean and variance hints at the possibility of as-yet-unknown mechanisms exerting some control over phage lysis timing. Indeed, researchers have long speculated that protein E may mediate lysis timing not just through its interactions with MraY, but also through interactions with the host protein SlyD (Young and Wang 2006). SlyD plays an important role in cell division and cell cycle regulation (Roof ) and is required for stability of wild-type E protein (Bernhardt ).

Strengths and limitations of our experimental and statistical approach

The strength of the experimental and statistical approach taken in this study lies in its power and applicability. It is powerful in that, more effectively than most other methods, it simultaneously captures information related to lysis time variance and burst size as a function of time. The approach also has potential application to a range of microbiological questions. The principles and tools underlying the assay and inferential apparatus could be applied usefully to other cases in which one seeks to make inferences about dynamic, Poisson-distributed processes that are sampled destructively over a time course (e.g., tracking the fates of individual tagged mutants in a dynamic population of cells or viruses).

However, we concede our approach has limitations. There is increasing evidence that host cells’ physiology plays a vital role in determining phage infection dynamics (Hadas ; Rabinovitch ; Storms ; Bull ). Yet, in synchronizing and rapidly sampling infections, our assay exposed the cells to temperature changes and nutritional regimes on a short timescale, which may have altered cells’ physiology in ways important for the measured phenotypes (van Elsas ). While these sorts of manipulations are not uncommon in the literature (e.g., Wang 2006; Cherwa Jr et al. 2009), they limit the generality of some conclusions drawn from this work. The unique conditions of the assay complicate comparisons to previous estimates of the clones’ burst size and lysis time phenotypes, and they make it impractical to validate our methods against traditional single-step growth experiments with much confidence.

By the same token, the lysis times and burst sizes inferred in this work may differ from those φX174 exhibits in nature. Environmental differences in factors like host density and diversity, as well as temperature changes and nutritional regimes, could possibly lead to differences between phenotypes measured in the assay and in the wild. That said, it is unclear how well any particular assay approximates the phage’s natural environment, which remains largely uncharacterized (Michel ). Further study of the ecological niche of φX174 would be necessary to clarify which environmental differences between our assay and the wild are most relevant to the phage’s life history phenotypes.

Finally, our assumption of a linear increase in burst size with time, though substantiated by significant past work with φX174 and other phages (Hutchison and Sinsheimer 1966; Wang ; Wang 2006), may be subject to debate (Storms ). We hope further research will help better evaluate the appropriateness of this assumption.

When it comes to achieving the core aims of this work, we believe the strengths of our experimental and statistical approach outweigh its limitations. It allowed us to demonstrate that, under at least one set of conditions held constant for eight similar but distinct φX174 clones, burst size and lysis time phenotypes (including lysis time variance) differed significantly. Furthermore, it is reasonable to assume that the degree of genetic variability among the clones is comparable to that which is the substrate for natural selection acting in the wild.

Evolutionary implications of heritable variation in lysis time variance

Our finding that variance in the lysis time of φX174 can be subject to at least some genetic control inspired us to ask, “What, then, could be the evolutionary implications of mutations affecting variance in phage lysis time?” Trivially, if burst size grows more slowly than exponentially once a cell has been infected, then a mutation increasing lysis time variance will be selectively favored, because early lysis events contribute more to population growth than late lysis events detract (Bull ; but see Singh and Dennehy 2014). Empirically, phage burst size is thought to increase linearly (Hutchison and Sinsheimer 1966; Wang ; Wang 2006; but see Storms ), immediately suggesting one mode of selection. Mutations increasing lysis time variance in a given environment likely would fix quickly under this mode of selection.

In order to control for the influence of this mode of selection, we recently have explored a model in which burst size increases exponentially in time since infection (D. M. Weinreich, C. Brown, and L. M. Wahl, unpublished data). This has allowed us to infer the intrinsic selective consequences of mutations that influence lysis time variance. As expected, in phage populations growing at their stable age-of-infection distribution in unlimited host cells (Bull ; Bull 2006), we found such mutations to be selectively neutral, since, by construction, the contribution of early bursts to population growth is now exactly balanced by (exponentially larger) late bursts.

However, the moment a mutation appears, the mutant subpopulation is by definition not at stable age distribution, and indeed we have found in simulations that mutations increasing variance in lysis time are transiently enriched in the population (D. M. Weinreich and C. Brown, unpublished data). This effect follows from the fact that the frequency of phages with high lysis time variance exceeds (sometimes considerably) the frequency of their low-variance brethren early in each lysis cycle. While the reverse is true late in the cycle, the effect is of smaller magnitude and much shorter lived. Finally, overlaying periodic demographic bottlenecks (e.g., as employed in lab conditions, but also surely in nature), we observed that mutations increasing lysis time variance can enjoy a dramatic increase in population frequency, despite their neutrality when they finally reach their stable age-of-infection distribution (D. M. Weinreich and C. Brown, unpublished). We therefore suggest a second mode of selection acting on mutations affecting variance in lysis time.

Future directions

The findings reported here suggest two complementary avenues for future investigation. One is microbiological, and it involves exploring φX174 lysis in relation to its genetic determinants and its host’s biology. It would be interesting to explore the mutants’ performance in a range of assays, such as those involving thermally regulated microscopy (Dennehy and Wang 2011), single-cell microfluidics-based techniques (Yin and Marshall 2012), or periodic filtration and sampling (Hutchison and Sinsheimer 1963). (The last and oldest of these techniques could prove especially fruitful, since it would allow direct access to both lysis time and burst size.) In addition, it would be highly informative to use methods like those designed by Storms to study how φX174 mutants’ burst sizes and lysis times vary as a direct function of the host’s cell cycle. These types of experiments could yield valuable insight into the regulation and host-dependence of φX174 life history phenotypes.

Theoretical evolutionary problems present another opportunity for future research. Does natural selection act on lysis time variance in nature? What is the fixation probability of a mutation that affects lysis time variance, and how do populations carrying such mutations evolve in the face of periodic bottlenecks? Existing computational models could be strengthened by more thoroughly exploring the models’ parameter sensitivity and investigating how mutations that pleiotropically influence adsorption, lysis time, and burst size affect predicted population dynamics. In addition, the models must be integrated with other theoretical work exploring the fixation probabilities of mechanistically distinct beneficial mutations in lytic phages (Wahl and DeHaan 2004; Hubbarde ; Patwa and Wahl 2008, 2009). Once elaborated and refined, models could be tested empirically through both “single-sampling” assays and longer-term evolution experiments. The former approach would assess the probability of loss of lysis time variants introduced at very low frequency in replicate phage populations, while the latter would explore whether real evolutionary outcomes match our theoretical predictions.