Proposals for the classification of human rhinovirus species A, B and C into genotypically assigned types

1Human rhinoviruses (HRVs) are highly prevalent respiratory pathogens traditionally associated with mild and self-limited upper respiratory tract infections (‘colds’). The replacement of relatively insensitive virus isolation methods with molecular methods for HRV detection has contributed to a recent reappraisal of their incidence of infection and genetic diversity. HRVs are now known to be intimately linked to the development and exacerbations of respiratory diseases such as asthma and responsible for more severe disease manifestations such as bronchiolitis in young children and in the immunosuppressed. Genetic characterization of HRVs detected by molecular methods has revealed much greater diversity than previously described; approximately one-third of HRV infections are now known to be caused by a third group of rhinoviruses (species C), uncultureable in vitro and therefore entirely missed by traditional virus isolation methods.

2HRVs are members of the genus Enterovirus of the family Picornaviridae, a large virus group that additionally contains human enteroviruses and a number of viruses infecting non-human primates, pigs, cattle, sheep and other mammals. The genus Enterovirus has recently been divided into a total of 12 species based on their genetic divergence, sites of replication and disease associations (http://www.picornaviridae.com). Among these species, three are designated Rhinovirus A, Rhinovirus B and Rhinovirus C, a nomenclature that reflects their primary tropism for the respiratory tract. Currently, humans are the only known natural host, although it is possible that primates could also be susceptible. Rhinovirus B and most Rhinovirus A variants use the intracellular adhesion molecule-1 (major) receptor for cell entry, whilst a subset of Rhinovirus A types use the low-density lipoprotein receptor. The receptor(s) used by HRV species C (HRV-C) types is currently unknown. HRVs are most commonly transmitted by the respiratory–salivary route, both by person-to-person contact and by airborne transmission. In temperate countries, infections occur primarily in two peaks, the first between April and May and the second between September and October.

3HRVs within each species show remarkable genetic and, where determined, antigenic heterogeneity. Species A and B isolates were originally characterized serologically, with 75 and 25 serotypes being classified by cross-neutralization assays by 1987, and prototype strains of each were lodged at the American Type Culture Collection (ATCC). However, as well as being laborious, time-consuming and requiring extensive panels of antisera specific to each serotype, classification by neutralization cannot be used for species C rhinoviruses because of their refractoriness to cell culture.

4Human enteroviruses (EVs) are classified into four species (Enterovirus A, Enterovirus B, Enterovirus C and Enterovirus D) and, like HRVs, isolates have also been traditionally characterized by neutralization assays. Some neutralization determinants are located on the external VP1 protein in the virus nucleocapsid, and it has been shown that sequence divergence in this region correlates with serotype designation of poliovirus isolates and other enteroviruses. Indeed, a threshold of 25 % nucleotide and 12 % amino acid sequence divergence generally corresponds to the previous serological division of EVs into members of the same or different serotypes. As well as providing a simpler method for EV type identification than the previously used serological assays, VP1 sequencing is now routinely used for identification and assignment of new EV types without serological confirmation of antigenic differences.

5In the current study, we investigated whether a comparable threshold of sequence diversity in the VP1 of HRVs would enable a similar genotypic division of HRV-A, -B and -C variants into types and additionally whether this correlated with sequence divergence in the VP4/VP2 region (routinely used for HRV typing). This investigation was prompted by the increasing use of molecular methods to characterize HRVs and the consequent lack of isolates or antisera with which to antigenically characterize them and, secondly, by the difficulty in culturing species C HRVs required for neutralization assays.

6This investigation follows on from previously published proposals for the classification of HRV-C based on a 13 % divergence threshold in VP1. These have been adopted by the Picornaviridae Study Group as the primary mechanism for defining HRV-C types. We have extended this classification system for species A and B rhinoviruses. The previously used term ‘serotype’ reflects the identification and classification of HRV variants by antigenic properties and therefore assignment by genetic comparison methods requires a different or combined term. Following the terminology now used for enteroviruses, we have adopted the term ‘type’ throughout the study to represent those HRV variants that have been identified and classified by either cross-neutralization or genetic comparisons.

7In order to determine whether a distinct threshold that divided pairwise p-distance comparisons in the VP1 region into intra- and inter-type values existed within each HRV species, distributions of pairwise nucleotide p-distances were constructed from a total of 435 HRV-A, 133 HRV-B and 206 HRV-C VP1 sequences (Fig. 1; sequences used for species C represent an expanded dataset from the one analysed previously; ). These showed a maximum within-species nucleotide p-distance of 39.4, 33.9 and 42.9 % for species A, B and C, respectively, and minimum values of between 10 and 14 % divergence. Detailed inspection allowed identification of minimum values in the distribution and assignment of a nucleotide divergence threshold for each species (Fig. 1, right-hand panels). A threshold of 13 % was identified for both HRV-A and -C sequences, whilst the less-divergent HRV-B isolates showed a minimum value that supported a 12 % threshold.

8The total datasets of available VP1 sequences contained unequal representations of different types within each species, ranging from one to a maximum of 13 for HRV-A (A49), 10 for HRV-B (B69) and 22 for HRV-C (C3). Particularly for pairwise distance comparisons, these may potentially distort the shape of the distribution of inter-type distances because comparisons with certain types are over-represented in the distance totals. An additional analysis was therefore carried out with only single sequences from each type (grey shading superimposed in Fig. 1, left graphs). For each, prototype strains of HRV-A and -B were those of the first full-genome sequence of a deposited ATCC isolate that was deposited in GenBank (highlighted in red in Table S1 in JGV Online). Prototype strains of HRV-C were those listed previously (; Table S1). This closely reproduced the analysis of the whole VP1 dataset in the shape and range of the inter-type distribution of distances (Fig. 1). Whilst distance distributions between VP1 sequences of HRV-A and -C showed a single peak, those of HRV-B showed a bimodal distribution. The peak between p-distance values of 0.19 and 0.26 corresponded to comparisons between HRV-B types within each of the four basally branching phylogenetic clusters (indicated by grey-outlined boxes in Fig. 2). The larger peak thus corresponded to comparisons between variants in different clusters.

9Phylogenetic trees were constructed of the VP1 region of HRV-A, -B and -C (Fig. 2). All HRV sequences fell into bootstrap-supported monophyletic groups that closely matched types assigned by sequence distances. These groups were numbered either by clustering with a known prototype strain (HRV-A/B 1–100), by order of submission of full-genome (HRV-A101–106, HRV-B100–104 and HRV-C1–11) or VP1 sequences (HRV-C12–51). The vast majority of sequences within each terminal clade showed intra-clade VP1 divergence substantially less than the proposed thresholds of 13 % in HRV-A and -C and 12 % in HRV-B (Table 1). Most of these also showed a minimum inter-clade VP1 divergence (with a different assigned type) that was much greater than the proposed threshold.

10Consistent with the thresholds identified for each species, three new types of HRV-A and three of HRV-B have been proposed previously; HRV-A101, HRV-A102, HRV-A103, HRV-B100, HRV-B101 and HRV-B102 (http://www.picornaviridae.com). In addition, the current analysis provided evidence for three further new HRV-A types (tentatively labelled HRV-A104–A106) and two novel HRV-B types (HRV-B103 and -B104), labelled in red in Fig. 2. The additional sequence data provided in this study and published elsewhere for HRV-C now create a phylogenetic tree that can be divided into 51 distinct clades. These have been assigned as types based on sequence divergence in the VP1 region (Table 1), adding a further 18 types to the total of 33 first defined in 2010 (http://www.picornaviridae.com).

11Although several of the putative new HRV-A and -B types included showed minimum VP1 distances far greater than the assignment thresholds shown in Fig. 1, several were close to or below the 13 and 12 % levels (Table 1). As examples of the former, proposed new HRV-B types B100, B101 and B102 showed minimum pairwise p-distances in VP1 of 0.214, 0.190 and 0.227 from their nearest confirmed HRV-B type (HRV-B35). Pairwise nucleotide p-distances between these three new HRV-B types ranged from 0.163 to 0.171. All proposed species C types were similarly clearly above the assignment threshold of 13 % for this species (Table 1). However, other variants proposed as new types showed minimum distances that were much closer to the threshold; for example, the six available sequences from the novel group putatively designated HRV-B103 showed VP1 distances from their nearest neighbour type (HRV-B6) ranging from 0.113 to 0.148. However, for these two types, only two pairwise comparisons fell below the 12 % threshold and both involved sequence F110_9318.

12The newly designated HRV-B104 represented a division of a previously single type, HRV-B52, based on the existence of two discrete phylogenetic clades with VP1 p-distances of 0.121–0.142 between them. In this case, the clade containing the original ATCC isolates, including strains FJ445188 (GenBank accession no.) and EF173424, has retained the HRV-B52 designation, whereas the clade containing field isolates FJ445137 and JF781506 has been reassigned as HRV-B104. Similar variability was observed in species A where two of the proposed new types showed distances close to the assignment threshold (HRV-A104 and -A105).

13Analysis of pairwise VP1 divergence between previously recognized HRV-A types also revealed several examples of VP1 distances below the 13 % threshold. For example, HRV-A8 and HRV-A95 showed an inter-type VP1 p-distance range of 0.015–0.019. We propose that all sequences of these two types are reassigned as HRV-A8, a proposal consistent with their previously documented serological cross-reactivity (Fig. 3a).

14Similarly, HRV-A29 and HRV-A44 have previously demonstrated cross-reactivity and prototype strains showed an inter-type VP1 p-distance of 0.0724 (Fig. 3a). We propose that these sequences be reassigned as HRV-A29. Interestingly, contemporary sequences that belonged to HRV-A29 by VP1 p-distance from the prototype strains (maximum p-distance of 0.127) showed a maximum divergence of 0.159 when compared with each other. All contemporary HRV-A29 sequences grouped separately on phylogenetic trees from both the HRV-A29 and HRV-A44 prototype sequences (Fig. 3a).

15HRV-A54 and HRV-A98 are a third type pair that are more closely related than expected, with a pairwise nucleotide VP1 p-distance of 0.107–0.124 (Fig. 3a). Although cross-neutralization properties of these two types have not been analysed in previous studies, we propose that all variants classified as these two types are reassigned as HRV-A54.

16In addition to those type pairs that have been combined into single type assignments, a number of other HRV types displayed intermediate divergence values that violated the proposed thresholds. These included HRV-A1A and -A1B, which are serologically cross-reactive but showed several pairwise p-distances that fell above the VP1 threshold (Fig. 3). We propose that members of these two clades retain their assignment as HRV-A1, although this should be reviewed as more sequences and confirmatory serological data become available. As a contrasting problem, the type pair HRV-A25 and -A62 showed VP1 pairwise distances ranging from 0.093 to 0.103 from each other (Fig. 3b), below the assignment threshold, despite showing no cross-neutralization. We propose leaving these assignments unchanged pending confirmation of the serological relationship between these two types. Similarly, prototype isolates of types HRV-A36 and -A89 (GenBank accession nos DQ473508 and M162488) showed VP1 p-distances of 0.1292, just below the assignment threshold (Fig. 3b). For these types, the HRV-A36 antiserum can neutralize HRV-A89 (as well as HRV-A50 and HRV-A58). However, the antiserum to HRV-A89 cannot neutralize HRV-A36 and almost all pairwise distances between more recent isolates of HRV-A36 (n = 9) and HRV-A89 (n = 5) are above the assignment threshold. We propose to retain their current designations as HRV-A36 and HRV-A89.

17Only one HRV-C type showed p-distance values that violated proposed thresholds. HRV-C8 strains were divided into two clades, one previously classified as HRV-Cpat28 and both containing contemporary sequences (Fig. 3b). Although the VP1 p-distance between the two clades ranged from 0.129 to 0.132, we propose that these two clades should remain classified as HRV-C8 on the basis of borderline VP1 divergence.

18Despite these examples, the number of pairwise distances in the overall distribution that were close to the assignment threshold (and therefore potentially ambiguously assigned) was small (Fig. 4). The majority represented comparisons between type pairs that should be combined (Fig. 4, green). Only six HRV-A types showed any between-type distances that were below the threshold (pink), and there were three with within-type distances that were above (Fig. 4, purple). As described above, their current classification can incorporate phylogenetic data and any available information on neutralization properties in addition to VP1 distances. Combined, these may potentially resolve any remaining uncertainties over type assignments.

19At the conclusion of this analysis of all available HRV VP1 sequences, totals of 77 HRV-A types, 29 HRV-B and 51 HRV-C types have been listed in full in Tables S1 and S2). Novel strains of HRV were classified based on date of submission to the Picornaviridae Study Group, rather than the earliest isolated strain.

20The well-conserved and readily amplified VP4/VP2 region of HRV is commonly used in studies of its epidemiology and clinical associations, with over 3900 HRV sequences from this region in GenBank, over 10 times the number of VP1 sequences (n = 384). To investigate the value of VP4/VP2 sequences for HRV type identification, we performed phylogenetic analysis on the entire dataset and compared groupings of sequences of known type to those observed in VP1 (Fig. 5).

21The majority of HRV strains grouped congruently between VP1 and VP4/VP2. Only one sequence showed evidence of phylogenetic incongruence; JX291115 is a member of HRV-C51 type in the VP1 region, with a pairwise nucleotide p-distance of 0.01 to JF317015 (HRV-C51 prototype strain). However, in VP4/VP2, this sequence showed a VP1 p-distance of 0.146 from JF317015 and lay closest to HRV-C15. Its status as a recombinant sequence was demonstrated by triplet analysis using the rdp software package. In contrast, comparable analysis of all other HRV-C and HRV-A/B variants with both full VP4/2 and VP1 sequences available by rdp and gard (164 types in species A, 70 in species B and 145 in species C) demonstrated no evidence of recombination between regions.

22We previously proposed a divergence threshold of 10 % for identifying different HRV-C types using VP4/VP2 sequences. Using the current much larger dataset of HRV-C sequences and a shorter fragment (positions 616–1004, based on typical coordinates of sequences available in GenBank), pairwise p-distance distributions were plotted, along with those for species A and B (Fig. 6). Detailed inspection of distributions of pairwise distances showed possible assignment thresholds of 10.5, 9.5 and 10.5 % for HRV-A, -B and -C, respectively, based on minimum values. The shorter length and greater sequence conservation than VP1 resulted in less-well-defined thresholds and a greater proportion of pairwise comparisons close to the assignment threshold, which may cause problems with type identification.

23For example, four type pairs showed VP4/VP2 divergence that was below the proposed threshold [HRV-A31/A47 (0.083), HRV-B52/B104 (0.093), HRV-C4/C49 (0.100) and HRV-C32/C39 (0.090)] even though they were previously assigned as different types based on VP1 divergence and separate clustering by phylogenetic analysis (Table 1 and Fig. 2). More problematic was the larger number of HRV variants of the same type with VP4/VP2 divergence above the assigned thresholds (22, 10 and six types of HRV-A, -B and -C, respectively; Table S3). Most of these types contained several contemporary sequences that had diversified substantially from their corresponding prototype sequence. In these cases, p-distances above the threshold occurred between contemporary isolates diversifying on separate evolutionary lineages, while retaining similarity in the intra-type range with their corresponding prototype type strain [VP4/VP2 distances of <10.5 % (species A) or 9.5 % (species B)]. This is similar to the pattern observed within the VP1 region of contemporary strains of HRV-A29.

24The original HRV-C classification included 28 PATs, which were temporarily assigned on the basis of VP4/VP2 sequence divergence alone until a matching VP1 sequence was obtained. Of these original 28 PATS, 16 have since been reassigned as new (confirmed) types based on VP1 divergence (new HRV-C type designations listed in Table S1 and S2), whilst one was combined into an existing type (HRV-C8/Cpat28; see above). VP1 sequences are unavailable for the remaining 11 PATs (Table S4).

25Eight further divergent VP4/VP2 groups were identified among species A (n = 4), and species B (n = 4) (Table S4). A further three provisionally assigned HRV-C types were also detected (designated HRV-Cpat29-pat31). All showed divergence above the VP4/VP2 threshold proposed for their species. These await sequence data from VP1 to confirm their assignment as different types.

26The distribution of HRV types (identified using the larger dataset of VP4/VP2 sequences and excluding PATs) was compared for the UK, the rest of Europe, the USA, Asia and other regions (Australia, the Middle East and Africa). Each geographical zone contained similar relative frequencies of HRV-A, -B and -C (Fig. 7). Individual types showed considerable variability in detection frequency, although those that were common showed similar detection frequencies in different geographical areas. Similarly, none of the prevalent HRV types showed any geographical restriction, findings that suggest that most or all HRV types circulate freely worldwide.

27The ability to perform genotyping directly from clinical specimens removes the need for cell culture passage and therefore greatly reduces the time required for type identification for all three HRV species. It additionally allows variants refractory to virus isolation, such as species C, to be identified. This parallels the transition from isolation and neutralization assay to VP1 sequencing for the identification of EV (sero)types. Because of its shorter length, greater sequence conservation and the possibility of using the same set of primers for amplification of all three HRV species, the VP4/VP2 region has been most commonly used for genetic characterization of HRV in clinical samples. Comparison of its predictive value compared with type assignments based on VP1 sequences provides evidence for its suitability of HRV type identification along with some restrictions. The latter include uncertainties in type identifications of variants showing sequence divergences from other HRV sequences close to the threshold assigned to that species. Examples of both inter-type distances below this threshold and intra-type distances above have been presented (Table S3), although patterns of phylogenetic clustering of variants provide further contributory information for type identification.

28Overall, however, the thresholds we have proposed divide HRV sequences into groups that correspond most consistently with phylogenetic groupings (with the small number of exceptions we have specifically itemised). Placing thresholds in any other position in the distance distributions would create groupings that conflict with those defined by phylogenetic analysis and by serology relationships. Placing the threshold at a minimum point in the distribution additionally minimizes the number of conflicting between- and within-group distances within an assigned type.

29The identification and definition of HRV types by using genetic data within the capsid region has the potential to revolutionize HRV typing. The ease with which this can be carried out, compared with traditional virus culture and neutralization assays, should make the endeavour of classifying all detected HRV sequences much more feasible and attractive. A more extensive global or sentinel system of type identification should allow large-scale investigations of HRV epidemiology, transmission and evolution; sequence data for comparisons are available worldwide almost instantly and will be of considerable value in surveillance for the emergence of new types of species with differing clinical presentations from those traditionally attributed to HRV.

30In addition to type identification, the definition of robust nucleotide divergence thresholds for the definition of new EV types (25 %) now provides the primary means for the assignment of new EV types. This has been used for several years and has created the assignment of 36 new EV types numbered from EV-73 onwards. We propose an equivalent system for the assignment of new HRV types in which criteria similar to those proposed for the assignment of species C rhinoviruses are extended to HRV-A and -B. HRV does, however, differ from EVs in showing lower divergence thresholds for type assignments (12–13 % nucleotide sequence divergence in VP1) and a lack of capsid recombination that enables VP4/VP2 sequences to be used at least for type identification. The criteria for HRV have been formulated as follows:

31An HRV type should be phylogenetically distinct from other HRV types.

32Analysis of divergence in the VP1 region should include at least 90 % of the full VP1 coding region for each variant.

33In the VP1 region, an HRV type should have at least 13 % (HRV-A), 12 % (HRV-B) or 13 % (HRV-C) nucleotide divergence from all other HRV types.

34In cases where minimum nucleotide divergence in VP1 of a putative new type lies close to the threshold for that species (±1 %), patterns of phylogenetic grouping with existing types may contribute to its assignment as would any available information on serological cross-reactivity.

35VP4/VP2 sequences are predictive of HRV species and type and can be used for type identification in epidemiological studies.

36HRV variants divergent from all known VP4/VP2 sequences (i.e. above type-assignment thresholds for the identified species in this region) but without an available VP1 sequence may be designated a PAT.

37Subsequent determination of VP1 sequences for PATs would allow their reassignment as new types if they show VP1 sequence divergence above the threshold for that species).

