Genetic Analysis of Influenza A/H1N1pdm Strains Isolated in Bangladesh in Early 2020

Influenza is one of the most common respiratory virus infections. We analyzed hemagglutinin (HA) and neuraminidase (NA) gene segments of viruses isolated from influenza patients who visited Evercare Hospital Dhaka, Bangladesh, in early 2020 immediately before the coronavirus disease 2019 (COVID-19) pandemic. All of them were influenza virus type A (IAV) H1N1pdm. Sequence analysis of the HA segments of the virus strains isolated from the clinical specimens and the subsequent phylogenic analyses of the obtained sequences revealed that all of the H1N1pdm recent subclades 6B.1A5A + 187V/A, 6B.1A5A + 156K, and 6B.1A5A + 156K with K209M were already present in Bangladesh in January 2020. Molecular clock analysis results suggested that the subclade 6B.1A5A + 156K emerged in Denmark, Australia, or the United States in July 2019, while subclades 6B.1A5A + 187V/A and 6B.1A5A + 156K with K209M emerged in East Asia in April and September 2019, respectively. On the other hand, sequence analysis of NA segments showed that the viruses lacked the H275Y mutation that confers oseltamivir resistance. Since the number of influenza cases in Bangladesh is usually small between November and January, these results indicated that the IAV H1N1pdm had spread extremely rapidly without acquiring oseltamivir resistance during a time of active international flow of people before the COVID-19 pandemic.


