Within-host evolutionary dynamics of seasonal and pandemic human influenza A viruses in young children

  1. Alvin X Han  Is a corresponding author
  2. Zandra C Felix Garza
  3. Matthijs RA Welkers
  4. René M Vigeveno
  5. Nhu Duong Tran
  6. Thi Quynh Mai Le
  7. Thai Pham Quang
  8. Dinh Thoang Dang
  9. Thi Ngoc Anh Tran
  10. Manh Tuan Ha
  11. Thanh Hung Nguyen
  12. Quoc Thinh Le
  13. Thanh Hai Le
  14. Thi Bich Ngoc Hoang
  15. Kulkanya Chokephaibulkit
  16. Pilaipan Puthavathana
  17. Van Vinh Chau Nguyen
  18. My Ngoc Nghiem
  19. Van Kinh Nguyen
  20. Tuyet Trinh Dao
  21. Tinh Hien Tran
  22. Heiman FL Wertheim
  23. Peter W Horby
  24. Annette Fox
  25. H Rogier van Doorn
  26. Dirk Eggink
  27. Menno D de Jong
  28. Colin A Russell  Is a corresponding author
  1. Department of Medical Microbiology & Infection Prevention, Amsterdam University Medical Center, Netherlands
  2. National Institute of Hygiene and Epidemiology, Viet Nam
  3. Ha Nam Centre for Disease Control, Viet Nam
  4. Children's Hospital 2, Viet Nam
  5. Children's Hospital 1, Viet Nam
  6. Vietnam National Children's Hospital, Viet Nam
  7. Siriraj Hospital, Mahidol University, Thailand
  8. Hospital for Tropical Diseases, Viet Nam
  9. National Hospital for Tropical Diseases, Viet Nam
  10. Oxford University Clinical Research Unit, Viet Nam
  11. Radboud Medical Centre, Radboud University, Netherlands
  12. Nuffield Department of Medicine, University of Oxford, United Kingdom
  13. Peter Doherty Institute for Infection and Immunity, University of Melbourne, Australia
  14. WHO Collaborating Centre for Reference and Research on Influenza, Australia
  15. Centre for Infectious Disease Control, National Institute for Public Health and the Environment, Netherlands
8 figures, 2 tables and 6 additional files

Figures

Figure 1 with 5 supplements
Genetic diversity of within-host influenza A virus populations.

Box plots summarizing the number of intra-host single-nucleotide variants (iSNVs; median, interquartile range [IQR], and whiskers extending within median ±1.5 × IQR) identified in samples with adequate breadth of coverage across the whole influenza virus genome in (A) seasonal A/H3N2 and (B) pandemic A/H1N1pdm09 virus samples, stratified by day(s) since symptom onset or qPCR cycle threshold (Ct) values. (C, D) Ct values as a function of day(s) since symptom onset for A/H3N2 viruses (C) and A/H1N1pdm09 viruses (D).

Figure 1—figure supplement 1
Sequence coverage across all influenza gene segments and samples.

Black line plots the mean coverage for a sliding window of 50 base pairs (stepsize = 25 base pairs). The interquartile range is shaded in dark pink while the full range is denoted in light pink. (A) H3N2; (B) H1N1pdm09.

Figure 1—figure supplement 2
Frequencies of nucleotide variants found in A/H3N2 viral reads sequenced from overlapping amplicons.

Each circle represents a nucleotide variant site (with frequency estimated between 0.02 and 0.98) found in reads attributed to at least two different amplicons (at least 100× coverage for each amplicon) and is colored by the cycle threshold (CT) value of the sample from which the variant was found. Scatter plot on the left panel compares the variant frequencies between any two amplicons while the plot on the right panel compares the variant frequencies of each amplicon to that when combining across all overlapping amplicons (i.e., the frequencies used for main analyses). The dashed line is the one-to-one expected value.

Figure 1—figure supplement 3
Maximum-likelihood phylogeny of putative majority (consensus) and minority whole-genome sequences (by concatenating all eight gene segments) of A/H3N2 virus samples.

