 Research
 Open access
 Published:
A fast and accurate enumerationbased algorithm for haplotyping a triploid individual
Algorithms for Molecular Biology volume 13, Article number: 10 (2018)
Abstract
Background
Haplotype assembly, reconstructing haplotypes from sequence data, is one of the major computational problems in bioinformatics. Most of the current methodologies for haplotype assembly are designed for diploid individuals. In recent years, genomes having more than two sets of homologous chromosomes have attracted many research groups that are interested in the genomics of disease, phylogenetics, botany and evolution. However, there is still a lack of methods for reconstructing polyploid haplotypes.
Results
In this work, the minimum error correction with genotype information (MEC/GI) model, an important combinatorial model for haplotyping a single individual, is used to study the triploid individual haplotype reconstruction problem. A fast and accurate enumerationbased algorithm enumeration haplotyping triploid with least difference (EHTLD) is proposed for solving the MEC/GI model. The EHTLD algorithm tries to reconstruct the three haplotypes according to the order of single nucleotide polymorphism (SNP) loci along them. When reconstructing a given SNP site, the EHTLD algorithm enumerates three kinds of SNP values in terms of the corresponding site’s genotype value, and chooses the one, which leads to the minimum difference between the reconstructed haplotypes and the sequenced fragments covering that SNP site, to fill the SNP loci being reconstructed.
Conclusion
Extensive experimental comparisons were performed between the EHTLD algorithm and the well known HapCompass and HapTree. Compared with algorithms HapCompass and HapTree, the EHTLD algorithm can reconstruct more accurate haplotypes, which were proven by a number of experiments.
Background
As a large number of sequencing data are available, the investigation of genetic variations has become one of the main topics in bioinformatics. Single nucleotide polymorphism (SNP), the most widespread form of variation, is believed to be the major genetic cause to phenotypic variability. A sequence of SNPs along a chromosome is referred to as a haplotype, which is more important for complete comprehending the complex genetic polymorphisms than isolated SNPs. Increasing evidence shows that haplotypes play a crucial role in studying the variations relating to diseases prediction and gene expression [1]. Therefore, computational methods to infer haplotypes are needed, for determining haplotypes is both time consuming and expensive by direct using biological experiments. In recent decade, the presented computational haplotyping algorithms generally fall into three categories [2]: (1) population haplotyping with genotype data [3, 4]; (2) population haplotyping with fragment data [5]; (3) individual haplotyping with fragment data [6]. In this paper, individual haplotyping problem is studied for a triploid individual.
The problem of individual haplotyping is also called as haplotype assembly problem or haplotype reconstruction problem. It has received extensive study in the recent decade. Most of the existing research results are regarding diploid individuals [1, 7, 8], and there is still a lack of research studies for reconstructing triploid ones. Several algorithms for assembling Kindividual haplotypes were proposed. Based on the minimum error correction (MEC) model and the minimum error correction with genotype information (MEC/GI) model, Wang et al. [9] and Qian et al. [10] respectively proposed a genetic algorithm and a particle swarm optimization algorithm to reconstruct diploid individual haplotypes, both of which can be adapted to reconstructing the Kindividual ones. The code length of the two algorithms is very long in practical applications, for it is equal to the number of sequencing SNP fragments. This brings huge solution space to these two algorithms and negatively affects the performance of them. Based on the minimum fragment removal (MFR) model [11], an exact exponential algorithm was introduced by Li et al. [11]. The time complexity of which is \(O(2^{2t}m^{2}n+2^{(K+1)t}m^{K+1})\), where m denotes the number of SNP fragments, n denotes the number of SNP sites, and t is the max number of holes covered by a fragment. The algorithm can not perform well with large m, n and t. In 2013, Aguiar et al. [12] introduced the HapCompass model and the minimum weighted edge removal (MWER) optimization for haplotyping polyploid genomes. Algorithm HapCompass aims to remove a minimal weighted set of edges from the compass graph such that a unique phasing may be constructed. The HapCompass algorithm performs on the spanningtree cycle basis of the compass graph to iteratively delete errors. However, in the same conflict cycle basis, there may be more than one edge having the same absolute value of weight. It may lead to the wrong SNP phasing to select the removed edge randomly. In 2014, Berger et al. [13] described a maximumlikelihood estimation framework HapTree for haplotyping a single polyploid. It can obtain better performance than the HapCompass algorithm [13]. In 2014, based on the MEC model, Wu et al. [14] presented a genetic algorithm GTIHR for reconstructing triploid haplotypes. Since the code length of algorithm GTIHR equals to the number of heterozygous sites in haplotype, the performance of the GTIHR algorithm is negatively affected by haplotype length and heterozygous rate. In this paper, the triploid individual haplotype assembly problem is studied based on the MEC/GI model. An enumerationbased algorithm enumeration haplotyping triploid with least difference (EHTLD) is proposed for solving it. Algorithm EHTLD reconstructs the three haplotypes according to the order of SNP loci along them. For reconstructing the three alleles of a given site, it enumerates three kinds of SNP values by using the site’s genotype, and chooses the kind of value resulting in the minimum difference between the reconstructed haplotypes and the sequenced fragments covering that SNP site. The experimental comparisons were performed between the EHTLD, the HapCompass and the HapTree algorithms. The results proved that the performance of algorithm EHTLD was superior to those of algorithms HapCompass and HapTree. The rest of this paper is arranged as follows. “Definitions and notations” section provides definitions and notations used later. “Algorithm EHTLD” section introduces the EHTLD algorithm. “Experimental results” section presents the experimental results of the EHTLD, the HapCompass and the HapTree algorithms. Some conclusions are drawn in the last section.
Definitions and notations
Triploid somatic cells contain three sets of chromosomes, i.e., a triploid organism has three copies of each chromosome. Since haplotype consists of the sequence of all SNPs along a chromosome, a triploid individual owns three haplotypes. It is commonly regarded that a SNP locus shows merely two possible alleles, hence the major allele can be represented as \(\text{`}0\text{'}\) and the minor one can be represented as \(\text{`}1\text{'}\). A haplotype can be encoded as a string over a 2letter alphabet \(\{0,1\}\) instead of four real bases {A,T,C,G}. A genotype is the conflation of three haplotypes on the homologous chromosomes. When three alleles at a SNP site have identical values, this SNP site is called a homozygous site, otherwise it is called a heterozygous site. For example, \((000)^{T}\) or \((111)^{T}\) represents the genotype value at a homozygous SNP site, while \((001)^{T}\) or \((011)^{T}\) represents the genotype value at a heterozygous SNP site. Suppose that m aligned SNP fragments, coming from three haplotypes of length n, are generated by DNA sequencing experiments. Let M denote an \(m\times n\) SNP matrix over the alphabet \(\{0,1,\}\) (− denotes the value is null). As shown in Fig. 1a, each row represents a SNP fragment, each column represents a SNP site, and each entry M[i, j] denotes the SNP allele of the ith fragment at the jth SNP site. Let \(G=(g_{1},g_{2},\ldots ,g_{n})\) denote the genotype matrix corresponding to M, where \(g_{j}=(g_{j1},g_{j2},g_{j3})^{T} (g_{jk}\in \{0,1\}, \quad k=1,2,3, \quad j=1,2,\ldots ,n)\) denotes the genotype value at the jth SNP site. Figure 1b shows an example of the genotype matrix.
Given a column \(M[,j]\,(j=1,2,\ldots ,n)\) of the matrix M, define r(j) as the set of fragments that cover the jth column. Given a row \(M[i,]\,(i=1,2,\ldots ,m)\) of the matrix M, let l(i) indicate the index of the leftmost SNP \(j\,(j=1,2,\ldots ,n)\) such that \(M[i,j]\ne \). Given two strings \(X=x_{1},x_{2},\ldots ,x_{n}\) and \(Y=y_{1},y_{2},\ldots ,y_{n}\), where \(x_{j}, y_{j}\in \{0,1,\}\, (j=1,2,\ldots ,n)\), the distance metric HD(X, Y, s, e) is defined as Formula (1). Take fragment \(f_{1}(10011)\) and fragment \(f_{2}(01010)\) in Fig. 1a for example. \(HD(f_{1}, f_{2}, 2, 5)=3\).
where
Let the strings X and Y be regarded as two SNP fragments, they are said to compatible if HD(X, Y, 1, n) = 0. The larger HD(X, Y, 1, n) is, the greater the probability of fragments X and Y coming from different chromosome copies or having sequencing errors is. If there are no errors in the data, the rows of M can be divided into three classes of compatible fragments. Three haplotypes can be reconstructed by assembling the fragments in the three classes. In this situation, the SNP matrix M is called feasible or errorfree. Given haplotype \(h=(h_{1}, h_{2}, h_{3})\,(h_{k}=(h_{k1}, h_{k2}, \ldots , h_{kn}), k=1,2,3)\) and genotype \(G=(g_{1}, g_{2}, \ldots , g_{n}) (g_{j}=(g_{j1}, g_{j2}, g_{j3})^{T}, j=1,2,\ldots ,n)\), if \(\sum \nolimits _{k=1}^3 h_{kj}=\sum \nolimits _{k=1}^3 g_{jk} (j=1,2,\ldots ,n), h\) and G are regarded as compatible.
Based on the above mentioned concepts for the haplotype reconstruction problem, the MEC/GI model can be described as follows [9]:
MEC/GI: Given a SNP matrix M and a genotype matrix G, correct the minimum number of entries in M (0 into 1 and vice versa) so that the resulting matrix is feasible, and the three reconstructed haplotypes are compatible with the genotype G.
Algorithm EHTLD
In this section, the EHTLD algorithm is described. The input consists of a SNP matrix M and a genotype matrix G. The output is three assembled haplotypes h = (\(h_{1}\), \(h_{2}\), \(h_{3}\)) of length n. In the first step of this algorithm, the matrices M and G are preprocessed by removing the homozygous SNPs, which do not play a role in assembling haplotypes. Subsequently, enumerates three kinds of values for the jth (j = 2, 3,…,n) SNP site in terms of its genotype, and chooses the one leading to the minimum difference between the reconstructed haplotypes and the fragments covering the jth site. After this iteration process is completed, three haplotypes \(h'=(h'_{1}, h'_{2}h'_{3}\)) having only heterozygous SNP sites are built, for only heterozygous SNPs are remained in the preprocessed matrices. Finally, \(h'\) is augmented by inserting the SNPs discarded in preprocessing step and the final haplotypes h is obtained. Some key steps of the EHTLD algorithm will be introduced in detail as follows.
Preprocessing
Since homozygous sites do not contribute to haplotype reconstruction, they are deleted from matrices M and G to improve the efficiency of assembly. Drop column \(j (j=1, 2,\ldots ,n)\) from G where \(g_{j1}=g_{j2}=g_{j3}\), and the corresponding column is also dropped from matrix M. The deleted column j (\(j=1, 2,\ldots ,n\)) is recorded as \(g_{j1}\). After dropping columns from matrix M, some SNP fragments with only—elements are also deleted, for they are also redundant information. The remained SNPs are all heterozygous sites. For convenience of description, the preprocessed matrices are still denoted by M and G. Sort the rows of M by their l(.) values in ascending order. For each column j (\(j=1, 2,\ldots ,n\)) of M, calculate set r(j) which contains the rows covering the jth column.
Enumerating and computing
The EHTLD algorithm iteratively reconstructs each heterozygous site of haplotypes \(h'\) = (\(h'_{1}\), \(h'_{2}\), \(h'_{3}\)). Each step concerns reconstructing the current empty site, starting from the left first site. Suppose that the first j − 1 sites of the three haplotypes \(h'\) have already been filled, i.e., (\(h'_{k1}\), \(h'_{k2}\),…, \(h'_{kj1}\)) (k = 1, 2, 3, j = 2, 3,…, n) has been assembled, and the jth site is under consideration. The calculating method comprises the following two steps.

