New statistical methods for estimation of recombination fractions in F2 population.

BACKGROUND: Dominant markers in an F2 population or a hybrid population have much less linkage information in repulsion phase than in coupling phase. Linkage analysis produces two separate complementary marker linkage maps that have little use in disease association analysis and breeding. There is a need to develop efficient statistical methods and computational algorithms to construct or merge a complete linkage dominant marker maps. The key for doing so is to efficiently estimate recombination fractions between dominant markers in repulsion phases. RESULT: We proposed an expectation least square (ELS) algorithm and binomial analysis of three-point gametes (BAT) for estimating gamete frequencies from F2 dominant and codominant marker data, respectively. The results obtained from simulated and real genotype datasets showed that the ELS algorithm was able to accurately estimate frequencies of gametes and outperformed the EM algorithm in estimating recombination fractions between dominant loci and recovering true linkage maps of 6 dominant loci in coupling and unknown linkage phases. Our BAT method also had smaller variances in estimation of two-point recombination fractions than the EM algorithm.
CONCLUSION: ELS is a powerful method for accurate estimation of gamete frequencies in dominant three-locus system in an F2 population and BAT is a computationally efficient and fast method for estimating frequencies of three-point codominant gametes.

Background

A great advance has been made in building genetic maps of various species due to the development of large-scale molecular marker technologies [1-7] and statistical methods [4, 8–18]. However, mapping of numerous molecular markers has been complicated by linkage phases of dominance [14–16, 19]. In two-point analysis, markers in repulsion phase provide quite less linkage information than in coupling phase [14, 15, 20, 21]. This is especially true for dominant markers in F2 population [14]. In practical mapping experiments, although the linkage phase for each dominant marker is random, a half of markers are derived from one of two coupling phases. The phase between couplings is repulsion [14, 15]. This situation results in two separate partner linkage maps for dominant markers: high linkage information content of markers in the coupling phase and low linkage information content of markers in the repulsion phase. Thus one has to build two complementary linkage maps [14, 15, 21, 22]. To date, there has not yet been an effective way to integrate both into a complete map. Mester et al. [15] attempted to use pairs of codominant and dominant (CD) markers to merge such two complementary maps because pairs of the CD markers in repulsion phase have much higher linkage information content than pairs of dominant-only markers in repulsion phase. However, this strategy demands that all dominant markers be paired with codominant markers, which is not a general case in mapping practice, otherwise, local and global disturbance will then violently affect the reliability of the integrated map.

The two-point analysis implemented by the expectation maximization (EM) algorithm [11–13, 23–25] is a highly powerful approach to estimate recombination fractions between codominant loci and between dominant loci in coupling phase, but the EM algorithm has very low power in estimation of recombination fractions between dominant loci in repulsion phase. This is because it is difficult for the EM algorithm to distinguish genotypes in coupling phase from those in repulsion phase for dominant markers.

Therefore, the key of developing a powerful method for mapping dominant loci in an intersection population is to overcome the difficulty of distinguishing coupling phase from repulsion phase. Since two-point analysis, as pointed out above, performs very poorly in the estimation of recombination fractions between dominant loci, three-point analysis is alternatively taken into account. However, few three-point EM algorithms can be applied to dominant markers because dominant markers are less informative for maximum likelihood estimation [26]. One effective way to carry out three-point analysis is to dissect three-point genotypes into various gamete components that are informative for distinction between coupling and repulsion phases, and then, to estimate their frequencies. With these estimated gamete frequencies, one can immediately estimate recombination fractions between dominant loci in couple and repulsion phases. A key to this strategy is to obtain estimate of gamete frequencies. On the basis of dissection of genotypes, Tan and Fu proposed a binomial analysis of three-point (BAT) to estimate frequencies of dominant gametes [19]. However, this binomial approach is limited to the frequency of the three-point recessive gamete abc. The accuracy of estimation is completely dependent on the observed frequency of its phenotype (aabbcc). We have developed a new method called “expectation least square” (ELS) to address this problem. ELS estimation, similarly to expectation maximum algorithm, is realized on the basis of Tan and Fu’s BAT method [19]. That is, the expectation of phenotype frequencies can be given by using Eqs. (1-9) in the BAT of Tan and Fu [19], and the difference between estimated and expected values of phenotype frequencies is given using least square. The expectation and least square steps are iterated so that the difference between estimated and expected values is less than tolerant value. In addition, we have also developed a fast binomial approach to estimate frequencies of codominant gametes.