Tip names are given in the format: ‘Patient ID_Days since symptom onset_putative consensus or minority sequence.’ The tree is rooted to the A/Brisbane/10/2007 virus (H3N2_Bris07; EPI_ISL_103644). Subject 1673 (green tips) and the D8 sample of subject 1878 (pink tips) might have arose from mixed infections or were contaminated by other strains.

Figure 1—figure supplement 4
Maximum-likelihood phylogeny of putative majority (consensus) and minority whole-genome sequences (by concatenating all eight gene segments) of H1N1pdm09 virus samples.

Tip names are given in the format: ‘Patient ID_Days since symptom onset_putative consensus or minority sequence.’ The tree is rooted to the A/California/04/2009 virus (H1N1pdm09_Cali09; EPI_ISL_376192). Encircled tips denote the consensus majority sequence of the sample.

Figure 1—figure supplement 5
Pearson’s correlation between the first day of oseltamivir treatment administered to patients and the last day on which viral samples with cycle threshold (CT) values ≤35 were collected.

Time is measured by number of days since symptom onset. Each point represents a patient included in this study who was treated with oseltamivir (Supplementary file 4).

Figure 2 with 4 supplements
Box plots (median, interquartile range [IQR], and whiskers extending within median ±1.5 × IQR) summarizing the empirical within-host evolutionary rates of (A) seasonal A/H3N2 viruses and (B) pandemic A/H1N1pdm09 viruses.

Top panel shows the evolutionary rate of individual gene segments over all timepoints (rg) while the bottom panel depicts the genome-wide evolutionary rate (rt) for each day since symptom onset. All rates are stratified by substitution type (synonymous – blue; nonsynonymous – red; gray – stop codon). Wilcoxon signed-rank tests were performed to assess if the paired synonymous and nonsynonymous evolutionary rates are significantly distinct per individual gene segment or timepoint (annotated with * if p<0.05). This was done for all sets of nonsynonymous and synonymous rate pairs except for those computed per day since symptom onset for A/H1N1pdm09 viruses due to the low number of data points available (median number of A/H1N1pdm09 virus samples collected per day since symptom onset = 2). Note that the scales of the y axes differ between (A) and (B) to better show rate trends.

Figure 2—figure supplement 1
Box plots (median, interquartile range [IQR], and whiskers extending within median ±1.5 × IQR) summarizing the empirical within-host evolutionary rates (rg,t) of different H3N2 viral gene segments.

All rates are stratified by substitution type (synonymous – blue; nonsynonymous – red; stop codon – gray). Wilcoxon signed-rank tests were performed to assess if the paired synonymous and nonsynonymous evolutionary rates are significantly distinct per timepoint (annotated with * if p<0.05).

Figure 2—figure supplement 2
Box plots (median, interquartile range [IQR], and whiskers extending within median ±1.5 × IQR) summarizing the empirical within-host evolutionary rates (rg,t) of different H1N1pdm09 viral gene segments.

All rates are stratified by substitution type (synonymous – blue; nonsynonymous – red; stop codon – gray). Wilcoxon signed-rank test was not performed here due to the low number of samples collected (i.e., median number of samples per day post-illness onset = 2).

Figure 2—figure supplement 3
Linear regression of within-host synonymous and nonsynonymous evolutionary rates of within-host A/H3N2 virus samples.

Each plotted line is the linearly regressed line to the evolutionary rates computed for each A/H3N2-infected individual. Based on our findings, we expect that synonymous rates correlate negatively with time while nonsynonymous rates have a positive temporal correlation. Colored lines represent those that fall within this expectation while dashed black lines represent those that did not.

Figure 2—figure supplement 4
Box plots (median, interquartile range [IQR], and whiskers extending within median ±1.5 × IQR) summarizing the empirical daily within-host evolutionary rates of seasonal A/H3N2 viruses.

