Variability analysis and inter-genotype comparison of human respiratory syncytial virus small hydrophobic gene

Human respiratory syncytial virus (HRSV) belongs to the Orthopneumovirus genus of the Pnemoviridae family [1] and has non-segmented RNA of negative polarity as its genome. The genome is approx. 15,200 nucleotides (nts) long, comprising 10 genes that encode for 11 proteins. The HRSV strains are divided into groups A and B based on genetic and antigenic differences [2].

There are three proteins on the surface of mature HRSV virions: glycoprotein (G), fusion protein (F) and small hydrophobic (SH) protein. SH protein is present at very low amounts [3, 4]. G and F are important for attachment and fusion with the target cell, respectively [5, 6] and they can elicit the production of neutralizing antibodies [7, 8].

SH protein is a short transmembrane protein (typically 64 amino acids for group A and 65 for group B) that is anchored by a hydrophobic signal-anchor sequence near the N-terminus, with the C- terminus oriented extracellularly. The SH ectodomain is only weakly immunogenic [4, 9]. In infected cells, SH protein accumulates mostly in lipid rafts of the Golgi and of the endoplasmic reticulum [3]. It is a viroporin, forming a transmembrane pentameric ring that functions as a channel for cations and small molecules [10, 11]. HRSV SH protein blocks or delays apoptosis through inhibition of the TNF-α signalling pathway [12] but the mechanism is still not fully elucidated.

In molecular epidemiology studies, HRSV surveillance and genotyping are based on sequences of the second, C-terminal hypervariable region of the G gene (HVR2). In frame duplications of 72 and 60 nts have occurred within this region, leading to emergence of group A strains belonging to genotype ON1 (which originated from NA1, and was firstly detected in 2010 [13]) and of group B strains belonging to BA genotypes, respectively.

To what extent do the strains identical in HVR2 differ elsewhere in their genomes is rarely investigated. In the analysis by Agoti et al. such strains differed by at least 1 and up to 9 nts across the complete genome [14].

Next to HVR2, the most diverse HRSV genomic region is the SH gene [15]. Unlike all other HRSV genes, the SH coding region constitutes less than 50% of the complete gene length (the coding regions of other HRSV genes constitute ca. 75–99% of the gene length). Following our previous analysis of HVR2 sequences of strains detected within a limited geographical area (the Zagreb region) and limited time frame (March 2011 to March 2014) [16], the goal of this study was to investigate whether virus variability and inter-genotype differences observed for HVR2 are also present in the SH gene. Within this pool of viruses, rapid and significant genetic changes have occurred in HVR2 [16]. The subset of strains analysed in this research belong to 3 HRSV A (NA1, ON1 and GA5) and 2 HRSV B genotypes (BA9 and BA10) and it included all the genotypes that were detected in Zagreb, 03/2011–03/2014. These are also the most common HRSV genotypes detected worldwide in the last decade.

Due to the fact that SH gene UTRs are relatively long (when compared to the coding region of the gene) and that the majority of mutations occur in them, HRSV SH gene is characterized by high variability while the SH protein remains highly conserved.

Besides the complete SH transcript, the genomic segment we analysed also included its two downstream regions, stopping immediately before the start codon of the G protein mRNA. The last two regions were included because they are known to play a role in the transcriptional regulation of both upstream and downstream genes [18–20].

We analysed whether a difference in the diversity of the SH/G junction (i.e. of the sequences between the stop codon of the upstream gene and the start codon of the downstream gene) among HRSV genotypes could be observed. A sequence characteristic only for GA5 strains was observed, but not all previously reported GA5 strains possess this sequence [20]. Whether it represents an evolutionary novelty, specific for recent GA5 strains requires further investigation (the oldest GA5 strain in our analysis was from Mexico, from 2004).

The overall nucleotide diversity of SHseg was investigated using two methods: conservation percentage that is not influenced by the number of identical sequences; and p-distance that is, and therefore, defines more precisely diversity within this particular HRSV population. The results of both methods were concordant. The level of genetic diversity observed in SHseg was not comparable to the one detected in HVR2, except in the 5’ UTR of the genotypes with prolonged G gene (ON1 and BA9) which gained striking dominance within the short timespan we analysed. On the contrary, the other SH UTR and SH-G intergenic region showed less variability in ON1 and BA9 strains, when compared to NA1 genotype.

