Genomewide Function Conservation and Phylogeny in the Herpesviridae

Identification of Homologous Protein Families and Function Assignment

Sequence homology among all proteins derived from complete herpesvirus genomes (Table 1) was determined and used as a basis to construct HPFs (Fig. 1). We identified 243 homologous families that contained two or more proteins, comprising 1496 proteins out of a total of 1692 predicted ORFs in the 19 genomes studied. We observed that ∼80% of the total herpesvirus proteins had homologs in a different herpesvirus, whereas 20% appeared to be unique to particular genomes, sometimes existing as multiple copies (paralogs). Three-dimensional structural information for a subset of herpesvirus proteins validated the homologous family groups. It was not possible to collapse the homologous families into smaller groups based on such structural information. We used GenBank header files to manually assign functions to the different HPFs, including those with only one protein member. Functions consisted of both a short definition, such as DNA polymerase, and a broad functional class, for example, replication. Uncharacterized proteins were assigned to the unknown class.

All HPFs that belong to different functional classes can be retrieved from http://www.biochem.ucl.ac.uk/bsm/virus_database. In addition, the HPFs can be searched using virus name, functional annotation, keywords, or GenBank protein entry numbers. Each HPF has been assigned a distinct family number (HPF 1, HPF 2, etc.).

Phylogenetic Distribution of Protein Homologous Families

We used the homologous families to build protein phylogenetic profiles (Pellegrini et al. 1999), in which for each homologous family the presence or absence in every genome was recorded in the form of a binary matrix, where 1 means presence of at least one protein from the genome and 0 means no protein. In this type of analysis, paralogous proteins, resulting from multiple copies of a gene in the same genome, will only be counted once. The profiles were used to determine the number of gene functions conserved in pairs of genomes and to construct phylogenetic trees. In addition, the profiles were used to determine the minimum number of functions conserved at the subfamily/lineage level and to study the degree of conservation with respect to the functional class of the gene.

The distribution of the number of shared functions, based on sequence homology, between any two genomes across the different Herpesviridae was in accordance with the main evolutionary herpesvirus lineages (subfamilies) and sublineages (individual viruses). A relatively high number of homologous families were conserved within subfamilies and a much lower number conserved between members of different subfamilies (Table 2, lower triangle). Using our sequence-comparison algorithm, the minimum number of shared homologous families was 26 conserved between the Alpha- and Betaherpesvirinae, and the maximum, 96, was found between the closely related HHV6-A and HHV6-B. The number of shared homologous families was more variable within subfamilies than between members of different subfamilies. For example, the number of shared homologous families within any two Alphaherpesvirinae viruses ranged from 52 to 77, but between members of this subfamily and the other two subfamilies, the range of conserved homologous families was much narrower, between 26 and 30. We also calculated the percentage of homologous families conserved between any two genomes, taken relative to the genome with a smaller number of different families (Table 2, upper triangle). At least one-third of the homologous families were conserved between any two genomes with respect to the smallest genome of the pair. Within subfamilies the percentage varied between 54% (HHV-5 vs. HHV-6) and 100% (HSV-1 vs. HSV-2).

Conservation of Function Within and Across Subfamilies

We next focused our attention at the number of homologous families, which formed the core set of proteins in the different herpesvirus subfamilies. We detected 26 different ORFs that were conserved across the Herpesviridae (Table 3), which is close to previous estimations of the minimal herpesvirus genome on the basis of clear sequence homology (Hannenhalli et al. 1995; McGeoch and Davison 1999a). Each of these ORFs formed a separate homologous family, except for the major and the minor capsid proteins that share a region of ∼66 amino acids and, therefore, are part of the same homologous family. Apart from this common set of genes, other ORFs were conserved in two subfamilies but were absent in the third. In particular, we found three homologous families that were specific for Alpha- and Gammaherpesviruses and 10 specific for Beta- and Gammaherpesviruses. We did not identify any homologous families present in all members of the Alpha- and Betaherpesviruses but not present in the Gammaherpesviruses. According to this, the Gamma and Beta lineages clearly share more genes with detectable sequence homology than either of the two with the Alphaherpesviruses. By computing the ORFs that were only conserved in all members of one subfamily but in no other herpesvirus, we determined the subfamily-specific homologous families. There were 22 such homologs for the Alphaherpesviruses, 23 for the Betaherpesviruses, and only 8 for the Gammaherpesviruses. By adding the homologous families conserved at the level of two or three subfamilies, we obtained 51 families totally conserved for Alphaherpesviruses, 59 for Betaherpesviruses, and 46 for Gammaherpesviruses.

