Cyanophages from a less virulent clade dominate over their sister clade in global oceans

Environmental virus communities are highly diverse. However, the infection physiology underlying the evolution of diverse phage lineages and their ecological consequences are largely unknown. T7-like cyanophages are abundant in nature and infect the marine unicellular cyanobacteria, Synechococcus and Prochlorococcus, important primary producers in the oceans. Viruses belonging to this genus are divided into two distinct phylogenetic clades: clade A and clade B. These viruses have narrow host-ranges with clade A phages primarily infecting Synechococcus genotypes, while clade B phages are more diverse and can infect either Synechococcus or Prochlorococcus genotypes. Here we investigated infection properties (life history traits) and environmental abundances of these two clades of T7-like cyanophages. We show that clade A cyanophages have more rapid infection dynamics, larger burst sizes and greater virulence than clade B cyanophages. However, clade B cyanophages were at least 10-fold more abundant in all seasons, and infected more cyanobacteria, than clade A cyanophages in the Red Sea. Models predicted that steady-state cyanophage abundances, infection frequency, and virus-induced mortality, peak at intermediate virulence values. Our findings indicate that differences in infection properties are reflected in virus phylogeny at the clade level. They further indicate that infection properties, together with differences in subclade diversity and host repertoire, have important ecological consequences with the less aggressive, more diverse virus clade having greater ecological impacts.


Supplementary Methods
Cyanobacteria growth measurements Synechococcus growth was measured by chlorophyll a autofluorescence of the cells using a Turner 10AU field fluorometer at 340-500 nm excitation and >665 nm emission or a SynergyMx plate reader at 440±20 nm excitation and 680±20 nm emission. Determination of Synechococcus concentrations in cultures was done by flow cytometry using the LSR-II (BD) for cultures of WH8109 and CC9605 or with the ImageStreamX (Amnis, Seattle, WA) for WH7803 cultures since the latter strain formed aggregates in 33% of the cases (75 out of 225 cells investigated in 4 different experiments). The number of aggregates and the number of cells per aggregate were taken into consideration when determining the per cell burst size and virulence for WH7803 (see below). Flow cytometry data were collected and analyzed as described below for field samples. ImageStreamX data were analyzed using the IDEAS6.2 software.

Virus growth curve experiments and adsorption kinetics
Virus growth curve experiments [62] were used to assess the dynamics of the infection cycle. The concentration of phages in the extracellular medium is determined with time after phage addition to the host culture. The time point with the minimum concentration of free viruses is considered the time of maximal adsorption. The time point prior to the release of new phage progeny into the extracellular medium is considered the end of the latent period.
Cultures were infected with cyanophages at the mid-log growth phase of hosts (~5x10 7 to 10 8 cellsꞏml -1 ) at an MOI of 2. Samples were collected at various time points, diluted 100-fold in ASW+NO3 medium [63,64] and filtered through 0.2 μm pore-sized Acrodisc Syringe Filters (Millex-GV 33 mm 0.22 µm PVDF). The filtrate was serially diluted and cyanophages enumerated using the plaque assay (see below). Three to eight independent experiments were performed for each host-phage interaction.
Adsorption kinetics of viruses to host cells were calculated by fitting a linear regression to the ln transformed phage abundances versus time for the period between phage addition and maximal adsorption, with the slope being the rate of adsorption. Due to the rapid adsorption of S-TIP28, its adsorption kinetics was determined from the difference of the ln transformed initial abundances and abundances at maximal adsorption divided by time. Data for these calculations were obtained from the growth curve experiments shown in Fig. 1 except for the Syn5 phage. Adsorption kinetics for Syn5 was calculated from specific adsorption experiments in which samples were collected every 3 min at an MOI of 0.01.
Comparison of phage growth curves was done using a multi-level modeling approach [65,66], in which we tested the combined effect of the type of the virus and time after infection with random effect of the replicate cultures. The model was performed using the 'lme4' package [67], the pvalues were calculated using the 'broom.mixed' package [68]. The figures were produced using the 'ggplot2' package [69].

