Bioinformatics analysis of SARS coronavirus genome polymorphism

Background We have compared 38 isolates of the SARS-CoV complete genome. The main goal was twofold: first, to analyze and compare nucleotide sequences and to identify positions of single nucleotide polymorphism (SNP), insertions and deletions, and second, to group them according to sequence similarity, eventually pointing to phylogeny of SARS-CoV isolates. The comparison is based on genome polymorphism such as insertions or deletions and the number and positions of SNPs. Results The nucleotide structure of all 38 isolates is presented. Based on insertions and deletions and dissimilarity due to SNPs, the dataset of all the isolates has been qualitatively classified into three groups each having their own subgroups. These are the A-group with "regular" isolates (no insertions / deletions except for 5' and 3' ends), the B-group of isolates with "long insertions", and the C-group of isolates with "many individual" insertions and deletions. The isolate with the smallest average number of SNPs, compared to other isolates, has been identified (TWH). The density distribution of SNPs, insertions and deletions for each group or subgroup, as well as cumulatively for all the isolates is also presented, along with the gene map for TWH. Since individual SNPs may have occurred at random, positions corresponding to multiple SNPs (occurring in two or more isolates) are identified and presented. This result revises some previous results of a similar type. Amino acid changes caused by multiple SNPs are also identified (for the annotated sequences, as well as presupposed amino acid changes for non-annotated ones). Exact SNP positions for the isolates in each group or subgroup are presented. Finally, a phylogenetic tree for the SARS-CoV isolates has been produced using the CLUSTALW program, showing high compatibility with former qualitative classification. Conclusions The comparative study of SARS-CoV isolates provides essential information for genome polymorphism, indication of strain differences and variants evolution. It may help with the development of effective treatment.


Background
Severe Acute Respiratory Syndrome (SARS) is a new infectious disease reported first in the autumn of 2002 and diagnosed for the first time in March 2003 [1]. It is still a serious threat to human health and SARS coronavirus (CoV) has been associated with the pathogenesis of SARS according to Koch's postulate [2].
Significant research efforts have been made into investigation of the SARS-CoV genome sequence, aimed at establishing its origin and evolution to help eventually in preventing or curing the disease it causes. Although the task is a hard one, it opens up the opportunity, amongst others, for comparative investigation of different SARS-CoV isolates aimed at identification of genome regions properties expressing different levels of sequence polymorphism [3][4][5][6][7][8].
The genome of SARS-CoV consists of a single positive RNA strand approximately 30 Kb in length, consisting of about 10 open reading frames (ORF), and about 10 intergenic regions (IGRs). The first two overlapping ORFs at the 5' end encompass two-thirds of the genome, while the rest of the ORFs at the 3' end account for the remaining third.  [9] in gbk and fasta formats ( Table 1). The main goal was twofold: first, to analyze and compare nucleotide sequences, to identify SNPs positions, insertions and deletions, and second, to group them according to sequence similarity, eventually pointing to phylogeny of SARS-CoV isolates.
According to the length of isolates (insertions and deletions) and the presence of SNPs, we classified them into three main groups with subgroups: "regular" isolates with no insertions or deletions (with different numbers of SNPs), isolates with "long insertions" and isolates with "many individual" insertions and deletions (with different positions of SNPs), which is close to phylogenetic analysis results.

