Characterizing genetic and antigenic divergence from vaccine strain of influenza A and B viruses circulating in Thailand, 2017–2020

We monitored the circulating strains and genetic variation among seasonal influenza A and B viruses in Thailand between July 2017 and March 2020. The hemagglutinin gene was amplified and sequenced. We identified amino acid (AA) changes and computed antigenic relatedness using the Pepitope model. Phylogenetic analyses revealed multiple clades/subclades of influenza A(H1N1)pdm09 and A(H3N2) were circulating simultaneously and evolved away from their vaccine strain, but not the influenza B virus. The predominant circulating strains of A(H1N1)pdm09 belonged to 6B.1A1 (2017–2018) and 6B.1A5 (2019–2020) with additional AA substitutions. Clade 3C.2a1b and 3C.2a2 viruses co-circulated in A(H3N2) and clade 3C.3a virus was found in 2020. The B/Victoria-like lineage predominated since 2019 with an additional three AA deletions. Antigenic drift was dominantly facilitated at epitopes Sa and Sb of A(H1N1)pdm09, epitopes A, B, D and E of A(H3N2), and the 120 loop and 190 helix of influenza B virus. Moderate computed antigenic relatedness was observed in A(H1N1)pdm09. The computed antigenic relatedness of A(H3N2) indicated a significant decline in 2019 (9.17%) and 2020 (− 18.94%) whereas the circulating influenza B virus was antigenically similar (94.81%) with its vaccine strain. Our findings offer insights into the genetic divergence from vaccine strains, which could aid vaccine updating.