Methods

Real data collection

Mouse genotype data: A RFLP dataset of 333 F2 mice was obtained from MAPMAKER/EXP (version 3.0b) [13].

Simulation

For dominant loci, we just took unknown phase into account in simulation and followed a point process model [27] and scheme of Tan and Fu [19] to perform simulations. In N F1 meioses, recombination events occurred at random between two adjacent loci. Here for the simplicity, we allowed for only independent crossovers during procedure of recombination occurrence between nonsister chromatides. We generated N F2 individuals with ratio = phenotype A: phenotype a = 3:1 at each dominant locus or A(homozygote): H(heterozygote): B(homozygote) = 1:2:1 at each codominant locus. We set three levels for sample size: N = 100, 200, and 300 F2 individuals and 100 iterations and used variance (equivalent to mean square error, MSE) that quantifies deviation of estimated recombination fraction between two adjacent loci from its true value to evaluate these estimators. Since the ELS and BAT estimators work in three-point system, three-point recombination fractions were incorporated to two-point recombination fractions by using Tan and Fu [19] method. Simulation of codominant and dominant F2 populations and the ELS and BAT estimations of gamete frequencies in F2population were implemented by our R functions (Additional file 1, source code).

Results

Estimation of the frequencies of three-locus gametes in an F2 population

Since our ELS method for accurate estimation of the frequencies of three-locus gametes in a population with random union of gametes is based on dissection of phenotypes, for convenience, we start by presenting the BAT method of Tan and Fu [19].

ELS estimation of frequencies of dominant marker gametes

Our study here is restricted to three biallelic dominant markers. We use A and a, B and b, C and c to represent two alleles at three loci where upper letters (A, B and C) stand for dominant alleles and lower letters (a, b and c) for recessive alleles. A triple-heterozygote individual via meiosis produces eight types of gametes at the three loci: ABC, ABc, Abc, AbC, aBC, abC, aBc and abc. Gametes ABC and abc are a pair of sister gametes on which two alleles at the all three loci are different and come from two different parents. Similarly, Abc and aBC, abC and ABc, AbC and aBc are also pairs of sister gametes. Two sister gametes theoretically have equal frequency in an F2population because no mutation, no migration, no gene conversion and no selection occur in such a random mating population. From the expectation that sister-gametes have equal frequencies, we have in an F2 population f(ABC) = f(abc) = q 1, f(ABC) = f(aBC) = q 2, f(ABc) = f(aBC) = q 3, f(AbC) = f(aBc) = q 4. These gamete frequencies are constrained by 2q 1 + 2q 2 + 2q 3 + 2q 4 = 1. The individuals in the population can be classified into four categories: category 0 in which all individuals possess 0 dominant locus, that is, all individuals have three recessive loci; categories 1, 2 and 3 in which all individuals have respectively only one, two and three homozygous or heterozygous dominant loci. To accurately estimate gamete frequencies, we dissect a phenotype into different zygote types (genotypes) in each category using sister gametes. In category 1, for example, aabbC_ has only locus c with one or two dominant alleles. Therefore it can be dissected into three zygote types:

Phenotypes aaB_cc and A_bbcc are dissected in a similar fashion. Category 2 also has three phenotypes and each of them can be dissected into four zygote types that are comprised of five pairs of sister gametes. For instance, phenotype type A_B_cc can be dissected into

Category 3 has only one phenotype. The phenotype is comprised of 8 zygote types (genotypes) and therefore it is not useful for estimate of gamete frequencies. We use Q 1, Q 2, Q 3, Q 4, Q 5, Q 6, and Q 7 to respectively represent the frequency expectations of phenotypes aabbcc, aabbC_, aaB_cc, A_bbcc, A_B_cc, A_bbC_, and aaB_C_ in a population. The frequency of phenotype aabbcc is

