A fast and accurate algorithm for single individual haplotyping
- Minzhu Xie^{1}Email author,
- Jianxin Wang^{2} and
- Tao Jiang^{3}
https://doi.org/10.1186/1752-0509-6-S2-S8
© Xie et al.; licensee BioMed Central Ltd. 2012
Published: 12 December 2012
Abstract
Background
Due to the difficulty in separating two (paternal and maternal) copies of a chromosome, most published human genome sequences only provide genotype information, i.e., the mixed information of the underlying two haplotypes. However, phased haplotype information is needed to completely understand complex genetic polymorphisms and to increase the power of genome-wide association studies for complex diseases. With the rapid development of DNA sequencing technologies, reconstructing a pair of haplotypes from an individual's aligned DNA fragments by computer algorithms (i.e., Single Individual Haplotyping) has become a practical haplotyping approach.
Results
In the paper, we combine two measures "errors corrected" and "fragments cut" and propose a new optimization model, called Balanced Optimal Partition (BOP), for single individual haplotyping. The model generalizes two existing models, Minimum Error Correction (MEC) and Maximum Fragments Cut (MFC), and could be made either model by using some extreme parameter values. To solve the model, we design a heuristic dynamic programming algorithm H-BOP. By limiting the number of intermediate solutions at each iteration to an appropriately chosen small integer k, H-BOP is able to solve the model efficiently.
Conclusions
Extensive experimental results on simulated and real data show that when k = 8, H-BOP is generally faster and more accurate than a recent state-of-art algorithm ReFHap in haplotype reconstruction. The running time of H-BOP is linearly dependent on some of the key parameters controlling the input size and H-BOP scales well to large input data. The code of H-BOP is available to the public for free upon request to the corresponding author.
Keywords
Background
Each human somatic cell contains 23 pairs of chromosomes, and there are about 0.5% differences between the DNA sequences of two copies of each chromosome [1]. The dominant DNA differences are single nucleotide polymorphisms (SNPs). Identification of the combination of alleles at the SNP loci on the same chromosome copy, i.e., haplotyping, is needed to fully understand the human genetic variation patterns and enhance the power of genome-wide association studies for complex diseases [2, 3]. Currently, it is expensive and labor-intensive to separate two copies of chromosomes by biological techniques [4], and most published human individuals' genomes contain only the mixed information, i.e., genotype information, of the underlying two copies of chromosomes [5]. Therefore, to reduce the cost, accurate and fast computational haplotyping methods are of urgent importance.
There have been many computational haplotyping models [6–8] and they can be grouped into two main classes: haplotype inference and haplotype assembly [6]. Haplotype inference is to phase the haplotypes of individuals in a pedigree or a population from their genotypes. Computer algorithms of haplotype inference have been used in the International HapMap Project [9] and the 1000 Genomes Project [5] to identify haplotypes. Haplotype assembly is also called Single Individual Haplotyping (SIH). SIH assembles a pair of haplotypes from an individual's aligned DNA sequence fragments. With the dramatically dropped cost of human whole genome sequencing, more and more human individual's DNAs have been sequenced. With mate-pairs sequencing and read length improvements of the next-generation sequencing technologies, and with the development of new sequencing technologies, SIH methods have been used to build haplotype-resolved genome of human beings [10, 11]. When there are enough DNA sequence fragments that cover two or more consecutive variant loci, SIH builds longer and more accurate haplotype blocks than haplotype inference does [12].
Since Lancia et al., first formalized the SIH problem [13], many optimization models and algorithms have been introduced to solve the problem [7, 14–25]. The main models are MEC (minimum error correction) [24], MFR (minimum fragment removal), and MSR (minimum SNP removal) [13]. Recently, Duitama et al., proposed a new model MFC (maximum fragments cut) [16]. Most of the models are NP-hard and APX-hard [16, 21, 26], and their exact algorithms run in time exponentially dependend on at least one input parameter [15, 17, 19–21]. Therefore, a large number of heuristic algorithms have been designed to deal with the problem [16, 22, 23, 25]. According to [16], one of the most accurate heuristic algorithms is HapCUT [25], while ReFHap [16] runs much faster than HapCUT without loss of accuracy. In this paper, we consider both quality measures "errors corrected" and "fragments cut", and propose a new optimization model, called Balanced Optimal Partition (BOP), for the SIH problem. The model generalizes the most popular model MEC and the recent model MFC. In fact, it could be made either model by setting some parameters to extreme values. To solve the model, we propose a dynamic programming algorithm H-BOP. By limiting the number of intermediate solutions at each iteration to an appropriately chosen small integer k, H-BOP is able to solve the model efficiently. The time complexity of H-BOP linearly depends on some of the key parameters controlling the input size and the algorithm scales well to large input data.
Results and discussion
We use a public available Java package SIH [27] to test the performance of H-BOP. The package contains a simulated data generator and implements algorithm ReFHap [12, 16]. The simulated data are generated according to five parameters: number of SNPs (haplotype length) n, number of fragments m, average fragment length l, sequencing error rate e, and gap rate g. In our experiments, since we only consider heterozygous SNPs, for each data set, a haplotype h_{1} containing n SNPs is generated randomly first and then the other haplotype h_{2} is obtained by flipping each allele of h_{1}. For each haplotype, m/ 2 fragments are randomly sampled from the haplotype and their lengths follow a normal distribution with mean l and variance 1. Finally for each fragment, every allele is flipped with probability e to introduce sequencing errors and, except at the first and last positions, every allele is deleted with probability g to introduce gaps. Given fragments generated as above, the average call coverage c is calculated by dividing the total number of alleles of the fragments by the haplotype length n. Please see [16] for more details.
Among many algorithms for the SIH problem, HapCUT and ReFHap are two of the most accurate heuristic algorithms [12, 16]. Since ReFHap is much faster than HapCUT, we only compare our algorithm with ReFHap. We implemented our algorithm H-BOP in Java and embedded it in the java package SIH and tested the accuracies, phased haplotype lengths and running time of H-BOP and ReFHap on simulated data and a real data set provided by [12]. All tests are carried out on a Windows 7 64 bit PC (3GHz CPU, 4GB RAM). To measure the haplotype reconstruction accuracy of an SIH algorithm, the hamming distance between the reconstructed haplotype pair and the real haplotype pair was previously used widely in the literature [14, 23]. However, a recent study [28] showed it over-penalizes simple switch errors. Therefore, as in [12, 16], we use switch errors to measure the accuracy of an algorithm. A switch error is an inconsistency between the reconstructed haplotype pair and the real haplotype pair over two contiguous SNPs. There may be some SNP sites where an algorithm is unable to determine the alleles of a haplotype. The phased haplotype length is defined as the number of SNP sites where the alleles of the reconstructed haplotype pair are determined. And the number of switch errors divided by the phased haplotype length is called switch error rate.
If the allele of a fragment f at a SNP site s is known, we say f covers s. When there are no fragments covering two consecutive loci, it is not possible to determine the haplotype containing these consecutive loci for all SIH models. Therefore, for each test we divide a haplotype into blocks according to the input fragments as in [16]. A block corresponds to a connected component of a graph G = (V, E) where V is the set of the SNP sites and there is a edge between two SNP sites s_{1} and s_{2} if and only if there is a fragment covers both s_{1} and s_{2}. The switch errors of an algorithm are the sum of the switch errors in all blocks. In the following simulation tests, the haplotype length n = 100, the gap rate g = 0.1 and each result is the average over 100 repeated experiments if there is no explicit specification.
Parameters of the algorithm H-BOP
There are two parameters w and k in H-BOP. The parameter w is a weighting factor. H-BOP tries to seek a solution with the minimum number of errors corrected when w = 0, and a solution with a maximum cut of the weighted conflict graph corresponding to the input fragments [16] when w is set sufficiently large. k is the maximum number of intermediate solutions that we will keep at each iteration of H-BOP. When k is large enough, H-BOP in fact becomes an exact algorithm. To choose appropriate values for w and k, we test H-BOP on different combinations of w and k.
Simulation results
Average phased haplotype lengths of ReFHap and H-BOP
l = 3 | l = 6 | l = 10 | |||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|
m = | 140 | 211 | 281 | 351 | 74 | 111 | 148 | 185 | 44 | 67 | 89 | 111 | |
e = .005 | ReFHap | 98.16 | 99.59 | 99.79 | 99.97 | 97.2 | 99.04 | 99.39 | 99.66 | 96.24 | 97.7 | 98.51 | 99.06 |
H-BOP | 98.17 | 99.61 | 99.8 | 99.97 | 97.2 | 99.04 | 99.39 | 99.66 | 96.24 | 97.7 | 98.51 | 99.06 | |
e = .01 | ReFHap | 97.99 | 99.6 | 99.68 | 99.9 | 97 | 98.87 | 99.32 | 99.59 | 95.68 | 97.81 | 98.73 | 99.13 |
H-BOP | 98.07 | 99.63 | 99.7 | 99.9 | 97 | 98.87 | 99.32 | 99.59 | 95.69 | 97.81 | 98.73 | 99.13 | |
e = .02 | ReFHap | 97.42 | 99.39 | 99.71 | 99.86 | 96.94 | 98.59 | 99.48 | 99.63 | 95.31 | 98.06 | 98.24 | 99.22 |
H-BOP | 97.57 | 99.48 | 99.74 | 99.86 | 96.95 | 98.59 | 99.48 | 99.63 | 95.31 | 98.06 | 98.24 | 99.22 | |
e = .05 | ReFHap | 96.6 | 98.88 | 99.49 | 99.77 | 95.61 | 98.31 | 99.14 | 99.5 | 94.25 | 97.56 | 98.62 | 98.87 |
H-BOP | 96.94 | 99.09 | 99.63 | 99.84 | 95.76 | 98.42 | 99.15 | 99.51 | 94.33 | 97.56 | 98.62 | 98.87 |
Results on real data
We downloaded a real data set from the SIH website [27], which contains the aligned sorted fosmid-based NGS DNA sequence fragments and gold-standard haplotypes of a HapMap trio child, NA12878 [12]. The total heterozygous SNP sites of the data are 1,704,166, the total fragments are 285,341, the average fragment length is 18.03, and the average call coverage is 3.02. Since the coverage is low, there are many consecutive heterozygous SNP site pairs not covered by any fragments, and hence the 23 pairs of chromosomes are divided into 17,839 haplotype blocks. Due to the low coverage, both H-BOP and ReFHap ran very fast and completed the reconstruction of haplotypes for all 23 pairs of chromosomes in about half a minute. The total phased haplotype lengths of H-BOP and ReFHap are 1,563,741 and 1,562,402 respectively. Compared with the gold-standard haplotypes, the total switch errors of the haplotypes built by H-BOP and ReFHap are 21,859 and 21,835 respectively. Though the switch errors of H-BOP is larger than that of ReFHap, the switch error rates of H-BOP and ReFHap are both 0.014.
To account for both completeness and quality, Duitama et al.[12] proposed an alternative measure QAN50 (quality adjusted N50). QAN50 is calculated as follows [12]:
(1) Break every haplotype block into the longest possible segments containing no switch errors.
(2) Calculate span (in reference base pairs) from the first phased SNP to the last phased SNP for every segment.
(3) Adjust each span by multiplying the span with phased SNPs ratio (the number of phased SNPs divided by the number of total SNPs) inside the segment (to correct for un-phased SNPs).
(4) Sort segments from the largest to the smallest adjusted span.
(5) Traverse the list and count the number of phased SNPs. When the count is more than a half of the total number of SNPs, the adjust span of the current segment is QAN50.
Clearly, algorithms with larger QAN50 values are more desirable. The QAN50 values of H-BOP and ReFHap on the above real data set are 114,261.09 and 113,831.67, respectively.
Conclusions
Haplotyping is regarded as one of the hardest challenges in personal genome sequencing. Though some single molecule sequencing technologies have been developed, they are still too expensive and labor-consuming. Computer algorithms are widely used to reconstruct haplotypes in personal genome sequencing. SIH uses computer algorithms to build a pair of haplotypes from an individual's aligned DNA sequence fragments. There are many different combinatorial optimization models for SIH, among which MEC is the most popular and MFC is the most recent introduction [16]. In this paper, we combine the two quality measures "errors corrected" and "fragments cut" used in MEC and MFC and introduce a new model BOP. We design a heuristic dynamic programming algorithm H-BOP to solve the model. By setting appropriate parameters of the algorithm, H-BOP is accurate and fast. Extensive simulation experiments show that H-BOP is generally faster and more accurate in assembling haplotypes than a recent state-of-art algorithm ReFHap. When the average fragment length is small and sequencing error rate is relatively high, H-BOP is significantly more accurate than ReFHap. The running time of H-BOP increases linearly as the average call coverage c increases, and H-BOP runs much faster than ReFHap when c is large. The test on a real data set also shows that H-BOP is superior to ReFHap considering both completeness and quality of the reconstructed haplotypes.
Methods
Formulation and problem
With the input of aligned DNA sequence fragments derived from a pair of chromosomes, SIH tries to reconstruct a pair of haplotypes of their underlying chromosomes. If there are no sequencing errors and for any two consecutive (but not necessarily adjacent) heterozygous SNP loci, there is at least one fragment covering both, we can easily determine the linkage relationship between two consecutive heterozygous SNPs and thus SIH is easy. However, sequencing errors are unavoidable which makes the problem complicated. In the following, we first introduce some notations and definitions similar to those in [15, 16, 20, 23], and then propose a new optimization model.
Since we only consider the alleles at SNP loci in the SIH problem, the input aligned fragments are encoded as an m × n SNP matrix M[15, 16, 25], where m is the number of fragments and n the number of SNP loci. An entry at the i th row and the j th column of M is denoted as M[i, j]. M[i, j] takes a value from {0, 1, −}, where '0' (or '1') encodes the major allele (or the minor allele, respectively) in the population and '−' represents an unknown allele. As in previous work [16, 25], we assume that all SNP loci are heterozygous and every fragment covers at least two heterozygous SNP loci. If this is not the case, a simple preprocessing as in [22] and [29] can be used to remove homozygous loci and fragments covering only one heterozygous SNP locus. In the remainder of the section, the i th row of M is equivalent to the i th fragment and the j column of M is equivalent to the j th SNP locus without explicit specification for briefness.
Given an m × n SNP matrix M, let the underlying haplotype pair be $\mathcal{H}=\left({H}_{1},{H}_{2}\right)$. The allele at the j th SNP locus of H_{1} (or H_{2}) is denoted as H_{1}[j] (or H_{2}[j], respectively). Notice that since all SNP loci of interest are heterozygous, for any j ∈ {1, ..., n}, H_{1}[j] ≠ H_{2}[j], i.e. when H_{1}[j] = 0, H_{2}[j] = 1 and when H_{1}[j] = 1, H_{2}[j] = 0. Let f_{ i } denote the i th fragment (i.e. the i th row of M). c(M[i_{1}, j], M[i_{2}, j]) = 1 means that ${f}_{{i}_{1}}$ and ${f}_{{i}_{2}}$conflict at SNP locus j (i.e. column j), and that if both fragments come from the same chromosome, either M[i_{1}, j] or M[i_{2}, j] is a sequencing error. Similarly, c(M[i, j], H_{ p } [j]) = 1 (p = 1 or 2) means that f_{ i } and H_{ p } conflict at SNP locus j, and that if f_{ i } comes from the chromosome with haplotype H_{ p } , M[i, j] must be a sequencing error.
If the real haplotype pair is $\mathcal{H}$, it is easy to verify that the number of sequencing errors in the input fragments is at least ${s}_{c}\left(\mathcal{H},\phantom{\rule{2.77695pt}{0ex}}M\right)$, which we call the errors corrected measure. Based on the principle of parsimony, it is a natural optimization objective that to minimize the number of errors and hence Minimum Error Correction [20, 23, 24] is the most popular model for SIH.
Minimum Error Correction (MEC): Given an m × n SNP matrix M, find a haplotype pair $\mathcal{H}$ such that ${s}_{c}\left(\mathcal{H},\phantom{\rule{2.77695pt}{0ex}}M\right)$ is minimized.
Based on the underlying haplotype pair $\mathcal{H}=\left({H}_{1},{H}_{2}\right)$, it is easy to partition the fragments into two groups G_{1}, G_{2} according to the following rule: For each fragment f_{ i } , if s_{ c } (H_{1}, f_{ i } ) < s_{ c } (H_{2}, f_{ i } ), add f_{ i } to G_{1}, otherwise add f_{ i } to G_{2}. Let ${\mathbb{P}}_{\mathcal{H}}$ denote the partition (G_{1}, G_{2}) obtained by the above rule. A partition $\mathcal{P}=\left({G}_{1},\phantom{\rule{2.77695pt}{0ex}}{G}_{2}\right)$ of a fragment set R is encoded as a map such that $\mathcal{P}\left(f\right)=0$ (or 1) if the fragment f in R belongs to G_{1} (or G_{2}, respectively).
Conversely, given a partition $\mathcal{P}=\left({G}_{1},\phantom{\rule{2.77695pt}{0ex}}{G}_{2}\right)$, a haplotype pair $\mathcal{H}=\left({H}_{1},{H}_{2}\right)$ can be constructed as follows. Let N_{ g,v } [j] be the number of fragments of G_{ g } whose allele at the j th locus is v for g = 1, 2 and v = 0, 1. For each SNP locus j, if N_{1,1}[j] + N_{2,0}[j] ≤ N_{1,0}[j] + N_{2,1}[j], H_{1}[j] = 0 and H_{2}[j] = 1; otherwise, H_{1}[j] = 1 and H_{2}[j] = 0. Let ${\mathbb{H}}_{\mathcal{P}}$ denote the haplotype pair obtained by the above method.
Theorem 1 Given an SNP matrix M and a haplotype pair $\mathcal{H},{s}_{c}\left(\mathcal{H},\phantom{\rule{2.77695pt}{0ex}}M\right)\ge {s}_{c}(\mathcal{H},\phantom{\rule{2.77695pt}{0ex}}M)$.
There is another formulation of MEC equivalent to the above one: Given an SNP matrix M, find a partition $\mathcal{P}$ of the rows in M such that ${s}_{c}\left({\mathbb{H}}_{\mathcal{P}},M\right)$ is minimized. In the following, ${s}_{c}\left({\mathbb{H}}_{\mathcal{P}},M\right)$ is called the errors corrected measure of$\mathcal{P}$. While MEC aims to find a partition such that the conflict between the fragments in the same group is minimized, Maximum Fragments Cut (MFC) [16] aims to find a partition such that the conflict between the fragments in G_{1} and the fragments in G_{2} is maximized.
Maximum Fragments Cut (MFC)[16]: Given an SNP matrix M, find a partition $\mathcal{P}$ of the rows of M such that ${s}_{d}\left(\mathcal{P},M\right)$ is maximized.
where w is a weight factor that is used to adjust the weight of the fragments cut measure. In the following we propose a new optimization model for the SIH problem.
Balanced Optimal Partition (BOP): Given an SNP matrix M, find a partition $\mathcal{P}$ of the rows in M such that ${s}_{p}\left(\mathcal{P},M\right)$ is minimized. A solution to BOP of M is denoted by BOP(M), i.e. a partition with the minimum partition score.
Note that when w = 0, BOP becomes MEC which has been proved NP-hard [24] and APX-hard [26]. Therefore, BOP is NP-hard and APX-hard.
Algorithm
Given an m × n SNP matrix M, there are 2^{m- 1}different partitions of m rows in M. Therefore, when m is large, it is impractical to enumerate all possible partitions and choose one with the minimum partition score. To solve the BOP model of M efficiently, we propose a dynamic programming algorithm in the subsection. We first consider the first row of M, then the first two rows and so on until we have considered all rows of M.
We need some definitions and notations. Let M [1..i, :] denote the SNP matrix consisting of only the first i rows of M. The first and last columns at which the i th row of M takes non '−' values are denoted by l(i) and r(i), respectively. For a column j, if l(i) ≤ j ≤ r(i), row i spans column j. In the following, we assume that all the rows of M are sorted such that if i_{1}< i_{2}, l(i_{1}) < l(i_{2}) or l(i_{1}) = l(i_{2}) ∧ r(i_{1}) ≤ r(i_{2}). Let R(i) denote the row set containing the rows in M[1..i, :] that span column l(i).
Let $\mathcal{P}$ be a partition of a row set R and ${\mathcal{P}}^{\prime}$ a partition of a subset ${R}^{\prime}$ of R. If for every row i ∈ R', $\mathcal{P}\left(i\right)={\mathcal{P}}^{\prime}\left(i\right)$, ${\mathcal{P}}^{\prime}$ is called the projection of $\mathcal{P}$ on R' and $\mathcal{P}$ is called an extension of ${P}^{\prime}$ on R. Fix a row i and let $\mathcal{P}$ be a partition of R(i). ${\mathcal{P}}^{\prime}$ is an optimal extension of $\mathcal{P}$, if the following conditions hold: (1) ${\mathcal{P}}^{\prime}$ is an extension of $\mathcal{P}$ on the row set R = {1, ..., i}; (2) for any possible extension ${\mathcal{P}}^{\u2033}$ of $\mathcal{P}$ on R, ${s}_{p}\left({\mathcal{P}}^{\prime},\phantom{\rule{2.77695pt}{0ex}}M\left[1..i,\phantom{\rule{2.77695pt}{0ex}}:\right]\right)\le {s}_{p}\left({\mathcal{P}}^{\u2033},\phantom{\rule{2.77695pt}{0ex}}M\left[1..i,\phantom{\rule{2.77695pt}{0ex}}:\right]\right)$.
Given a partition $\mathcal{P}$ of R(i), let ${\epsilon}_{i}\left(\mathcal{P}\right)$ denote an optimal extension of $\mathcal{P}$. We call ${s}_{p}\left({\mathcal{E}}_{i}\left(\mathcal{P}\right),\phantom{\rule{2.77695pt}{0ex}}M\left[\mathsf{\text{l}}..i,\phantom{\rule{2.77695pt}{0ex}}:\right]\right)$ the partition score of $\mathcal{P}$ and denote it as ${s}_{p}^{i}\left(\mathcal{P}\right)$ for briefness.
Theorem 2 For an m × n SNP matrix M, let$\mathcal{P}$be a partition of$R\left(m\right).\phantom{\rule{2.77695pt}{0ex}}{\mathcal{E}}_{m}\left(\mathcal{P}\right)$ is a solution to BOP of M if the following condition holds: for any possible partition ${\mathcal{P}}^{\prime}$of R(i), ${s}_{p}^{m}\left(\mathcal{P}\right)\le {s}_{p}^{m}\left({\mathcal{P}}^{\prime}\right)$.
In Equation (10) (or 11), ${\delta}_{c}\left({{\mathcal{P}}_{1}}^{\prime}\right)\phantom{\rule{2.77695pt}{0ex}}\left(\mathsf{\text{or}}{\delta}_{c}\left({{\mathcal{P}}_{2}}^{\prime}\right)\right)$ is the difference between the errors corrected measures of ${\mathcal{E}}_{i+1}\left({{\mathcal{P}}_{1}}^{\prime}\right)\phantom{\rule{2.77695pt}{0ex}}\left(\mathsf{\text{or}}{\mathcal{E}}_{i+1}\left({{\mathcal{P}}_{2}}^{\prime}\right)\right)$ and ${\mathcal{E}}_{i}\left({\mathcal{P}}^{\u2033}\right)$; and ${\delta}_{d}\left({{\mathcal{P}}_{1}}^{\prime}\right)\phantom{\rule{2.77695pt}{0ex}}\left(\mathsf{\text{or}}{\delta}_{d}\left({{\mathcal{P}}_{2}}^{\prime}\right)\right)$ is the difference between the fragments cut measures of ${\mathcal{E}}_{i+1}\left({{\mathcal{P}}_{1}}^{\prime}\right)\phantom{\rule{2.77695pt}{0ex}}\left(\mathsf{\text{or}}{\mathcal{E}}_{i+1}\left({{\mathcal{P}}_{2}}^{\prime}\right)\right)$ and ${\mathcal{E}}_{i}\left({\mathcal{P}}^{\u2033}\right)$. The values ${\delta}_{c}\left(\mathcal{P}\right)$ and ${\delta}_{d}\left(\mathcal{P}\right)$ are calculated by calling the following functions $DeltaEC\left(i,\mathcal{P}\right)$ and $DeltaFC\left(i,\mathcal{P}\right)$, respectively.
$DeltaEC\left(i,\mathcal{P}\right)$
{ δ = 0;
for j = l(i + 1),..., r(i + 1) do
{ if M[i + 1, j] = = '-' then continue;
N_{1,0} = N_{2,0} = N_{1,1} = N_{2,1} = 0;
for each row l ∈ R_{ c } (i, i + 1) do
{ if M[l, j] = = '-' then continue;
v = M[l, j], g = $\mathcal{P}\left(l\right)$,N_{ g,v } + +; }
if N_{1,1} + N_{2,0} ≤ N_{1,0} + N_{2,1}then
{ δ = δ - (N_{1,1} + N_{2,0)}; }
else { δ = δ - (N_{1,0} + N_{2,1)}; }
v = M[i + 1, j], g = $\mathcal{P}\left(i+1\right)$, N_{ g,v }+ +;
if N_{1,1} + N_{2,0} ≤ N_{1,0} + N_{2,1}then
{ δ = δ + (N_{1,1} + N_{2,0}); }
else { δ = δ + (N_{1,0} + N_{2,1}); }
}
return δ; }
$DeltaFC\left(i,\mathcal{P}\right)$
{ δ = 0, g_{0} = $\mathcal{P}\left(i+1\right)$;
for j = l(i + 1), ..., r(i + 1) do
{ if M[i + 1, j] = = '-' then continue;
for each row l ∈ R_{ c } (i, i + 1) do
{ if M[l, j] = = '-' then continue;
v = M[l, j], g = $\mathcal{P}\left(l\right)$;
if g == g_{0}then continue;
if v == M[i + 1, j] then δ - -;
else δ + +; }
}
return δ; }
The time complexity of H-BOP is O(mkk_{1}k_{2}), and the space complexity is O(mk_{1} + mk), where ${k}_{1}=\underset{i=1,\dots ,m}{\text{max}}\left(r\left(i\right)-l\left(i\right)+1\right)$ and ${k}_{2}=\underset{i=1,\dots ,m}{\text{max}}\left(\left|R\left(i\right)\right|\right)$.
Declarations
Acknowledgements
This research was supported by the National Natural Science Foundation of China under Grant Nos. 61070145 and 61232001.
This article has been published as part of BMC Systems Biology Volume 6 Supplement 2, 2012: Proceedings of the 23rd International Conference on Genome Informatics (GIW 2012). The full contents of the supplement are available online at http://www.biomedcentral.com/bmcsystbiol/supplements/6/S2.
Authors’ Affiliations
References
- Levy S, et al.: The diploid genome sequence of an individual human. PLoS Biology. 2007, 5 (10): e254-e254. 10.1371/journal.pbio.0050254.PubMed CentralView ArticlePubMedGoogle Scholar
- Tewhey R, Bansal V, Torkamani A, Topol EJ, Schork NJ: The importance of phase information for human genomics. Nat Rev Genet. 2011, 12 (3): 215-23. 10.1038/nrg2950.PubMed CentralView ArticlePubMedGoogle Scholar
- Casey JP, et al.: A novel approach of homozygous haplotype sharing identifies candidate genes in autism spectrum disorder. Hum Genet. 2012, 131 (4): 565-79. 10.1007/s00439-011-1094-6.PubMed CentralView ArticlePubMedGoogle Scholar
- Fan HC, Wang J, Potanina A, Quake SR: Whole-genome molecular haplotyping of single cells. Nat Biotechnol. 2011, 29: 51-7. 10.1038/nbt.1739.PubMed CentralView ArticlePubMedGoogle Scholar
- The 1000 Genomes Project Consortium: A map of human genome variation from population-scale sequencing. Nature. 2010, 467 (7319): 1061-73. 10.1038/nature09534.PubMed CentralView ArticleGoogle Scholar
- Zhang XS, Wang RS, Wu LY, Chen L: Models and algorithms for haplotyping problem. Current Bioinformatics. 2006, 1: 105-114. 10.2174/157489306775330570.View ArticleGoogle Scholar
- Xie M, Wang J, Chen J, Wu J, Liu X: Computational models and algorithms for the single individual haplotyping problem. Current Bioinformatics. 2010, 5: 18-28. 10.2174/157489310790596411.View ArticleGoogle Scholar
- Browning SR, Browning BL: Haplotype phasing: existing methods and new developments. Nat Rev Genet. 2011, 12 (10): 703-14. 10.1038/nrg3054.PubMed CentralView ArticlePubMedGoogle Scholar
- The International HapMap Consortium: A haplotype map of the human genome. Nature. 2005, 437 (7063): 1299-1320. 10.1038/nature04226.PubMed CentralView ArticleGoogle Scholar
- Kitzman JO, et al.: Haplotype-resolved genome sequencing of a Gujarati Indian individual. Nat Biotechnol. 2011, 29: 59-63. 10.1038/nbt.1740.PubMed CentralView ArticlePubMedGoogle Scholar
- Suk EK, et al.: A comprehensively molecular haplotype-resolved genome of a European individual. Genome Res. 2011, 21 (10): 1672-85. 10.1101/gr.125047.111.PubMed CentralView ArticlePubMedGoogle Scholar
- Duitama J, et al.: Fosmid-based whole genome haplotyping of a HapMap trio child: evaluation of Single Individual Haplotyping techniques. Nucleic Acids Res. 2012, 40 (5): 2041-53. 10.1093/nar/gkr1042.PubMed CentralView ArticlePubMedGoogle Scholar
- Lancia G, Bafna V, Istrail S, Lippert R, Schwartz R: SNPs problems, complexity and algorithms. Proc. Ann. European Symp. on Algorithms (ESA), Volume 2161 of Lecture Notes in Computer Science. Edited by: auf der Heide FM. 2001, Berlin/Heidelberg: Springer, 182-193.Google Scholar
- Geraci F: A comparison of several algorithms for the single individual SNP haplotyping reconstruction problem. Bioinformatics. 2010, 26 (18): 2217-25. 10.1093/bioinformatics/btq411.PubMed CentralView ArticlePubMedGoogle Scholar
- Xie M, Wang J, Chen J: A model of higher accuracy for the individual haplotyping problem based on weighted SNP fragments and genotype with errors. Bioinformatics. 2008, 24 (13): i105-13. 10.1093/bioinformatics/btn147.PubMed CentralView ArticlePubMedGoogle Scholar
- Duitama J, Huebsch T, McEwen G, Suk EK, Hoehe MR: ReFHap: a reliable and fast algorithm for single individual haplotyping. Proceedings of the First ACM international Conference on Bioinformatics and Computational Biology. 2010, Niagara Falls, New York: ACM, 160-169.View ArticleGoogle Scholar
- He D, Choi A, Pipatsrisawat K, Darwiche A, Eskin E: Optimal algorithms for haplotype assembly from whole-genome sequence data. Bioinformatics. 2010, 26 (12): i183-90. 10.1093/bioinformatics/btq215.PubMed CentralView ArticlePubMedGoogle Scholar
- Bansal V, Bafna V: HapCUT: an efficient and accurate algorithm for the haplotype assembly problem. Bioinformatics. 2008, 24 (16): i153-9. 10.1093/bioinformatics/btn298.View ArticlePubMedGoogle Scholar
- Wang J, Xie M, Chen J: A practical exact algorithm for the individual haplotyping problem MEC/GI. Algorithmica. 2010, 56 (3): 283-296. 10.1007/s00453-009-9288-1.View ArticleGoogle Scholar
- Xie M, Wang J, Chen J: A practical parameterised algorithm for the individual haplotyping problem MLF. Mathematical Structures in Computer Science. 2010, 20 (5): 851-863. 10.1017/S096012951000023X.View ArticleGoogle Scholar
- Bafna V, Istrail S, Lancia G, Rizzi R: Polynomial and APX-hard cases of the individual haplotyping problem. Theoretical Computer Science. 2005, 335: 109-125. 10.1016/j.tcs.2004.12.017.View ArticleGoogle Scholar
- Panconesi A, Sozio M: Fast hare: a fast heuristic for single individual SNP haplotype reconstruction. Proc. WABI, Volume 3240 of Lecture Notes in Computer Science. Edited by: Jonassen I, Kim J. 2004, Berlin/Heidelberg: Springer, 266-277.Google Scholar
- Wang RS, Wu LY, Li ZP, Zhang XS: Haplotype reconstruction from SNP fragments by minimum error correction. Bioinformatics. 2005, 21 (10): 2456-2462. 10.1093/bioinformatics/bti352.View ArticlePubMedGoogle Scholar
- Lippert R, Schwartz R, Istrail S: Algorithmic strategies for the single nucleotide polymorphism haplotype assembly problem. Brief. Bioinform. 2002, 3: 1-9.View ArticleGoogle Scholar
- Genovese LM, Geraci F, Pellegrini M: SpeedHap: an accurate heuristic for the single individual SNP haplotyping problem with many gaps, high reading error rate and low coverage. IEEE/ACM Trans Comput Biol Bioinform. 2008, 5 (4): 492-502.View ArticlePubMedGoogle Scholar
- Cilibrasi R, van Iersel L, Kelk S, Tromp J: The complexity of the single individual SNP haplotyping problem. Algorithmica. 2007, 49: 13-36. 10.1007/s00453-007-0029-z.View ArticleGoogle Scholar
- SIH. [http://www.molgen.mpg.de/~genetic-variation/SIH/]
- Lo C, Bashir A, Bansal V, Bafna V: Strobe sequence design for haplotype assembly. BMC Bioinformatics. 2011, 12 (Suppl 1): S24-10.1186/1471-2105-12-S1-S24.PubMed CentralView ArticlePubMedGoogle Scholar
- Xie M, Wang J: An improved (and practical) parameterized algorithm for the individual haplotyping problem MFR with mate-pairs. Algorithmica. 2008, 52 (2): 250-266. 10.1007/s00453-007-9150-2.View ArticleGoogle Scholar
Copyright
This article is published under license to BioMed Central Ltd. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.