Genetic characterization. Phylogenetic trees for each influenza virus were generated using the nucleotide sequences obtained from this study and additional HA sequences from other Thai strains previously identified and available from the National Center for Biotechnology Information (www.ncbi.nlm.nih.gov) and the Global Initiative for Sharing All Influenza Data (GISAID) (http://platf orm.gisai d.org) database. A total of 183 nucleotide sequences of influenza A(H1N1)pdm09, 232 nucleotide sequences of influenza A(H3N2), and 177 nucleotide sequences of the influenza B virus were aligned with the reference and vaccine strains for each influenza virus using MUSCLE. Nucleotide substitution models for the HA gene of influenza A(H1N1)pdm09 (TN93 + G), influenza A(H3N2) (HKY + G + I) and influenza B (HKY + G) were implemented in MEGAX 29 . Trees were constructed using the maximum-likelihood method and bootstrapping of 1,000 replicates. Deduced amino acid sequences were compared to the reference and vaccine strains to identify substitutions, which were noted at each branch node on the phylogenetic trees. Potential N-linked glycosylation sites on the HA gene was analyzed using the NetNGlyc 1.0 server 30  Determination of selection pressure. The ratio of non-synonymous/synonymous substitutions (dN/ dS) was considered when evaluating codon under selective pressure. dN/dS was analyzed using the mixedeffects model of evolution (MEME) and the fixed effects likelihood (FEL) methods. Both algorithms were in the HYPHY software implemented in the Datamonkey webserver (https ://www.datam onkey .org/) 31 . Positively selected residue was considered significant at P = 0.1.

Estimation of computed antigenic relatedness.
The antigenic relatedness of influenza virus was computed using the P epitope model, which took into account the distinct antigenic sites between circulating strains and the vaccine strain by considering the epitope sites. The P epitope model was calculated by the fraction of number of amino-acid substitutions in the dominant epitope and the total number of amino acids in that dominant epitope. The association between P epitope model and antigenic distance measured by hemagglutinin inhibition assay or antigenic relatedness (efficacy, E) was determined by a mathematical formula. For A(H1N1) pdm09, E = − 1.19 × P epitope + 0.53 in which efficacy is 53% when the P epitope = 0 32 . For A(H3N2), the association between the antigenic relatedness and P epitope is given by E = − 2.47 × P epitope + 0.47 in which efficacy is 47% when P epitope = 0 33 . For the influenza B virus, E = − 0.864 × P epitope + 0.6824 in which efficacy is 68.24% when P epitope = 0 34 . The trend of computed antigenic relatedness was identified using R v3.6.0 (R Foundation for Statistical Computing, Vienna, Austria; https ://www.r-proje ct.org). The difference among annually computed antigenic relatedness was calculated using one-way ANOVA (P < 0.05 was considered statistically significant). Statistical analysis was done using Prism 8.0 (GraphPad, San Diego, CA, USA; https ://www.graph pad.com).
A(H3N2). Phylogenetic analysis of 232 HA sequences of influenza A(H3N2) comprising the Thai strains identified in this study and previously elsewhere in Thailand, along with the vaccine and reference strains, showed that 98.7% (229/232) belonged to clade 3C.2a (the remaining three Thai strains were 3C.3a, all identified in 2020) (Fig. 3). These Thai strains diverged into subclades 3C.2a1 (defined by N121K, R142G, and N171K on HA1, and I77V and G155E on HA2) and 3C.2a2 (defined by T131K, R142K, and R261Q on HA1) (Fig. S1 , we identified N156K, L161I and S164T substitutions on the Sa antigenic site (Fig. 5A), which often drives A(H1N1)pdm09 evolution, particularly at residue 156 13 . We found that T185I, D187A, and Q189E mapped to the Sb antigenic site and overlapped the RBS, while S74R is on epitope Cb. No residue changes appeared on epitopes Ca1 nor Ca2. Using A/Hong Kong/4801/2014 as a representative A(H3N2) vaccine strain for clade 3C.2a, five residue changes mapped to antigenic epitope A (T131K, T135K, S137F, A138S, R142G/K), three residues mapped to epitope B (T128A, T160K, and L194P), one residue mapped to epitope C (H311Q), four residues mapped to epitope D (N96S, N121K, N171K, and A212T), and three residues mapped to epitope E (E62G, K92R, and R261Q) (Fig. 5B). Using A/Switzerland/9,715,293/2013 as a representative of A(H3N2) vaccine strain for subclade 3C.3a, which was the vaccine strain in 2015, four residue changes identified among the Thai strains were mapped onto the HA structure (S144K on site A, F193S on site B, N246E on site D, and S91N on site E) ( On the HA structure of B/Brisbane/60/2008, circulating B/Victoria strains possessed substitutions at residues I117V, N129D, and K136E, which aggregate in the vicinity of the 120 loop region (Fig. 5D). Meanwhile, D197N which is a reversion of the substitution in the vaccine strain acquired as a result of egg-adaption was mapped to the 190 helix region. In contrast, none of the residue changes identified among the circulating B/Yamagata strains occurred on the antigenic or the RBS (Fig. 5E).
Determination of the selection pressure on influenza virus. Amino acid changes on the HA often result from selective pressure being exerted on the virus by the host immunity during infection. We therefore evaluated the potential positive selection of these residues by examining their rate of change (dN/dS) ( Table S2)

Discussion
Year-round influenza activity in Thailand can give rise to antigenically drifted influenza virus strains, which are sufficiently different to be categorized into emerging new subclades. These strains sometimes differ significantly from the vaccine strains and are able to escape the host immunity elicited by the annual influenza vaccine. www.nature.com/scientificreports/ Consequently, the effectiveness of influenza immunization is substantially reduced 3,35 . Our study aimed to characterize the predominantly circulating A(H1N1)pdm09 and A(H3N2), as well as both lineages of the influenza B virus and how residue changes can affect the antigenic relatedness. Analysis of their HA gene sequence and deduced amino acids identified genetic heterogeneity, which were different from those found in the southern hemisphere vaccine strains for their respective years. Among the influenza A(H1N1)pdm09 virus strains in this study, most amino acid substitutions contributing to the phylogenetic cluster transitions were accumulated at the antigenic sites Sa and Sb, while A(H3N2) demonstrated greater diversity with most mutations dominantly located on epitopes A, B, D and E. Meanwhile, antigenic drifts of the influenza B virus occurred at the 120 loop and the 190 helix. Observed changes on the antigenic and the RBS can be important since even one residue substitution at these sites could potentially drive new antigenic variants 12,36 .
The The computed antigenic relatedness for A(H3N2) during this study period, on the other hand, performed worse than A(H1N1)pdm09. This was not surprising given that the evolutionary rate for A(H3N2) in the HA1 domain is considerably greater than that for A(H1N1)pdm09 43 . Among several circulating A(H3N2), clade 3C.2a predominated since 2015 27 and has now accumulated sufficient changes to be classified as 3C.2a1 and 3C.2a2 17 . Since then, 3C.2a1b was established, which has gained additional substitutions (either T131K or T135K with T128A) in HA1 and these changes are located in the antigenic and the RBS region 7 . The 3C.2a1b virus was also reportedly the predominant subclade in Europe for 2017-2018 winter month 37 . Due to multiple clades/subclades circulation, our study found the mismatch of the 2019 southern hemisphere vaccine (A/Switzerland/8060/2017 (subclade 3C.2a2)) and circulating strains in 2019 (3C.2a1b + T131K ). Furthermore, although A/South Australia/34/2019 which is a 3C.2a1b + T131K strain was announced for the 2020 southern hemisphere vaccine strain 22 , the predominate strains in the first three months of 2020 were subclade 3C.2a1b with additional T135K and T128A which slightly far away from their vaccine strain. T135K and T128A, which are located at antigenic sites A and B respectively, are predicted to cause a loss of glycosylation that might alter the HA antigenic properties and affect antibody recognition 6 . Further selective pressure analysis on the HA1 of A(H3N2) suggests that the residues changes at positions 131, 135, 144 and 193 have dN/dS ≥ 1 which indicated the positive selection due to immune escape mechanisms. Detection of clade 3C.3a strains in 2020 has been less frequent so far. A previous report showed that antibodies generated against 3C.2a less neutralize 3C.3a virus, particularly when 3C.3a possesses an additional F193S. It also does not neutralize 3C.2a1 well compared to 3C.2a2 44 .
The A(H3N2) vaccine strain was not effective during the 2018-19 influenza season in Europe, and low VE (− 58 to 57%) continues during 2019-20 season 39,42 . Our study also found that the computed antigenic relatedness against A(H3N2) was < 50% between 2019 and 2020 due to the circulation of multiple clades/subclades. These data suggested that the genetic diversity in A(H3N2) might hamper identification of a well-matched virus as a vaccine component in the 2014-2015 influenza season 45 A meta-analysis indicated reduced protection and substantial variation of VE against A(H3N2). At the same time, the influenza vaccine provided moderate-tohigh protection against A(H1N1)pdm09 and influenza B viruses 35 . Alternatives to egg-based manufacturing should proceed since a negative impact of egg-induced mutations in the H3N2 vaccine strain has been found 46 .
Newly emerging influenza B/Victoria strains are antigenically distinct and possess either double deletion at residues 162-163 or triple deletion at residues 162-164 within the HA1 domain. In this study, the triple deletion strains were more commonly found than the ancestor B/Brisbane-like strain during the 2019-2020 season. New B/Victoria lineages are actively circulating and warrant changing vaccine strain recommendations for this lineage since 2018 8,17 . Interim VE estimates of the 2019-2020 influenza season in the U.S. against the predominant B/Victoria lineage was 39-59% which showed a higher VE than influenza A virus 41 . This is consistent with a European study showing the VE of 62-83% against the influenza B virus for all ages during 2019-2020 season 42 . Our study also revealed a high computed antigenic relatedness against the new B/Victoria virus.
Control measures to mitigate coronavirus disease 2019 (COVID-19) pandemic have resulted in a decrease of influenza activity in 2020 thus far compared to the corresponding period of the previous year [47][48][49] . Monitoring influenza virus genetic drift remains critical in tracking novel residue changes, which could potentially affect virulence and evasion of host immunity.

Scientific Reports
| (2021) 11:735 | https://doi.org/10.1038/s41598-020-80895-w www.nature.com/scientificreports/ Our study had several limitations. First, there was no information on the vaccination status of the individuals with the samples sent to us for influenza virus testing. Second, the results for computed antigenic relatedness, which was assessed using the accumulated substitutions on antigenic sites, ideally would require additional confirmatory antigenic characterization such as hemagglutinin inhibition or virus neutralization assay in order to complement and strengthen our existing data. Finally, we did not investigate the genetic changes in the neuraminidase gene, which encodes a surface glycoprotein that is also immunogenic.
By monitoring the genetic and antigenic changes of the influenza virus circulating in the tropics, data from this study suggests that new clades and subclades for both influenza A(H1N1)pmd09 and A(H3N2) viruses increasingly differ from those of the chosen vaccine strain, but not influenza B virus. These findings highlight the need for improved vaccine strain match particularly for A(H3N2).

Data availability
All data generated during this study are contained within this manuscript and its Supplementary Information files.