Skip to main content
  • Methodology article
  • Open access
  • Published:

Mapping quantitative trait loci in line cross with repeat records

Abstract

Background

Phenotypes with repeat records from one individual or multiple individuals were often encountered in practices of mapping QTL in linecross. The current genetic mapping method for a trait with repeat records is adopted by simply replacing the phenotype by the average value of the repeat records. This simple treatment has not sufficiently utilized the information from the replication and ignored the impacts of the permanent environmental effects on the accuracy of the estimated QTL.

Results

We propose to map QTL by using the repeatability model to directly analyze the repeat records rather than simply analyze the mean phenotype, improving the efficiency of QTL detecting because of adequately utilizing the information from data and allowing for the permanent environmental effects. A maximum likelihood method implemented via the expectation-maximization (EM) algorithm is applied to perform the parameter estimation of the repeatability model. The superiority of the mapping method based on the repeatability model over simple analysis using the mean phenotype was demonstrated by a series of simulations.

Conclusion

Our results suggest that the proposed method can serve as a powerful alternative to existing methods. By mean of the repeatability model, utilizing the repeat records on individual may improve the efficiency of QTL detecting in line cross.

Background

Replication is the fundamental of the experimental design, the important advantages of which are that it allows for an estimate of experimental error and increases the reliability of information obtained at each experimental point [1, 2]. Replication denotes sampling or measuring multiple times under the same experimental condition (within one treatment), where the experimental unit may be either one individual or multiple individuals with the identical genetic background.

Often plants or animals are observed more than once for a particular trait. For examples, fleece weight of sheep in different years, blood pressure and pulse of a human over time, litter size of sows over time, antler size of deer in different seasons, racing results of horses from several races, exam scores of students during university and so on. These records observed belong to replicate ones if they are not influenced by the measuring environments, such as the years, seasons, parities, races.

In classical quantitative genetics, a trait with repeat records is generally analysed by means of the repeatability model [3, 4], in which, there is an additional permanent environmental effect besides an individual's additive genetic value for a trait. The permanent environmental effect as a measure of the differences among experimental units, is a non-genetic effect common to all observations on the same individual [5]. Such environmental effects are usually accounted for in the model to ensure accurate prediction of breeding values [4]. However, the repeatability model has not been paid adequate attention to mapping QTL by using data with repeat records.

The current genetic mapping method for a trait with repeat records is adopted by simply replacing the phenotype by the average value of the repeat records [6, 7]. This simple treatment has not sufficiently utilized the information from the replication and ignored the impacts of the permanent environmental effects on the accuracy of the estimated QTL, although it enables to improve the power of detecting QTL with a certain extent.

In this study, we apply the repeatability model to mapping quantitative trait loci with repeat records and demonstrate the higher efficiency of this model by the simulations.

Theory and methods

Mapping QTL based on the mean phenotype

Take a simple F2 population of size n derived from two homozygous lines as an example. There are the three possible genotypes denoted by Q1Q1, Q1Q2, and Q2Q2, respectively, at a quantitative trait locus Q. The phenotypic value of an individual i is usually described by the following linear model,

y i = μ + z i a + w i d + e i ,

Where μ is the population mean, a and d are additive and dominant effects of the QTL, e i is the residual error with a N(0, σ2) distribution, and

