Introduction

Highly Pathogenic Avian Influenza (HPAI) viruses are a global concern for both animal and human health [1]. HPAI H5N1 viruses of the A/goose/Guangdong/1/1996 lineage were first discovered in 1996 in China and since have continued to circulate in poultry and wild birds across Asia, the Middle East, Europe and Africa [2, 3]. In total, 68 countries have been affected and millions of birds have succumbed to the disease or have been culled to prevent further spread of the disease [4]. HPAI H5N1 viruses infect humans sporadically and may cause severe disease with a high case fatality rate among confirmed hospitalized patients. To date, there are 861 confirmed human cases, of which 455 have died [5]. HPAI H5N1 viruses continue to evolve through genetic drift and reassortment events with other avian influenza A viruses, resulting in multiple genetic clades and subtypes [6, 7].

HPAI H5N1 viruses were first identified in Indonesia from poultry outbreaks on the Java Island in 2003 and had since spread to other parts of the country [8, 9]. During subsequent years, clade 2.1 viruses became enzootic in Indonesia [10]. However, a new HPAI H5N1 clade 2.3 virus was detected during poultry outbreaks since 2012. 200 human infections with HPAI H5N1 viruses have been reported in Indonesia so far with higher case fatality rates among reported cases (84%) compared to the rest of the world afflicted by the virus. We previously showed that high nasopharyngeal viral load was associated with more severe outcome of human H5N1 infections in Indonesia and that the virus was more commonly detected in blood relative to other geographical regions affected by HPAI H5N1 [11, 12]. Strikingly, although the number of detections in humans peaked in 2005–2006 in Indonesia and subsequently declined in the following years, the case fatality rate in Indonesia increased over time from 65% in 2005 to 100% since 2012. This increase was associated with higher viral load prior to treatment and the presence of mutations in the matrix protein that confers adamantane resistance [11]. However, reasons for the higher viral load and case fatality rate are still unclear. More detailed sequence analyses are warranted to investigate the presence of known virulence markers and substitutions related to possible human adaptation, which can help explain the higher viral loads and increased mortality.

In response to the outbreaks, the Indonesian government implemented a strategy to reduce the incidence of HPAI H5N1 virus infections in poultry including stamping out of infected poultry, culling of contiguous flocks and poultry vaccination [13]. Several vaccines were developed and implemented to match the circulating strain in the poultry and pandemic preparedness over the time [14, 15]. The initial vaccine used was based on the A/Chicken/Legok/2003 isolate, a clade 2.1.1 virus [13, 16]. By 2010, the vaccine strain was updated to subclade 2.1.2 and 2.1.3 viruses, based on isolates A/Chicken/West Java/30/07 and A/Chicken/Nagrak/30/07, respectively [17]. To date, a new subclade of H5N1 has emerged (2.3.2.1) and a new vaccine was developed based on isolate A/duck/Sukoharjo/BBVW-1428-9/2012 [10]. However, as a consequence of the large-scale vaccination, antigenic drift was induced in poultry and consequently the vaccines became less effective [18].

Despite vaccination efforts, the number of poultry outbreaks remained high and the epidemic in poultry continued to spread among 32 out of 34 Indonesian provinces with over 11,000 reported poultry outbreaks since 2007 [19]. The large number of poultry outbreaks continues to pose a threat for future zoonotic infections in humans, antigenic drift and possible host adaptations that could increase the pandemic risk of circulating viruses [20, 21].

Improved insights into the genetic and antigenic characteristics of HPAI H5N1 viruses from Indonesia provide a better understanding of its epidemiology, the high case fatality rate and for a pandemic risk assessment [22]. Sequencing data contain valuable information about viral genetic characteristics, including presence of known human adaptive markers, resistance against available antiviral drugs or other changes that can explain the high and rising mortality, while antigenic characterization will help assess the potential protection of pre-pandemic vaccines.

