ABSTRACT
Influenza B virus (IBV) undergoes seasonal antigenic drift more slowly than influenza A virus, but the reasons for this difference are unclear. While the evolutionary dynamics of influenza viruses play out globally, they are fundamentally driven by mutation, reassortment, drift, and selection at the level of individual hosts. These processes have recently been described for influenza A virus, but little is known about the evolutionary dynamics of IBV during individual infections and transmission events. Here, we define the within-host evolutionary dynamics of IBV by sequencing virus populations from naturally infected individuals enrolled in a prospective, community-based cohort over 8,176 person-seasons of observation. Through analysis of high depth-of-coverage sequencing data from samples from 91 individuals with influenza B, we find that IBV accumulates lower genetic diversity than previously observed for influenza A virus during acute infections. Consistent with studies of influenza A viruses, the within-host evolution of IBVs is characterized by purifying selection and the general absence of widespread positive selection of within-host variants. Analysis of shared genetic diversity across 15 sequence-validated transmission pairs suggests that IBV experiences a tight transmission bottleneck similar to that of influenza A virus. These patterns of local-scale evolution are consistent with the lower global evolutionary rate of IBV.
IMPORTANCE The evolution of influenza virus is a significant public health problem and necessitates the annual evaluation of influenza vaccine formulation to keep pace with viral escape from herd immunity. Influenza B virus is a serious health concern for children, in particular, yet remains understudied compared to influenza A virus. Influenza B virus evolves more slowly than influenza A virus, but the factors underlying this are not completely understood. We studied how the within-host diversity of influenza B virus relates to its global evolution by sequencing viruses from a community-based cohort. We found that influenza B virus populations have lower within-host genetic diversity than influenza A virus and experience a tight genetic bottleneck during transmission. Our work provides insights into the varying dynamics of influenza viruses in human infection.
INTRODUCTION
Influenza viruses rapidly mutate and evolve through selection, genetic drift, and reassortment (1). On a global scale, influenza A virus (IAV) and influenza B virus (IBV) evolve under strong positive selection driven by pressure for escape from preexisting population immunity (2, 3). Selection of new antigenic variants contributes to reduced effectiveness of seasonal influenza vaccines, necessitating annual updates of vaccine strains (4). IAV and IBV both undergo seasonal antigenic drift and share a similar genomic architecture, but their ecology and evolution differ in important ways (5). While IBV accounts for roughly one-third of influenza’s burden of morbidity and mortality (6, 7), it circulates only in humans and seals and is considered to be a lower pandemic risk than IAV due to its limited animal reservoirs. Like IAV, there are cocirculating, antigenically distinct lineages of IBV that are included in the quadrivalent influenza vaccine. Two lineages of IBV diverged in the 1980s, B/Victoria/2/87-like and B/Yamagata/16/88-like, here referred to as B/Victoria and B/Yamagata, respectively (8).
IBV evolves more slowly than IAV on a global scale and has a lower rate of antigenic drift, but the reasons for this are poorly understood (5, 9). Similar evolutionary forces are involved in the antigenic evolution of IAV and IBV, generally characterized by nonsynonymous substitutions at antigenic sites in the surface hemagglutinin (HA) protein (10, 11) and reassortment within and between lineages (12–14). The IBV polymerase has a lower mutation rate than that of IAV (15). However, it is unclear whether the slower global evolution of IBV is driven by its lower mutation rate or other differences in selection at the global scale.
All new seasonal influenza variants are ultimately derived from de novo mutations within individual hosts (16). Therefore, understanding how new variants arise within individuals and transmit between them is essential to defining how novel viruses spread in host populations. For example, if the relative mutation rate is a major factor underlying the global evolutionary differences across IAV and IBV, we might also expect to see differences in their within-host dynamics. We and others have used next-generation sequencing (NGS) to investigate the within- and between-host evolutionary dynamics of IAV in humans (16–21). We have found that there is little accumulation of intrahost variants during acute infections of immunocompetent individuals (18, 19), and we have not found evidence of changes in intrahost diversity by vaccination status or other proxies for immunological history (19, 20, 22). The IAV transmission bottleneck is stringent (19), which generally means that few variants that arise within hosts are able to transmit. Together, these studies suggest that positive selection of novel variants is an inefficient process in IAV-infected hosts, contrasting with its patterns of significant positive selection at the global level (23). Despite the importance of intrahost processes to influenza virus evolution, these dynamics have not been systematically investigated in IBV.
Here, we use NGS to define the within-host diversity of IBV populations from individuals enrolled in the Household Influenza Vaccine Evaluation (HIVE) study, a community-based household cohort initiated in 2010. We apply a previously validated sequencing approach and bioinformatic pipeline (19, 20, 24) to identify intrahost single-nucleotide variants (iSNV) arising during infection with B/Victoria and B/Yamagata viruses. We find that IBV has significantly lower intrahost diversity than IAV, consistent with its lower mutation rate and lower rate of evolution. We analyze shared iSNV across 15 genetically validated household transmission pairs and find that like IAV, IBV is also subject to a tight genetic bottleneck at transmission. These data provide the first systematic evaluation of the genetic architecture of IBV populations during natural human infection and provide insights into the comparative epidemiology and evolution of influenza viruses.
(This article was submitted to an online preprint archive [25].)
RESULTS
We used high depth-of-coverage sequencing to define the intrahost genetic diversity in IBV-positive samples collected from individuals in the HIVE, a prospective, household cohort in southeastern Michigan that follows 200 to 350 households annually (Table 1). This cohort provides an opportunity to investigate natural infections and transmission events in a community context. Individuals that meet symptom-based criteria for an upper respiratory illness during the surveillance period undergo collection of nasal and throat swabs for molecular detection of respiratory viruses by reverse transcription-PCR (RT-PCR). Starting in 2014-2015, individuals also provided a sample collected at home prior to subsequent collection of a second specimen at the on-site clinic.
Influenza B viruses over seven seasons in a household cohort
Over seven seasons (2010-2011 through 2016-2017) and 8,176 person-seasons of observation, we identified 111 individuals infected with B/Yamagata and 67 infected with B/Victoria (Table 1). Several households had clusters of infections of two or three IBV-positive individuals within 7 days of each other, suggestive of within-household transmission. Because variant identification is sensitive to input viral titer (24), we first measured viral loads of all available IBV-positive samples by RT-quantitative PCR (RT-qPCR) (Fig. 1A). Any samples with a viral load below 103 copies/μl were not submitted for sequencing. For samples with a viral load in the range of 103 to 105 copies/μl, we performed two independent RT-PCRs and sequenced replicate libraries on separate sequencing runs. We sequenced samples with viral loads above 105 copies/μl of transport medium in a single replicate. From the available IBV-positive samples, we were able to obtain sequence data on 106 samples from 91 individuals, consisting of 35 individuals infected with B/Victoria and 56 infected with B/Yamagata (Table 1).
Viral load and sequencing coverage. (A) Boxplot of viral load (genome copies per microliter of swab transport medium, y axis) by day of sampling relative to symptom onset (x axis). The boxes display median and 25th and 75th percentiles, with whiskers extending to the most extreme point within the range of the median ± 1.5 times the interquartile range. (B) Sequencing coverage is plotted with read depth on the y axis and location within a concatenated influenza B virus genome on the x axis. The mean coverage for each sample was calculated over a sliding window of size 200 and a step size of 100. The data are displayed for all samples at each window as a boxplot, showing the median and 25th and 75th percentiles, with whiskers extending to the most extreme point within the range of the median ± 1.5 times the interquartile range; all values outside this range are shown as individual points.
We identified intrahost single nucleotide variants (iSNV) using our previously validated bioinformatic pipeline. As in our previous work, we report iSNV at frequencies of 2% or above, for which we have well-defined sensitivity and specificity (19). We consider sites with >98% frequency to be essentially fixed, setting the frequency at those sites to 100% (see Materials and Methods). We achieved a mean coverage of 10,000× per sample across most genome segments, with generally lower coverage on segments encoding nucleoprotein (NP) and nonstructural protein (NS) (Fig. 1B). We restricted our analysis of iSNV to samples with an average genome coverage of greater than 1,000×, which includes 99 of the original 106 sequenced samples.
Within-host genetic diversity of IBV in natural infections.All samples exhibited low genetic diversity. The vast majority had no iSNV above the 2% cutoff. Of the 99 samples with high-quality NGS data, 70 had no minority iSNV, 17 had one iSNV, 7 had two iSNV, and 3 samples had 3 iSNV (median, 0; interquartile range [IQR], 0 to 1) (Table 2). Two outliers had a large number of iSNV, with 8 and 20 iSNV. These two samples came from the same individual, with one collected at home and the second at the study clinic 2 days later. Most of the iSNV in these two samples were present at similar frequencies, 3 to 5% in the home sample and 17 to 23% in the clinic sample (Table 3), both of which were sequenced in duplicate on separate Illumina runs. The high number of mutations present at similar frequencies is suggestive of a mixed infection with distinct haplotypes or strains as opposed to de novo mutations arising on a single genetic background. The iSNV in the home-collected sample are all found in the subsequent clinic-collected sample, each with a similar change in frequency across the two samples. This further supports the conclusion that these mutations are on the same genome in a mixed infection with two distinct strains.
Identified iSNV, excluding samples from one putative mixed infection
Identified iSNV in one vaccinated individuala with a putative mixed infection during the 2014–2015 season
We examined how within-host diversity changes by day of sampling during IBV infections, as the virus population rapidly expands and contracts. As we have previously shown that specimen viral load can affect the sensitivity and specificity of variant identification (24), we sought to control for this variable in our analysis. Although viral load generally decreased with time after symptom onset (Fig. 1A), we found that within-host diversity as measured by number of identified minority iSNV did not vary with viral load (Fig. 2A) (P = 0.2996, adjusted r2 = 0.0009) or with day of infection (Fig. 2B) (P = 0.62, analysis of variance [ANOVA]). The frequencies of the identified iSNV were consistent across replicate libraries from the same samples, indicating that our measurements of iSNV frequency are precise (Fig. 2C).
Intrahost minority SNV by day post-symptom onset and viral load. (A) The number of minority iSNV per sample is plotted on the y axis by day post-symptom onset on the x axis. Data are displayed as boxplots representing the median and 25th and 75th percentiles, with whiskers extending to the most extreme point within the range of the median ± 1.5 times the interquartile range. The raw data points are shown in black overlaid on top of the boxplots; points from mixed infection samples are shown in orange. (B) Scatterplot relating the number of minority iSNV per sample on the y axis to the log10 of viral load, in genome copies per microliter, on the x axis. Data points from the mixed infection are shown in orange. (C) Frequency of minority iSNV in samples sequenced in duplicate. Orange dots represent variants identified in samples with a viral load of 103 to 104 genome copies per microliter, and blue dots represent variants in samples with a viral load of 104 to 105 genome copies per microliter.
We detected minority iSNV across all eight genome segments (Fig. 3). The ratio of nonsynonymous to synonymous iSNV was 0.74, which given the excess of nonsynonymous sites across the genome suggests significant purifying selection. There was only one minority iSNV present in more than one individual; we identified a variant encoding a synonymous mutation in PB1 in two individuals from separate households infected with B/Yamagata in the 2016–2017 season. We did not identify any nonsynonymous minority iSNV in the known antigenic sites of IBV hemagglutinin, which suggests that positive selective pressure for variants that escape antibody-mediated immunity is not particularly strong within hosts. We found that there is no difference in distribution of the number of iSNV per sample between vaccinated and nonvaccinated individuals (Fig. 4A). During the first few seasons of the study, some individuals received trivalent vaccines, which contain only one of the two IBV lineages. We therefore repeated this analysis, excluding three individuals for whom we had no information about specific vaccine product and reclassifying six individuals who received trivalent vaccines and were infected with a lineage not included in that season’s trivalent formulation as “unvaccinated.” We again found no difference in the number of iSNV between groups (Mann-Whitney U [MWU] test, P = 0.9103). Together, these data indicate that vaccine-induced immunity is not a major diversifying force for IBV within hosts in our study population. This is consistent with our previous work on IAV in the HIVE as well as a randomized-controlled trial of vaccine efficacy (FLU-VACS), both of which showed no difference in intrahost diversity based on same-season vaccination status (19, 20). Intrahost diversity was similar between B/Victoria and B/Yamagata virus populations (Fig. 4B), consistent with our previous comparison of subtype A/H3N2 and A/H1N1 viruses (19).
Intrahost SNV frequency by genome position and mutation type. All minority (<50%) iSNV from 99 samples are displayed with their frequency on the y axis and their position within a concatenated influenza B virus genome on the x axis. Synonymous mutations are shown in orange, and nonsynonymous mutations are shown in blue.
Intrahost SNV by vaccination status and IBV lineage. (A) Numbers of minority iSNV per sample across all 99 samples are shown (y axis) by current-season vaccination status of the host (x axis). Samples from the mixed infection are shown in orange. (B) Numbers of minority iSNV per sample are shown (y axis) by IBV lineage (x axis). Samples from the mixed infection are shown in orange. (C) Pairwise nucleotide diversity (π, y axis) by influenza virus type (x axis), stratified by iSNV frequency cutoff (top). Medians are shown as red lines. Data for influenza A virus are from 243 samples described by McCrone et al. (18). Data on influenza B virus are from 97 high-quality samples in the present study. Samples from mixed infections in both studies are excluded. (D) Numbers of minority iSNV in 43 of the 99 high-quality samples (y axis), consisting of B/Yamagata from the 2014–2015 season, B/Victoria from the 2015–2016 season, and B/Yamagata from the 2016–2017 season based on alignments to the original references from the 2012–2013 season versus season-matched reference genomes (x axis).
We compared the within-host genetic diversity of IBV to our previously published data on IAV from the HIVE cohort (19). Here, IBV exhibits lower within-host diversity than IAV (Table 4) (P < 0.001, MWU test). IBV samples had a lower median number of minority iSNV (median, 0; IQR, 0 to 1) than IAV samples (median, 2; IQR, 1 to 3); 71% of IBV samples contained no minority iSNV, compared to 20% of IAV samples. The difference in within-host diversity was robust to relaxation of the iSNV frequency cutoff to 1% and 0.5%. IBV also exhibited a lower within-host diversity, as measured by nucleotide diversity (π), which takes both the number and the frequency of iSNV into account (Fig. 4C and Table 4). To ensure that our results were not an artifact of overly stringent quality thresholds, we also identified minority iSNV with less conservative read mapping quality (MapQ) and base quality (Phred) scores. We identified the same set of minority iSNV with a MapQ cutoff of 20 as with the original cutoff of 30. Similarly, reduction of the Phred base-quality cutoff to >25 in addition to a MapQ score cutoff of >20 resulted in only 20 more minority iSNV, 8 of which were found in the individual with a mixed infection. The other additional 12 minority iSNV were dispersed across specimens and did not significantly change the overall distribution of within-host diversity. We also examined whether our results were biased by use of a single B/Yamagata and B/Victoria reference for alignment and variant calling, which were both drawn from the 2012-2013 season (see Materials and Methods). We realigned sequence data from 43 of the original 99 samples to season-specific reference genomes isolated in southeastern Michigan. We found that the overall alignment rate for any given specimen was similar between the original reference and the new season-matched reference. Variant identification based on the new references and the original quality thresholds resulted in the same distribution of within-host diversity, although the identities of some iSNV were different (Fig. 4D).
Within-host diversity of IAV versus IBV
Together, these results indicate that our measurements of within-host diversity are robust to several technical aspects of variant identification, which are unlikely to account for the lower observed diversity of IBV. Because these data are from the same cohort and were generated using the same sequencing approach and analytic pipeline as our previous IAV data sets, the observed differences likely reflect true biological differences between IAV and IBV.
Identification of household transmission pairs.We compared viral diversity across samples from individuals in the same household to investigate the genetic bottleneck that influenza B viruses experience during natural transmission. Over the seven influenza seasons, 39 households in the HIVE cohort had two or more individuals positive for the same IBV lineage within a 7-day interval (Table 1). This epidemiologic linkage is suggestive of transmission events but does not rule out coincident community-acquired infection (19). We identified 16 putative transmission pairs for which we sequenced at least one sample from each individual. In one of these pairs, the putative recipient was the individual with a mixed infection. The donor did not have evidence of a mixed infection based on number of iSNV, which would imply that the recipient may have been infected twice or that the second virus was lost from the donor by the time of sampling. This pair was excluded from the between-host analysis, leaving 15 putative transmission pairs for which we have high-quality sequencing data on both donor and recipient influenza populations.
We used our sequencing data to determine which of these epidemiologically linked household pairs were actual IBV transmission pairs. We generated maximum likelihood phylogenetic trees for samples from the two IBV lineages using the concatenated coding consensus sequences. Phylogenetic analysis provided genetic evidence that the 15 epidemiologically linked pairs were indeed true transmission pairs, as epidemiologically linked pairs were found nearest each other in each tree (Fig. 5A and B; vertical bars with household identifiers [ID]). We also validated these transmission pairs by analyzing the genetic distance across viral populations. True transmission pairs should have genetically similar populations exhibiting low genetic distance, while individuals with coincident community acquisition are more likely to have populations with a higher genetic distance. We compared the genetic distance between epidemiologically linked household pairs and random community pairs from the same season and infected with the same IBV lineage, using L1-norm as the measurement of genetic distance (Fig. 5C). The distribution of random community pairs functions as a null model of genetic distances among locally circulating strains. All of the 15 putative transmission pairs fell on the tail of this distribution, below the 5th percentile of the community pair L1-norm distribution, indicating that they are true transmission pairs (Fig. 5C). While the L1-norm is a function of both the consensus sequence and the iSNV, this signal was driven predominantly by consensus differences, as reflected in the phylogenetic analysis.
Identification of household transmission pairs. Maximum likelihood phylogenetic tree of all B/Victoria (A) and B/Yamagata (B) samples from this study. Concatenated consensus coding sequences were aligned with MUSCLE, and phylogenetic trees were constructed with RAxML. Tip labels are denoted as enrollee ID, household ID, season, and lineage, separated by underscores; tip labels are colored by season. (C) Histogram of genetic distance, as measured by L1-norm, between household pairs and random community pairs from the same season and lineage. The bar heights for each group are normalized to the maximum for each group for comparison. Community pairs are shown in orange, and household pairs are shown in blue. The dotted red line indicates the 5th percentile of the community pair distribution.
Comparison of viral diversity across transmission pairs.Transmission bottlenecks restrict the genetic diversity that is passed between hosts. With a loose transmission bottleneck, many unique genomes will be passed from donor to recipient. Because this will allow two variants at a given site to be transmitted, sites that are polymorphic in the donor are more likely to be polymorphic in the recipient. However, in the case of a tight or stringent bottleneck, sites that are polymorphic in the donor will likely be either fixed or absent in the recipient. We have previously demonstrated that influenza A virus experiences a tight transmission bottleneck of 1 to 2 unique genomes (19). Across our 15 IBV transmission pairs, we found no sites that were polymorphic in the donor and the recipient (Fig. 6). iSNV present in the donor were either fixed (100%) or absent (0%) in the recipient. These data suggest a stringent transmission bottleneck for influenza B virus, similar to that of influenza A virus. As there were fewer samples, transmission pairs, and iSNV in our IBV data set, we were unable to obtain a robust and precise estimate of bottleneck size.
Shared diversity across household transmission pairs with influenza B virus. iSNV for 15 validated transmission pairs using samples closest to the time of transmission (inferred based on day of symptom onset) are shown. Each iSNV is plotted as a point, with its frequency in the recipient (y axis) versus its frequency in the donor (x axis).
DISCUSSION
Here, we define the within-host genetic diversity of IBV in natural infections by sequencing 106 samples collected over 8,176 person-seasons of observation in a household cohort. Because the HIVE study prospectively identifies individuals with acute respiratory illness regardless of severity, these samples capture IBV dynamics in a natural setting, reflective of infections occurring in the community. We show that within-host diversity of IBV is remarkably low, with most samples displaying no intrahost variants above our level of detection. We also find that IBV experiences a tight transmission bottleneck, limiting the diversity that is passed between hosts. IBV exhibits significantly lower within-host diversity than IAV. These findings reflect the lower relative evolutionary rate of IBV compared to IAV.
Our findings are largely consistent with what has been observed in IAV infections in humans (17–20). We found that only a minority of samples contain iSNV, the majority of which encode synonymous changes, consistent with a predominance of purifying selection within hosts. If immune-driven selective pressures were sufficiently strong to drive positive selection of antigenic variants at the individual level, we would expect to see enrichment of variants in antigenic regions. However, variants were no more common in the antigenic proteins hemagglutinin and neuraminidase, and we found no intrahost variants in known antigenic regions of hemagglutinin. We also found that the extent of within-host diversity did not vary with current-season vaccination status, further suggesting that immune selection is not particularly strong within hosts (19, 20, 22). Our data suggest that selective sweeps occur infrequently at the individual level, with selection evident only over a broader scale of time and space (19, 26). We recognize, however, that it is possible for individual level selective pressure to vary in magnitude by age, locale, influenza infection history, or immune status (27).
We do find that there are important differences in the within-host evolution of IAV and IBV. IBV displays significantly lower within-host diversity than IAV. Since measurements of within-host diversity can vary based on host population, sequencing approach, and variant calling algorithm (28), a strength of our study is that our comparison is based on samples from the same cohort with the same sequencing approach and analytic pipeline. In both of our studies, we have sequenced swab samples directly without prior culture, accounted for the confounding effect of viral load, and used a standardized, empirically validated analytic pipeline for variant identification (24). The only difference in methodology between these two studies is the multiplex amplification primers. In both viruses, these primers target the highly conserved ends of influenza virus genome segments, making it unlikely that this factor would drastically alter the amplification efficiency of within-host variants. Our analytic pipeline includes rigorous quality criteria to reduce false positives that can be introduced by amplification and Illumina sequencing. Importantly, these empirical quality criteria did not mask diversity actually present in these samples, strengthening the conclusion that IBV exhibits lower within-host diversity than IAV.
The most likely biological explanation for IBV’s lower within-host diversity is its de novo mutation rate, which is thought to be at least 2-fold lower than that of IAV (15). Viral mutation rates are critical to the diversification of rapidly evolving viruses within hosts. Under a neutral model, the number and frequency of minority variants are dependent on the mutation rate and demographics of the population (16). In such a model, the expected number of variants is highly sensitive to variation in the mutation rate across the range commonly estimated in RNA viruses. In light of our results, a more thorough comparison of mutation rates across influenza viruses is needed.
Another possible factor underlying IBV’s reduced diversity is the mutational robustness of the IBV genome relative to that of IAV. If IBV were less robust to mutation, stronger negative selection on multiple genes in IBV could result in more limited within-host diversity, perhaps located to certain regions of the genome. However, we found that the distributions of iSNV across IAV and IBV genomes are relatively similar. Furthermore, we have previously shown that the distribution of mutational fitness effects in influenza A/WSN/33/H1N1 matches that of other RNA and single-stranded DNA (ssDNA) viruses (29). Given that viruses across families with vastly different genomic architecture have similar mutational robustness, this is unlikely to account for the differences in within-host diversity between IAV and IBV.
We find that IBV experiences a stringent genetic bottleneck between hosts. A stringent transmission bottleneck places a constraint on the rate of adaptation of viral populations within and between individual hosts. Population bottlenecks reduce the effective population size, which increases random genetic drift and decreases the efficiency of selection (30). This results in a reduced ability of selection to fix beneficial mutations and to remove deleterious ones, which can decrease population fitness. However, there are potential evolutionary advantages to stringent bottlenecks, including removal of defective interfering particles (31, 32). While we were not able to estimate the size of the transmission bottleneck for IBV as precisely as for IAV, it is likely that the bottleneck sizes are comparable across the two viruses given the similarities in their transmission routes and ecology in the human population. Data from many more transmission pairs are necessary for a more robust estimate.
Together, our results are consistent with the lower rate of global evolution observed in IBV lineages than that in both seasonal A/H1N1 and A/H3N2 (10, 12, 14, 33). We suggest that a lower intrinsic mutation rate leads to reduced within-host diversity. With a comparably tight bottleneck, fewer de novo variants will rise to a level where they can be transmitted and spread through host populations. Combined with a lower incidence of IBV than of IAV, this would result in fewer variants that eventually spread and influence global dynamics. However, further investigation in larger populations are required to evaluate the within-host dynamics of both types of seasonal influenza viruses and how they contribute to larger-scale evolutionary patterns.
MATERIALS AND METHODS
Description of the HIVE cohort.The HIVE study is a prospective, community-based household cohort in southeastern Michigan based at the University of Michigan School of Public Health (34–39). The cohort was initiated in 2010, with enrollment of households with children occurring on an annual basis and an active surveillance period lasting from October through May. In 2014, active surveillance was expanded to take place year-round. Participating adults provided informed consent for themselves and their children, and children ages 7 to 17 years provided oral assent. Individuals in each household were followed prospectively for acute respiratory illness, defined as two or more of the following: cough, fever or feverishness, nasal congestion, chills, headache, body aches, or sore throat. Study participants meeting the criteria for acute respiratory illness attended a study research clinic at the University of Michigan School of Public Health, where a combined throat and nasal swab, or a nasal swab only for children less than 3 years old, was collected by the study team. Beginning in the 2014-2015 season, study participants with acute respiratory illnesses took an additional nasal swab at home at the time of illness onset, collected either by themselves or by a parent. The study was approved by the Institutional Review Board of the University of Michigan Medical School.
Viral detection, lineage typing, and viral load quantification.We processed upper respiratory specimens (combined nasal and throat swab or nasal swab) for confirmation of influenza virus infection by reverse transcription-PCR (RT-PCR). We extracted viral RNA with either QIAamp viral RNA mini kits (Qiagen) or PureLink Pro 96 viral RNA/DNA purification kits (Invitrogen) and tested samples using the SuperScript III Platinum one-step quantitative RT-PCR system with ROX (Invitrogen) and primers and probes for universal detection of influenza A and B viruses (CDC protocol, 28 April 2009). Specimens positive for influenza virus were tested using subtype/lineage primer and probe sets, which are designed to detect influenza A (H3N2), A (H1N1)pdm09, B (Yamagata), and B (Victoria) viruses. An RNAseP primer/probe set was run for each specimen to confirm specimen quality and successful RNA extraction.
We quantified the viral load in each sample by RT-qPCR using primers specific for the open reading frame of segment 8 (NS1/NEP): forward primer 5′-TCCTCAACTCACTCTTCGAGCG-3′, reverse primer 5′-CGGTGCTCTTGACCAAATTGG-3′, and probe 5′-(FAM)-CCAATTCGAGCAGCTGAAACTGCGGTG-(BHQ1)-3′. Each reaction mixture contained 5.4 μl of nuclease-free water, 0.5 μl of each primer at 50 μM, 0.1 μl of ROX dye, 0.5 μl SuperScript III RT/Platinum Taq enzyme mix, 0.5 μl of 10 μM probe, 12.5 μl of 2× PCR buffer master mix, and 5 μl of extracted viral RNA. To relate genome copy number to threshold cycle (CT) value, we used a standard curve based on serial dilutions of a plasmid control, run in duplicate on the same plate.
Amplification, library preparation, and sequencing.We amplified viral cDNA from all eight genomic segments using the SuperScript III one-step RT-PCR Platinum Taq high-fidelity kit (Invitrogen). Each reaction mixture contained 5 μl of extracted viral RNA, 12.5 μl of 2× PCR buffer, 2 μl of primer cocktail, 0.5 μl of enzyme mix, and 5 μl of nuclease-free water. The primer cocktail was a mixture of B-PBs-UniF, B-PBs-UniR, B-PA-UniF, B-PA-UniR, B-HANA-UniF, B-HANA-UniR, B-NP-UniF, B-NP-UniR, B-M-Uni3F, B-Mg-Uni3F, B-M-Uni3R, B-NS-Uni3F, and B-NS-Uni3R (sequences and proportions are listed in reference 40). The thermocycler protocol was as follows: 45°C for 60 min, 55°C for 30 min, and 94°C for 2 min, followed by 5 cycles of 94°C for 20 s, 40°C for 30 s, and 68°C for 3 min 30 s, followed by 40 cycles of 94°C for 20 s, 58°C for 30 s, and 68°C for 3 min 30 s, and a final extension of 68°C for 10 min. We confirmed IBV genome amplification by gel electrophoresis. We sheared amplified cDNA (100 to 500 ng) on a Covaris ultrasonicator with the following settings: time, 80 s; duty cycle, 10%; intensity, 4; cycles per burst, 200. We prepared sequencing libraries with NEBNext Ultra DNA library prep kits (NEB) and sequenced them on an Illumina NextSeq with 2 × 150 paired end reads (mid-output run, v2 chemistry). To increase the specificity of variant identification, samples with a viral load between 103 and 105 genome copies/μl of transport medium were amplified and sequenced in duplicate. Samples amplified from B/Victoria and B/Yamagata plasmid clones were included on each sequencing run to account for sequencing errors. The plasmids used in the control reactions were generated by segment-specific RT-PCR from clinical samples of B/Victoria and B/Yamagata strains from the 2012-2013 season, followed by gel extraction and TOPO-TA cloning (Invitrogen). The sequence of each plasmid was determined by Sanger sequencing. We generated the plasmid control amplicons included on each Illumina sequencing run by using the same multiplex amplification protocol but with cloned plasmid DNA as the template.
Identification of iSNV.Intrahost single-nucleotide variants (iSNV) were identified using a previously described analytic pipeline (24). We identified iSNV in samples that had an average genome coverage greater than 1,000× and a viral load greater than 103 genome copies/μl of transport medium in the original sample. Sequencing adapters were removed with Cutadapt (41), and reads were aligned to the sequences derived from the B/Victoria and B/Yamagata plasmid controls with Bowtie2 (42). Duplicate reads were marked and removed with Picard and samtools (43). Putative variants were identified with the R package deepSNV using data from the clonal plasmid controls of each sequencing run (44). Minority iSNV (<50% frequency) were identified using the following empirically derived criteria: deepSNV P value, <0.01; average mapping quality, >30; average Phred score, >35; and average read position, in the middle 50% (positions 37 and 113 for 150-bp reads). For samples processed in duplicate, we used only variants that were present in both replicates; the frequency of the variant in the replicate with greater coverage at that site was used. Lastly, variants with a frequency of <2%, which have a higher false-positive rate from RT-PCR and/or sequencing errors, were not included in downstream analyses.
In our previous work on IAV, we found that there were multiple sites with mutations that were essentially fixed (>0.95) relative to those of the plasmid control and in which the base in the plasmid control was therefore identified as a minority variant in the sample (19). At these sites, deepSNV is unable to estimate the base-specific error rate and cannot distinguish true minority iSNV; however, we found that we could accurately identify minority variants at these sites at a frequency of 2% or above (19). This frequency threshold was incorporated into the pipeline for iSNV identification at these sites. Therefore, we identify intrahost variants with frequencies between 2 and 98%. Minority iSNV are the subset of these variants with a frequency between 2 and 50% relative to the sample consensus, which we use as a metric of within-host diversity. For each minority iSNV, we identify the majority iSNV present at a frequency of 50 to 98%. Any sites that were monomorphic after application of quality filters were assigned a frequency of 100%. Nucleotide diversity (π) was calculated using identified iSNV in each sample by the formula described by Zhao and Illingworth (45).
Data and code availability.Raw sequence data, with human content filtered out, are available at the NCBI Sequence Read Archive under BioProject accession number PRJNA561158. Code for the variant identification pipeline is available at http://github.com/lauringlab/variant_pipeline. Analysis code is available at http://github.com/lauringlab/Host_level_IBV_evolution.
ACKNOWLEDGMENTS
We acknowledge the individuals enrolled in the HIVE study for their participation. We thank Maria Virgilio for technical assistance.
This work was supported by grants R01 AI 118886 and NIH R21 AI 141832 from the NIH and a Burroughs Wellcome Fund Investigator in the Pathogenesis of Infectious Diseases award to A.S.L. The HIVE cohort was supported by grant U01 IP 001034 from the CDC to A.S.M. and E.T.M. and by grant U01 IP000474 from the CDC and grant R01 097150 from the NIH to A.S.M. A.L.V. was supported by grant T32 GM 007863 from the NIH.
A.S.M. reports personal fees from Sanofi and personal fees from Seqirus. E.T.M. reports consulting fees from Pfizer and research support from Merck. A.S.L. reports consulting fees from Sanofi and Roche. In all cases, these are outside the submitted work.
FOOTNOTES
- Received 6 October 2019.
- Accepted 29 November 2019.
- Accepted manuscript posted online 4 December 2019.
- Copyright © 2020 American Society for Microbiology.