The other 6 phenotypes have their frequencies:

Using Q = 2 (q 2 q 3 + q 2 q 4 + q 3 q 4), Eq. (4) is simplified as

Estimates of q 1 ,  …  , q 4 can be obtained from the above sets of equations by replacing Q with their observed frequencies where k = 1, 2,…,7 for 7 phenotypes. Theoretically, eqs. (1) and (3) are sufficient to make solutions for the frequencies of four types of gametes. However, Eq. (5) can be used to further minimize noise in the observed frequencies. That is, Q 2, Q 3,and Q 4 can be alternatively estimated aswhere Q = Q 5 + Q 6 + Q 7 + Q 1 − 0.25 [19]. It implicates that Q 2, Q 3, and Q 4 can also be estimated from the estimated frequencies of Q 1, Q 5, Q 6, and Q 7. Thus, we can combine the two sets of estimates of Q 2, Q 3, and Q 4 into one set:where a and b are weights of and , respectively, where k = 2, 3, and 4. and are respectively estimates of Q and . In general case, a = b (see Additional file 3: Appendix B). An alternative method for weighting is and b  = 1 − a . When the sample is small, it is likely that 0 or 0. In such a case, one can set a = 1 and b = 0 for 0, or a = 0 and b = 1 for 0 and = 0. Since , q 3 can be given by

Similarly,

Q 1, Q 2, Q 3, and Q 4 are respectively estimated by , , , , therefore q 3, q 2, q 4, and q 1 are respectively estimated by

In Eq. (9), accurate estimation of q 1 is a key contribution to accurate estimations of q 2, q 3, and q 4. Equations (3) and (4) show that Q 2 ~ Q 7 can also provide information of solution to q 1. But it is impossible to directly obtain a solution for q 1 from Q 2 ~ Q 7. To estimate q 1 from Q 1 ~ Q 7, we here proposed a seeking method, named “expectation least square” (ELS) method.

Similar to the EM method [11, 25, 28, 29], the ELS method also consists of two steps. The first step is the expectation step, denoted by E-step, and the second step is the least-square step, denoted by LS-step. q 1 is initialized to be . We use to estimate q 2, q 3, and q 4 and get , , and from Eqs. (9). Then, we calculate the expected values of Q 2 ~ Q 7 from Eqs. (3) ~ (4) with , , and . At iteration j, we realize E-step and LS-step to get , , and :

E-step:

Calculate the expected values ~ of Q 2 ~ Q 7 by replacing , , , and into Eqs. (3) ~ (4) where , , and are obtained by wherewhere i = 2 , …, 4 and where .

LS-step:

Calculate square value using

Note that is a value we want to seek for, therefore, Eq. (10) does not contain . As it is very difficult to directly get solutions for these four q-values from the derivative approach, we use an iteration approach to minimize square value:

Use to calculate , , and where j is the jth iteration, j = 1 , …, and Δ is specified with a very small value. Here our algorithm to realize LS-step is

If , then

if , then ,

otherwise,

else if , then

if , then ,

otherwise, .

Note that there are not and in this algorithm. The iteration will stop at where t is a given tolerant value. Once the final estimate () of q 1 is found at a given tolerant value where j = f, the final estimates of q 2, q 3, and q 4 are obtained. Then we let , , , and .

BAT for estimation of the frequencies of codominant marker gametes in F2 population

To avoid confusing notations in codominant loci with those in dominant loci, we let 0 and 1 code for homozygote from two parents, respectively, and 2 code for heterozygote at a locus. Since homozygote and heterozygote at three loci can be recognized, most of zygotes are informative for estimation of the frequencies of four pairs of sister gametes. We still assume that the sister-gametes have equal frequencies, that is, q 1 = f(111) = f(000), q 2 = f(100) = f(011), q 3 = f(110) = f(001), q 4 = f(101) =f(010) in F2 population. Here these complementary zygote type pairs are listed as follows:

Let P 1, P 2, P 3 and P 4 represent the frequencies of complementary homozygote types (111/000), (100/011), (110/001), and (101/010) in each of which all three loci are homozygous; let P 12, P 13, P 14, P 23, P 24, and P 34 be the frequencies of complementary two-locus homozygote types (200/211), (002/112), (121/020), (021/120), (102/012), and (201/210) in each of which only one locus are heterozygous and let P 1234, P 1324, P 1423 be the frequencies of complementary one-locus homozygote types (122/022), (221/220) and (212/202) in each of which two loci are heterozygous. Then, , , , , P 12 = 4q 1 q 2, P 13 = 4q 1 q 3, P 14 = 4q 1 q 4, P 23 = 4q 2 q 3, P 24 = 4q 2 q 4, P 34 = 4q 3 q 4, P 1234 = 4q 1 q 2 + 4q 3 q 4, P 1324 = 4q 1 q 3 + 4q 2 q 4, P 1423 = 4q 1 q 4 + 4q 2 q 3. From the zygote type pair list above, we find that the frequencies of these 12 pairs of zygote types can constitute two sets of 6 binomial equations:

We use arithmetic mean to get frequencies of these zygote types in F2 population:where and b  = 1 − a . +b (q  + q )2= (a  + b ) (q  + q )2 = (q  + q )2 where i and j are gamete types i and j (i = 1, 2, 3 and j = 2, 3, 4 and i ≠ j). Thus, the frequencies of four types of non-sister gametes in a codominant three-locus system in an F2population are easily and fast estimated by where and are respective estimates of Q and P in F2 population where k = 1,…,4 denote gamete types 1, …, 4.

A modified BAT method (BAT II) for estimating the frequencies of eight gamete types without assumption that the sister gametes have equal frequencies in any generation population is given in Additional file 2, Appendix A.

Estimation of recombination fractions

Since these four qs are estimated separately, sum of them does not always satisfy a constraint of . For this reason, we normalize our estimates as

For three linked loci, the frequencies of the four gamete pairs can be used to find the double crossover types by distinguishing coupling phase from repulsion phase between loci. For example, for an order a-b-c of the three loci a, b and c, p 4 is determined to be the frequency of double crossover types if its value is the smallest and/or p 1 is the largest, which are produced at three coupling loci or p 1 is found to be the frequency of double crossover types if its value is the smallest and/or p 4 is the largest, which are formed at loci a and c in coupling phase and locus b in repulsion phase. In a similar way, we can also define p 3 or p 2 as the frequency of double crossover types.

If p 4 is frequency of double crossover types, then the recombination fractions between loci a and b, between loci b and c, and between loci a and c can be estimated by

For the linkage orders a-c-b and b-a-c, the recombination fractions between loci are also estimated in a similar way.

In the repulsion phase, the linkage a-b-c order of three loci determines p 1 to be the frequency of double crossover types, so estimates of recombination fractions between loci a and b, between loci b and c, and between loci a and c are

For the linkage orders b-a-c and a-c-b, the recombination fractions between three loci in the repulsion phase can be estimated in this way.

r  , r  , and r are simple notations of three recombination fractions in a triple. However, when n markers on a chromosome or a fragment are genotyped, it is difficult to use these notations of three recombination frequencies to denote recombination fractions in n(n − 1)(n − 2)/6 triples. To notate recombination fractions in multiple triples, we let r  = r where c is referred to as a reference marker for recombination fraction between markers a and b, r  = r where b as reference marker for that between loci a and c, and r  = r where a as reference locus for that between markers b and c, in a three-locus system consisting of markers a, b, and c [19]. In more general fashion, we denote i for the first marker, j for the second maker, and k for the last marker. Thus, the rest n − 2 markers are combined with loci i and j into n − 2 three-points, therefore, there are n − 2 estimates of the recombination fraction between markers i and j. Hence estimate of recombination fraction between loci i and j is given by Tan and Fu’s method [19]:

Practical examples