Variants that could potentially be PCR artifacts were removed (i.e., those found under the 75th percentile [6.3%] of frequency range of variants located in overlapping amplicons but were only detected in one amplicon, see Figure 1—figure supplement 2). All rates are stratified by substitution type (synonymous – blue; nonsynonymous – red; gray – stop codon). Wilcoxon signed-rank tests were performed to assess if the paired synonymous and nonsynonymous evolutionary rates are significantly distinct per individual gene segment or timepoint (annotated with * if p<0.05). This was done for all sets of nonsynonymous and synonymous rate pairs computed between days 3 and 9 since symptom onset.

Figure 3 with 2 supplements
Histogram of the mean number of minority intra-host single-nucleotide variants (iSNVs) identified per sample across all.

(A) A/H3N2 and (B) A/H1N1pdm09 virus specimens, sorted by frequency bins of 5% and substitution type (synonymous – blue; nonsynonymous – red; stop-codon – gray).

Figure 3—figure supplement 1
Histogram of the mean number of minority intra-host single-nucleotide variants (iSNVs) identified per sample across all H3N2 viral gene segments across all samples sorted by frequency bins of 5% and substitution type (synonymous – blue; nonsynonymous – red; gray – stop codon).
Figure 3—figure supplement 2
Histogram of the mean number of minority intra-host single-nucleotide variants (iSNVs) identified across all H1N1pdm09 viral gene segments across all samples sorted by frequency bins of 5% and substitution type (synonymous – blue; nonsynonymous – red; stop-codon – gray).
Figure 4 with 9 supplements
Intra-host single-nucleotide variants in within-host IAV populations.

(A) Breakdown of intra-host single-nucleotide variants (iSNVs) identified in seasonal A/H3N2 virus samples. The top panels plot the nucleotide positions where iSNVs were found in at least two subjects. The bottom panels show the frequencies at which iSNVs were identified. For sites with iSNVs that were found in two or more subjects, the interquartile ranges of variant frequencies are plotted as vertical lines and the median frequencies are marked with a dash. If the iSNV was only found in one subject, its corresponding frequency is plotted as a circle. All iSNVs are stratified to either synonymous (blue), nonsynonymous (red), or stop codon (gray) variants. Only the nonsynonymous variants are plotted if both types of variants are found in a site. Positions of antigenic sites of the hemagglutinin (HA)_ gene segment (Igarashi et al., 2010) are marked in green on the top panels. (B) Box plots of the frequencies of synonymous and nonsynonymous variants between antigenic and non-antigenic sites of seasonal A/H3N2 HA gene segment. (C) Similar plots to (A) for iSNVs found in pandemic A/H1N1pdm09 virus samples. (D) Similar plots to (B) for HA iSNVs identified in the pandemic A/H1N1pdm09 virus samples.

Figure 4—figure supplement 1
Plots of intra-host hemagglutinin (HA) amino acid variants in A/H3N2-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Different canonical antigenic sites of the HA protein are colored (HA numbering based on H3 numbering without signal peptide). Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced. The variant frequencies of all putative antigenic sites are also plotted.

Figure 4—figure supplement 2
Plots of intra-host neuraminidase (NA) amino acid variants in A/H3N2-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced. The first day of oseltamivir treatment for individuals with resistance mutations is annotated below the x-axis.

Figure 4—figure supplement 3
Plots of intra-host nucleoprotein (NP) amino acid variants in A/H3N2-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced.

Figure 4—figure supplement 4
Plots of M2 protein intra-host amino acid variants in A/H3N2-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced.

Figure 4—figure supplement 5
Plots of hemagglutinin (HA) intra-host amino acid variants in A/H1N1pdm09-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Different canonical antigenic sites of the HA protein are colored (HA numbering based on H3 numbering without signal peptide). Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced. The variant frequencies of all putative antigenic sites are also plotted. The frequencies of the five putative HA antigenic variants of A/H1N1pdm09 viruses are marked by arrows for better clarity.