38New HRV type names should be numbered sequentially and separately for each species – the next assignments for the three rhinovirus species would be HRV-A107, HRV-B105 and HRV-C52.

39The ‘prototype strain’ of each new HRV type should be the first reported full-genome sequence or, if this is unavailable, the first reported full VP1 sequence.

40The Picornaviridae Study Group should continue to oversee assignment of new HRV types, as the group includes members with specific expertise in EV type assignment and a number of scientists currently active in HRV research.

41The use of only capsid-coding regions in type-assignment criteria should not detract from the importance of continued investigation of other genomic regions, in particular where these may contribute to the phenotype and disease associations of HRV. In addition, the use of the VP1 region for definition of new HRV types should not discourage the widespread continued use of 5′-untranslated region and VP4/VP2 screening protocols. Screening with these relatively conserved regions allows a much greater opportunity for discovery of previously unknown HRV types through further sequence characterization of VP1 and other genome regions. Although this represents, to our knowledge, the most comprehensive survey of HRV sequence data currently possible, the guidelines should be subject to further review as additional data becomes available.

42In all three species, the phylogeny of the full VP1 region reliably separated isolates of known differing type into distinct bootstrap-supported clades. This is consistent with the presence of a large number of the neutralizing immunogenic sites in HRV-A and -B that underlie their serological properties. The difficulties with culturing species C rhinoviruses have to date prevented their serological properties from being directly investigated. It is possible that the sinus mucosal organ culture system used to culture HRV-C15 might prove more broadly applicable for other species C variants. Limited serological reinvestigation of a number of previous HRV-A and -B type assignments would be useful for isolates where cross-reactivity is uncertain (such as HRV-A54/A98) or conflicts with their sequence divergence in VP1.