Here we used RFLP (restriction fragment length polymorphisms) data of 333 F2 mice from MAPMAKER/EXP (version 3.0b), LANDER et al. [13] to illustrate performances of our ELS and BAT methods to estimate recombination fractions between dominant and codominant loci. RFLP markers are codominant markers. In genotype data of 333 F2 mice, “A” stands for homozygote A (two alleles from parent A), “H” for heterozygote H (an allele from parent A and the other from parent B), and “B” for homozygote B (two alleles from parent B). We arbitrarily selected 6 codominant markers from the original genotype data. To evaluate our ELS algorithm, we converted the codominant genotype data into dominant genotype data by changing B to H. For convenience, we used arabic digits (1, 2,…,6) to label these six markers: marker 1, marker 2, …, marker 6. Sometime we also used locus 1, locus 2, …, locus 6 to mark these six marker loci. The frequencies of 20 non-sister gametes were estimated by respectively performing ELS on the dominant data and BAT on the codominant genotype data, normalized by using Eq. (16) and the results are summarized in Tables 1 and 2. For the ELS estimation, three non-sister gametes containing loci 4 and 6 (146, 246 and 346) fitted well the ratio of 1:1:1:1 (Chi-square test p-value >0.084, Table 1), indicating that loci 4 and 6 are unlinked to loci 1, 2 and 3. In addition, the frequencies of gametes 256, 356, and 345 also fitted the ratio of 1:1:1:1 with p-value ≥ 0.063 (Chi-square test, Table 1), but gametes 156, 245 and 145 had the ratios significantly deviating against 1:1:1:1 (Chi-square test p-value <0.0212, Table 1), we could infer that locus 5 was linked to loci 1 but independent of locus 3 and unascertained at locus 2. Thus, we definitely excluded loci 4 and 6 in the linkage. By using eqs. (17) – (19), the recombination fractions in four triples (123), (125), (135), and (235) were calculated by following the five given steps: the first step is to determine the linkage order of three loci in triple. For example, in triple (123), p 1 = f(abc) = 0.208668 is the largest value while p 2 = f(Abc) = 0.086162 is the smallest one, that is to say, gamete Abc is double crossover type and abc is parental type, so their order is 2(b)-1(a)-3(c). Step2 is to determine linkage phase: since gamete bac is parental type and bAc is double crossover type, gamete BAC or bac is couple phase. At step 3, we abstracted frequencies of gametes 123, 125, 135, 235 (Table 3) from Table 1. At step 4, recombination fractions between loci in a triple were estimated as

Table 1

The ELS estimated frequencies of four nonsister gametes in 20 triplets of 6 dominant loci in 333 F2 micea

locus			frequency of non-sister gamete				Chi-square test
a	b	c	p1 = f(abc)	p2 = f(Abc)	p4 = f(aBc)	p3 = f(abC)	p-value	ratio
1	2	3	0.208668	0.086162	0.094698	0.110472	0.000339
1	2	4	0.14751	0.087958	0.160866	0.103665	0.028
1	2	5	0.200976	0.080676	0.108494	0.109854	0.00092
1	2	6	0.192229	0.084408	0.126033	0.09733	0.0023
1	3	4	0.140566	0.079252	0.181795	0.098387	0.0038
1	3	5	0.209093	0.065783	0.12237	0.102753	0.00012
1	3	6	0.16895	0.063323	0.165648	0.102079	0.0011
1	4	5	0.173539	0.100396	0.079665	0.1464	0.0069
1	4	6	0.16482	0.102771	0.113819	0.118591	0.0837	1:1:1:1
1	5	6	0.173447	0.07932	0.148645	0.098588	0.0059
2	3	4	0.141958	0.088919	0.172784	0.096339	0.012
2	3	5	0.202566	0.085775	0.112098	0.099561	0.00084
2	3	6	0.139672	0.086032	0.16843	0.105866	0.0212
2	4	5	0.173634	0.117152	0.085607	0.123607	0.0212
2	4	6	0.10113	0.133668	0.133668	0.131535	0.3134	1:1:1:1
2	5	6	0.156859	0.096648	0.142758	0.103735	0.0634	1:1:1:1
3	4	5	0.143582	0.120278	0.098137	0.138002	0.1954	1:1:1:1
3	4	6	0.105025	0.139707	0.129181	0.126086	0.3581	1:1:1:1
3	5	6	0.140841	0.109125	0.152351	0.097682	0.1028	1:1:1:1
4	5	6	0.146847	0.096392	0.156905	0.099856	0.0476

