Coronavirus phylogeny based on 2D graphical representation of DNA sequence

Abstract A novel coronavirus has been identified as the cause of the outbreak of severe acute respiratory syndrome (SARS). Previous phylogenetic analyses based on sequence alignments show that SARS‐CoVs form a new group distantly related to the other three groups of previously characterized coronaviruses. In this aritcle, a new approach based on the 2D graphical representation of the whole genome sequence is proposed to analyze the phylogenetic relationships of coronaviruses. The evolutionary distances are obtained through measuring the differences among the two‐dimensional curves. © 2006 Wiley Periodicals, Inc. J Comput Chem 27: 1196–1202, 2006


Introduction
The outbreak of atypical pneumonia, referred to as severe acute respiratory syndrome (SARS) was first identified in Guangdong Province, China, and spread to several countries later. A novel coronavirus was isolated and found to be the cause of SARS. The SARS-coronavirus is a new member of the order Nidovirales, family Coronaviridae, and genus Coronavirus. Some researchers have considered the mutation analysis and phylogenetic analysis. [1][2][3][4][5][6] Phylogenetic analysis using biological sequences can be divided into two groups. The algorithms in the first group calculate a matrix representing the distance between each pair of sequences and then transform this matrix into a tree. In the second type of approaches, instead of building a tree, the tree that can best explain the observed sequences under the evolutionary assumption is found by evaluating the fitness of different topologies. For example, Jukes and Cantor, 7 Kimura, 8 Barry and Hartigan, 9 Kishino and Hasegawa, 10 and Lake 11 proposed various distance measures. Camin and Sokal, 12 Eck and Dayhoff, 13 Cavalli-Sforza and Edwards, 14 and Fitch 15 gave parsimony methods. Felsenstein et al. [16][17][18] proposed maximum likelihood methods.
But, all of these methods require a multiple alignment of the sequences and assume some sort of an evolutionary model. In addition to problems in multiple alignment (computational complexity and inherent ambiguity of the alignment cost criteria), these methods become insufficient for phylogenies using complete genomes. Multiple alignment become misleading due to gene rearrangement, inversion, transposition, and translocation at the substring level, unequal length of sequences, etc, and statistical evolutionary models are yet to be suggested for complete genomes. On the other hand, whole genome-based phylogenic analyses are appearing because single gene sequences generally do not possess enough information to construct an evolutionary history of organisms. Factors such as different rates of evolution and horizontal gene transfer make phylogenetic analysis of species using single gene sequences difficult.
Mathematical analysis of the large volume genomic DNA sequence data is one of the challenges for bioscientists. Graphical representation of DNA sequence provides a simple way of viewing, sorting, and comparing various gene structures. In recent years several authors outlined different graphical representation of DNA sequences based on 2D, 3D, or 4D. [19][20][21][22][23][24][25][26][27][28][29][30][31][32] Graphical techniques have emerged as a very powerful tool for the visualization and analysis of long DNA sequences. These techniques provide useful insights into local and global characteristics and the occurrences, variations, and repetition of the nucleotides along a sequence that are not as easily obtainable by other methods. 29,33 Based on these graphical representation several authors outlined some approaches to make comparison of DNA sequences [34][35][36][37][38] . Recently, we present a new two-dimensional graphical representation of DNA sequences, which has no circuit or degeneracy. 19 Here, a new approach based on the 2D graphical representation of the whole genome sequence is proposed to analyze the phylogenetic relationships of genomes. The evolutionary distances are obtained through measuring the differences among the 2D curves. The examination of the phylogenetic relationships of coronaviruses illustrates the utility of our approach.