The current HRSV genotyping and molecular surveillance system is based on HVR2 [21, 22], a quite short genomic segment. It partitions HRSV strains in over 15 group A and over 20 group B genotypes. Regarding phylogeny, SHseg in not as informative as HVR2. E.g. in analyses based on HVR2, ON1 and NA1 strains are placed on different branches, even when one copy of the duplicated segment is excluded from the analysis [23]. Still, basing molecular epidemiology studies only on HVR2 may be quite misleading. In our pool of HRSV strains, over 20% of them would be incorrectly assumed to be identical to another strain if only the HVR2 region was analysed. An alternative HRSV genotyping scheme has recently been proposed [24], based on the complete G protein ectodomainwhich is twice longer than HVR2. Although this would increase the discrimination between strains, it is still focused solely on one protein gene. The G gene is non-representative of the HRSV genome, as it is marked by elevated substitution rates (attributed to relaxed selective constraints [25]) and a flip-flop substitution pattern [26]. Broadening the target genome sequence to SH (whose length in our analysis was 401–411 nts in group A, 412 or 413 in group B) or only to SH 5’ UTR (in our analysis, this was a region of 83 or 84 nts in group A, 85 in group B) would provide better insight into the interrelation between molecular epidemiology and molecular and evolutionary dynamics.

After the translation of our genomic sequences, only 15 different protein sequences were obtained. Although group A and group B strains diverged approximately 350 years ago [27], the majority of amino acids which differ from the consensus in group A have the same identity as corresponding amino acids in group B strains. This indicates that the variability of the SH protein is not only low, but is also restricted. Based on the fact that HVR2 has the highest evolutionary rates, we concentrated on strains possessing the same HVR2 nucleotide sequences, but still differing in their SH protein sequences. Those are the strains most likely to be highly similar and therefore we hypothesized that the observed differences could be a result of selective pressure. Out of 4 sets of strains with differences in the SH ectodomain, in 3 sets the differences were in SH protein region 49–51. In fact, besides S35F the only amino acid differences observed in group A strains which differed from the group B consensus, were located within 49–51 region. In our viral pool, which was very limited both geographically and temporally, 6.7% of group A strains and 51.3% of group B strains possessed mutation in region 49–51.

Our hypothesis that SH region 49–51 is an antigenic site still needs confirmation, but the same amino acids were shown to be highly variable by Lima et al. [28]. Their analysis was performed on strains belonging to HRSV genotypes different from ours and detected in Brazil in 2004–2005. Besides finding T49I variability (as observed by us), P and A at position 49 were also reported among group B strains. Chen et al. [29] detected variability in the same region in American HRSV strains from 1998 to 1999 [29]. In group A, the most common amino acid at position 49 was V (the same as in our analysis), but I was also observed. Among group B strains, at position 51, besides the most common amino acid C and F (the same as in our analysis), Chen et al. also identified Y [29].

A deletion of 9 nts coding for amino acids 49–51 was found in the SH gene sequence of an NA1 strain from January 2013. The deletion was found only in this strain that caused a bronchiolitis in a 2-month-old child. The deletion was confirmed by 3 independent RNA isolations and preparations of DNA for sequencing. During the same HRSV epidemic, strain RSV3234–12 (from November 2012) was detected in a sample from a 4-month-old child with bronchiolitis. The two strains are identical in the SH, G an F genes, except for the fact that RSV3234–12 does not possess this deletion in the SH gene. Other genomic regions were not sequenced. In group A strains, H51 (which is missing from RSV5–13) is a crucial residue in the regulation of SH ion channel activity [30]. Together with H22, H51 plays a key role in the opening and closing mechanism of the SH pentameric pore, since they are located in strategic places within the chains, close to the N- and C-termini [30, 31]. An alternative amino acid at position 51 (Y, a non-basic amino acid) has been reported before, in a Brazilian GA5 strain from 2004 [28], but so far there were no reports of the deletion observed in RSV5–13. In vitro experiments regarding the functionality of the SH protein from RSV5–13 or modelling studies that would compare its structure to other HRSV SH proteins have not been performed.

Basing HRSV molecular epidemiology studies solely on HVR2 largely underestimates the complexity of circulating virus populations. In strain identification, broadening of the genomic target sequence to SH gene would provide a more comprehensive insight into viral pool versatility.

Unlike its gene, SH protein is characterized by intra-group conservation and by highly restricted inter-group variability except for amino acids 49–51, implicating that this protein region may be relevant for the antigenicity of the virus.