Analysis of Different Functional Classes

The number of homologous families with known function identified in viruses from different lineages was variable. The Alphaherpesviruses were the best characterized, with between 60% and 80% of the proteins of any virus having an assigned function. This percentage was between 55% and 70% for the Gammaherpesviruses. The Betaherpesviruses contained the largest number of uncharacterized proteins among the different lineages. Only about half of the predicted HHV-6 and HHV-7 ORFs and about one-fourth of the HHV-5 (human cytomegalovirus, HCMV) ORFs have a documented function.

Next we compared the degree of conservation of the different homologous families across the whole Herpesviridae. To do this, we analyzed separately the phylogenetic profiles of homologous families that belonged to different functional classes. The analysis is shown for the structural class (Fig. 2). The size distribution of the homologous family, taken as the number of different viruses represented, was markedly different for the seven functional classes (Fig. 3). Genes involved in nucleotide metabolism and DNA repair were the most conserved, with most of them being in large groups containing viruses from two or three subfamilies. Structural proteins, including capsid and tegument proteins, were also well conserved, as were proteins from the replication functional class. However, glycoproteins showed a much lower conservation and most of them belonged to families with a size of 1–3 viruses, clearly below the size of a herpesvirus subfamily. Proteins identified as being involved in transcription, as well as proteins in the others group, which included genes involved in virus-host interactions, were also poorly conserved. Finally, the majority of homologous families with an as-yet-unknown function (unknown class) fell into the 1–3 viruses size range.

Interestingly, the three proteins that have been conserved in all Alpha- and Gammaherpesviruses but not in Betaherpesviruses belonged to the same functional group, nucleotide metabolism/DNA repair, namely ribonucleotide reductase small subunit, dUTPase, and thymidine kinase (HPF 28, 29, and 31, respectively). In contrast, homologous families that are exclusively conserved between the Beta- and Gammaherpesviruses were structural or of unknown function. One HPF, 9, contained the DNA origin-binding protein from the Alphaherpesviruses and the Betaherpesviruses HHV-6/HHV-7. However, this protein showed no homology to any Gammaherpesvirus protein or to proteins from the Betaherpesvirus HHV-5 (human cytomegalovirus). Homologous families that appeared to be exclusive to particular herpesvirus subfamilies occurred across different functional classes, although 12 homologous families corresponded to structural proteins in Alphaherpesviruses and 12 to genes of unknown function in Betaherpesviruses.

Phylogenetic Reconstruction Based on Function Conservation