Introduction
Seasonal influenza is one of the most common respiratory virus infections. It is estimated to induce about 3 to 5 million cases of severe illness and about 290,000 to 650,000 respiratory deaths per year prior to the coronavirus disease 2019 (COVID- 19) pandemic (WHO fact sheets. Influenza (Seasonal). https://www.who.int/news-room/ fact-sheets/detail/influenza-(seasonal); accessed on 19 July 2021). Seasonal influenza is caused by influenza virus, which belongs to the family Orthomyxoviridae and is classified into types A, B, and C (International Committee on Taxonomy of Viruses (ICTV). https: //talk.ictvonline.org/taxonomy/; accessed on 19 June 2021). Influenza type A (IAV) is the primary type that causes seasonal influenza. It possesses eight negative-sense, singlestranded RNA segments that encode hemagglutinin (HA), neuraminidase (NA), RNA polymerase subunit PB2 (PB2), RNA polymerase subunit PB1 (PB1) and PB1-F2, RNA

Ethics Statement
This study proposal was approved by the Research and Ethics Committee of Evercare Hospital Dhaka (ERC 25/2020-1).

Clinical Specimens
Clinical specimens of nasal or nasopharyngeal fluid were collected at Evercare Hospital Dhaka between January and February 2020 from 94 patients with influenza-like illness as defined by the World Health Organization (WHO) (Global Epidemiological Surveillance Standards for Influenza; https://www.who.int/publications/i/item/9789241506601; accessed on 20 July 2021). Two swabs were taken from each patient: one specimen was tested using the SD BIOLINE influenza antigen test kit (Standard Diagnostics, Inc., Yongin-si, (Applied Biosystems). Obtained sequences were deposited in the EpiFlu™ database of the GISAID (www.gisaid.org) (accessed on 21 February 2022) with the accession number EPI1986952, EPI1986953, EPI1986977-EPI1986982.

Sequence Characterization
Phylogenetic analyses were conducted on the HA segment of IAV H1N1pdm09. The IAV H1N1pdm09 sequences were downloaded from the EpiFlu TM database of the Global Initiative on Sharing All Influenza Data (GISAID) database (https://www.gisaid.org/; accessed on 13 September 2021)) [14]. DNA and amino acid sequence analyses were performed with GENETYX ® software (GENETYX, Tokyo, Japan) and MEGA X: Molecular Evolutionary Genetics Analysis across computing platforms [15].

Phylogenetic Analysis
The nucleotide sequence dataset was aligned using Muscle in AliView v1.26 [16]. The substitution model was selected, and the maximum likelihood (ML) tree was constructed in IQ-TREE with 1000 ultrafast bootstrap replicates [17]. The time-scaled tree was constructed using the BEAST package with a Bayesian Markov chain Monte Carlo approach [18]. The input dataset was inspected for a positive temporal signal in Tempest v1.5.3 [19], and the nucleotide substitution model was determined by using ModelFinder [20]. An exponential population growth model and the relaxed molecular clock were employed as previously described [21]. Triplicate runs of Markov chain Monte Carlo chain lengths of 30,000,000 generations with sampling every 3000 generations were performed, with individually obtained effective sample sizes over 200 traced in Tracer and combined in LogCombiner. A maximum clade credibility tree was constructed in TreeAnnotator and visualized in Figtree v1.4.4 (http://tree.bio.ed.ac.uk/software/figtree; accessed on 16 October 2021).

Background of the Clinical Specimens, RT-PCR Analysis, and Virus Isolation
Of 94 patients with influenza-like symptoms, 21 patients were IVA-positive using a rapid diagnosis test kit. There were no influenza B-positive samples. The IVA-positive specimens were collected from 14 males and seven females with mean age of 24.7 years (range, 0.5-71 years), and there was an average of 2.48 ± 1.97 days between the onset of clinical signs and sample collection (Supplementary Table S2). Subsequently, the IVA-positive specimens were subjected to RT-PCR for subtyping and to SYBR green real-time RT-PCR for the quantitation of influenza virus RNA. The subtyping RT-PCR results indicated that 17 of 21 specimens contained H1N1pdm, while the remaining four were all negative for H1N1pdm, seasonal H1N1 before 2009, and seasonal H3N2. Quantitative real-time RT-PCR results indicated that the virus quantity ranged from 0.3 FFU/mL to 4.0 × 10 3 FFU/mL. The four samples in which we failed to determine subtype showed the lowest, second lowest, third lowest, and fourth lowest virus quantity (range, 0.3-0.5 FFU/mL) on quantitative real-time PCR (Supplementary Table S2). Next, we attempted to isolate the virus from the 17 clinical specimens that were H1N1pdm-positive by RT-PCR subtyping; viruses were successfully isolated from four clinical specimens (No. 10, 15, 46, and 59; Supplementary  Table S2). All these isolates were obtained after the second passages.

Genetic Characterization of the HA Segments Obtained in the Present Study
The four viruses isolated in the present study were subjected to sequencing and phylogenetic analyses. The amino acid sequences of the HA segments were characterized and classified into subclades according to the classification of the February 2019 WHO vaccine composition meeting; we identified seven subclades (6B.1A1-6B.1A7) defined by amino acid substitutions [22]. The sequences of the four isolates were all classified into subclade 6B.1A, as they possessed amino acid substitutions S74R, S164T, and I295V in HA1 as compared to the prototype A/Michigan/45/2015. They were further classified into subclade 6B.1A5A by the presence of the additional amino acid substitutions of S183P, N129D, T185I, and N260D (T185I was in the Sb antigenic site). They also carried K130N, N156K, L161I, V250A, and E506D (N156K and L161I were in the Sa antigenic site), which are specific to subclade 6B.1A5A + 156K (Figure 1, Supplementary Table S3). These results indicated that all the four strains obtained in the present study belonged to this subclade 6B.1A5A + 156K. classified into subclades according to the classification of the February 2019 WHO vaccine composition meeting; we identified seven subclades (6B.1A1-6B.1A7) defined by amino acid substitutions [22]. The sequences of the four isolates were all classified into subclade 6B.1A, as they possessed amino acid substitutions S74R, S164T, and I295V in HA1 as compared to the prototype A/Michigan/45/2015. They were further classified into subclade 6B.1A5A by the presence of the additional amino acid substitutions of S183P, N129D, T185I, and N260D (T185I was in the Sb antigenic site). They also carried K130N, N156K, L161I, V250A, and E506D (N156K and L161I were in the Sa antigenic site), which are specific to subclade 6B.1A5A+156K (Figure 1, Supplementary Table S3). These results indicated that all the four strains obtained in the present study belonged to this subclade 6B.1A5A+156K.

Human IAV H1N1pdm09 Detected in Bangladesh from 2018 to 2020
To investigate the subclades that have circulated in Bangladesh in the past, the HA sequences of Bangladesh strains of human IAV H1N1pdm registered in the EpiFlu GISAID database between January 2018 and March 2021 were collected, and a phylogenetic analysis was performed. The sequence dataset consisted of 376 strains (124, 212, and 18 strains registered in 2018, 2019, and 2020, respectively), the four sequences obtained in the present study, and 20 reference strains of 6B.1A and 6B.1A1-A7. Results showed that

Human IAV H1N1pdm09 Detected in Bangladesh from 2018 to 2020
To investigate the subclades that have circulated in Bangladesh in the past, the HA sequences of Bangladesh strains of human IAV H1N1pdm registered in the EpiFlu GI-SAID database between January 2018 and March 2021 were collected, and a phylogenetic analysis was performed. The sequence dataset consisted of 376 strains (124, 212, and 18 strains registered in 2018, 2019, and 2020, respectively), the four sequences obtained in the present study, and 20 reference strains of 6B.1A and 6B.1A1-A7. Results showed that 122 of 2018 strains were clustered with Ireland/84630/2018, the reference strain of the subclade 6B.1A6 defined by the amino acid substitutions T120A and S183P (Figure 2 122 of 2018 strains were clustered with Ireland/84630/2018, the reference strain of the subclade 6B.1A6 defined by the amino acid substitutions T120A and S183P (Figure 2). The remaining two 2018 strains were classified as either 6B.1A2 or 6B.1A5A (Figure 2, Supplementary Table S3). On the other hand, all of the 2019 strains belonged to subclade 6B.1A5A. Remarkably, the subclade 6B.1A6 that dominated in 2018 had completely disappeared in 2019. Eighteen strains registered in 2020 and the four strains of the present study belonged to subclade 6B.1A5A, which was further phylogenetically classified into three genetic groups: (1) two strains (9.1%) were related to the parental subclade 6B.1A5A strain; (2)     To examine the origin, divergence time, and related strains of these three genetic groups of Bangladesh 2020 strains, a new dataset was prepared for Bayesian phylogenetic analysis; the dataset included the most related strains obtained from a BLAST analysis conducted in EpiFlu GISAID, the earlier strains of 6B.1A5A + 187V/A and 6B.1A5A + 156K in influenza genomic epidemiology provided by Nextstrain in GISAID [24], and the reference strains of 6B.1A and 6B.1A1-A7. In this new dataset, the correlation coefficient between the collection date and root-to-tip divergence was 0.725 ( Figure 3A). The phylogenetic tree ( Figure 3B Table S4). These results indicated that the human IAV H1N1pdm09 subclades that dominated in Bangladesh in 2020 also dominated all over the world.

Genetic Characterization of NA Segments Obtained in the Present Study
NA segment sequences of Bangladesh strains of human IAV H1N1pdm registered in the EpiFlu GISAID database between January 2018 and March 2021 were collected, and a phylogenetic analysis was performed. The sequence dataset consisted of 356 strains (129, 209, and 18 strains registered in 2018, 2019, and 2020, respectively) and the four sequences obtained in the present study (Supplementary Table S4). The results showed almost the same pattern as HA segments (Supplementary Figure S1). Nearly all the strains in 2018 (126 out of 129 strains) formed a distinct cluster from the rest of strains. This result indicated that NA sequences that dominated in 2018 totally disappeared from Bangladesh in 2019 and 2020 as observed in HA sequences. Furthermore, the strains collected in 2020 formed three distinct clusters, just like HA genes. NA segment analysis indicated that all four sequences from the present study lacked the H275Y amino acid substitution (Supplementary Table S5) that confers oseltamivir resistance [25]. The H275Y mutation was found only one strain (A/Bangladesh/4005 2019|EPI ISL 395161 NA) in the Bangladesh strains registered in the GISAID database between January 2018 and March 2021 (Supplementary  Table S5); however, it was found in seven (five 6B.1A5A + 156K, one 6B.1A5A + 187V/A, and one non-6B clade) of the 737 strains registered from all over the world between March 2020 and March 2021 (Supplementary Table S6). These results indicated that oseltamivir resistance was still rare in these human IAV H1N1pdm09 subclades.

Discussion
In the present study, we determined the nucleotide sequences of HA and NA genome segments of IAV H1N1pdm specimens collected in Bangladesh between January and February 2020. The HA gene sequences were then compared to the sequences in EpiFlu GISAID database. The viruses obtained in the present study belonged to subclade 6B.1A5A + 156K. Our molecular clock analysis suggested that another subclade, 6B.1A5A + 187V/A, which was also present in Bangladesh in 2020, Subsequently, these two subclades spread worldwide, particularly in America, East Asia, South Asia, and Southeast Asia. Notably, the more recent clade 6B.1A5A + 156K with K209M was also observed, and its time of emergence was inferred to be September 2019 (August to November 2019). This virus was detected for the first time in Hong Kong and South Korea in late November and early December 2019, and it later became dominant in Bangladesh, especially in March 2020. Therefore, the present study showed that all of the recent subclades of H1N1pdm, 6B.1A5A + 187V/A, 6B.1A5A + 156K, and 6B.1A5A + 156K with K209M, all of which first emerged between April and September 2019 in countries other than Bangladesh, were also present in Bangladesh in January 2020. Since the number of influenza cases in Bangladesh is usually small between November and January, [26], these results demonstrated that the IAV had spread extremely rapidly.
The amino acid substitution of N156K in subclade 6B.1A5A + 156K is reported to have a significant impact on vaccine effectiveness. Only 2 years after the 2009 H1N1pdm pandemic, Strengell et al. showed that the N156K mutation had reduced the effectiveness of the influenza vaccine [27]. Then, in 2013, Guarnaccia et al. showed that the N156K mutation altered the binding efficiency of HA-specific antibodies [28]. However, the prevalence of the N156K mutation was as low as 0.15% among the global isolates in registered databases in 2009 to 2012 [28] and was 8.2% among Vietnamese isolates in 2010 to 2013 [29]. On the other hand, the results of the present study showed that in Bangladesh, 0%, 1.9%, and 60.9% of the isolates had this mutation in 2018, 2019, and 2020, respectively, and after March 2020, 47.7% of the isolates reported worldwide had this mutation. These results indicated that the prevalence of this mutation has rapidly increased since 2020. Tenford et al. reported that the sera from individuals vaccinated with the 2019-2020 season vaccine recognized the 6B.1A5A + 156K strain poorly, and the vaccine failed to protect against infection by this virus [30]. Therefore, the WHO has changed their recommendation on the composition of the influenza vaccines for the 2021 southern hemisphere influenza season and the 2021-2022 northern hemisphere influenza season to vaccines using the A/Victoria/2570/2019 strain or A/Wisconsin/588/2019 strain as reference strains for H1N1 [23,31]. Vaccines using these strains showed better recognition of 6B.1A5A + 156K but poor recognition of other strains, such as 6B.1A5A + 187V/A, which lacks the N156K mutation [23,30]. However, 6B.1A5A + 187V/A strains have also been detected in many cases since 2020; therefore, monitoring of these strains may be important for assessing the effectiveness of future vaccines.
COVID-19 is a respiratory infection similar to influenza. The preventative measures that have been adopted against the spread of COVID-19 are thought to have also reduced the number of influenza cases [23]. In addition, asymptomatic patients play an important role in influenza transmission [32,33], and a decrease in the total number of influenza patients also decreases the number of asymptomatic patients, which may result in a reduction in the herd immunity to influenza viruses. Due to the shift of human and equipment resources to COVID-19 testing, there has been a reduction of influenza surveillance and reporting activities, and the number of influenza virus strains being registered in databases recently has greatly decreased when compared to previous years [23]. Although it remains unclear how SARS-CoV-2 infection affects influenza virus infection, several cases of coinfection with both of these viruses have been reported [34]. Therefore, it is expected that influenza surveillance will become even more important in the future [35].

Conclusions
We determined nucleotide sequences of HA and NA genes of IAV H1N1pdm specimens collected in Bangladesh between January and February 2020. Phylogenic analyses revealed that all of three H1N1pdm recent subclades were already present in Bangladesh in January 2020. Molecular clock analyses suggested that these subclade emerged in April-September 2019 in countries other than Bangladesh. These results indicated that the IAV H1N1pdm had spread extremely rapidly during a time of active international flow of people before the COVID-19 pandemic.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/tropicalmed7030038/s1, Table S1: Primers used for RT-PCR subtyping, real-time PCR, cDNA synthesis, and sequencing; Table S2: Data on the clinical specimens and the immunochromatographic test, RT-PCR subtyping, quantitative real-time RT-PCR, and virus isolation results; Table S3: Amino acid sequences of HA of the Bangladeshi IAVs registered from 2018 to 2020; Table S4: Amino acid sequences of HA of the global IAVs registered from March 2020 to March 2021; Table S5: Amino acid sequences of NA of the Bangladeshi IAVs registered from 2018 to 2020; Table S6: Amino acid sequences of NA of the global IAVs registered from March 2020 to March 2021; Figure S1: The H1N1pdm09 NA genes in Bangladesh in 2018 to 2020. A maximum-likelihood tree of H1N1pdm09 was generated based on the NA gene with 1000 ultrafast bootstrap replicates using IQ-TREE. Nucleotide sequences of the NA genes retrieved from EpiFlu GISAID are labeled as follows: isolate name|isolate ID. The Bangladesh 2018, 2019, and 2020 sequences are shown as brown circles, green squares, and blue rhombuses, respectively. The four sequences obtained in the present study are shown as red rhombuses. The H1N1pdm09 reference sequences are shown as black triangles. The strain with H275Y amino acid substitution is shown as a right green square. Funding: This research was funded by the Japan Agency for Medical Research and Development under grant numbers JP19fm0108003, 20wm0225010h0101, and 21wm0225010h0102, with support from the Osaka University ASEAN campus project.

Institutional Review Board Statement:
The study was conducted in accordance with the Declaration of Helsinki and approved by the Research and Ethics Committee of Evercare Hospital Dhaka (ERC 25/2020-1).
Informed Consent Statement: Patient consent was waived since we used left-over specimens in the study.

Data Availability Statement:
The data presented in this study are available in this article and Supplementary Materials.