z i = { + 1 f o r  Q 1 Q 1 0 f o r  Q 1 Q 2 − 1 f o r  Q 2 Q 2 and w i = { − 1 f o r  Q 1 Q 1 + 1 f o r  Q 1 Q 2 − 1 f o r  Q 2 Q 2 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqabeqadaaabaGaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGH9aqpdaGabeqaauaabeqadiaaaeaacqGHRaWkcqaIXaqmaeaacqWGMbGzcqWGVbWBcqWGYbGCcqqGGaaicqqGrbqudaWgaaWcbaGaeGymaedabeaakiabbgfarnaaBaaaleaacqaIXaqmaeqaaaGcbaGaeGimaadabaGaemOzayMaem4Ba8MaemOCaiNaeeiiaaIaeeyuae1aaSbaaSqaaiabigdaXaqabaGccqqGrbqudaWgaaWcbaGaeGOmaidabeaaaOqaaiabgkHiTiabigdaXaqaaiabdAgaMjabd+gaVjabdkhaYjabbccaGiabbgfarnaaBaaaleaacqaIYaGmaeqaaOGaeeyuae1aaSbaaSqaaiabikdaYaqabaaaaaGccaGL7baaaeaacqqGHbqycqqGUbGBcqqGKbazaeaacqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabg2da9maaceqabaqbaeqabmGaaaqaaiabgkHiTiabigdaXaqaaiabdAgaMjabd+gaVjabdkhaYjabbccaGiabbgfarnaaBaaaleaacqaIXaqmaeqaaOGaeeyuae1aaSbaaSqaaiabigdaXaqabaaakeaacqGHRaWkcqaIXaqmaeaacqWGMbGzcqWGVbWBcqWGYbGCcqqGGaaicqqGrbqudaWgaaWcbaGaeGymaedabeaakiabbgfarnaaBaaaleaacqaIYaGmaeqaaaGcbaGaeyOeI0IaeGymaedabaGaemOzayMaem4Ba8MaemOCaiNaeeiiaaIaeeyuae1aaSbaaSqaaiabikdaYaqabaGccqqGrbqudaWgaaWcbaGaeGOmaidabeaaaaaakiaawUhaaaaacqGGUaGlaaa@7FA1@

If m i records are repeatedly sampled from each individual and the phenotypic value of an individual i is measured by the average of m i records, the model is modified as

y ¯ i = μ + z i a + w i d + e ¯ i , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaOGaeyypa0dcciGae8hVd0Maey4kaSIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqWGHbqycqGHRaWkcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabdsgaKjabgUcaRiqbdwgaLzaaraWaaSbaaSqaaiabdMgaPbqabaGccqGGSaalaaa@41C7@
(2)

where

y ¯ i = 1 m i ∑ j = 1 m i y i j , e ¯ i = 1 m i ∑ j = 1 m i e i j , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaOGaeyypa0ZaaSaaaeaacqaIXaqmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaaGcdaaeWbqaaiabdMha5naaBaaaleaacqWGPbqAcqWGQbGAaeqaaaqaaiabdQgaQjabg2da9iabigdaXaqaaiabd2gaTnaaBaaameaacqWGPbqAaeqaaaqdcqGHris5aOGaeiilaWIafmyzauMbaebadaWgaaWcbaGaemyAaKgabeaakiabg2da9maalaaabaGaeGymaedabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaaaaOWaaabCaeaacqWGLbqzdaWgaaWcbaGaemyAaKMaemOAaOgabeaaaeaacqWGQbGAcqGH9aqpcqaIXaqmaeaacqWGTbqBdaWgaaadbaGaemyAaKgabeaaa0GaeyyeIuoakiabcYcaSaaa@5814@

and the variable with additional subscript j indicates the corresponding variable for the j th record of the i th F2 individual. The residual error now follows a N(0, σ2/m i ) distribution, given that e ij ~ N(0, σ2).

Let

f ( y ¯ i | θ , z i , w i ) = m i 2 π σ exp [ − m i 2 σ 2 ( y ¯ i − μ − z i a − w i d ) 2 ]
(3)

be the conditional density of y ¯ i MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaaaa@2FC6@ , where θ = [μ a d σ2]Tare the parameters; the log likelihood function defined under the missing variables z i and w i is

L ( θ ) = ∑ i = 1 n ln { E [ f ( y ¯ i | θ , z i , w i ) ] } MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGmbatcqGGOaakiiGacqWF4oqCcqGGPaqkcqGH9aqpdaaeWbqaaiGbcYgaSjabc6gaUbWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdGcdaGadeqaaiabdweafjabcUfaBjabdAgaMjabcIcaOiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaGccqGG8baFcqWF4oqCcqGGSaalcqWG6bGEdaWgaaWcbaGaemyAaKgabeaakiabcYcaSiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaeiykaKIaeiyxa0facaGL7bGaayzFaaaaaa@5301@
(4)

The expectation-maximization (EM) algorithm [8] can be used to obtain the MLE, as shown below,

[ ∑ i = 1 n m i ∑ i = 1 n m i E ( z i ) ∑ i = 1 n m i E ( w i ) ∑ i = 1 n m i E ( z i ) ∑ i = 1 n m i E ( z i 2 ) ∑ i = 1 n m i E ( z i w i ) ∑ i = 1 n m i E ( w i ) ∑ i = 1 n m i E ( z i w i ) ∑ i = 1 n m i E ( w i 2 ) ] [ μ a d ] = [ ∑ i = 1 n m i y ¯ i ∑ i = 1 n m i E ( z i ) y ¯ i ∑ i = 1 n m i E ( w i ) y ¯ i ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaadaWadaqaauaabeqadmaaaeaadaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaaqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGcbaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiabdweafjabcIcaOiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaeiykaKcaleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOqaamaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqWGfbqrcqGGOaakcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPaWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdaakeaadaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOGaemyrauKaeiikaGIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGcbaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiabdweafjabcIcaOiabdQha6naaDaaaleaacqWGPbqAaeaacqaIYaGmaaGccqGGPaqkaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGcbaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiabdweafjabcIcaOiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaem4DaC3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGcbaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiabdweafjabcIcaOiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaeiykaKcaleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOqaamaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqWGfbqrcqGGOaakcqWG6bGEdaWgaaWcbaGaemyAaKgabeaakiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaeiykaKcaleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOqaamaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqWGfbqrcqGGOaakcqWG3bWDdaqhaaWcbaGaemyAaKgabaGaeGOmaidaaOGaeiykaKcaleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaaaakiaawUfacaGLDbaadaWadaqaauaabeqadeaaaeaaiiGacqWF8oqBaeaacqWGHbqyaeaacqWGKbazaaaacaGLBbGaayzxaaGaeyypa0ZaamWaaeaafaqabeWabaaabaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaaabaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdaakeaadaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOGaemyrauKaeiikaGIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaaqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGcbaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiabdweafjabcIcaOiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaeiykaKIafmyEaKNbaebadaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaaaakiaawUfacaGLDbaaaaa@F7D2@
(5)

and

σ 2 = [ ∑ i = 1 n m i ] − 1 ∑ i = 1 n m i E [ ( y ¯ i − μ − z i a − w i d ) 2 ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaahaaWcbeqaaiabikdaYaaakiabg2da9maadmaabaWaaabCaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOGaay5waiaaw2faamaaCaaaleqabaGaeyOeI0IaeGymaedaaOWaaabCaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiabdweafbWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdGcdaWadaqaaiabcIcaOiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWF8oqBcqGHsislcqWG6bGEdaWgaaWcbaGaemyAaKgabeaakiabdggaHjabgkHiTiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaemizaqMaeiykaKYaaWbaaSqabeaacqaIYaGmaaaakiaawUfacaGLDbaaaaa@5EB0@
(6)

The expectation shown in Equation 6 can be further expressed as

E [ y ¯ i − μ − z i a − w i d ) 2 ] = ( y ¯ i − μ ) 2 + a 2 E ( z i 2 ) + d 2 E ( w i 2 ) − 2 ( y ¯ i − μ ) [ a E ( z i ) + d E ( w i ) ] + 2 a d E ( z i w i ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGabaaabaGaemyrauKaei4waSLafmyEaKNbaebadaWgaaWcbaGaemyAaKgabeaakiabgkHiTGGaciab=X7aTjabgkHiTiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaemyyaeMaeyOeI0Iaem4DaC3aaSbaaSqaaiabdMgaPbqabaGccqWGKbazcqGGPaqkdaahaaWcbeqaaiabikdaYaaakiabc2faDbqacmaaWdaaOfaamgGaaCzcaiabg2da9iabcIcaOiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWF8oqBcqGGPaqkdaahaaWcbeqaaiabikdaYaaakiabgUcaRiabdggaHnaaCaaaleqabaGaeGOmaidaaOGaemyrauKaeiikaGIaemOEaO3aa0baaSqaaiabdMgaPbqaaiabikdaYaaakiabcMcaPiabgUcaRiabdsgaKnaaCaaaleqabaGaeGOmaidaaOGaemyrauKaeiikaGIaem4DaC3aa0baaSqaaiabdMgaPbqaaiabikdaYaaakiabcMcaPiabgkHiTiabikdaYiabcIcaOiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWF8oqBcqGGPaqkcqGGBbWwcqWGHbqycqWGfbqrcqGGOaakcqWG6bGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabgUcaRiabdsgaKjabdweafjabcIcaOiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaeiykaKIaeiyxa0Laey4kaSIaeGOmaiJaemyyaeMaemizaqMaemyrauKaeiikaGIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPaaaaaa@8B4C@

Define the posterior probabilities of the three QTL genotypes for j th individual as

p i k ∗ = p i k f ( y ¯ i | θ , z i , w i ) ∑ k = 1 3 p i k f ( y ¯ i | θ , z i , w i ) ( k = 1 , 2 , 3 ) , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqabeqacaaabaGaemiCaa3aa0baaSqaaiabdMgaPjabdUgaRbqaaiabgEHiQaaakiabg2da9maalaaabaGaemiCaa3aaSbaaSqaaiabdMgaPjabdUgaRbqabaGccqWGMbGzcqGGOaakcuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaOGaeiiFaWhcciGae8hUdeNaeiilaWIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGGSaalcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPaqaamaaqadabaGaemiCaa3aaSbaaSqaaiabdMgaPjabdUgaRbqabaaabaGaem4AaSMaeyypa0JaeGymaedabaGaeG4mamdaniabggHiLdGccqWGMbGzcqGGOaakcuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaOGaeiiFaWNae8hUdeNaeiilaWIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGGSaalcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPaaaaeaacqGGOaakcqWGRbWAcqGH9aqpcqaIXaqmcqGGSaalcqaIYaGmcqGGSaalcqaIZaWmcqGGPaqkcqGGSaalaaaaaa@6E07@
(7)

where p ik are the conditional probabilities inferred by marker information, then

E ( z i ) = p i 1 ∗ − p i 3 ∗ , E ( z i 2 ) = p i 1 ∗ + p i 3 ∗ , E ( w i ) = p i 2 ∗ − p i 1 ∗ − p i 3 ∗ , E ( w i 2 ) = 1  and  E ( z i w i ) = p i 3 ∗ − p i 1 ∗ . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGabaaabaGaemyrauKaeiikaGIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcqGH9aqpcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcYcaSiabdweafjabcIcaOiabdQha6naaDaaaleaacqWGPbqAaeaacqaIYaGmaaGccqGGPaqkcqGH9aqpcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaey4kaSIaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcYcaSaqaaiabdweafjabcIcaOiabdEha3naaBaaaleaacqWGPbqAaeqaaOGaeiykaKIaeyypa0JaemiCaa3aa0baaSqaaiabdMgaPjabikdaYaqaaiabgEHiQaaakiabgkHiTiabdchaWnaaDaaaleaacqWGPbqAcqaIXaqmaeaacqGHxiIkaaGccqGHsislcqWGWbaCdaqhaaWcbaGaemyAaKMaeG4mamdabaGaey4fIOcaaOGaeiilaWIaemyrauKaeiikaGIaem4DaC3aa0baaSqaaiabdMgaPbqaaiabikdaYaaakiabcMcaPiabg2da9iabigdaXiabbccaGiabbggaHjabb6gaUjabbsgaKjabbccaGiabdweafjabcIcaOiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaem4DaC3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcqGH9aqpcqWGWbaCdaqhaaWcbaGaemyAaKMaeG4mamdabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabigdaXaqaaiabgEHiQaaakiabc6caUaaaaaa@8E26@
(8)

Because p i k ∗ MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGWbaCdaqhaaWcbaGaemyAaKMaem4AaSgabaGaey4fIOcaaaaa@31EB@ is a function of the unknown parameters, iterations are required for EM algorithm. The iterations are described as

Step 0: Set up initials for θ(0).

Step 1: Calculate the posterior probabilities p i k ∗ MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGWbaCdaqhaaWcbaGaemyAaKMaem4AaSgabaGaey4fIOcaaaaa@31EB@ with equation (7).

Step 2: Substituting (8) into equation (5), estimate

μ ( 1 ) = [ ∑ i = 1 n m i ] − 1 [ ∑ i = 1 n m i y ¯ i − a ( 0 ) ∑ i = 1 n m i ( p i 1 ∗ − p i 3 ∗ ) − d ( 0 ) ∑ i = 1 n m i ( p i 2 ∗ − p i 1 ∗ − p i 3 ∗ ) ] a ( 1 ) = [ ∑ i = 1 n m i ( p i 1 ∗ + p i 3 ∗ ) ] − 1 [ ∑ i = 1 n m i ( p i 1 ∗ − p i 3 ∗ ) y ¯ i − μ ( 1 ) ∑ i = 1 n m i ( p i 1 ∗ − p i 3 ∗ ) − d ( 0 ) ∑ i = 1 n m i ( p i 3 ∗ − p i 1 ∗ ) ] d ( 1 ) = [ ∑ i = 1 n m i ] − 1 [ ∑ i = 1 n m i ( p i 2 ∗ − p i 1 ∗ − p i 3 ∗ ) y ¯ i − μ ( 1 ) ∑ i = 1 n m i ( p i 2 ∗ − p i 1 ∗ − p i 3 ∗ ) − a ( 1 ) ∑ i = 1 n m i ( p i 3 ∗ − p i 1 ∗ ) ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeWabaaabaacciGae8hVd02aaWbaaSqabeaacqGGOaakcqaIXaqmcqGGPaqkaaGccqGH9aqpdaWadaqaamaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaaabaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdaakiaawUfacaGLDbaadaahaaWcbeqaaiabgkHiTiabigdaXaaakmaadmaabaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaakiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaaabaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdGccqGHsislcqWGHbqydaahaaWcbeqaaiabcIcaOiabicdaWiabcMcaPaaakmaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqGGOaakcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcMcaPiabgkHiTiabdsgaKnaaCaaaleqabaGaeiikaGIaeGimaaJaeiykaKcaaOWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoakiabcIcaOiabdchaWnaaDaaaleaacqWGPbqAcqaIYaGmaeaacqGHxiIkaaGccqGHsislcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcMcaPaWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdaakiaawUfacaGLDbaaaeaacqWGHbqydaahaaWcbeqaaiabcIcaOiabigdaXiabcMcaPaaakiabg2da9maadmaabaWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoakiabcIcaOiabdchaWnaaDaaaleaacqWGPbqAcqaIXaqmaeaacqGHxiIkaaGccqGHRaWkcqWGWbaCdaqhaaWcbaGaemyAaKMaeG4mamdabaGaey4fIOcaaOGaeiykaKcacaGLBbGaayzxaaWaaWbaaSqabeaacqGHsislcqaIXaqmaaGcdaWadaqaamaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqGGOaakcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcMcaPiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaaabaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdGccqGHsislcqWF8oqBdaahaaWcbeqaaiabcIcaOiabigdaXiabcMcaPaaakmaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqGGOaakcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcMcaPiabgkHiTiabdsgaKnaaCaaaleqabaGaeiikaGIaeGimaaJaeiykaKcaaOWaaabmaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoakiabcIcaOiabdchaWnaaDaaaleaacqWGPbqAcqaIZaWmaeaacqGHxiIkaaGccqGHsislcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeiykaKcaleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOGaay5waiaaw2faaaqaaiabdsgaKnaaCaaaleqabaGaeiikaGIaeGymaeJaeiykaKcaaOGaeyypa0ZaamWaaeaadaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaaqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGccaGLBbGaayzxaaWaaWbaaSqabeaacqGHsislcqaIXaqmaaGcdaWadaqaamaaqadabaGaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqGGOaakcqWGWbaCdaqhaaWcbaGaemyAaKMaeGOmaidabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabigdaXaqaaiabgEHiQaaakiabgkHiTiabdchaWnaaDaaaleaacqWGPbqAcqaIZaWmaeaacqGHxiIkaaGccqGGPaqkcuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaaqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aOGaeyOeI0Iae8hVd02aaWbaaSqabeaacqGGOaakcqaIXaqmcqGGPaqkaaGcdaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOGaeiikaGIaemiCaa3aa0baaSqaaiabdMgaPjabikdaYaqaaiabgEHiQaaakiabgkHiTiabdchaWnaaDaaaleaacqWGPbqAcqaIXaqmaeaacqGHxiIkaaGccqGHsislcqWGWbaCdaqhaaWcbaGaemyAaKMaeG4mamdabaGaey4fIOcaaOGaeiykaKIaeyOeI0Iaemyyae2aaWbaaSqabeaacqGGOaakcqaIXaqmcqGGPaqkaaGcdaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaaqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aOGaeiikaGIaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabgkHiTiabdchaWnaaDaaaleaacqWGPbqAcqaIXaqmaeaacqGHxiIkaaGccqGGPaqkaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaGccaGLBbGaayzxaaaaaaaa@6B66@
(9)

Step 3: Substituting (8) into equation (6), estimate

σ 2 ( 1 ) = [ ∑ i = 1 n m i ] − 1 { ∑ i = 1 n m i [ p i 1 ∗ ( y ¯ i − μ ( 1 ) − a ( 1 ) + d ( 1 ) ) 2 + p i 2 ∗ ( y ¯ i − μ ( 1 ) − d ( 1 ) ) 2 + p i 3 ∗ ( y ¯ i − μ ( 1 ) + a ( 1 ) + d ( 1 ) ) 2 ] } MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaahaaWcbeqaaiabikdaYiabcIcaOiabigdaXiabcMcaPaaakiabg2da9maadmaabaWaaabCaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOGaay5waiaaw2faamaaCaaaleqabaGaeyOeI0IaeGymaedaaOWaaiWabeaadaaeWbqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOWaamWaaeaacqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOWaaeWaaeaacuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaOGaeyOeI0IamGjG=X7aTnaaCaaaleqabaGaeiikaGIaeGymaeJaeiykaKcaaOGaeyOeI0Iaemyyae2aaWbaaSqabeaacqGGOaakcqaIXaqmcqGGPaqkaaGccqGHRaWkcqWGKbazdaahaaWcbeqaaiabcIcaOiabigdaXiabcMcaPaaaaOGaayjkaiaawMcaamaaCaaaleqabaGaeGOmaidaaOGaey4kaSIaemiCaa3aa0baaSqaaiabdMgaPjabikdaYaqaaiabgEHiQaaakmaabmaabaGafmyEaKNbaebadaWgaaWcbaGaemyAaKgabeaakiabgkHiTiab=X7aTnaaCaaaleqabaGaeiikaGIaeGymaeJaeiykaKcaaOGaeyOeI0Iaemizaq2aaWbaaSqabeaacqGGOaakcqaIXaqmcqGGPaqkaaaakiaawIcacaGLPaaadaahaaWcbeqaaiabikdaYaaakiabgUcaRiabdchaWnaaDaaaleaacqWGPbqAcqaIZaWmaeaacqGHxiIkaaGcdaqadaqaaiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWF8oqBdaahaaWcbeqaaiabcIcaOiabigdaXiabcMcaPaaakiabgUcaRiabdggaHnaaCaaaleqabaGaeiikaGIaeGymaeJaeiykaKcaaOGaey4kaSIaemizaq2aaWbaaSqabeaacqGGOaakcqaIXaqmcqGGPaqkaaaakiaawIcacaGLPaaadaahaaWcbeqaaiabikdaYaaaaOGaay5waiaaw2faaaWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdaakiaawUhacaGL9baaaaa@9BED@

Step 4: Go to step 1, which complete one round of iteration.

Mapping QTL based on the repeatability model

Partitioning residual error e i in model (1) into an individual-specific permanent environmental effect ζ i and random environmental effect ε ij , the j th phenotypic value of an individual i is represented as

y ij = μ + z i a + w i d + ζ i + ε ij

This is a mixed effects model, also called repeatability model, with a and d being treated as the fixed effects and p i as the random effect. i.i.d. N(0, σ ζ 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaidaaaaa@314D@ ) distribution and ε ij i.i.d. N(0, σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ ) distribution.

We use an m i × 1 vector y i = [yi 1yi 2… y im ]T, for n = 1, 2, …, n to denote the array of phenotypic values for the i th individual and define ϕ i = [1 1 … 1]Tas a vector of dimension m i . In matrix notation, model (9) can be written as

y i = ϕ i μ + z i ϕ i a + w i ϕ i d + ϕ i ζ i + ε i (11)

where ε i = [εi 1εi 2 … ε im ]Tis an m i × 1 vector for the random environmental effects which follows N(0, I i , σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ ) with I i being an (m i × 1) × (m i × 1) identity matrix. The conditional expectation of model (11) given the fixed effects is

E(y i ) = M i = ϕ i μ + z i ϕ i a + w i ϕ i d (12)

and the variance-covariance matrix is

V a r ( y i ) = V i = φ i φ i T σ ζ 2 + I i σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGwbGvcqWGHbqycqWGYbGCcqGGOaakcqWG5bqEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabg2da9iabdAfawnaaBaaaleaacqWGPbqAaeqaaOGaeyypa0dcciGae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWFgpGzdaqhaaWcbaGaemyAaKgabaGaemivaqfaaOGae83Wdm3aa0baaSqaaiab=z7a6bqaaiabikdaYaaakiabgUcaRiabdMeajnaaBaaaleaacqWGPbqAaeqaaOGae83Wdm3aa0baaSqaaiab=v7aLbqaaiabikdaYaaaaaa@4E9D@
(13)

which applies to all i = 1, 2, …, n.

The conditional density of y i based on M i and V i is

f ( y i | θ , z i , w i ) = 1 2 π | V i | 1 / 2 exp [ − 1 2 ( y i − M i ) T V i − 1 ( y i − M i ) ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGMbGzcqGGOaakcqWG5bqEdaWgaaWcbaGaemyAaKgabeaakiabcYha8HGaciab=H7aXjabcYcaSiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaeiilaWIaem4DaC3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcqGH9aqpdaWcaaqaaiabigdaXaqaamaakaaabaGaeGOmaiJae8hWdahaleqaaOWaaqWaaeaacqWGwbGvdaWgaaWcbaGaemyAaKgabeaaaOGaay5bSlaawIa7amaaCaaaleqabaGaeGymaeJaei4la8IaeGOmaidaaaaakiGbcwgaLjabcIha4jabcchaWnaadmaabaGaeyOeI0YaaSaaaeaacqaIXaqmaeaacqaIYaGmaaGaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWGnbqtdaWgaaWcbaGaemyAaKgabeaakiabcMcaPmaaCaaaleqabaGaemivaqfaaOGaemOvay1aa0baaSqaaiabdMgaPbqaaiabgkHiTiabigdaXaaakiabcIcaOiabdMha5naaBaaaleaacqWGPbqAaeqaaOGaeyOeI0Iaemyta00aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkaiaawUfacaGLDbaaaaa@6B5D@
(14)

where θ = [μ a d σ p 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGaemiCaahabaGaeGOmaidaaaaa@30FE@ σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ ]. Corresponding log-likelihood function defined is

L ( θ ) = ∑ i = 1 n ln { E [ f ( y i | θ , z i , w i ) ] } MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGmbatcqGGOaakiiGacqWF4oqCcqGGPaqkcqGH9aqpdaaeWbqaaiGbcYgaSjabc6gaUbWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdGcdaGadeqaaiabdweafjabcUfaBjabdAgaMjabcIcaOiabdMha5naaBaaaleaacqWGPbqAaeqaaOGaeiiFaWNae8hUdeNaeiilaWIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqGGSaalcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabc2faDbGaay5Eaiaaw2haaaaa@52E9@
(15)

With derivative for μ, a and d, we can obtain

[ ∑ i − 1 n φ i T V i − 1 φ i ∑ i = 1 n E ( z i ) φ i T V i − 1 φ i ∑ i = 1 n E ( w i ) φ i T V i − 1 φ i ∑ i = 1 n E ( z i ) φ i T V i − 1 φ i ∑ i = 1 n E ( z i 2 ) φ i T V i − 1 φ i ∑ i = 1 n E ( z i w i ) φ i T V i − 1 φ i ∑ i = 1 n E ( w i ) φ i T V i − 1 φ i ∑ i = 1 n E ( z i w i ) φ i T V i − 1 φ i ∑ i = 1 n E ( w i 2 ) φ i T V i − 1 φ i ] [ μ a d ] = [ ∑ i = 1 n φ i T V i − 1 y i ∑ i = 1 n E ( z i ) φ i T V i − 1 y i ∑ i = 1 n E ( w i ) φ i T V i − 1 y i ] ,
(16)

but the explicit equations for σ ζ 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaidaaaaa@314D@ and σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ can not be derived in the same way. Instead of above likelihood function, we construct the following likelihood function by using joint conditional density of y i ,

L ( θ ) = ∑ i = 1 n ln { E [ f ( y i | θ 1 , z i , w i ) g ( ζ i | σ ζ 2 ) ] } MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGmbatcqGGOaakiiGacqWF4oqCcqGGPaqkcqGH9aqpdaaeWbqaaiGbcYgaSjabc6gaUnaacmqabaGaemyrauKaei4waSLaemOzayMaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGG8baFcqWF4oqCdaWgaaWcbaacbaGae4xmaedabeaakiabcYcaSiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaeiilaWIaem4DaC3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcqWGNbWzcqGGOaakcqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcYha8jab=n8aZnaaDaaaleaacqWF2oGEaeaacqaIYaGmaaGccqGGPaqkcqGGDbqxaiaawUhacaGL9baaaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaaa@6075@
(17)

Where θ1 = [μ a d ζ i σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ ]

f ( y i | θ 1 , z i , w i ) = 1 2 π σ ε exp [ − 1 2 σ ε 2 ( y i − M i − φ i ζ i ) T ( y i − M i − φ i ζ i ) ] g ( ζ i | σ ζ 2 ) = 1 2 π σ ζ exp ( − ζ i 2 2 σ ζ 2 ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqabeGabaaabaGaemOzayMaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGG8baFiiGacqWF4oqCdaWgaaWcbaGaeGymaedabeaakiabcYcaSiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaeiilaWIaem4DaC3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcqGH9aqpdaWcaaqaaiabigdaXaqaamaakaaabaGaeGOmaiJae8hWdahaleqaaOGae83Wdm3aaSbaaSqaaiab=v7aLbqabaaaaOGagiyzauMaeiiEaGNaeiiCaa3aamWaaeaacqGHsisldaWcaaqaaiabigdaXaqaaiabikdaYiab=n8aZnaaDaaaleaacqWF1oqzaeaacqaIYaGmaaaaaOGaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWGnbqtdaWgaaWcbaGaemyAaKgabeaakiabgkHiTiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGae8NTdO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkdaahaaWcbeqaaiabdsfaubaakiabcIcaOiabdMha5naaBaaaleaacqWGPbqAaeqaaOGaeyOeI0Iaemyta00aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=z7a6naaBaaaleaacqWGPbqAaeqaaOGaeiykaKcacaGLBbGaayzxaaaabaGaem4zaCMaeiikaGIae8NTdO3aaSbaaSqaaiabdMgaPbqabaGccqGG8baFcqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaidaaOGaeiykaKIaeyypa0ZaaSaaaeaacqaIXaqmaeaadaGcaaqaaiabikdaYiab=b8aWbWcbeaakiab=n8aZnaaBaaaleaacqWF2oGEaeqaaaaakiGbcwgaLjabcIha4jabcchaWnaabmaabaGaeyOeI0YaaSaaaeaacqWF2oGEdaqhaaWcbaGaemyAaKgabaGaeGOmaidaaaGcbaGaeGOmaiJae83Wdm3aa0baaSqaaiab=z7a6bqaaiabikdaYaaaaaaakiaawIcacaGLPaaaaaaaaa@9B93@

With derivative for θ1, we obtain

[ ∑ i = 1 n m i ∑ i = 1 n m i E ( z i ) ∑ i = 1 n m i E ( w i ) m 1 ⋯ m n ∑ i = 1 n m i E ( z i ) ∑ i = 1 n m i E ( z i 2 ) ∑ i = 1 n m i E ( z i w i ) m 1 E ( z 1 ) ⋯ m n E ( z i ) ∑ i = 1 n m i E ( w i ) ∑ i = 1 n m i E ( z i w i ) ∑ i = 1 n m i E ( w i 2 ) m 1 E ( w 1 ) ⋯ m n E ( w i ) m 1 m 1 E ( z 1 ) m 1 E ( w 1 ) m 1 + σ ε 2 σ ζ 2 ⋯ 0 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ m n m n E ( z i ) m n E ( w i ) 0 ⋯ m n + σ ε 2 σ ζ 2 ] [ μ a d ζ 1 ⋯ ζ n ] = [ ∑ i = 1 n m i y ¯ i ∑ i = 1 n m i E ( z i ) y ¯ i ∑ i = 1 n m i E ( w i ) y ¯ i m 1 y ¯ 1 ⋯ m n y ¯ n ]
(18)

and

σ ε 2 = 1 ∑ i = 1 n m i ∑ i = 1 n E [ ( y i − φ i μ − z i φ i a − w i φ i d − φ i ζ i ) T ( y ˜ i − φ i μ − z i φ i a − w i φ i d − φ i ζ i ) ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaOGaeyypa0ZaaSaaaeaacqaIXaqmaeaadaaeWaqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaaqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaaakmaaqahabaGaemyrauKaei4waSLaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=X7aTjabgkHiTiabdQha6naaBaaaleaacqWGPbqAaeqaaOGae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWGHbqycqGHsislcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGaemizaqMaeyOeI0Iae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPmaaCaaaleqabaGaemivaqfaaOGaeiikaGIafmyEaKNbaGaadaWgaaWcbaGaemyAaKgabeaakiabgkHiTiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGae8hVd0MaeyOeI0IaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiabdggaHjabgkHiTiabdEha3naaBaaaleaacqWGPbqAaeqaaOGae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWGKbazcqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=z7a6naaBaaaleaacqWGPbqAaeqaaOGaeiykaKIaeiyxa0faleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaaa@8F4E@
(19)

Where

E [ ( y i − φ i μ − z i φ i a − w i φ i d − φ i ζ i ) T ( y i − φ i μ − z i φ i a − w i φ i d − φ i ζ i ) ] = ( y i − φ i μ ) T ( y i − φ i μ ) 2 + m i a 2 E ( z i 2 ) + m i d 2 E ( w i 2 ) − 2 ( y i − φ i μ ) T [ φ i a E ( z i ) + φ i d E ( w i ) ] + 2 m i a d E ( z i w i ) + m i E ( ζ i 2 ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGabaaabaGaemyrauKaei4waSLaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsisliiGacqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=X7aTjabgkHiTiabdQha6naaBaaaleaacqWGPbqAaeqaaOGae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWGHbqycqGHsislcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGaemizaqMaeyOeI0Iae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPmaaCaaaleqabaGaemivaqfaaOGaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=X7aTjabgkHiTiabdQha6naaBaaaleaacqWGPbqAaeqaaOGae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWGHbqycqGHsislcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGaemizaqMaeyOeI0Iae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabc2faDbqaceaa8gGaaCzcaiabg2da9iabcIcaOiabdMha5naaBaaaleaacqWGPbqAaeqaaOGaeyOeI0Iae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWF8oqBcqGGPaqkdaahaaWcbeqaaiabdsfaubaakiabcIcaOiabdMha5naaBaaaleaacqWGPbqAaeqaaOGaeyOeI0Iae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWF8oqBcqGGPaqkdaahaaWcbeqaaiabikdaYaaakiabgUcaRiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOGaemyyae2aaWbaaSqabeaacqaIYaGmaaGccqWGfbqrcqGGOaakcqWG6bGEdaqhaaWcbaGaemyAaKgabaGaeGOmaidaaOGaeiykaKIaey4kaSIaemyBa02aaSbaaSqaaiabdMgaPbqabaGccqWGKbazdaahaaWcbeqaaiabikdaYaaakiabdweafjabcIcaOiabdEha3naaDaaaleaacqWGPbqAaeaacqaIYaGmaaGccqGGPaqkcqGHsislcqaIYaGmcqGGOaakcqWG5bqEdaWgaaWcbaGaemyAaKgabeaakiabgkHiTiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGae8hVd0MaeiykaKYaaWbaaSqabeaacqWGubavaaGccqGGBbWwcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiabdggaHjabdweafjabcIcaOiabdQha6naaBaaaleaacqWGPbqAaeqaaOGaeiykaKIaey4kaSIae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWGKbazcqWGfbqrcqGGOaakcqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabc2faDjabgUcaRiabikdaYiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOGaemyyaeMaemizaqMaemyrauKaeiikaGIaemOEaO3aaSbaaSqaaiabdMgaPbqabaGccqWG3bWDdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabgUcaRiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOGaemyrauKaeiikaGIae8NTdO3aa0baaSqaaiabdMgaPbqaaiabikdaYaaakiabcMcaPaaaaaa@EDA4@
(20)

so, we can simply utilize existing mixed model EM algorithm to find the MLE of parameters [9]. Followings are the EM steps for the mixed model analysis.

Step 0: Initialize all parameters with values in their legal domain, denoted by θ(0).

Step 1: Compute the posterior probabilities of the three genotypes for each individual

p i k ∗ = p i k f ( y i | θ 1 , z i , w i ) ∑ k = 1 3 p i k f ( y i | θ 1 , z i , w i ) for  k = 1 , 2 , 3
(21)

Step 2: Compute all the expectations involved in the following maximization steps (same with the equation (8)).

Step 3: Find the posterior distribution of the random effect p i from equation (18). This posterior distribution turns out to be a mixture of three normal distributions with a mean

E ( ζ i ) = σ p 2 φ V − 1 [ y i − φ μ − ( p i 1 ∗ − p i 3 ∗ ) φ a − ( p i 2 ∗ − p i 1 ∗ − p i 3 ∗ ) φ d ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGfbqrcqGGOaakiiGacqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabg2da9iab=n8aZnaaDaaaleaacqWGWbaCaeaacqaIYaGmaaGccqWFgpGzcqWGwbGvdaahaaWcbeqaaiabgkHiTiabigdaXaaakmaadmaabaGaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWFgpGzcqWF8oqBcqGHsislcqGGOaakcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcMcaPiab=z8aMjabdggaHjabgkHiTiabcIcaOiabdchaWnaaDaaaleaacqWGPbqAcqaIYaGmaeaacqGHxiIkaaGccqGHsislcqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeyOeI0IaemiCaa3aa0baaSqaaiabdMgaPjabiodaZaqaaiabgEHiQaaakiabcMcaPiab=z8aMjabdsgaKbGaay5waiaaw2faaaaa@6CBE@
(22)

and a variance

V a r ( ζ i ) = σ ζ 2 − φ i V − 1 φ i T σ ζ 4 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGwbGvcqWGHbqycqWGYbGCcqGGOaakiiGacqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabg2da9iab=n8aZnaaDaaaleaacqWF2oGEaeaacqaIYaGmaaGccqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiabdAfawnaaCaaaleqabaGaeyOeI0IaeGymaedaaOGae8NXdy2aa0baaSqaaiabdMgaPbqaaiabdsfaubaakiab=n8aZnaaDaaaleaacqWF2oGEaeaacqaI0aanaaaaaa@4BD0@
(23)

Step 4: Update the population mean, additive effect and dominance effect by equation (16). The resulting equations are equivalent to equations (9) replacing m i with φ T V i − 1 φ MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGacaGaaiaabeqaaeqabiWaaaGcbaGaeqOXdy2aa0baaSqaaaqaaiabdsfaubaakiabdAfawnaaDaaaleaacqWGPbqAaeaacqGHsislcqaIXaqmaaGccqaHgpGzdaWgaaWcbaaabeaaaaa@35CB@ .

Step 5: Update the covariance matrix of the random effect

σ ζ 2 ( 1 ) = 1 n ∑ i = 1 n E ( ζ i 2 ) = 1 n ∑ i = 1 n [ V a r ( ζ i ) + E 2 ( ζ i ) ] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaiJaeiikaGIaeGymaeJaeiykaKcaaOGaeyypa0ZaaSaaaeaacqaIXaqmaeaacqWGUbGBaaWaaabmaeaacqWGfbqrcqGGOaakcqWF2oGEdaqhaaWcbaGaemyAaKgabaGaeGOmaidaaOGaeiykaKcaleaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoakiabg2da9maalaaabaGaeGymaedabaGaemOBa4gaamaaqadabaWaamWaaeaacqWGwbGvcqWGHbqycqWGYbGCcqGGOaakcqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPiabgUcaRiabdweafnaaCaaaleqabaGaeGOmaidaaOGaeiikaGIae8NTdO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkaiaawUfacaGLDbaaaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabd6gaUbqdcqGHris5aaaa@623D@
(24)

Step 6: Update the residual variance by equation (19)

σ 2 = [ ∑ i = 1 n m i ] − 1 { ∑ i = 1 n m i [ p i 1 ∗ ( y i − μ − a + d − φ i ζ i ) T ( y i − μ − a + d − φ i ζ i ) + + p i 2 ∗ ( y ¯ i − μ − d − φ i ζ i ) T ( y ¯ i − μ − d − φ i ζ i ) + p i 3 ∗ ( y ¯ i − μ + a + d − φ i ζ i ) T ( y ¯ i − μ + a + d − φ i ζ i ) ] } MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeGacaa4baaIcuaabmqadeaaaeaaiiGacqWFdpWCdaahaaWcbeqaaiabikdaYaaakiabg2da9maadmaabaWaaabCaeaacqWGTbqBdaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoaaOGaay5waiaaw2faamaaCaaaleqabaGaeyOeI0IaeGymaedaaOWaaiqabeaadaaeWbqaaiabd2gaTnaaBaaaleaacqWGPbqAaeqaaOWaamqaaeaacqWGWbaCdaqhaaWcbaGaemyAaKMaeGymaedabaGaey4fIOcaaOGaeiikaGIaemyEaK3aaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWF8oqBcqGHsislcqWGHbqycqGHRaWkcqWGKbazcqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=z7a6naaBaaaleaacqWGPbqAaeqaaOGaeiykaKYaaWbaaSqabeaacqWGubavaaGccqGGOaakcqWG5bqEdaWgaaWcbaGaemyAaKgabeaakiabgkHiTiab=X7aTjabgkHiTiabdggaHjabgUcaRiabdsgaKjabgkHiTiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGae8NTdO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkcqGHRaWkaiaawUfaaaWcbaGaemyAaKMaeyypa0JaeGymaedabaGaemOBa4ganiabggHiLdaakiaawUhaaaqaciaazbaanOGaaCzcaiaaxMaacqGHRaWkcqWGWbaCdaqhaaWcbaGaemyAaKMaeGOmaidabaGaey4fIOcaaOGaeiikaGIafmyEaKNbaebadaWgaaWcbaGaemyAaKgabeaakiabgkHiTiabeY7aTjabgkHiTiabdsgaKjabgkHiTiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGae8NTdO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkdaahaaWcbeqaaiabdsfaubaakiabcIcaOiqbdMha5zaaraWaaSbaaSqaaiabdMgaPbqabaGccqGHsislcqWF8oqBcqGHsislcqWGKbazcqGHsislcqWFgpGzdaWgaaWcbaGaemyAaKgabeaakiab=z7a6naaBaaaleaacqWGPbqAaeqaaOGaeiykaKcabiqbaKcaaeGcaKHcauIcaeCdcaWLjaGaaCzcaiaaxMaacaWLjaGaey4kaSYaaiGabeaadaWacaqaaiabdchaWnaaDaaaleaacqWGPbqAcqaIZaWmaeaacqGHxiIkaaGccqGGOaakcuWG5bqEgaqeamaaBaaaleaacqWGPbqAaeqaaOGaeyOeI0Iae8hVd0Maey4kaSIaemyyaeMaey4kaSIaemizaqMaeyOeI0Iae8NXdy2aaSbaaSqaaiabdMgaPbqabaGccqWF2oGEdaWgaaWcbaGaemyAaKgabeaakiabcMcaPmaaCaaaleqabaGaemivaqfaaOGaeiikaGIafmyEaKNbaebadaWgaaWcbaGaemyAaKgabeaakiabgkHiTiab=X7aTjabgUcaRiabdggaHjabgUcaRiabdsgaKjabgkHiTiab=z8aMnaaBaaaleaacqWGPbqAaeqaaOGae8NTdO3aaSbaaSqaaiabdMgaPbqabaGccqGGPaqkaiaaw2faaaGaayzFaaaaaaaa@D8C2@

Step7: Repeat from step 1 to step 6 until a certain convergence criterion is reached.

MLE of parameters in both model (2) and (10) are iteratively solved at specific location on chromosomes using EM algorithm and the QTL position and effects are determined by means of likelihood ratio statistics in chromosome or genome scanning.

Simulation studies

A series of simulation experiments were used to compare the efficiency and behaviour of two mapping methods based on the repeatability model with simple analysis using the mean phenotype for a trait with repeat records. We simulated a single chromosome of 100 cM long with 11 evenly spaced codominant markers for an F2 population with sample size n = 100 and a single QTL was put at position 25 cM (between markers 3 and 4). Under the null model, the QTL was assigned a value of zero for both the additive and dominance effects. The empirical critical values of likelihood ratio statistics for testing the presence of the QTL were obtained by simulating 1000 replicates. Under the alternative model, nonzero and equal additive and dominance effects were simulated. The simulations were replicated 100 times. Empirical power was calculated by counting the number of runs in which test statistics were greater than the critical values.

Factor considered include the QTL size, measured as the proportion of the phenotypic variance explained by the QTL (also called the QTL heritability), the number of replicates and σ ζ 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaidaaaaa@314D@ : σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ i.e the variance ratio of permanent environmental effect to random environmental effect. The QTL size was set at three levels: a = d = 0.265, 0.577, 0.943 correspond to the three levels of h2 = 0.05, 0.10, 0.20 respectively. The number of replicates was examined at five levels: m = 1, 3, 5, 10, 15, and σ ζ 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaidaaaaa@314D@ : σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ = 1:4, 2:3, 2.5:2.5, 3:2, 4:1, remaining σ ζ 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8NTdOhabaGaeGOmaidaaaaa@314D@ + σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ = 5.0.

The j th phenotypic value of individual i was simulated by using the repeatability model:

y ij = μ + z i a + w i d + ξ i σ ζ + η ij σ ε (25)

Where both ξ i and η ij are the random numbers from standard normal distribution.

The results of all simulations consistently show that under the same experimental condition, (1) using the repeatability model can significantly increase the statistical power of QTL detecting compared with simple analysis using the mean phenotype, (2) the position and effects of QTL, especially the proportion of phenotypic variance contributed by QTL were more accuracy estimated by using the repeatability model than using the genetic mapping model without permanent environmental effects to analyze mean phenotype. The superiority of the repeatability model over the simple analysis using the mean phenotype performs in evidence under the condition of the low QTL heritability.

The effects of number of replications on the efficiency and behaviour of the two methods were investigated only at variance ratio of permanent environmental effect to random environmental effect of 1:1. The results of simulations were listed in Table 1 and 2, respectively, by different mapping method. Notices that the simulated results at m = 1 (no replication) only correspond to the mapping method based on the mean phenotype for no solution by using the repeatability model. As expected, the statistical power of QTL detecting with replication is higher than no replication, based on either the mean phenotype or the repeatability model. The estimation of QTL parameters show a general tendency to improve as the number of replications increases.

Table 1 Effects of the number of replications on the mapping analysis based on the repeatability model
Table 2 Effects of the number of replications on the simple analysis using the mean phenotype

We have also investigated the impact of the variance ratio of permanent environmental effect to random environmental effect on differences in mapping performance between the two methods. The results of simulations fixing five replications were listed in Table 3. The difference in variance between permanent environmental effect and random environmental effect is greater under fixing total variance of random effects, the superiority of the mapping method based on the repeatability model over the mean phenotype is clearer in the statistical power of QTL detecting. The possible reasons are that either the large variance of random environmental effect made reliability of the individual's mean phenotype value low or the variance of residual error in model (2) increases with the variance of permanent environmental effect increased.

Table 3 Comparisons of the mapping analysis based on the repeatability model with the simple analysis using the mean phenotype under the conditions of different the variance ratios of permanent environmental effects to random environmental effects

Discussion

For a trait with repeat records, we proposed use of the repeatability model to map QTL, which distinguishes from simple analysis using the mean phenotype not only in the data analyzed but essentially in the model adopted. Simple analysis using the mean phenotype was based on regular genetic model for mapping QTL in linecross, which excluded the permanent environmental effects. The excluded permanent environmental effects were deposited to the residual error, decreasing the accuracy of estimation for QTL parameters, which was strictly proved in the relevant books to statistic models [e.g., [10, 11]]. Of course, the loss of data information has also influenced the performance of mapping QTL based on the mean phenotype.

Replication required either the experimental conditions must be the same when multiple records were observed only from one individual or the genetic backgrounds must be the identical for each individual while those records were from multiple individuals. If the former was not satisfied, then such "repeat" records observed became longitudinal data, such as test-day records of milk production and body weight in cattle, were genetically analysed using the random regresion model which is essentially the repeatability model nested submodels of time [12–14]. Besides cloned individuals and progencies from each plant in RIL, the later was hard to be satisfied. For example, there were incompletely same genetic backgrounds among individuals within a family and F3 progenies from one F2 individual. To improve the efficiency of detecting QTL using such data, the genetic backgrounds should be at least taken into account in the analysis [7], furthermore, the repeatability model may be a good choice for directly analyzing such "repeat" records.

Although we demonstrate the statistical method of QTL mapping using a F2 population as an example, other more simple or complex designs, such as backcross population and full-sib family can also be extended. Assuming only one QTL in the model considered here is to conveniently investigate efficiency of presented method based on various estimates. If a trait is controlled by multiple loci, the composite interval mapping [15, 16] or Bayesian mapping [e.g., [17, 18]] will be proposed for mapping those QTLs by incorporating marker-cofactors outside the scanning interval or all the QTLs into the model (9).

References

  1. Fisher RA: The design of experiments. 1971, New York, Hafner Publishing Company, 9

    Google Scholar 

  2. Steel RGD, Torrie JH: Principles and procedures of statistics: a biometrical approach. 1980, Tokyo, McGraw-Hill Kogakusha, 2

    Google Scholar 

  3. Henderson CR: Applications of Linear Models in Animal Breeding. 1984, Guelph ON Univ of Guelph

    Google Scholar 

  4. Mrode RA: Linear Models for the Prediction of Animal Breeding Values. 1996, UK, CAB International

    Google Scholar 

  5. Falconer DS: Introduction to Quantitative Genetics. 1960, London,Oliver & Boyd

    Google Scholar 

  6. Zhang TY, Yuan J, Yu W, Guo Z, Kohel RJ: Molecular tagging of a major QTL for fiber strong in upland cotton and its marker-assisted selection. Theor Appl Genet. 2003, 106: 262-268.

    CAS  PubMed  Google Scholar 

  7. Zhang Y, Xu S: Mapping Quantitative Trait Loci in F2 Incorporating Phenotypes of F3 Progeny. Genetics. 2004, 166: 1981-1993. 10.1534/genetics.166.4.1981.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  8. Dempster AP, Laird NM, Rubin DB: Maximum likelihood from incomplete data via EM algorithm. J R Stat Soc Ser B. 1977, 39: 1-38.

    Google Scholar 

  9. Henderson CR: Recent developments in variance and covariance estimation. J Anim Sci. 1986, 63: 208-216.

    Google Scholar 

  10. Zar JH: Biostatistical Analysis. 1996, Prentice Hall, 3

    Google Scholar 

  11. Neter J, Kutner MH, Nachtsheim CJ, Wasserman W: Applied Linear Statistical Models. 1996, RD Irwin, Homewood, IL, 4

    Google Scholar 

  12. Henderson CR: Analysis of covariance in the mixed model: Higher level, no homogenous, and random regressions. Biometrics. 1982, 38: 623-640. 10.2307/2530044.

    Article  PubMed  Google Scholar 

  13. Schaeffer LR: Application of random regression model in animal breeding. Livest Prod Sci. 2004, 86: 35-45. 10.1016/S0301-6226(03)00151-9.

    Article  Google Scholar 

  14. Macgregor S, Knott SA, White I, Visscher PM: Quantitative trait locus analysis of longitudinal quantitative trait data in complex pedigrees. Genetics. 2005, 171: 1365-1376. 10.1534/genetics.105.043828.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  15. Jansen RC: Controlling the type I and type II errors in mapping quantitative trait loci. Genetics. 1994, 138: 871-881.

    PubMed Central  CAS  PubMed  Google Scholar 

  16. Zeng ZB: Precision mapping of quantitative trait loci. Genetics. 1994, 136: 1457-1468.

    PubMed Central  CAS  PubMed  Google Scholar 

  17. Satagopan JM, Yandell BS, Newton MA, Osborn TC: A Bayesian approach to detect quantitative trait loci using Markov chain Monte Carlo. Genetics. 1996, 144: 805-816.

    PubMed Central  CAS  PubMed  Google Scholar 

  18. Yi N, Xu S: Bayesian mapping of quantitative trait loci under complicated mating designs. Genetics. 2001, 157: 1759-1771.

    PubMed Central  CAS  PubMed  Google Scholar 

Download references

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Runqing Yang.

Additional information

Authors' contributions

RQY coordinated the study, developed the foundational principle of the method and wrote the computing program and the paper. FM was responsible for the simulation experiment and carried out the analysis of results.

Rights and permissions

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.

Reprints and permissions

About this article

Cite this article

Yang, R., Fang, M. Mapping quantitative trait loci in line cross with repeat records. BMC Genet 8, 47 (2007). https://doi.org/10.1186/1471-2156-8-47

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/1471-2156-8-47

Keywords