Here, we conducted a study to characterize the viral genetics of HPAI H5N1 viruses isolated from patients in Indonesia between 2008 and 2015 that could explain the virulence leading to the high case fatality rate. We investigated the presence of known molecular determinants of virulence, receptor binding properties and antiviral susceptibility using whole genome sequencing of the HPAI H5N1 viruses. In addition, we investigated the antigenic properties of these human virus isolates and compared them to previous antigenic changes in HPAI H5N1 viruses from poultry, in order to assess the usefulness of and protection by current available pre-pandemic virus vaccines.

Materials and methods

Specimen and data collection

As part of the national procedure for avian influenza case investigation in Indonesia, respiratory specimens were collected from suspected H5N1 cases admitted to hospitals throughout Indonesia and sent to the national reference laboratory for influenza at the National Institute of Health Research and Development (NIHRD) in Jakarta. Suspected cases were defined according to World Health Organization (WHO) criteria [23]. The NIHRD is the reference laboratory under the Indonesian Ministry of Health responsible for laboratory testing and event-based surveillance of emerging infectious diseases in humans, including avian influenza A/H5N1 virus. Because Indonesian clinical specimens are obtained from suspected H5N1 cases as part of the national outbreak procedure for HPAI H5N1 case investigations, requirement for informed consent has been waived by the Indonesian Ministry of Health.

The specimens and data were collected from January 2008 to December 2015, according to the national outbreak investigation protocol following circulation of HPAI H5N1 viruses in South East Asia [24]. All of the specimens collected were stored and analyzed at the NIHRD. Laboratory identification and confirmation was determined using real-time reverse transcriptase-polymerase chain reaction (RT-PCR) typing and subtyping assay according to the Centers for Disease Control (CDC) (Atlanta, United States) protocol [25]. Specimens for all laboratory-confirmed cases were selected for subsequent virus isolation and genetic analyses, based on the specimen with the lowest cycle threshold (Ct) value according to the real-time RT-PCR, available for each patient.

Specimen isolation and characterization

The 71 selected specimens positive for influenza A(H5N1) virus with a Ct value below 35 were grown in 9- to 10-day-old specific pathogen-free (SPF) embryonated chicken eggs in a Biosafety Level 3 (BSL3) facility [26]. After incubation at 37 °C for 30 h, the egg allantoic fluid was harvested and hemagglutination titers were determined by hemagglutination assay. A total of 35 positive cultures were obtained due to variable specimen quality and limited availability of specimen volumes. The viral RNA was extracted from 200 µl of influenza virus positive allantoic fluids using High Pure RNA isolation kit (Roche) with on-column DNase treatment according to the manufacturer’s instructions. The RNA was reverse transcribed into cDNA using Uni12M primer (AGCRAAAGCAGG) [27] using SuperScript III Reverse Transcriptase (Invitrogen, Carlsbad, CA, USA) according to the manufacturer’s protocol. PCR amplification was performed using gene specific whole genome degenerative primer sets (primer sequences available upon request) [28,29,30] using Platinum Taq DNA Polymerase High Fidelity (Invitrogen). The PCR products were then purified with the ExoSAP-IT purification kit (Affimetrix, Inc, Santa Clara, CA) according to the manufacturer’s protocol.

The complete coding sequences were sequenced using the Big Dye Terminator v3.1 cycle sequencing kit (Applied Biosystem, Foster City, CA, USA). The products of the sequencing reactions were cleaned using Big Dye X Terminator Kit (Applied Biosystem, Foster City, CA, USA) according to manufacturer’s instructions and sequenced in a 16-capillary 3130xl Genetic Analyzer (Applied Biosystem, Foster City, CA, USA).

All nucleotide sequences obtained from this study have been deposited in the GISAID database (see Supplemental Table S2).

Phylogenetic analyses