a: The data came from MAPMAKER/EXP(3.0b) [27]

Table 2

The BAT estimated frequencies of nonsister gametes in 20 triplets of 6 codominant loci in 333 F2 micea

locus			frequency of non-sister gamete				Chi-square test
a	b	c	p1 = f(000)	p2 = f(100)	p3 = f(001)	p4 = f(010)	p-value	ratio
1	2	3	0.242929	0.066568	0.080146	0.110358	3.348e-07
1	2	4	0.145838	0.08977	0.162134	0.102258	0.0291
1	2	5	0.196094	0.091387	0.121051	0.091467	0.0017
1	2	6	0.17224	0.104184	0.143308	0.080268	0.0098
1	3	4	0.165099	0.068697	0.191983	0.074221	7.1334e-05
1	3	5	0.222931	0.079297	0.147828	0.049944	1.0943e-06
1	3	6	0.177615	0.065713	0.187929	0.068743	2.3780e-05
1	4	5	0.158699	0.103969	0.114759	0.122573	0.1336	1:1:1:1
1	4	6	0.165089	0.113628	0.139128	0.082155	0.0249
1	5	6	0.155874	0.091722	0.16263	0.089774	0.0121
2	3	4	0.142565	0.093943	0.179432	0.08406	0.0050
2	3	5	0.216411	0.069533	0.134853	0.079203	2.0404e-05
2	3	6	0.160337	0.092614	0.172787	0.074262	0.0020
2	4	5	0.172459	0.100044	0.105018	0.12248	0.0365
2	4	6	0.154284	0.140079	0.121173	0.084464	0.0559	1:1:1:1
2	5	6	0.167782	0.072156	0.153072	0.10699	0.0057
3	4	5	0.154314	0.118649	0.10895	0.118086	0.2051	1:1:1:1
3	4	6	0.144358	0.108635	0.131647	0.115359	0.3080	1:1:1:1
3	5	6	0.16738	0.053399	0.176627	0.102594	0.0002
4	5	6	0.153613	0.092081	0.124467	0.129838	0.1092	1:1:1:1

a: The data came from MAPMAKER/EXP(3.0b) [27]

Table 3

The ELS estimated frequencies of nonsister gametes in triplets of dominant loci 1, 2, 3 and 5 in 333 F2 mice

locus			frequency of gamete
a	b	c	p1 = f(abc)	p2 = f(Abc)	p3 = f(abC)	p4 = f(aBc)
1	2	3	0.208668	0.086162	0.094698	0.110472
1	2	5	0.200976	0.080676	0.108494	0.109854
1	3	5	0.209093	0.065783	0.12237	0.102753
2	3	5	0.202566	0.085775	0.112098	0.099561

Similarly, we also estimated the recombination fractions in triples (125), (135), and (235) (Table 3). Finally, the three-point estimates of the recombination fractions were incorporated into two-point estimates by applying Eq. (19) to the data in Table 4:

Table 4

The estimated recombination fractions between dominant loci in four triples

triple			Recombination fraction between loci
a	b	c	b-a	a-c	b-c
1	2	3	0.39326	0.36172	0.41034
1	2	5	0.38106	0.37834	0.43669
1	3	5	0.33707	0.37631	0.45024
2	3	5	0.37067	0.39575	0.42332

Table 2 displays frequencies of codominant gametes estimated by our BAT method. It is clear to see that frequencies of gametes 145, 246, 345, 346, and 456 fitted well ratio of 1:1:1:1 with p-value ≥ 0.0559 (Chi-square test), however, the frequencies of gametes 156, 256 and 356 did not fit the ratio of 1:1:1:1 with p-value < 0.0121 (Chi-square test, Table 2), inferring that loci 4 and 6 are unlinked to loci 1, 2 and 3 but locus 5 could not be inferred to linked to them. Again, in codominant genotype data, locus 5 was still unascertained. Following the steps above, we obtained estimates of recombination fractions between these four loci (Table 5). Both ELS estimates of recombination fractions between dominant loci and BAT estimates between codominant loci show that locus 5 could not be tightly linked to any one of loci 1, 2 and 3. Loci 1, 2 and 3 could be determined to have linkage order of 2-1-3. Simulation data also showed that the codominant estimator had higher precision than the dominant estimator (see Simulation data section), suggesting that codominant markers indeed contain higher linkage information than dominant ones.