Figure 4—figure supplement 6
Plots of neuraminidase (NA) intra-host amino acid variants in A/H1N1pdm09-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced. The first day of oseltamivir treatment for individuals with resistance mutations is annotated below the x-axis.

Figure 4—figure supplement 7
Plots of M2 protein intra-host amino acid variants in A/H1N1pdm09-infected individuals.

Top panel shows the number of subjects where nonsynonymous variants were found in the respective protein site. Bottom panel plots selected as well as parallel amino acid mutations found in multiple patients against days since illness onset. Filled circles represent days on which samples were collected and sequenced.

Figure 4—figure supplement 8
Frequency distributions of intra-host single-nucleotide variants (iSNVs) below the 2% variant calling threshold found in nucleotide positions NP-1150 and M-917 that encode for amino acid sites NP-384 and M2-77, respectively.

All A/H3N2 virus samples collected from all patients with site coverage above the 100× are included. The distributions were compared to that of neighboring sites, ±10 nucleotide positions adjacent to NP1150 and M-917.

Figure 4—figure supplement 9
Plots of within-host recurring A/H3N2 amino acid variants NP-G384R and M2-R77* based on variant calls and frequencies after remapping sample reads to their respective sample consensus sequence.
The trimeric and monomeric crystal structures of nucleoprotein (PDB: 3ZDP) (Chenavas et al., 2013) of influenza A viruses.

Amino acid sites with potentially linked epistatic amino acid variants as tabulated in Table 1 are separately colored, with their corresponding positions annotated on the monomeric structure.

Figure 6 with 1 supplement
Evolutionary rates computed from forward-time Monte Carlo within-host simulations for different mean deleterious effects (sd-) of nonsynonymous mutations.

We assumed that synonymous mutations are neutral for all simulations. (A) Neutral expectation where all nonsynonymous mutations are neutral (fneu,NS=100%). We tested our hypotheses where the majority of nonsynonymous mutations are non-neutral (fneu,NS=1%) and they are either (B) strongly (sd-=-0.1) or (C) weakly (sd-=-0.01) deleterious.

Figure 6—figure supplement 1
Number of virions (N) against replicative generation (t) based on a target cell-limited within-host model.

Blue line with markers denotes the population size computed from the model. When N>107 virions, we assumed that N remained constant at 107 (pink dashed line) to reduce computational costs of forward-time simulations.

Appendix 1—figure 1
Genetic diversity of within-host influenza A virus populations as estimated by nucleotide diversity π statistic.

Box plots summarizing the π statistic (intra-host single-nucleotide variants [iSNVs]; median, interquartile range [IQR], and whiskers extending within median ±1.5 × IQR) computed for samples with adequate breadth of coverage across the whole influenza genome. (A) Seasonal A/H3N2 and (B) pandemic A/H1N1pdm09 viruses. All box plots are either stratified by day(s) since symptom onset or qPCR cycle threshold (Ct) values.

Appendix 1—figure 2
A/H1N1pdm09 virus household transmission pairs.

(A) Schematic of A/H1N1pdm09 virus household transmission pairs identified by epidemiological linkage and plotted based on timing of sample collection (left panel) Intra-host single-nucleotide variant (iSNV) frequencies found in the donors and recipients of the five transmission pairs (right panel). (B) Violin plots of L1-norm pairwise genetic distance per site between different A/H1N1pdm09 virus sample pairs (each circle = 1 pair of virus samples). Transmission pairs are those represented in Figure 5A. Longitudinal pairs are made up of sample pairs collected from the same individual on the first and any other later timepoints on which the patient was sampled. These pairs are stratified by whether they were collected from households located in the same community (i.e., household) or combined with the rest of the analyzed A/H1N1pdm09 virus samples collected from hospitals (i.e., all). For the same aforementioned categories, we also plotted the distribution of L1-norm distances for pairs of viruses collected from different individuals. All circles are colored by the difference in days between which the sample pairs were collected (ΔT). All p-values reported are based on Mann–Whitney U tests, which were used to determine if the L1-norm genetic distance distributions of the two categories marked by the ends of the horizontal line above are statistically distinct.