Phylogenetic profiles were used to construct phylogenetic trees based on whole-genome homologous family content. In this type of tree, the distances between the different viruses are based on the degree of conservation of gene functions. Therefore, the topology of the tree will be affected by gene loss, gene capture (typically from the host genome in herpesviruses), and extensive sequence divergence beyond the recognition by the sequence comparison methods used here. The phylogenetic profiles were bootstrapped 100 times before constructing the trees. To build neighbor-joining trees, we explored the use of two types of intergenomic distance, the fraction of nonshared functions, and the fraction of dissimilar functions (Fig. 4A,B, respectively). The branching order of the two trees was the same for the two approaches and the main differences were in the branch lengths. As expected, the distance method that used the total of dissimilar functions, which was not standardized to the size of the smaller genome, reflected the difference in the number of genes per genome much better (Fig. 4B). For example, the branch length for HHV-5, which has approximately twice as many genes than any other herpesvirus genome, was longer than branches in other parts of the tree. Bootstrap supports, in general, were very high, with the exception of the split of the two ungulate herpesviruses, alcelaphine herpesvirus 1 (AHV-1) and equine herpesvirus 2 (EHV-2), with bootstrap values of 44% and 37%, respectively. In addition to neighbor-joining trees, we built up a maximum parsimony tree from the same set of data (Fig. 4C). Again the branching pattern was the same, except for the independent split of AHV-1 and EHV-2 as sisters, although, again, the bootstrap value was relatively low (64%).

The trees based on the phylogenetic profile clearly resolved the splits between herpesvirus subfamilies and sublineages (Table 1). In addition, our data regarding the number of shared functions between different subfamilies supported previous observations of an early split of the Beta- and Gammaherpesviruses from the Alphaherpesviruses. This branching pattern was observed when we simulated a root by using an artificial outgroup genome that had none of the homologous proteins, that is, a row of 0s in the phylogenetic profile. To compare these trees with a sequence-comparison-based tree, we constructed an alignment of all conserved domains in the 26 ORFs identified as clear homologs in all herpesviruses. These genes have been preserved throughout herpesvirus evolution and are present in one copy per genome. The alignment contained 8900 positions and, using 100 bootstrapped data sets, neighbor-joining, UPGMA, and maximum parsimony trees were constructed. All trees showed the same topology and the neighbor-joining tree is shown in Figure 4D. The trees were representative and agreed with previous phylogenetic trees produced using a smaller set of highly conserved herpesvirus proteins (McGeoch and Davison 1999a).

There was complete consistency between the trees based on either function conservation or on sequence alignment, for the Alpha- and Betaherpesviruses, with all trees producing the same branching pattern. Among the Gammaherpesviruses, branch differences occurred in the positions of the ungulate herpesviruses (AHV-1 and EHV-2) and the murine herpesvirus MHV-68 when comparing the various trees. The position of the MHV-8, previously unresolved (McGeoch and Davison 1999b), appeared basal to the rhadinoviruses (all Gammaherpesviruses except HHV-4 in this study) in the alignment-based tree with a bootstrap value of 99%. Instead MHV-68 clustered together with the human and primate viruses in the other trees (bootstrap values of 87% and 69%). AHV-1 and EHV-2 formed a cluster in the neighbor-joining trees based on homologous family conservation. This association is in accordance with the hypothesis that herpesviruses have coevolved with their hosts (McGeoch and Cook 1994; McGeoch and Davison 1999b). However, the bootstrap values were low and the cluster was not observed in the other two trees. Therefore, the result is suggestive but requires further investigation.

Identification of Homologous Families

A total of 19 complete genomes representative of viruses in the Herpesviridae were retrieved from GenBank (see Table 1). Protein sequences from all identified ORFs were extracted and used to build up a protein-sequence dataset containing a total of 1692 proteins. XDOM (Gouzy et al. 1997) was used to identify homology between the proteins and to identify regions of sequence similarity that were common to related proteins. XDOM is based on BLASTP (Altschul et al. 1990) and had previously been used to identify regions of protein-sequence similarity in different complete genomes from bacteria, archaea, and eukarya (Gouzy et al. 1999). Initially, we empirically tested several parameters of the program so as to maximize sensitivity without compromising accuracy. After the initial observations, XDOM was used with the parameters SCORE = 75 and SCORE2 = 40 instead of the default values (90 and 50, respectively). We found that these parameters increased sensitivity although they still prevented the appearance of spurious matches between functionally unrelated proteins. A C++ program, PSC BUILDER, was written to cluster protein sequence domains together into HPFs. We clustered all proteins that shared at least one sequence domain, so that in each HPF there is at least one conserved region that is present in all proteins (Fig. 1). The method used identifies all proteins that share sequence similarity. Therefore, orthologous and paralogous sequences, derived from recent gene duplications, may be found in the same HPF. Proteins that did not share sequence homology to any other protein were treated as single-protein families. In these cases, the equivalent of the HPF-conserved sequence region will be the complete protein sequence.