Table 5

Comparison between two estimators of recombination fractions between markers

two loci		the ELS estimate in dominant genotype data	the BAT estimate in codominant genotype data
1	2	0.387164	0.271317
1	3	0.349396	0.275955
1	5	0.377323	0.439563
2	3	0.390506	0.339240
2	5	0.416221	0.426574
3	5	0.436782	0.402158

Simulation data

We performed simulation study to compare the two estimators of recombination fractions. We followed the simulation scheme of Tan and Fu [19]. Briefly, we set two linkage maps comprised of 6 dominant loci and 6 codominant loci, respectively. Five possible map distances 10, 15, 20, 25, and 30 cM (1 cM = 1%) were randomly assigned to the five adjacent intervals on these two linkage models with equal probability (see Methods for detail). The point process model [27] was used to generate F2 population. We did not consider recombination interference and linkage disequilibrium. Recombination fractions between adjacent loci in an unknown linkage phase (or say random phase) were estimated by the two-point EM [14, 23] and ELS estimators in 100 repeated samples of 100, 200, and 300 individuals drawn from the simulated F2 population. These two estimators were rated by the variance that quantifies deviation of estimated recombination fraction between two adjacent loci from its true value and is equivalent to mean squared error (MSE). For dominant markers, simulation shows that the ELS algorithm had much smaller variances in estimation of true recombination fractions between adjacent loci in samples of 100, 200 and 300 F2 individuals than two-point EM algorithm (Fig. 1). In Table 6, one can find that ELS had slightly higher probability of recovering true linkage maps of 6 loci than EM [14, 23] and BAT in the case of coupling phase and samples of 100 and 200 F2 individuals. When sample size reached 300 individuals, both ELS and EM recovered true coupling linkage maps with 100% probability and BAT also had 97.9% recovery rate. However, in unknown phase, ELS recovered true linkage maps of 6 loci with 23.4% probability in sample of 100 F2 individuals and reached 85% recovery rate in sample of 300 F2 individuals. By contrast, EM had very low recovery rate (23.4%) even when sample size was 300. Therefore, ELS performed much better than two-point EM algorithm in all given scenarios. An inexact comparison can be done between ELS and three-point EM algorithm of Lu et al. [30], Table 4 in Lu et al. showed that their three-point EM algorithms had 98.5% probability of finding the correct linkage map of three dominant markers in coupling phase from a sample of full-sib 100 individuals (corresponding to 100 F2 individuals), our ELS had 96.7% probability of recovering true linkage map of 6 dominant markers in coupling phase in 100 F2 individuals (Table 6). The probability to find a given linkage map will remarkably decrease as number of markers increases. So we can predict that the three-point EM algorithm would not have over 96.7% of the probability to find a given linkage map of 6 dominant markers. For the repulsion phase (or trans × trans), Lu et al.’s three-point EM algorithm had 99.5% probability of finding a correct linkage map of three markers in 100 full-sib individuals, which is higher than 98.6% in coupling phase. In theory, any EM algorithm should have much lower probability to find a given linkage order in repulsion phase than in coupling phase because the repulsion phase has much less linkage information content than the coupling phase [14, 26]. So, this result may be required to be confirmed in more simulations. Since Lu et al. did not implement simulation of random phase case and the repulsion phase is not random phase, the comparison cannot be made between the three-point EM and ELS algorithms in the random phase. For codominant markers, the BAT method performed with smaller variances than the two-point EM algorithm in the most cases. The results provided strong evidence for the conclusion that a method or algorithm based on three-point gametes can mitigate effect of low linkage information of repulsion phase on estimation of recombination fractions. Compared to the simulated results in Table 3 in [19], one can find that the ELS algorithm is better than the Tan and Fu’s BAT method. Table 3 in [19] showed that in case of unknown phase, the BAT method outperformed two-point EM.