2D Graphical Representation of DNA Sequences
As shown in Figure 1, which is similar with Yan's 34 method, we construct a pyrimidine-purine graph on two quadrants of the cartesian coordinate system, with pyrimidines(T and C) in the first quadrant and purines(A and G) in the fourth quadrant. The unit vectors representing four nucleotides A,G,C, and T are as follows: where m is a real number, n is a positive real number but not a perfect square number. Using this representation, we will reduce a DNA sequence into a series of nodes P 0 , P 1 , P 2 , . . . , P N , whose coordinates where N is the length of the DNA sequence being studied) satisfy where a i , c i , g i and t i are the cumulative occurrence numbers of A, C, G, and T, respectively, in the subsequence from the first base to the i-th base in the sequence. We define a 0 = c 0 = g 0 = t 0 = 0. We called the corresponding plot set a characteristic plot set. The curve connecting all plots of the characteristic plot set, in turn, is called the characteristic curve, which is determined by m,n, that satisfy the above mentioned condition. In Figure 2, we show the chimpanzee corresponding curves with different parameters n and  m. Observing Figure 2, we find that chimpanzees have similar curves despite corresponding different parameters of n and m. They have the same tendency despite different lengths. In Figure 3, we present the 2D curves for 24 complete coronavirus genomes (see Table 1) with parameters n = 1/2 and m = 3/4 chosen initially by Yan et al. 34 Observing Figure 3, we find that the curves of BCoV, BCoVL, BCoVM, and BCoVQ have some similar tendencies. The curves of MHV2, MHV, MHVM, and MHVP have some similar tendencies. The curves of BJ01, CUHK-Su10, CUHK-W1, SIN2679, SIN2748, SIN2774, HKU-39849, SIN2500, SIN2677, TW1, Urbani, and TOR2 have some similar tendencies.

Phylogenetic Tree of Coronaviruses
For any sequence, we have a set of points (x i , y i ), i = 1, 2, 3, . . . , N, where N is the length of the sequence. The coordinates of the geometrical center of the points, denoted by x 0 and y 0 , may be calculated as follows 29 The element of the covariance matrix CM of the points are defined: The above four numbers give a quantitative description of a set of point (x i , y i ), i = 1, 2, . . . , N, scattering in a 2D space. Obviously,    the matrix is a real symmetric 2 × 2 one. The eigenvectors and their associated eigenvalues are defined as follows: Corresponding to each eigenvalue λ k , there's an eigenvector EV k . Corresponding to λ 1 < λ 2 , the two eigenvectors are denoted by EV λ1 , EV λ2 , respectively. In Table 2, we list the (x 0 , y 0 ) and eigenvectors belonging to 24 species with parameters m = 1 2 , n = 3 4 . To facilitate the quantitative comparison of different species in terms of their collective parameters, we introduce a distance scale and an angle scale as defined below. Suppose that there are two species i and j, the parameters are is the geometrical center of the curve belonging to species i. λ i 1 , λ i 2 are the two eigenvalues of matrix CM i corresponding to species i. The distance d ij between the two points is. 39 where d ij denotes the distance between the geometric centers of the ith and the jth genomes, and M is the total number of all genomes (M = 24, here). Then we obtain a real M × M symmetric matrix whose elements are d ij .
To reflect the differences between the trends of every two 2D curves, the angles between the corresponding eigenvectors of every two genomes are used. The 2D vectors are denoted as follows: The angle between the two vectors is denoted as follows: The sum of θ k ij over k for given i, j can be used to reflect the trend information of the eigenvectors involved Consequently, two sets of parameters are obtained. The first reflects the difference of center positions represented by the Euclidean distance between the geometric centers. The second indicates the difference of the trends of the 2D curves represented by the related eigenvectors. The overall distance D ij between the species i and j is defined by Accordingly, a real symmetric M × M matrix D ij is obtained and used to reflect the evolutionary distance between the species i and j. The clustering tree is constructed using the UPGMA method in PHYLIP package (http://evolution.genetics.washington.edu/ phylip.html). The final phylogenetic tree is drawn using the DRAWGRAM program in the PHYLIP package. In Figure 4, we present the phylogenetic tree belonging to 24 species.

Conclusion
Most existing approaches for phylogenetic inference use multiple alignment of sequences and assume some sort of an evolutionary model. The multiple alignment strategy does not work for all types of data, for example, whole genome phylogeny, and the evolutionary models may not always be correct. Our representation provides a direct plotting method to denote DNA sequences without degeneracy. From the DNA graph, the A, T, G, and C usage as well as the original DNA sequence can be recaptured mathematically without loss of textual information. The current 2D graphical representation of DNA sequences provides different approaches for constructing the phylogenetic tree. Unlike most existing phylogeny construction methods, the proposed method does not require multiple alignment. Also, both computational scientists and molecular biologists can use it to analysis DNA sequences efficiently with different parameters of n and m.