Plaque assay and plaque size estimate
The plaque assay was used to quantify the number of infective viruses in a lysate or experimental sample. Serial dilutions of the lysate were mixed with the assay culture and combined with 0.28% ultraPure low melting point (LMP) agarose (Invitrogen) in ASW+NO3 medium and poured into petri dishes with 1 mM sodium sulfite [25,55].
Plaque sizes were compared for each pair of clade A and clade B phages after growth on lawns of the same host cultures for a period of 4 or 5 days (SI Appendix, Table S3). Petri dishes were photographed and the diameter of the plaques were measured using the ImageJ Java image processing program (https://imagej.nih.gov/ij/).

Virulence and burst size measurements
Virulence and burst size was evaluated using a single cell approach [30]. Cells at 5-10x10 7 cellsꞏml -1 were infected at an MOI of 2 so that all cells would contact an infective phage. At maximum adsorption, cells were sorted by flow cytometry into individual wells in 96-well plates in PBS (137 mM NaCl, 2.7 mm KCl, 10 mM Na2HPO4, 1.8 mM KH2PO4 at pH 7.4; Dulbecco's Phosphate Buffered Saline, Biological Industries) using the FACSAria-IIIu (BD) cell sorter with an 85 µm nozzle at a pressure of 45 psi.
For the virulence assay, cells were sorted into culture-containing 96-well plates and incubated under growth conditions for up to two weeks. If the sorted cell lyses, this initiates infection and lysis of the culture in the well. The percentage of wells with a significant decrease in cyanobacterial color relative to uninfected control plates was used to score the number of lysed cells in the original sorted culture. Results are shown as % of lysed cells for phages infecting WH8109 and CC9605, and as % lysed cells per sorted event for phages infecting WH7803 since 33% of cells produced aggregates of more than 1 cell (see above). Six to 21 experiments were carried out for each hostphage interaction, with each experiment consisting of 120 cells sorted into internal wells of two 96-well plates.
For the burst size, single cells were sorted into individual wells of medium-containing 96-well plates. The plates were incubated in host growth conditions overnight, a period longer than the length of the lytic cycle, to allow lysis. The number of viruses produced by each cell was determined by plating the contents of each well on separate petri dishes using the plaque assay. Plates with single plaques were excluded from the analysis since it is not possible to differentiate between cases where a cell produced a single plaque or lysis was delayed and we plated prior to cell lysis. The results of the assay are shown as the number of viruses produced per lysed cell. Since 33% of the WH7803 cells formed aggregates, we calculated the burst size for phages on this strain by dividing the number of phages produced per sorted event by 2, the combined average number of cells in over 50 sorted events per experiment, including singlets and aggregates.

Decay measurements
The decay rate was determined for freshly lysates. Triplicate cultures were infected at a MOI of ~0.1 and phages were allowed to replicate until the host culture declined, but prior to complete lysis of the cultures. Lysates were collected after filtering the cultures through 0.2 µm PVDF filters and were incubated in glass tubes under host growth conditions: a 14/10 h light/dark cycle at a light intensity of 20 µmol photonsꞏm -2 ꞏs -1 and a temperature of 21ºC. The loss of infectivity was quantified over the period of a month by the plaque assay. Decay rates were calculated by fitting linear regressions to the ln transformed plaque abundances versus time with the slope of the regression giving the rate of decay.

Field study site
The Gulf of Aqaba is the northwest basin of the Red Sea. The deepest point reaches 1800 m. On average the Gulf of Aqaba is 14 km wide and is an evaporation basin surrounded by deserts, with almost no entry of water runoff. The only influx of surface water is through the Straits of Tiran, which are characterized by a shallow sill at a depth of 252 m. The water temperature fluctuates from ~27 °C at the surface during summer to 20.5 °C in winter, the latter of which is the water temperature below the stratified layer throughout the year. The lack of cold deep water results in annual deep winter mixing [20,35,70]. Winter mixing is followed by temperature-induced stratification, formed by a combination of increased solar radiation and pulses of warm advected water that increases the temperature in the upper 100-150 m by 1-2°C in a matter of days [70].

Physical environmental conditions
Water column conditions including temperature, salinity and pressure, were measured in-situ by the CTD instrument (SBE 19plus V2 SeaCAT Profiler). The σ density anomaly was calculated using the Ocean Data View software algorithm [71]. Here we define the maximal mixing depth as the depth down to which a uniform distribution of σ density was observed, with a difference of less than 0.01 kgꞏm -3 per meter. We present environmental data from the beginning of the calendar year in 2013 to provide information on water column conditions during mixing prior to the beginning of our sampling period in the spring. and 60 m at ~09:30 AM (UTC) which was from 10:00-12:30 local time. Samples for chlorophyll a and flow cytometry counts were filtered through a 60 µm plankton mesh while samples for virus enumeration were prefiltered through a 20 µm mesh into darkened bottles and kept at room temperature while being processed for the different procedures (see below).
Chlorophyll a concentrations were determined using a cold acetone extraction procedure [60]. Water (300 ml) was filtered onto a GF/F filter (Whatman GF/F, 25 mm), placed in 90% acetone and incubated for 24 h at 4°C in the dark. Extracted chlorophyll a concentrations were determined using a non-acidification fluorometric method with a pre-calibrated Turner TD700 fluorometer.

Enumeration of photosynthetic cells and bacteria
Samples for determining photosynthetic cells (Synechococcus, Prochlorococcus, eukaryotic phytoplankton) and bacterial abundances were fixed in 0.125% glutaraldehyde, frozen in liquid nitrogen and stored at -80°C until analysis. The samples were analyzed with the LSR-II flow cytometer (BD Biosciences) equipped with a 488 nm laser from 200 µl subsamples in two or three technical replicates. As an internal reference for size and fluorescence intensity, 1 μm diameter yellow-green microspheres (Fluoresbrite) were added. The data from the LSR-II was analyzed using FCSexpress 5.
Abundances of Synechococcus, Prochlorococcus, and photosynthetic picoeukarya were determined based on their characteristic autofluorescence and size, determined from forward scatter. Prochlorococcus and photosynthetic picoeukaryotic cells were detected by red fluorescence of chlorophyll a (emission at 692/640nm) while Synechococcus was detected by orange fluorescence of phycoerythrin (emission at 580/30nm) [72]. The photosynthetic picoeukaryotes and Prochlorococcus were distinguished from each other by their difference in forward scatter and the intensity of their red fluorescence. Quantification of Prochlorococcus is underestimated in the upper 40 m during stratification due to their very low chlorophyll a fluorescence and insufficient sensitivity of the LSR-II for their detection. Due to this underestimation of Prochlorococcus we show correlation results for depths from 60 to 140 m in the main text. We note here that a significant correlation was also found when the 0 to 40 m samples were included and was ρ=0.74 (S=28362, p<2.2e-16, n=76) for the annual comparison.
Enumeration of total bacteria was done after DNA staining with SYBR Green I (Invitrogen) at a dilution of 10 -4 of the commercial stock solution [72] and a 15 min dark incubation. Detection of heterotrophic prokaryotes was based on the green fluorescence of SYBR Green I, excited with the 488 nm laser and detected by a 530/20 nm filter. Heterotrophic bacteria were quantified as the difference between the abundances of total bacteria and cyanobacteria. Therefore, during the stratification period they will be overestimated by the same degree that Prochlorococcus is underestimated.

Quantification of virus-like-particles (VLPs)
Bulk virioplankton was determined from VLP abundances following Patel [61]. Samples were filtered on-board over a 0.2 µm syringe filter (Millex-GV 33 mm 0.22 µm PVDF) into two 50 ml falcon tubes, fixed at a final concentration of 1.6% formaldehyde (Bio Lab) (that had been 0.02 µm filtered). Samples were incubated for 20 minutes in the dark and frozen and stored until analysis at -80°C. Samples were thawed at 25°C and 2-4 ml was immediately filtered onto 0.02 µm Anodisc aluminum oxide filters (Whatman, Kent, UK), stained with SYBR Green I (Molecular Probes Inc., Eugene, OR, USA), and enumerated by epifluorescence microscopy with a Cell Observer microscope (Zeiss) equipped with an AxioCam HS video camera, filters for blue excitation (470/40nm) and green emission (525/50nm) a 1000X magnification. Signals with a radius of 3-12 pixels were considered virus-like-particles, based on a comparison with manual counts. Images of 20 fields from each sample were analyzed using OpenCFU 3.9.0 software [73].

Quantification of T7-like cyanophages
The quantification of T7-like cyanophages was done using the polony method as described previously [37]. Samples collected from Niskin bottles were filtered over a 0.2 µm syringe filter as above. The <0.2 µm filtrate, containing the viral fraction, was frozen in liquid nitrogen after 24 hours and stored at -80°C until analysis. In the lab 8.3 µl of the <0.2 µm filtrate was added to 27.5 µl gel mix containing 10% acrylamide, 20 µM acrydite primer Acr-534-Rd, 0.2% BSA, 0.1% APS, 0.1% TEMED and poured into a well in custom-made microscope slides (Thermo Fisher Scientific) and polymerized in an argon chamber for 30 min. The acrydite primer is used to anchor the amplicon to the acrylamide gel. Degenerate primers targeting the DNA polymerase gene were used to amplify T7-like cyanophages at the genus level: forward primer 341Fd-15-NNN: NNNCCNAAYYTNGSNCAR, and reverse primer Acr-534Rd: [5Acr]TGNWRYTCRTCRTGNAYRAA. The final volume of the gel on the slide is 11.6 µl. The slides were washed with MilliQ and 0.025% Tween-20 to remove the excess gel and dried in air for 20-30 minutes.
After polymerization the other PCR reagents were diffused into the gel in a 20 µl volume. These are 1X PCR buffer for Jumpstart Taq polymerase (Sigma-Aldrich), 0.25mM dNTPs, 10 µM unmodified primer 341-Fd-15-NNN, 0.1% Tween-20, 0.2% BSA, 13.4 units of Jumpstart Taq polymerase (Sigma-Aldrich)). The slides were sealed with a Secure-Seal hybridization chamber (Grace Biolabs) and filled with mineral oil. Thermocycling was carried out in a twin-tower slide thermocycler (DNA Engine with dual block slide chamber, Bio-Rad). Prior to cycling, an initial denaturation step at 94°C for 5 min is used to make the viral DNA in the capsids accessible. A total of 50 cycles are performed with denaturation at 94°C for 45 s, annealing at 50°C for 45 s and elongation at 72°C for 2 min. A final 6 min elongation step was performed at 72°C after which the slides were stored at 10 °C overnight.
Reliable quantification is achieved with at least 10 polonies per gel with a limit of accurate quantification of 1×10 4 phagesꞏml -1 . When fewer polonies are observed, sample concentration is required for accurate quantification. A 40-50-fold concentration for clade A phage quantification was done for 4 representative depth profiles using iron flocculation [37,74], decreasing the lower threshold of accurate quantification to 200 phagesꞏml -1 . FeCl3 was added to the samples at a final concentration of 0.18 mM and incubated at room temperature for at least 1 h for samples collected in April, July and August 2013 and January 2014. After incubation the samples were centrifuged for 10 min at 17,000Xg at 10°C and then resuspended in a buffer containing 0.125 M Tris, 0.1 M Na2EDTA and 0.125 M oxalic acid, pH 6.

Quantification of cyanobacterial infection by T7-like cyanophages
The quantification of percent infection by clade A and clade B T7-like cyanophages was done using the iPolony method [38]. Samples collected from Niskin bottles were filtered through a 20 µm Nylon mesh. Cells in the filtrate were fixed with 0.1% glutaraldehyde (final concentration) in the dark for 30 minutes and then frozen in liquid nitrogen and stored at -80°C until analysis. Sorting was performed on a BD Influx flow cytometer equipped with a 488 nm laser and small particle detector on one-drop purity mode. Synechococcus and Prochlorococcus were gated based on their autofluorescence properties and size as described above. Thousands of sorted Synechococcus and Prochlorococcus cells from each sample were screened for the presence of intracellular T7-like cyanophage DNA using the polony method as described above. Percent infection was calculated by dividing the number of cyanophage amplicons by the number of cells input in each reaction. Corrections for co-sorted free cyanophages and differential detection across the latent period were then applied as described previously [38].

Host-virus population modeling
In our first model we investigated the dynamics of host-virus interactions when one host genotype (HA) is infected by a clade A phage (VA) and another host genotype (HB) is infected by a clade B phage (VB) (Fig. 6a). Phages have burst sizes (β), virulence ( ), and latent periods (1/η) that are equal to the mean measured in this study for the clade A and clade B phages. For clade A these are βA = 103.8 phage•cell -1 , A = 0.64, 1/ηA = 3.3 h, while for clade B these are βB = 35 phage•cell -1 , B = 0.30, 1/ηB = 7.7 h (Table S4). Contact rate (φ) and decay coefficient (m) were assumed to be equal for both phage clades and are based on literature values. The contact rate was 10 -7 ml·virus -1 ·h -1 assuming 1 µm and 40 nm diameter for cyanobacteria and T7-like cyanophage, respectively [75] and a decay coefficient of 0.02 h -1 was used [76]. The carrying capacity (K) was set to 10 6 cells·ml -1 , a common maximal abundance of cyanobacteria in subtropical waters [20, 77] and a growth rate (r) of 0.06 h -1 was used for the cyanobacteria based on Synechococcus [78,79]. See Table S4 for definitions of all parameters and their initial values.

Eqs. 2, and 5 describe the dynamics of cells infected by clade A and clade B phages. The term
H is the production of newly infected cells and is the product of the number of infected cells, I, lysed at rates . As modeled here, the length of an infection is exponentially distributed with mean infection periods of 1/ , respectively, where is the burst size.
Eqs. 3 and 6 describe the dynamics of the free phage. The term describe the production of new phages due to the lysis of infected cells where is the burst size and is the loss of phage due to decay.
is the combined loss of viruses due to decay and infection of a new host. The terms are the loss of phage due to adsorption to uninfected cells.
The solution of the equations at steady-state for both host-virus interactions are: The ratio of clade A and clade B abundances was determined from Eqs. 9 and 12, while the relative infection levels were calculated from Eq. 8 and Eq. 11. We also used our model to predict the mortality rates of cyanobacteria caused by clade A versus clade B phages. At steady-state, the growth of the cyanobacterial host is balanced by the mortality caused by viral lysis. Consequently, the ratio of total mortality rates is the ratio of the steady-state host densities * / * , as shown below in Eq 13, because the same growth rate and carrying capacity is assumed for cyanobacterial hosts of both clade A and clade B phages.
(13) * * * * The second model experiment considers a situation whereby a single cyanobacterial genotype (HAB) is infected by both a clade A (VA) and a clade B (VB) phage (Fig. 6b).
The solution of the equations, for > , at steady-state for both host-virus interactions are: Our models are necessarily simplifications of the complexity of interactions found in nature [58]. They do not explicitly account for the influence of other mortality agents on dynamics, such as grazing and viral lysis via other cyanophage families including those with broad host-ranges, which would stabilize population oscillations around steady state solutions (Fig. S6). Other factors, such as the rate of phage decay, host growth rates, and carrying capacity also influence steadystate outcomes and oscillations of these interactions (Fig. S6). As such, we use these population dynamic models to understand the directionality of the relationship between infection physiology of cyanophages and their ecological impacts. All numerical solutions to the differential equations were calculated using the deSolve [80]. The code for these models is available at https://github.com/lindelllab/Maidanik-et-al-2021.git.        Table S4). The values changed or parameter added relative to this are shown in bold in the legend. Dynamics oscillate less with increased host growth rate (h -1 ), reduced host carrying capacity (cells ml -1 ), increased phage decay rates (h -1 ) and when the grazing rate (h -1 ) is included. Note that a clade B phage with parameter values of r = 0.06; k = 5e5; m = 0.5; g = 0.02 declines to become extinct and is not shown. The addition of the grazing terms of gH (grazing of susceptible hosts) and gI (grazing of infected hosts) changes Eq. 1, 2, 4 and 5 (see Methods) to: Host abundances (cells•ml -1 )