Fig. 1

Variances of estimated recombination fractions between adjacent dominant loci in unknown linkage phase deviated from their respective true values Variance of estimated recombination fraction between adjacent dominant loci is given by simulating 100 estimates around true recombination fraction between adjacent loci. The variance here is equivalent to mean square error (MSE)

Table 6

Efficiencies of estimators of recombination fractions in recovering the true linkage maps of 6 dominant loci in the case of random distance

Estimator		Linkage phase	Sample size
			100	200	300
Two-point EM		CP	92.3	97.8	100.0
		UP	15.7	22.9	23.4
	ELS	CP	96.7	100.0	100.0
		UP	50.5	77.0	85.1
	BAT	CP	82.1	95.9	97.9
		UP	26.0	40.9	42.3

CP: Coupling phase and UP: unknown phase

Discussion

Accurate estimation of recombination fractions is a key for mapping multiple markers. Therefore, powerful method for estimating recombination fractions is required. For dominant loci, the EM and ML methods have been verified to have low power to estimate frequencies of recombination between loci in repulsion phase [14, 19]. This is because the EM method cannot distinguish dominant homozygous genotypes from dominant heterozygous genotypes.

Compared to the EM algorithm, the ELS algorithm based on Tan and Fu’s method [19] has small bias for estimating recombination fractions between dominant loci on a chromosome in a larger F2 population due to the following reasons: (a) gamete analysis can effectively distinguish marker linkage phases; (b) accurately estimate q 1, and (c) average of estimates of recombination fraction between two loci over all reference loci [Eq. (19)] effectively balances sampling error. Estimation of q 1 is restriction of the Tan and Fu’s method. We here proposed iteration expectation-least square algorithm (ELS) to seek for accurate q 1 estimation. This new algorithm is similar to expectation maximum algorithm and its statistical properties will be given by more simulation comparisons in elsewhere. In addition, importance for high efficiency of recombination fraction estimation is . ELS had much higher recovery rate by using than by using in both coupling and unknown phase (Table 6). Correlation analysis also indicated that indeed has the linkage behavior similar to (Additional file 3, Appendix B). Furthermore, we found that obtained from a data set of 100 simulated samples of 100 F2 individuals has remarkably smaller variance than (data not shown). To fully confirm that is the optimal choice in our ELS method, was taken into account where = if >0, otherwise, = . The simulated result showed that ~31% of linkage maps recovered true order of 6 dominant loci in samples of 100 F2 individuals, which is apparently lower than that by using . For this reason, we chose in our ELS algorithm. Besides the ELS algorithm, average of recombination fraction between two loci over all reference loci greatly reduces noise of recombination fractions.

BATII given in Additional file 2, Appendix A, can be used to estimate frequencies of 8 codominant gamete types in any nature population because it does not require the assumption that the sister gametes have equal frequencies in a population. However, its estimation accuracy is not higher than the first BAT method in F2 population because sister-gametes really have equal frequencies and two-locus heterozygote types are not useful in the BATII. In a natural population, for example, human population, the frequencies of these gametes are not purely derived from recombination events but may be due to selection, genetic drift, migration and mutation. If, however, sister gametes are found to be equal in statistics, then these frequencies can still be used to inference recombination fractions between loci and recombination inference.

Conclusions

Accurate estimation of recombination fractions between loci is given by methodologies developed for accurate estimation of gamete frequencies in a population. Analyses of simulated and real dominant and codominant data show that the ELS method proposed here is a powerful algorithm for accurate estimation of frequencies of gametes with unknown phase in dominant three-locus system in F2 population and BAT is a computationally efficient and powerful method for estimating frequencies of non-sister three-point codominant gametes.

Source code: three R functions: BAT.R, ELS.R, simulatF2.R. (ZIP 4 kb)

Appendix A. Binomial analysis of three-point method (BATII) is described in detail. BATII is used to estimate frequencies of sister gametes at codominant loci in natural populations. (DOCX 184 kb)

Appendix B. A proof of a proposition that equal weights of two datasets combined into a dataset have maximum linkage information and minimum error for linkage analysis is given. (DOCX 41 kb)