Function Identification

Protein function, if known, was extracted for each herpesvirus protein from the original sequence-entry annotations. As no major disagreements were found in the annotated function of different proteins in the same homologous family, we considered that a function could be used to define most herpesvirus HPFs. Functions were simple definitions such as DNA polymerase or capsid protein. All protein functions were classified into seven major pathways or functional classes: replication, nucleotide metabolism and DNA repair, transcription, structural (including capsid, tegument, and virus assembly proteins), glycoproteins, others (including proteins involved in host-virus interactions such as immune modulation proteins), and unknowns.

Phylogenetic Profiles of the Homologous Families

Phylogenetic profiles can be defined from the presence or absence of a HPF in each virus genome (Pellegrini et al. 1999). A matrix was constructed, which for each homologous family, the presence of proteins from each given genome was expressed as 1 (presence) or 0 (absence). The matrix consisted of 439 columns for the total of homologous families, including those with only one protein, and 19 rows for the number of herpesvirus genomes. The presence of more than one protein from the same genome in the same homologous family (presumably due to paralogous genes) was not taken into account for the purpose of matrix construction. For the separate analysis of functional class conservation, the complete matrix was split into class submatrices. The number of shared gene functions across all genomes was determined as a whole number, representing all homologous families in which both genomes were present and also as a percentage of the number of shared functions.

Phylogenetic Analysis of Herpesvirus Genomes on a Functional Basis

The phylogenetic profiles were used to conduct phylogenetic analysis of the different viruses. The different protein families can be considered as molecular function characters for which the different viruses are positive (1) or negative (0). The data was bootstrapped 100 times using our own scripts and maximum parsimony, and distance methods (neighbor-joining) were applied.

For the distance methods, two distance measures were used: (1) Fraction of nonshared functions dx,y = 1−[(positive in X and in Y)/(minimum between total positives in X and total positives in Y)] and (2) fraction of dissimilar functions dx,y = [(positive in X but not in Y) + (positive in Y but not in X)]/total of homologous families.

In both cases, a positive refers to a 1 in the matrix (presence of a gene from the homologous family in that genome).

The first measure was previously used to build trees from gene content in unicellular organisms (Snel et al. 1999); the second was chosen because it may better satisfy the property of additivity of distance (Rzhetsky and Nei 1993). We used the programs NEIGHBOR and DNAPARS from the PHYLI8P package (Felsenstein 1993) for neighbor-joining and maximum parsimony methods, respectively. Consensus trees were derived using CONSENSE from the same package. The final trees were drawn with TREEVIEW (Page 1996).

Phylogenetic Analysis Based on Protein Sequence Alignments

We used the 26 ORFs identified as homologous in all Herpesviridae to construct a phylogeny based on sequence similarity. Alignments from a total of 28 conserved domains from the 26 ORFs and derived with MKDOM (Gouzy et al. 1997) were concatenated to form a single alignment of 8900 amino acids, including gaps. The alignment was bootstrapped 100 times and distances were computed with CLUSTALX default metric based on the Gonnet matrices (Benner et al. 1994) and corrected for multiple substitutions. Neighbor-joining trees were constructed using CLUSTALX (Thompson et al. 1997); UPGMA and maximum parsimony trees were constructed using NEIGHBOR and PROTPARS, respectively, from the PHYLIP package (Felsenstein 1993). Consensus trees were obtained with CONSENSE from PHYLIP and trees visualized with TREEVIEW (Page 1996).