The assembly and editing process of sequences from all eight gene segments was performed using Codon Code software (Gene Codes, USA). All sequences were aligned using ClustalW as available within BioEdit software version 7.0.8.0 [31]. To infer the evolutionary relationships between the viruses, maximum likelihood (ML) phylogenetic trees were constructed using RAxML 8.2.12 with the GTRGAMMA nucleotide substitution model [32, 33]. A ML phylogenetic tree was constructed using the combined nucleotide alignment of hemagglutinin (HA) sequences from the newly sampled viruses and reference sequences used to defined the H5 nomenclature system (https://www.who.int/influenza/gisrs_laboratory/201101_h5smalltreealignment.txt; Fig. 1) [34, 35]. Sequence data of human and avian H5N1 viruses from Indonesia with all eight influenza virus gene segments (200 viral isolates as of January 2020) was downloaded from the (GISAID) EpiFlu Database [36]. Individual ML trees were reconstructed for each gene segment to compare the genetic diversity of the newly sampled viruses against those previously collected from Indonesia (Fig. S1). Tanglegrams were visualized using the Baltic toolkit (https://github.com/evogytis/baltic).

Fig. 1
figure 1

Maximum-likelihood phylogenetic tree of HA sequences of the newly sampled human HPAI H5N1 viruses. New virus isolates are indicated with encircled tips and colored by their respective year of sample collection. WHO reference strains are used to define the H5 nomenclature system [34, 35]

Full size image

Residue and molecular analysis

Amino acid sequences were analyzed to identify substitutions potentially linked to human adaptation, virulence, antiviral resistance and antigenic properties as listed in the CDC H5N1 Genetic Change Inventory [37]. In addition to this inventory, we also used FluSurver to identify potentially relevant substitutions present in our sequence dataset (https://www.gisaid.org, https://flusurver.bii.a-star.edu.sg). FluSurver is a web-based tool to rapidly screen the sequences for potential mutations based on the curated and published literature.

Antigenic assays

Virus titers were determined by hemagglutination assay and antigenic characterization was performed by hemagglutination inhibition (HI) assays according to WHO protocols [38, 39]. The ferret antisera specifically reactive to defined H5 hemagglutinin clades were raised as described previously [40]. All antisera were pretreated overnight at 37 °C with receptor destroying enzyme (RDE Vibrio cholerae neuraminidase), followed by inactivation for 1 h at 56 °C. The HI assays were performed using the following procedures: twofold serial dilutions of 50 µl antisera starting at a 1:20 were mixed with 25 µl of a virus containing 4 hemagglutinating units (HAU) and were incubated at 37 °C for 30 min. Then, 25 µl of 1% turkey erythrocytes was added and incubated at 4 °C for 1 h. The HI titer is determined as the reciprocal value of the highest serum dilution that completely inhibited the hemagglutination of the turkey erythrocytes. Antigenic properties were determined for 25 representative novel isolates. Selection was based upon available HA titer of virus stocks and availability of at least two independent replicate experiments, measuring HI titers for all available ferret sera.

Analysis of antigenic properties was conducted using antigenic cartography methods as described previously [40, 41]. Briefly, the HI titers are converted to a distance matrix in which the distance between one antigen and one antiserum corresponds to the difference between the log2 value of the maximum observed titer to the antiserum from any antigen and the titer of the antigen to the antiserum. This distance matrix is used as input for multidimensional scaling algorithms, which arrange the antiserum and antigen points in space to best satisfy the target distances specified by the HI data by minimizing the error. Therefore, the distances between the points in an antigenic map represent antigenic distance as measured by the HI assay, in which the distances between antigens (virus isolates) and antisera are inversely related to the log2 HI titer. Although only distances between antigens and antisera are measured in the HI assay, antigenic maps allow the indirect measure of antigenic distances between two viruses.

Result

Epidemiological characteristic of confirmed cases

During the course of this study between 2008 and 2015, over 8000 poultry outbreaks of HPAI H5N1 viruses were reported [42] and 82 cases of laboratory-confirmed human HPAI H5N1 virus infection were collected. Of these 82 cases, we successfully cultured 35 virus isolates for genetic and antigenic characterization. Table 1 shows a summary of the epidemiological and other data of these 35 cases. Among the 35 patients, the median age was 21 (range 2–40), 14 (60%) were male, and 21 (40%) were female, and 18 (51%) received Oseltamivir treatment. Specimens were collected at median 8 days post onset of symptoms (range 0–17). There were a total of 6 specimens collected in 2008, 6 specimens in 2009, 5 specimens in 2010, 5 specimens in 2011, 5 specimens in 2012, 5 specimens in 2013, 2 specimens in 2014 and 2 specimens in 2015. These samples were collected from regions with high incidence of poultry H5N1 outbreaks [16, 21], including West Java (31%), followed by Jakarta (26%) and the Banten province (14%).

Table 1 Overview of demographic and clinical characteristics of H5N1-infected patients January 2008–December 2015 included in this study
Full size table

Phylogenetic analyses show close relationship between human and avian HPAI H5N1 viruses

A ML phylogenetic tree based on hemagglutinin (HA) sequences was constructed to infer the evolutionary relationships between the newly isolated viruses and HPAI H5N1 viruses circulating globally (see “Materials and methods”). The ML tree of the HA gene showed that all novel viruses belonged to the monophyletic lineage of clade 2.1.3.x viruses known to have circulated in Indonesia between 2008 and 2015, including clade 2.1.3.2 (10 viruses), 2.1.3.2a. (8 viruses), 2.1.3.2b (16 viruses), and 2.1.3.3 (1 virus) (Fig. 1).

To further elucidate the phylogenetic relationships between the novel H5N1 viruses and those collected from Indonesia previously, ML phylogenetic trees were constructed for each individual influenza virus gene segment (Fig. S1). The phylogenetic groupings of individual gene segments were largely similar to the clade nomenclature assigned to the HA gene segment for most of the novel viruses. Some of the novel viruses are reassortants derived from gene segments attributed to different subclades of Indonesian H5N1 viruses. The NS, NP, MP, NA, PB2 genes of A/Indonesia/NIHRD14122/2014 and A/Indonesia/NIHRD13157/2013, both of which HA genes belonged to clade 2.1.3.2a, were found in the monophyletic lineage of clade 2.1.3.2b viruses in the respective gene phylogenies. Additionally, A/Indonesia/NIHRD7781/2008 (HA-2.1.3.3) had PB1, NP, NA and MP gene segments found in viruses belonging to clade 2.1.3.2 (Fig. S2). There were also no clear distinct phylogenetic groupings between human and avian viruses in any of the gene segment analyzed, indicating that viruses infecting both host types in Indonesia were genetically similar (Fig. S1D).

Residue and molecular analyses show limited presence of known human adaptive markers

Next, we compared the amino acid substitutions found in these newly isolated viruses against molecular markers known to alter viral phenotypes such as virulence, drug resistance and human host adaptation (Table S1) [43,44,45,46,47].

Hemagglutinin

The HA protein can affect the virulence and host range of HPAI H5N1 viruses due to (1) the presence of a multibasic cleavage site, as well as changes to (2) host cell receptor specificity, (3) N-linked glycosylation patterns and (4) HA stability.

The pathogenicity of avian influenza viruses is determined by the cleavability of the HA glycoprotein. The presence of multiple basic amino acid residues at the cleavage site of HA allows the glycoprotein to be cleaved into mature subunits HA1 and HA2 by furin-like proteases, which are ubiquitously expressed. To the contrary, HA containing a single basic residue are cleaved by trypsin-like proteases, predominantly expressed in the respiratory and intestinal tract of birds and the respiratory tract of humans. All of the new isolates analyzed in this study are highly pathogenic avian influenza viruses that encode a multibasic HA cleavage site [48, 49]. The cleavage site motif PQRESRRKKR↓G was found in 30 of the 35 newly isolated viruses while other variations (i.e., PQREGRRKKR↓G, PQRESKRKKR↓G, PQRESRRRKR↓G and PQRESRRKRR↓G) were observed in the remaining isolates.

Another key feature of HA related to both virulence and human adaptation is the receptor specificity. Conserved residues within the receptor binding site (RBS) of HA are required for binding to sialic acid receptors (SIA), while several other residues in domains surrounding the RBS are key determinants of receptor specificity. Residues in these domains, the 130-loop, 190-helix and the 220-loop, determine the specificity for either the avian-type receptor or human-type receptor, α2-3-linked SIA or α2-6-linked SIA, respectively. Several key residues within these domains have been identified at positions including 186, 190, 193, 224, 226 and 228 (H3 numbering). High conservation of amino acid sequences was found at the receptor binding site (RBS). All isolates possessed a conserved residue at position N186, E190, N224, Q226 and G228, the most apparent residues involved in receptor binding specificity (as reviewed in [50]), indicating preferential binding of the viruses to avian like α2-3-linked SIA [51]. Interestingly, polymorphism was observed at position 193 for which a methionine or isoleucine was observed, instead of the more common arginine or lysine. However, the exact role in receptor specificity for this residue needs to be determined [52, 53].

N-linked glycosylation of the influenza virus HA protein plays important roles in protein folding and modulates virus pathogenicity and evasion of neutralizing antibodies [54, 55]. In addition, glycans within the vicinity of the RBS region may alter receptor binding affinity and/or specificity. Like other clade 2.1 viruses, the H5N1 viruses from Indonesia contain seven potential N-linked glycosylation sites. N-linked glycosylation at positions 14, 15, 27, 290 and 488 are highly conserved among many HA subtypes [56]. In addition, H5N1 viruses can contain N-linked glycosylation sites at positions 158 and 169. The absence of glycosylation site 158 was linked to human receptor specificity and affinity, and aerosol transmission in ferrets, an animal model representative for aerosol transmission between humans [57, 58]. However, no substitution (i.e., N158X or T160X) removing this glycosylation site was observed in the new Indonesian samples.

Besides changes to receptor specificity and glycosylation patterns, the protein stability of HA is also important for human host adaptation, transmission and possibly virulence [46]. Nonetheless, we did not find in any of the 35 virus isolates any HA substitutions (i.e., H103Y, T315I [57, 58] and Y351H, H352Q, an K387I (H3 numbering) [59, 60] that are known to increase replication and virulence of avian influenza virus H5N1 or H7N9 in mammalian animal models by mediating HA protein stability. However, it is expected that other positions and substitutions within the HA trimer could affect stability and therefore be involved in human adaptation and transmissibility of H5N1 viruses, which would require further research to identify.

Neuraminidase

All of the novel HPAI H5N1 isolates were found to contain the deletion of amino acids between positions 49 and 68 in the stalk region of its NA glycoprotein. This shorter stalk length of NA was previously linked to increase virulence of H5N1 viruses in mammals [61,62,63].

Furthermore, the neuraminidase (NA) protein serves as a target for NA inhibitors (NAI) such as oseltamivir, zanamivir, peramivir and laninamivir, which block the NA enzyme active site to limit influenza virus egress. Eighteen of 35 patients were treated with oseltamivir in our study. However, none of the newly isolated viruses encoded known NAI resistance mutations (i.e., V116A, I117V, E119V, G136K, V149A, R156K, D198N, S246N, H275Y, R293K, N295S (N1 numbering)) [64]. This corresponds with our earlier study showing that acquisition of NAI resistance is extremely rare in H5N1-infected individuals in Indonesia, be it before or during treatment [11]. Of note, Q136H was observed in 14 of the 35 isolates. Although Q136L is associated with reduced sensitivity to zanamivir and oseltamivir, Q136H had no effect on sensitivity to NAI when tested in H1N1pdm2009 or H3N2 [65].

M2

We previously showed that substitutions related to amantadine resistance are common in H5N1 viruses in Indonesia [11] even though amantadine treatment is not administered anymore. The prevalence of amantadine resistance-related substitutions increased over time from 37.5% in 2005, to 95% in 2009 and 100% during subsequent years [11]. Various amantadine resistance substitutions in the M2 protein were also found in all 35 isolates, including V27A (34 viruses), V27T (1 virus), S31G (1 virus), and S31N (5 viruses). Interestingly, isolates in this study collected in more recent years often encode resistance mutations in both positions 27 and 31. These results indicate that Indonesian H5N1 viruses are sensitive to NA inhibitors but resistant to M2 inhibitor, despite the absence of amantadine treatment [11].

Polymerase complex

Besides receptor specificity, polymerase activity is known to be a hallmark for host adaptation and virulence. The polymerase complex is a heterotrimer that consists of the PB2, PB1 and PA subunits. The PB2 protein is an important determinant of virulence and host range. PB2 substitutions such as E627K and D701N [12] can dramatically increase polymerase activity of avian influenza viruses in mammalian cells. In particular, PB2-E627K is a key molecular determinant of host range [66] and a virulence factor during human infection with HPAI H5N1 [67].

Both PB2-E627K (5 of 35 viruses) and PB2-D701N (1 virus) substitutions were observed in a limited number of the novel viruses presented in this study. This is in strong contrast with human H5N1 viruses collected in other geographic regions where PB2-E627K substitution is common [68]. We also checked if there are other known substitutions found in the polymerase complex and nucleoprotein that enhance polymerase activity of avian influenza viruses in human cells reported in previous studies [69]. While there are some genetic variations present in some of these positions, there were no obvious markers of human adaptation or virulence that could be linked to the high case fatality rate of H5N1 infections in humans in Indonesia.

Immune regulatory proteins NS1 and PB1-F2

Both PB1-F2 and NS1 have immune regulatory roles for influenza virus. The full-length PB1-F2 protein (90 aa) inhibits type I interferon response mediated by the mitochondrial antiviral signaling protein [70, 71]. However, the open reading frame (ORF) of the auxiliary protein which occurs in the second ORF of the PB1 gene segment can be truncated or lost [70, 72]. All of the Indonesian H5N1 viruses were found to encode the full-length PB1-F2 protein. Furthermore, the N66S substitution in PB1-F2 known to increase virulence [70] was not found in any of the novel viruses.

On the other hand, the four-amino-acid sequence motif (ESEV) at the carboxyl terminus of NS1 facilitates the nonstructural protein to bind to cellular PDZ-containing proteins that are involved in host cellular signaling pathways [73]. The ESEV motif was found in all of the new Indonesian isolates.

NS1 mutations such as P42S, D87E, L98F and I101M were also found to modulate the virulence of H5N1 viruses [73]. Additionally, substitutions in NS1 (N200S, G205R) as well as NS2 (T47A, M51I) proteins may result in decreased antiviral responses in the host [74]. However, none of these substitutions were found in the 35 H5N1 viruses.

Antigenic analyses show two antigenically distinct virus populations which co-circulate

Besides the presence of virulence and human adaptive markers, it is important to monitor and understand the antigenic properties of circulating influenza viruses. Vaccination is a primary measure to control or prevent H5N1 outbreaks in poultry and could be used to protect humans from H5N1 infections, should these viruses become pandemic. Influenza viruses can easily escape from available vaccines by substitutions in the major antigenic sites on the globular head domain of HA [75]. Here we investigated the antigenic properties of these human HPAI H5N1 viruses.

To characterize the antigenic diversity of Indonesian human influenza H5N1 viruses by HI test, we selected a panel of ferret antisera able to detect antigenic variation between representative viruses [40]. We included ferret antisera against clade 2.1.3.1 (A/Chicken/South Sulawesi/157/2011), different clade 2.1.3.2 viruses (A/Indonesia/5/2005, A/Chicken/Central Java/51/2009, A/Chicken/East Java/121/2010, A/Chicken/West Java/30/2007 and A/Chicken/West Java/119/2010) and against clade 2.1.3.3 (A/Chicken/North Sumatra/72/2010). We determined the antigenic properties of 25 isolates analyzed in HI assays using this panel of ferret antisera. Antigenic cartography was used to visualize the antigenic relatedness of the HPAI H5N1 isolates in a 2D space (Fig. 2). The antigenic map showed that the human HPAI H5N1 viruses from Indonesia clustered into two antigenic groups. The first group of viruses clustered around two representative antisera of clade 2.1.3.2 (A/Chicken/East Java/121/2010, A/Chicken/West Java/119/2010) and contained most of the Indonesian H5N1 viruses (21 of 25 included here) of genetically clade 2.1.3.2, 2.1.3.2a and 2.1.3.2b. Although these viruses do cluster around the A/Chicken/East Java/121/2010 sera, quite some antigenic differences are apparent within this group. The second group (4 of 25) was antigenically closer to representative antisera of clade 2.1.3.1 (A/Chicken/South Sulawesi/157/2011), clade 2.1.3.2 (A/Indonesia/5/2005, A/Chicken/Central Java/51/2009), clade 2.1.3.3 (A/Indonesia/North Sumatra/72/2010). This finding is similar with a previous study, describing different antigenic clusters within avian H5N1 viruses in Indonesia isolated from poultry [40]. This study identified a small number of residues immediately adjacent to the RBS within the antigenic sites in the globular head domain of HA, which are primarily responsible for antigenic changes. We investigated genetic diversity at these positions. Amino acid differences were identified at these six antigenically important positions located near the receptor binding sites [40]: 129, 133, 151, 183, 185 and 189 (H5 numbering as there is no equivalent for position 129 in H3N2 viruses; Table 2). All viruses antigenically clustering into the first antigenic group possessed residues S129, S133, I151, N183, A185 and M189, except for isolates 11,046 and 12,452 that contained T183 and I189, respectively. The N183T substitution does not seem to have an antigenic effect, although the M189I could have a small antigenic effect as indicated by the placement on the outside of the cluster. The virus isolates of cluster 2 contained residues S129, S133, I151, D183, A185 and R189, typical for A/Indonesia/5/2005 antigenic-like viruses. The study by Koel et al. has previously shown that substitutions D183N and R189M are indeed responsible for the antigenic differences between these two clusters. Isolate 10,364 contained A133; however, this substitution does not seem to have major antigenic effects as this isolate clusters with other viruses. It was previously shown that a combination of substitutions at position 133 and 185 is necessary to result in antigenic changes [40].

Fig. 2
figure 2

Antigenic map of Indonesia human H5N1 viruses. Circles and open squares indicate the positions of viruses and antisera, respectively. Both axes represent antigenic distance: one square on the antigenic map represents a distance of one antigenic unit, corresponding to a twofold difference in the HI assay. The antigenic map was generated using antigenic cartography, a method that uses multidimensional scaling algorithms to place virus and antiserum points in a 2D space such that their relative position in the map reflects the HI titers with minimal error. The distance between a virus-and-antiserum pair is inversely related to the HI titer of the virus to that antiserum. The color coding of the human HPAI H5N1 isolates is based on their year of isolation as depicted in Fig. 1. Virus isolate names and antisera are abbreviated to isolate number/year

Full size image
Table 2 HA protein sequence diversity of important antigenic sites
Full size table

These data showed that viruses belonging to two distinct antigenic groups infected humans in Indonesia. Interestingly, the presence of viruses from different years of isolation in both clusters (as indicated by the color coding of the viruses in Fig. 2) indicates that these different antigenic variants were co-circulating. Full protection in humans would therefore have likely required a multivalent vaccine, including at least A/Indonesia/5/2005- and A/chicken/East Java/121/2010-like viruses, which are currently under development or approved for human and/or poultry use [17, 76].

Discussion

Indonesia has suffered numerous HPAI H5N1 virus outbreaks in poultry farms, live bird markets and backyard poultry, which have resulted in 200 reported human cases with a case fatality rate of over 80%. This high case fatality rate is in sharp contrast with lower case fatality rates in countries affected by other genetic clades of HPAI H5N1 viruses. Despite the large number of cases, limited genetic and antigenic information is available for HPAI H5N1 isolates from human cases in Indonesia. Here we present 35 HPAI H5N1 virus isolates collected from infected humans in Indonesia as part of outbreak investigations between 2008 and 2015. These novel HPAI H5N1 viruses were genetically typed as clade 2.1.x viruses (i.e., 2.1.1, 2.1.2, 2.1.3) using hemagglutinin sequences in previous analyses [34]. Clade 2.1.1 viruses were mostly identified from poultry outbreaks between 2003 and 2004, while the viruses from clade 2.1.2 and 2.1.3 circulated at the beginning of 2004 in birds and were also isolated from humans [77]. As frequency of clade 2.1.2 viruses identified in Indonesia started to decline in 2005, with none detected after 2008 [78], clade 2.1.3 viruses continued to circulate. Phylogenetic analyses of the HA gene showed that all of the new 35 isolates cluster with other 2.1.3.x viruses, similar to earlier published human and avian sequences from Indonesia. These new isolates are closely related to previously sequenced avian HPAI H5N1 viruses circulating in Indonesia.

We recently showed that the high case fatality rate of human Indonesian H5N1 cases correlated with viral load prior to treatment and increased from 65% in 2005 to 100% since 2012 [11]. We found that this high case fatality rate coincided with the high prevalence of amantadine resistance-conferring M2 substitutions; however, no mechanistic explanation for the role of such substitutions in virulence of HPAI H5N1 is known yet. The aforementioned study did not include any further sequencing data looking into specific virulence and human adaptive markers. Our analyses of the full genomic sequences for the 35 isolates did not indicate any potential genetic changes that might explain the increase in case fatality rate and virulence over time. We did not find any known genetic markers associated with human adaptation or virulence in the HA or polymerase complex genes. There were no changes in the RBS region of HA that was indicative of a switch towards the human-type receptor. Although some genetic diversity was observed in the polymerase genes, well-known substitutions such as PB2-E627K and PB2-D701N, which are often selected upon infection of humans and affects the virulence of avian influenza viruses such as H5N1 [12, 66, 67, 79], were not commonly found in the new samples. However, there could also be other currently unknown adaptive substitutions present in the new human samples. Confirmatory investigations into the effects of these substitutions on the activity of the polymerase complexes of both human and avian HPAI H5N1 viruses should be done in future studies.

Further sampling and research, involving the collection of more full genome sequences from avian and human viruses as well as both in vitro and in vivo characterization of virus replication and pathogenicity, are warranted to determine if Indonesian HPAI H5N1 viruses are indeed more virulent than H5N1 viruses of other genetic clades circulating in other geographic areas. This should also address whether there is a specific selection for more virulent viruses in humans only or whether Indonesia HPAI A/H5N1 viruses are more virulent in general, also in poultry in Indonesia, resulting in the consequence that zoonotic events happen with more virulent viruses, resulting in higher case fatality rates. Further characterization of virus isolates from different periods of time will have to show if more recent viruses are indeed more virulent and whether this could be contributed to specific molecular markers.

A possible reason why human H5N1 cases have declined is the implementation of large-scale vaccination of poultry against H5N1. Most countries affected by H5N1 virus outbreaks, including Indonesia, have implemented poultry vaccination as a key strategy for the control of H5N1 infections. Currently, the H5N1 vaccine for poultry in Indonesia is an inactivated bivalent vaccine containing H5N1 viruses belonging to clades 2.3.2 and 2.1.3. As for the selection of human seasonal vaccine strains, antigenic analyses are required to understand and predict vaccine effectiveness. From the current study, antigenic analyses of human H5N1 viruses from 2008 till 2015 identified two antigenic groups of human clade 2.1.3 viruses that co-circulate in Indonesia. Based on a previous study by Koel et al., the antigenic differences could be explained by alternative amino acids present at several key residues at the rim of the receptor binding site.

No evidence of antigenic change over time was observed or association with geographical location. Therefore, a combination of two of the current available pre-pandemic human H5N1 vaccines, A/Indonesia/5/20 05 and A/Indonesia/NIHRD11771/2011, would have been required to optimize protection against the two different antigenic groups.

In summary, we performed genetic and antigenic analyses of H5N1 influenza viruses isolated from humans between 2008 and 2015. We observed low levels of genetic diversity and only sporadically prevalence of known substitutions associated with human adaptation and virulence (e.g., PB2-627K). However, the analysis only captured the majority variants and did not include the presence of minority variants present during infection. Additionally, we have limited our genetic analyses to known substitutions only. To ascertain and better understand the high mortality associated with human HPAI H5N1 virus infections in Indonesia, it is essential to perform more in-depth analysis of genetic diversity during human infections with HPAI H5N1 virus and to functionally characterize the observed substitutions. Furthermore, our data showed that two antigenic groups co-circulated in Indonesia, with no evidence of antigenic change over time. A combination of available pre-pandemic vaccines was required to be protective against circulating viruses of study period.