1.
Enumerating three kinds of possible values according to \(g_{j}\):

a.
if \(\sum \nolimits _{k=1}^3 g_{jk}\) = 1, the three kinds of values are (\(h'_{1j}=0, h'_{2j}=0, h'_{3j}=1\)), (\(h'_{1j}=0, h'_{2j}=1, h'_{3j}=0\)) and (\(h'_{1j}=1, h'_{2j}=0, h'_{3j}=0\)).

b.
if \(\sum \nolimits _{k=1}^3 g_{jk}\) = 2, the three kinds of values are (\(h'_{1j}\) = 0, \(h'_{2j}\) = 1, \(h'_{3j}\) = 1), (\(h'_{1j}\) = 1, \(h'_{2j}\) = 0, \(h'_{3j}\) = 1) and (\(h'_{1j}\) = 1, \(h'_{2j}\) = 1, \(h'_{3j}\) = 0).

a.

2.
Given the jth site value (\(h'_{1j}\), \(h'_{2j}\),\(h'_{3j}\)), let D(\(h'_{1j}\), \(h'_{2j}\), \(h'_{3j}\)) measure the difference between the reconstructed haplotypes and the fragments covering the jth site, as defined in Formula (3). From the three kinds of values enumerated in step (1), choose the one with the minimum D(.) value.
$$\begin{aligned} D(h'_{1j},h'_{2j},h'_{3j})=\sum \limits _{i\in r(j)} min\{HD(h'_{k},M[i,],l(i),j)k=1,2,3\} \end{aligned}$$(3)
In the following, we give an example for enumerating and computing by using the matrices in Fig. 1. As shown in Fig. 2, assume that the first three sites of the three haplotypes \(h'\) = (\(h'_{1}\), \(h'_{2}\), \(h'_{3}\)) have been reconstructed, i.e., \(h'\) = (\(h_{1}'\)(011), \(h_{2}'\)(010), \(h_{3}'\)(100)), and the fourth site is under reconstruction. The genotype of the fourth SNP site is (101)\(^{T}\), hence haplotypes \(h'\) = (\(h'_{1}\), \(h'_{2}\), \(h'_{3}\)) have the following three kinds of possible values on the fourth SNP site: (\(h_{14}'\) = 0, \(h_{24}'\) = 1, \(h_{34}'\) = 1), (\(h_{14}'\) = 1, \(h_{24}'\) = 0, \(h_{34}'\) = 1) and (\(h_{14}'\) = 1, \(h_{24}'\) = 1, \(h_{34}'\) = 0). The values of D(0,1,1), D(1,0,1) and D(1,1,0) are computed respectively according to the fragments in Fig. 1a and the three haplotypes \(h'\) = (\(h'_{1}\), \(h'_{2}\), \(h'_{3}\)). D(0,1,1) = 3, D(1,0,1) = 2, D(1,1,0) = 1. Because D(1,1,0) is the smallest, (\(h_{14}'\) = 1, \(h_{24}'\) = 1, \(h_{34}'\) = 0) is chosen, and \(h'\) = (\(h_{1}'\)(0111), \(h_{2}'\)(0101), \(h_{3}'\)(1000)) are reconstructed.
Augmenting
The homozygous SNPs that are deleted by preprocessing must be reinserted. The reconstructed haplotypes \(h'\) = (\(h'_{1}\), \(h'_{2}\), \(h'_{3}\)) are augmented by the bits of the columns removed, and h = (\(h_{1}\), \(h_{2}\), \(h_{3}\)) are built. For a given position j, haplotypes \(h'_{1}\), \(h'_{2}\) and \(h'_{3}\) are inserted with \(g_{j1}\) when the discarded column j is recorded as \(g_{j1}\). Based on the above mentioned steps, the EHTLD algorithm for assembling triploid haplotypes is depicted in Fig. 3.
Now the time complexity of the EHTLD algorithm is discussed. In preprocessing, dropping redundant information and calculating set \(r(\cdot )\) take time \(O(m\times n)\), sorting the rows of M takes time \(O(m\times logm)\). During enumerating and computing, three haplotypes with only heterozygous SNP sites are reconstructed, which takes time \(O(c \times n\times len)\), here c denotes the fragments coverage, and len represents the average length of fragments. In augmenting, the discarded columns can be reinserted by scanning the columns only once, which takes time O(n). In summary, the time complexity of the algorithm is \(O(m\times n+m \times logm+ c \times n \times len)\).
Experimental results
In this section, the EHTLD algorithm is compared with two stateoftheart algorithms, i.e., the HapCompass [12] and the HapTree [13] algorithms. Algorithms EHTLD and HapCompass were implemented on a Windows 7 and the compiler was Microsoft Visual C# 2012. The Python program HapTree (v0.1), downloaded from http://groups.csail.mit.edu/cb/haptree/, was implemented on a Linux system. All the tests below are conducted on a 64 bit PC with Intel Core i5 2.50GHz CPU and 6GB RAM. One hundred data sets were generated for each parameter setting. The average over 100 runs at each parameter setting was calculated and presented.
The vector error (VE) [13], the reconstruction rate (RR) [1, 9, 15] and the minimum error correction (MEC) score [12] were used to measure the performance of the algorithms. The vector error, generalized from switch error, is a special kind of measurement for evaluating the accuracy of polyploid phasing. Given three reconstructed haplotypes, the vector error is equal to the minimum number of segments on them for which a switch must occur to correspond with the three true haplotypes, i.e., the minimum number of segments a reconstructed phase and the true phase have in common [13].
The reconstruction rate (RR), which measures the similarity degree between the pair of true haplotypes and the pair of reconstructed ones, is a widely adopted index to evaluate diploid phasing [1, 9, 15]. For triploid phasing, we generalized it to calculate the similarity degree between the three true haplotypes and the three reconstructed haplotypes. Assuming that h = (\(h_{1}\), \(h_{2}\), \(h_{3}\)) are the original haplotypes, and \(\hat{h}\) = (\(\hat{h}_{1}\), \(\hat{h}_{2}\), \(\hat{h}_{3}\)) are the reconstructed haplotypes. RR is defined as the proportion of nucleotides that are reconstructed correctly, as shown in Formula (4).
where \(r_{i_{k}j_{k}}\) = \(HD(h_{i_{k}},\hat{h}_{j_{k}},1,n)\).
The minimum error correction (MEC) score measures the minimum number of mismatches between the reconstructed haplotypes \(\hat{h}\) = (\(\hat{h}_{1}\), \(\hat{h}_{2}\), \(\hat{h}_{3}\)) and the SNP matrix M, as shown in Formula (5).
To the best of our knowledge, the real triploid haplotype data are not available in the public domain, Aguiar et al. [12] and Berger et al. [13] used computergenerated simulated data. Therefore, simulated data were also used in our experiments. Three simulation haplotypes h = (\(h_{1}\), \(h_{2}\), \(h_{3}\)) of length n were created by using the following method. \(h_{1}\) was generated at random firstly. \(h_{2}\) was generated by flipping each bit of \(h_{1}\) randomly so that the hamming distance between \(h_{1}\) and \(h_{2}\) was equal to a given parameter d. \(h_{3}\) also had the same length and \(h_{3j}\) was set to \(h_{1j}\) or \(h_{2j}\,(j=1,2,\ldots ,n)\) with uniform probability. With regard to fragment data, two kinds of sequencing simulators, CELSIM [16] and MetaSim [17], were adopted to generate simulation fragments, and the testing datasets were called as CELSIM instances and MetaSim instances, respectively.
CELSIM instances
In this section, the evaluation of the EHTLD, the HapCompass and the HapTree algorithms is described by using CELSIM instances. CELSIM was invoked to simulate shotgun sequencing platform. m_{1} single SNP fragments and m_{2} matepair SNP fragments were generated. A single fragment had a length ranging from \(f_{min}\) to \(f_{max}\), and a matepair fragment had a length of n/10. The coverage was c/2 for both kinds of fragments, and the total coverage was c. Reading errors were planted into the fragments with probability \(p_{s}\). In practical applications of shotgun sequencing, the values of \(f_{min}\) and \(f_{max}\) are 3 and 7, respectively, c ranges from 5 to 10, and \(p_{s}\) ranges between 2 and 5% [2, 18]. In the following tables, algorithms EHTLD, HapCompass and HapTree are abbreviated to EH, HC and HT, respectively.
In Table 1, 12 sets of parameters were set in dealing with error rate \(p_{s}\), where c = 10, \(f_{min}\) = 3, \(f_{max}\) = 7, n = 100 and d = 0.3. It can be seen from this table that algorithm EHTLD can achieve much higher reconstruction rates, smaller vector errors and MEC scores than the HapCompass and the HapTree algorithms in every \(p_{s}\) setting. When \(p_{s}\) = 0, algorithm EHTLD achieves reconstruction rate of 0.97, which is higher than both HapCompass and HapTree algorithms by about 9.0%, and vector error of 3, which is less than them by 8 times or so. In particular, the MEC score obtained by algorithm EHTLD reaches zero, while those of the other two algorithms are 126 and 57. Although the increase of \(p_{s}\) plays stronger negative effect on algorithm EHTLD than on algorithms HapCompass and HapTree, the EHTLD algorithm still obtains better performance than algorithms HapCompass and HapTree with high error rate. When \(p_{s}\) = 0.2, the RRs of algorithms EHTLD, HapCompass and HapTree are 0.92, 0.89 and 0.88, the vector errors of them are 14, 31 and 26, and the MEC scores of them are 335, 407 and 364, respectively. The three algorithms all execute very efficiently when \(p_{s}\) ranges from 0 to 0.2.
In Table 2, nine sets of parameters were set in dealing with coverage c, where n = 100, \(f_{min}\) = 3, \(f_{max}\) = 7, \(p_{s}\) = 0.05 and d = 0.3. From Table 2 we observe that algorithm EHTLD still obtains the highest reconstruction rate and the smallest vector error and MEC score under different coverage settings. When the coverage is 2, the RRs of algorithms EHTLD, HapCompass and HapTree are 0.94, 0.89 and 0.86, the vector errors of them are 10, 30 and 29, and the MEC scores of them are 16, 40 and 19. When the coverage increases, the RR of algorithm EHTLD increases gradually, while that of algorithm HapCompass fluctuates between 0.89 and 0.90, and that of algorithm HapTree varies between 0.85 and 0.91. Generally, the increase of coverage plays a positive role in the improvement of algorithm performance, for much more original fragment information can be utilized. However, it is not apparent for algorithms HapCompass and HapTree.
Table 3 compares the performance of the three algorithms with different haplotype lengths n, where c = 10, \(f_{min}\) = 3, \(f_{max}\) = 7, \(p_{s}\) = 0.05 and d = 0.3. As can be seen from this table, algorithm EHTLD still obtains superior results to the other two algorithms under each parameter setting. With the increase of haplotype length, the three algorithms experience a gradual degradation in the performance. When n is 100, the RR of algorithm EHTLD is 0.97, which is higher than both HapCompass and HapTree algorithms by about 7.8%, the vector error of algorithm EHTLD is 4, which is less than algorithms HapCompass and HapTree by about 86 and 85%, the MEC score of algorithm EHTLD is 85, which is less than algorithms HapCompass and HapTree by about 56 and 21%, respectively. When n is 1000, the RRs of them decrease to 0.92, 0.88 and 0.86, the vector errors of them increase to 136, 322 and 314, and the MEC scores of them go up to 479, 1029 and 595, respectively. The running time of the three algorithms increases significantly with the increase of n. When n = 100, the running time of algorithms EHTLD, HapCompass and HapTree is 0.01, 0.03 and 0.01 s, respectively, while n = 1000, it increases to 4.36, 4.67 and 4.51 s, respectively.
In Table 4, three groups of parameters were set in dealing with single fragment length range \([f_{min}, f_{max}]\), where c = 10, \(p_{s}\) = 0.05, n = 100 and d = 0.3. As shown in Table 4, algorithm EHTLD still performs the best under different parameter settings. When \([f_{min}, f_{max}]\,=\,[3,7]\), the RRs of algorithms EHTLD, HapCompass and HapTree are 0.97, 0.90 and 0.89, the vector errors of them are 4, 29 and 26, and the MEC scores of them are 85, 194 and 117, respectively. With the decrease of the length of single fragment, the decline of fragments overlap might be disadvantageous for haplotype reconstruction. When \([f_{min}, f_{max}]\,=\,[1,2]\), the RRs decrease to 0.94, 0.90 and 0.88, the vector errors increase to 14, 30 and 28, and the MEC scores drop to 65, 86 and 76, respectively. The decrease of the MEC scores explain the shorter the fragments are, the more probability the fragments agree with the reconstructed haplotypes. The change of single fragment length plays little effect on the running time of the three algorithms.
Table 5 compares the three algorithms with different hamming distances d, where \(f_{min}\) = 3, \(f_{max}\) = 7, c = 10, \(p_{s}\) = 0.05 and n = 100. It can be seen from this table that the performance of algorithm EHTLD remains relatively stable under different d, while that of algorithms HapCompass and HapTree suffers strong negative influence with the increase of hamming distance. For example, when d varies from 0.1 to 1.0, the RR of the EHTLD algorithm fluctuates between 0.97 and 1.0, while those of the HapCompass and the HapTree algorithms achieve decrease rate up to 35 and 20%, respectively.
MetaSim instances
MetaSim was used to simulate 454 sequencing platform. m SNP fragments, including \({m_1}\,=\,(1p_{m})\times m\) single ones and \({m_2}\) = \({p_m} \times m\) matepair ones, were generated, where \(p_{m}\) denoted the probability of matepair fragments and was set to 0.25 in the experiments. A single fragment had an expected length of \(f\_len\), and a matepair fragment had a length of \(3 \times f\_len\). Since each matepair fragment consists of two single fragments of the same haplotype, the coverage c equals to \([({m_1}+2\times {m_2})\times f\_len]/3\times n\).
Table 6 gives the comparisons with coverage ranging from 5 to 50, where n = 100, \(f\_len\) = 5, and d = 0.3. In Table 7, six sets of experimental results under different haplotype length settings are displayed, where c = 20, \(f\_len\) = 5, and d = 0.3. In Table 8, three instances were generated in dealing with single fragment length \(f\_len\), where n = 100, c = 20, and d = 0.3. In Table 9, the test results under different parameter d are shown, where n = 100, c = 20, and \(f\_len\) = 5. The experimental results obtained from MetaSim instances indicate that algorithm EHTLD still obtain much higher reconstruction rates, smaller vector errors and MEC scores than the HapCompass and the HapTree algorithms under different \(c, n, f\_len\) and d settings.
Conclusion
The minimum error correction with genotype information (MEC/GI) model is one of the important computational models for solving single individual SNP haplotyping problem. In this paper, an enumerationbased algorithm EHTLD is presented for haplotyping a triploid single individual by using this model. Algorithm EHTLD reconstructs the three haplotypes according to the order of SNP loci along them. For a SNP site being reconstructed, the EHTLD algorithm enumerates three possible values in terms of the site’s genotype, and chooses the one leading to the minimum difference between the reconstructed haplotypes and the fragments covering that SNP site. The reconstructed alleles of a SNP site mainly depend on the fragments which cover the site, and are little affected by other former reconstructed alleles. Therefore, the former wrongly reconstructed SNP alleles would not affect the latter reconstructed SNP alleles, i.e., reconstructed errors on the former SNP alleles would not spread to the latter ones. The kind of enumeration strategy can also be apply to reconstruct haplotypes of other ploidy, which will be studied in the future.
Compared with algorithms HapCompass and HapTree by using two kinds of simulated sequencing data, the EHTLD algorithm can get the highest reconstruction rates, the smallest vector errors and MEC scores, which was tested by a number of experiments. In addition, algorithm EHTLD still achieves satisfying performance even with high error rate, low fragment coverage, or long haplotype length. All of these advantages may contribute to the practical application of the EHTLD algorithm when haplotyping a triploid single individual.
References
Geraci F. A comparison of several algorithms for the single individual SNP haplotyping reconstruction problem. Bioinformatics. 2010;26(18):2217–25.
Wu JL, Liang BB. A fast and accurate algorithm for diploid individual haplotype reconstruction. J Bioinform Comput Biol. 2013;11(4):1350010.
Clark AG. Inference of haplotypes from PCRamplified samples of diploid populations. Mol Biol Evol. 1990;7:111–22.
Gusfield D. Inference of haplotypes from samples of diploid populations: complexity and algorithms. J Comput Biol. 2001;8:305–23.
O’Neil ST, Emrich SJ. Haplotype and minimumchimerism consensus determination using short sequence data. BMC Genomics. 2012;13(Suppl 2):S4.
Lancia G, Bafna V, Istrail S, Lippert R, Schwartz R. SNPs problems, complexity and algorithms. In: Heide FM, editor. Proceeding on 9th European symposium on algorithms, Aarhus, Denmark; 2001. p. 182–93.
Chen ZZ, Deng F, Wang LS. Exact algorithms for haplotype assembly from wholegenome sequence data. Bioinformatics. 2013;29(16):1938–45.
Mazrouee S, Wang W. FastHap: fast and accurate single individual haplotype reconstruction using fuzzy conflict graphs. Bioinformatics. 2014;30(ECCB):i371–8.
Wang RS, Wu LY, Li ZP, Zhang XS. Haplotype reconstruction from SNP fragments by minimum error correction. Bioinformatics. 2005;21(10):2456–62.
Qian WY, Yang YJ, Yang NN, Chun L. Particle swarm optimization for snp haplotype reconstruction problem. Appl Math Comput. 2008;196:266–72.
Li ZP, Wu LY, Zhao YY, Zhang XS. A dynamic programming algorithm for the khaplotyping problem. Acta Math Sinic English Series. 2006;22:405–12.
Aguiar D, Istrail S. Haplotype assembly in polyploid genomes and identical by descent shared tracts. Bioinformatics. 2013;29(13):i352–60.
Berger E, Yorukoglu D, Peng J, Berger B. HapTree: a novel Bayesian framework for single individual polyplotyping using NGS data. PLoS Comput Biol. 2014;10(3):e1003502.
Wu JL, Chen XX, Li XC. Haplotyping a single triploid individual based on genetic algorithm. BioMed Mater Eng. 2014;24(6):3753–62.
Xie MZ, Wang JX, Chen JE. A high accurate model of the individual haplotyping problem based on weighted SNP fragments and genotype with errors. Bioinformatics. 2008;24(ISMB):i105–13.
Myers G. A dataset generator for whole genome shotgun sequencing. In: Proceedings of the 7th international conference on intelligent systems for molecular biology; 1999. p. 202–10.
Richter DC, Ott F, Auch AF, Schmid R, Huson DH. MetaSim—a sequencing simulator for genomics and metagenomics. PLOS ONE. 2008;3(10):e3373.
Panconesi A, Sozio M. Fast hare: a fast heuristic for single individual SNP haplotype reconstruction. In: Proceedings of 4th workshop on algorithms in bioinformatics; 2004. p. 266–77.
Authors’ contributions
JW designed and implemented the algorithms and methods, QZ contributed on experimental design and data processing. JW wrote the most part of the manuscript. QZ helped in data preparing and modifying the manuscript. All the work was guided by JW in the whole process. Both authors read and approved the final manuscript.
Acknowledgements
The authors are grateful to anonymous referees for their helpful comments and to Professor Gene Myers for his kindly providing the source codes of CELSIM. This research is supported by Guangxi Collaborative Innovation Center of Multisource Information Integration and Intelligent Processing.
Competing interests
The authors declare that they have no competing interests.
Availability of data and materials
Not applicable.
Consent for publication
Not applicable.
Ethics approval and consent to participate
Not applicable.
Funding
This research is supported by the National Natural Science Foundation of China under Grant Nos. 61363035, 61762015 and 61502111, Guangxi Natural Science Foundation under Grant No. 2015GXNSFAA139288, Research Fund of Guangxi Key Lab of Multisource Information Mining & Security Nos. 14A0302 and 15A0302, “Bagui Scholar” Project Special Funds.
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Author information
Authors and Affiliations
Corresponding author
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Wu, J., Zhang, Q. A fast and accurate enumerationbased algorithm for haplotyping a triploid individual. Algorithms Mol Biol 13, 10 (2018). https://doi.org/10.1186/s1301501801290
Received:
Accepted:
Published:
DOI: https://doi.org/10.1186/s1301501801290