43There is a possibility that a smaller fragment of VP1 could be suitable for genotyping in HRV. reported that a 450 nt segment at the 3′ end of VP1 was effective for EV type identification and the results had a 100 % correlation with neutralization results. A second study showed full correlation of a 303 nt stretch of VP1 with neutralization data, but this was only applied to 59 strains. In addition, a simple, fast and effective method of typing EVs by pyrosequencing has recently been developed. This may become extremely useful as pyrosequencing becomes more widely available worldwide. Whether a shorter region than the whole VP1 would be similarly predictive of HRV types and could be used for future type assignments requires further investigation.

44VP4/VP2 sequences generated from several studies of HRV molecular epidemiology in the UK (Simmonds et al., unpublished data) (431 HRV-A and 113 HRV-B) were divided into groups based on phylogenetic clustering. Any isolate that did not show a close grouping with a prototype full genome or represented a distinct phylogenetic cluster of HRV-C was selected for amplification of the VP1 region. The procedure for amplification of the VP1 region of HRV-A, -B and -C and initial descriptions of sequences generated have been described previously.

45All available HRV VP4/VP2 and VP1 sequences were downloaded from GenBank in November 2012. The genome regions analysed in this study included the full VP1 region and a fragment of VP4/VP2 (nt 616–1004 numbered according to GenBank accession no. EF582385 for HRV-C4) that is commonly used in studies of HRV epidemiology. Any sequences that were <90 % complete across these regions or wereclassified as non-functional were excluded. All sequences were aligned using the sse v1.0 sequence editor package and pairwise nucleotide p-distances were calculated using the program Sequence Distances within the package. Phylogenetic trees were constructed in the mega v5.0 software package by the neighbour-joining method from 100 bootstrap resampled sequence alignments of maximum composite likelihood distances, with pairwise deletion for missing data. All phylogenetic analysis was undertaken using the Enterovirus D sequence (GenBank accession no. D00820 for EV-D70) as an outgroup.

46Phylogenetic trees were inspected visually for evidence of incongruence between the VP4/VP2 and VP1 regions. Datasets containing these two non-consecutive regions were concatenated to give one continuous sequence, which was then additionally analysed for the occurrence of recombination. rdp v4.0 analysis was used to screen each HRV species dataset for the occurrence of recombination within the capsid-coding region. Each sequence set was examined using a combination of algorithms implemented in the software package, specifically rdp, GeneConv, MaxChi, Chimaera, SiScan and Bootscan. A threshold P value of <0.05 and detection by at least two algorithms was required for a recombination event to be detected. Each HRV species dataset was further analysed using the gard program as part of the HyPhy package. Any predicted recombination breakpoints were analysed by construction of neighbour-joining maximum composite likelihood trees for sequence fragments 5′ and 3′ to the putative breakpoint.