Analysis of genomic polymorphism of the isolates resulted in the following facts
II) Similarity analysis showed that a significant number of isolates have the same length (29727 bases), the same beginning and ending subsequences (that seem to be exact starts and ends of the complete SARS-CoV genome up to the poly-"a" at the 3' end), thus forming a kind of referent group; these are the isolates TWH, TWC3, TWK, TWS, TWY, Urbani, Frankfurt 1 ( Figure 1). The fully sequenced isolate TWH then has been chosen as the referent isolate for sequence comparisons since its average number of SNPs compared to other isolates is the smallest. For example, TWH and Urbani have an average number of SNPs 15.7 and 17.6 respectively for all the isolates, and 5.7 and 10.5 respectively for the referent group.
For SNPs see the tables in the additional files 1 and 2.
IV) There is a group of isolates that have significant length insertions (29 nucleotides) inside the sequence. These are the isolates GD01, SZ3, SZ16. A significant number of individual insertions have been identified in ZJ01 and ZMY 1 isolates ( Figure 1, additional files 3,4,5).
Among the SNP contents of isolates, there is a significant difference in the number of SNPs for different pairs of isolates. For TWH as the referent isolate, this number varies from 1 to 80 SNPs. Isolates may be classified into three groups based on the number of SNPs with TWH ( Figure  2): 1. with less than 15 (TWC3, TWK, TWS, TWY, Urbani,  TWJ, TWC, TW1 Finally, besides the number, there are differences in positions of SNPs (potential mutation sites). In order to avoid nucleotide changes that probably arose during propagation of the virus in cell culture and sequencing, Figure 3 represents positions (on the relative scale of all isolates and on TWH scale) where two or more SNPs occurred, not taking into consideration isolates with long insertions (GD01, SZ3 and SZ16). The positions of multiple SNPs of these three isolates, similar as far as these three are concerned, are highly different from all the others and are represented in Figure 4. These results coincide with those published in Marra et al's paper [4] for Urbani and Tor2 isolates, but differ from those published in Ruan's paper [8] for the 14 isolates therein analyzed (Sin-group, BJgroup, Tor2, Urbani, CUHK-W1, HKU-39849, GD01), which were obviously based on different revisions of the PubMed NCBI Entrez database [9];

Grouping of isolates
The isolates from the dataset considered may be classified according to their sequence polymorphism and SNP contents properties just described. At first, properties (III, IV) may result in three different groups ( Figure 2): A. "regular isolates" whose nucleotide structure is close to the referent group (  Structural tree for SARS-CoV isolates Figure 2 Structural tree for SARS-CoV isolates. The tree is based on qualitative analysis of sequence variation of 36 isolates.   protein  1ab  1ab  1ab  1ab  1ab  1ab  1ab  1ab  1ab  1ab  1ab  1ab  1ab  S  S  S  hyp  hyp  E  M  M  hyp  hyp  hyp  N   Relative  scale  2562  3858  7930  8585  9417  9867  11461  11506  17590  18991  19090  19111  19865  21749  22250  24962  25329  26080  26233  26507  26630  27273  27843  27858  28328   TWH scale  2557  3852  7919  8572  9404  9854  11448  11493  17564  18965  19064  19084  19838  21721  22222  24933  25299  26050  26203  26477  26600  27243  27812  27827 28268 Positions with two or more SNPs in B group with amino acid changes Figure 4 Positions with two or more SNPs in B group with amino acid changes. Only SNPs in B group isolates, regarding TWH, have been counted. The same notation is applied as in Figure 3.   (Figure 7b) Finally, the positions of SNPs will move CUHK-W1 from A1 into A2 group (more than 50% of common SNP positions) while WHU will move from A2 into A1 (less than 30% of common SNP positions), giving the final grouping of isolates presented as a structural tree (Figure 2 Although qualitative in nature, the structural tree turns out to be close to the quantitative grouping which is a basis for (computational) phylogenetic classification. Tables in additional files 1  additional files 2 for isolates of A2 group and TWH, on the TWH scale, additional files 3 for isolates of group B with TWH, and additional files 4,5 for isolates of C1, C2 groups, respectively). Figures 5,6,7 represent density distribution of SNPs, insertions and deletions on the TWH scale, for the same groups of isolates. Figure 8 represents the overall density distribution of SNPs, insertions and deletions for all the 36 isolates, along with the gene map for TWH (which is quite similar to gene maps of other isolates). Density distributions do not show regularities yet (with respect to the number of available sequences) that could provide for precise statistical characterization. Still, they exhibit crowding regions close to the 3' end which is also characterized by the presence of a number of proteins of unknown function.

SARS-CoV isolates
It can also be noted that the proposed grouping of 36 isolates, based on different criteria, still conserves the previous classification T-T-T-T / C-G-C-C [8]. All the isolates from groups A1 and C have T-T-T-T configuration, while all the isolates from groups A2 and B have C-G-C-C configuration, except for GZ50, BJ04 being TGCC (Figure 2, The fact that the enlarged number of isolates exhibits the same properties relating to the four loci supports the assumption that the mutations could not have arisen by chance base substitution during propagation in cell culture and the sequencing procedure [8].

Changes in amino acids
We analyzed amino acid changes in proteins for the annotated isolates (19 out of 36), and presumed proteins in non-annotated ones for multiple SNPs in all the isolates.
resulted by SNPs. The only significant change in amino acid properties is in S protein Gly→Asp (A2, B groups, i.e., in CUHK-W1, GZ50, BJ01-BJ03, GD01, SZ3 and SZ16 isolates) and hypothetical protein Cys→Arg (the same isolates, BJ04 in addition). The only addition in nonannotated sequences is in hypothetical protein following S protein in TWH, exhibiting silent change, and in nonannotated BJ02 and BJ03, corresponding to the hypothetical protein, Gly→Glu. Similar analysis can be done for amino acid changes corresponding to SNPs at positions specific for B group isolates (Figure 4). Taking into account the only annotated isolate GD01, there are five amino acid changes in polyprotein 1ab, two amino acid properties changes in S protein (Ser→Leu and Tyr→Asp, the second being within the epitope region), one amino acid change in M protein and one amino acid property change (Cis→Arg) in BGI-PUP.

Phylogenetic analysis
The SARS-CoV isolates have been multialigned using the CLUSTALW program [10] as the very first step in obtaining a phylogenetic tree. The aligned sequences have been submitted then to CLUSTALW for bootstrapping and phylogenetic tree production. Enlargement of the sequence set resulted in the refinement of the phylogenetic tree produced, as compared to previous results such as Ruan [8] and Zhang&Zheng [14], obtained for 14 and 16 isolates, respectively. The phylogenetic tree obtained, drawn using the PhyloDraw program [15], is represented in Figure 9. It is similar to our structural tree based on qualitative analysis of the isolates (Figure 2).
The results of the analysis of dissimilarities, described in previous paragraphs, are in accordance with the alignment obtained by CLUSTALW, but regrouped and formatted in a way that facilitates further interpretation and application.

Conclusion
Comparative analysis of genome sequence variations of 38 SARS-CoV isolates resulted in some conclusions that might be of interest in further investigation of the SARS-CoV genome: 1. All of the SARS-CoV isolates are highly homologous (more than 99% pairwise). Most of them have similar nucleotide structure, with the same 5' and 3' ends and poly-"a" at the 3' end of different length (0-24), some of them with a single short deletion close to the 3' end of the sequence; out of 312 SNPs in total, only two are in IGRs.
2. Three of the 38 isolates have long insertions within the sequence; 3. Two of the isolates have a large number of individual insertions / deletions, exhibiting different SNP positions; 4. All the isolates may be grouped according to sequence polymorphism into three groups (with up to two subgroups), reflecting their similarities / dissimilarities. Since the isolate sequences have a high degree of homology, different properties of groups are represented in a more transparent way in the classification tree obtained by such a qualitative analysis, than in a bootstrapped phylogenetic tree obtained from multialigned sequences using the CLUSTALW program [10].
5. The total number of amino acid changes caused by multiple SNPs is 15 (in isolates of A, C groups) and 34 in isolates of B group. The total number of silent mutations is 10 (for A, C groups) and 7 (for B group).
6. Since S protein is of special interest regarding its receptor affinity and antigenecity, it is interesting to notice that all amino acid properties' changes are located in its outer membrane region, one for A, C groups and two for B group.
7. The results obtained may be useful in further investigation aiming at identification of SARS-CoV genome regions responsible for its infectious nature.

Dataset
We investigated the complete genomes of 38 SARS-CoV isolates. Nucleotide sequences are taken from the PubMed NCBI Entrez database [9] in gbk and fasta formats ( Table  1).
The coverage included all the isolates published by October 31 st 2003 (with updated revisions). The identifiers, accession numbers, genomic size (in nucleotides), revision dates and country or source of the isolates considered are included in the table, together with labels as referred in this paper. The fully sequenced isolate TWH has been chosen as the referent isolate, since its average number of SNPs was the lowest as compared to all other isolates.

Methods for similarity analysis
For similarity analysis of isolates, the following procedure has been applied consisting of two steps: 1. identification of structurally identical parts of isolates, i.e., insertion and deletion sites 2. identification of SNPs in structurally identical parts.
Step 1 has been carried out by a function performing similarity analysis of subsequences of a given length (e.g., 100 bps), and identifying significantly non-matching strings as being inserted in the corresponding sequence (i.e. deleted from the other). Since significant number of isolates have the same length (29727 bases) and starting and ending subsequences (that seem to be the exact starts and ends of the complete SARS-CoV genome up to the poly-"a" at the 3' end), they may be considered as forming a representative group. The nucleotide structure of all other isolates was analyzed with respect to this representative group. For each pair of isolates (x,y) (x from the representative group), a file InsDelx-y has been produced containing positions and lengths of each of the insertions or deletions in the isolate y.
Step 2 has been carried out by comparing structurally identical parts (of the same length) of pairs of isolates. The starting and ending positions of those parts have been taken from the file InsDelx-y (for comparison of x and y), produced in step 1. The procedure returns results in a file with SNPs in the two sequences (files Mismx-y).
We also used the CLUSTALW program [10] for multialignment as a control process, as well as for phylogenetic investigations.

Methods for phylogenetic investigation
In order to use similarity analysis results for drawing any phylogenetic conclusions about the SARS-CoV genome dataset, a CLUSTALW [10] multialigned output has been generated and a bootstrapped phylogenetic tree has been produced and drawn using the PhyloDraw program [15].