Tables

Table 1
FoldX stability predictions of likely linked nonsynonymous minority variants found in A/H3N2 nucleoprotein.

The mean ΔΔG and standard deviation (SD) values reported are based on the results of five distinct simulations. Variants with mean ΔΔG<-0.46 kcal/mol are deemed to be stabilizing while destabilizing mutants were estimated to yield ΔΔG>0.46 kcal/mol.

ΔΔG(kcal/mol)
VariantsMeanSD
G384R−3.840.06
M426I2.610.05
G384R,M426I−0.420.06
G102R4.870.00
G384R,G102R0.760.09
A493T11.960.30
G384R,A493T5.560.19
V197I−3.110.02
S353Y−1.970.68
V197I,S353Y−4.480.14
Table 2
Parameter values used in the within-host model.
ParameterMeaningValue (units)Source
-Number of hours per replicative generation6 hrAssumption
rAverage number of virions produced by an infected cell100 virionsFrensing et al., 2016
C0Initial target cell population size4 × 108 virionsHadjichrysanthou et al., 2016
dPer capita decay rateTwo per-generationAssumption
R0Within-host basic reproduction number5Hadjichrysanthou et al., 2016
µPer-site, per-generation mutation rate3 × 10−5 per-site, per-generationMcCrone et al., 2020

Additional files

Supplementary file 1

Mean number of nonsynonymous (NS), synonymous (S), and stop codon (Stop) variants per sample for each gene segment as well as the corresponding NS/S ratio.

https://cdn.elifesciences.org/articles/68917/elife-68917-supp1-v2.xlsx
Supplementary file 2

Potentially linked nonsynonymous variants in within-host A/H1N1pdm09 and A/H3N2 virus samples.

Sample names are given in the format of 'Patient ID_Days since symptom onset.' Both linkage disequilibrium (LD) and the normalized LD' measures are tabulated alongside the inferred maximum-likelihood haplotype frequencies (q10 and q01 are the haplotype frequencies with variant i or ii only while q11 is the frequency of haplotypes encoding both variants).

https://cdn.elifesciences.org/articles/68917/elife-68917-supp2-v2.xlsx
Supplementary file 3

A/H3N2 segment-specific primers.

https://cdn.elifesciences.org/articles/68917/elife-68917-supp3-v2.xlsx
Supplementary file 4

Patients metadata (provided as an Excel file).

https://cdn.elifesciences.org/articles/68917/elife-68917-supp4-v2.xlsx
Supplementary file 5

Acknowledgment table of reference sequences downloaded from GISAID.

https://cdn.elifesciences.org/articles/68917/elife-68917-supp5-v2.xlsx
Transparent reporting form
https://cdn.elifesciences.org/articles/68917/elife-68917-transrepform-v2.docx

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Alvin X Han
  2. Zandra C Felix Garza
  3. Matthijs RA Welkers
  4. René M Vigeveno
  5. Nhu Duong Tran
  6. Thi Quynh Mai Le
  7. Thai Pham Quang
  8. Dinh Thoang Dang
  9. Thi Ngoc Anh Tran
  10. Manh Tuan Ha
  11. Thanh Hung Nguyen
  12. Quoc Thinh Le
  13. Thanh Hai Le
  14. Thi Bich Ngoc Hoang
  15. Kulkanya Chokephaibulkit
  16. Pilaipan Puthavathana
  17. Van Vinh Chau Nguyen
  18. My Ngoc Nghiem
  19. Van Kinh Nguyen
  20. Tuyet Trinh Dao
  21. Tinh Hien Tran
  22. Heiman FL Wertheim
  23. Peter W Horby
  24. Annette Fox
  25. H Rogier van Doorn
  26. Dirk Eggink
  27. Menno D de Jong
  28. Colin A Russell
(2021)
Within-host evolutionary dynamics of seasonal and pandemic human influenza A viruses in young children
eLife 10:e68917.
https://doi.org/10.7554/eLife.68917