Multiple Interval Mapping (MIM) for QTL detection in a F2 Mimulus population

• In a previous exercise, we obtained the following linkage map for a Mimulus complex F2 population, derived from an interspecific cross between Mimulus guttatus and Mimulus nasutus. It is detailed in this page and in this vídeo.

Multiple Interval Mapping (MIM) for the phenotype corolla width

Based on the map above, a Multiple Interval Mapping (MIM) approach was performed to find QTLs, as follows.

• Exportation of the list with the linkage groups (markers and positions) from package OneMap

write_map(maps_list, "mimulus_onemap.map")

• Installation and loading the R/qtl package (more information about it in this page)

install.packages("qtl")

library("qtl")

• Directory where we can find the raw_file with genotypic and phenotypic information

raw_file <- paste(system.file("extdata", package = "onemap"),"m_feb06.raw", sep = "/")

• Data reading

f2_qtlmimulus <- read.cross("mm", file = raw_file, mapfile = "mimulus_onemap.map")

##  --Read the following data:
##  Type of cross:          f2 
##  Number of individuals:  287 
##  Number of markers:      418 
##  Number of phenotypes:   16 
##  --Cross type: f2

f2_qtlmimulus

##   This is an object of class "cross".
##   It is too complex to print, so we provide just this summary.
##     F2 intercross
## 
##     No. individuals:    287 
## 
##     No. phenotypes:     16 
##     Percent phenotyped: 96.2 93.4 96.2 93 95.8 87.8 96.2 95.8 96.2 96.2 
##                         90.2 96.2 95.8 96.2 95.8 96.2 
## 
##     No. chromosomes:    14 
##         Autosomes:      1 2 3 4 5 6 7 8 9 10 11 12 13 14 
## 
##     Total markers:      390 
##     No. markers:        20 22 26 35 21 46 32 26 25 21 35 24 19 38 
##     Percent genotyped:  89.6 
##     Genotypes (%):      AA:16.4  AB:24.8  BB:22.4  not BB:19.3  not AA:17.1

• Marker information plots

plotMap(f2_qtlmimulus,main=substitute(paste("Genetic linkage map - F2 ")(italic('Mimulus guttatus x Mimulus nasutus'))))

plotMissing(f2_qtlmimulus, main="Missing genotypes")

Corolla width (mm)

#plot_pheno<-plot.pheno(f2_qtlmimulus, pheno.col=16)
plot(plot_pheno, col="yellow", xlab="Corola width (mm)", main="ww")

• Jitter map

mimulus_jm <- jittermap(f2_qtlmimulus, amount=1e-6)
summary(mimulus_jm)

##     F2 intercross
## 
##     No. individuals:    287 
## 
##     No. phenotypes:     16 
##     Percent phenotyped: 96.2 93.4 96.2 93 95.8 87.8 96.2 95.8 96.2 96.2 
##                         90.2 96.2 95.8 96.2 95.8 96.2 
## 
##     No. chromosomes:    14 
##         Autosomes:      1 2 3 4 5 6 7 8 9 10 11 12 13 14 
## 
##     Total markers:      390 
##     No. markers:        20 22 26 35 21 46 32 26 25 21 35 24 19 38 
##     Percent genotyped:  89.6 
##     Genotypes (%):      AA:16.4  AB:24.8  BB:22.4  not BB:19.3  not AA:17.1

• Genotype simulation given observed marker data

mimulus_imp <- sim.geno(mimulus_jm, n.draws=16, step=5, off.end=0, error.prob=0.001, map.function="kosambi")

• Conditional genotype probability in each genome position

mimulus_prob <- calc.genoprob(mimulus_imp, step = 5, off.end=0, error.prob=0.001, map.function="kosambi", stepwidth="fixed")

• Using the function scantwo to perform a two-dimensional scan for corolla width phenotype.

out_mimulus<-scantwo(mimulus_prob, method="hk",pheno.col=16)

##  --Running scanone
##  --Running scantwo
##  (1,1)
##  (1,2)
##  (1,3)
##  (1,4)
##  (1,5)
##  (1,6)
##  (1,7)
##  (1,8)
##  (1,9)
##  (1,10)
##  (1,11)
##  (1,12)
##  (1,13)
##  (1,14)
##  (2,2)
##  (2,3)
##  (2,4)
##  (2,5)
##  (2,6)
##  (2,7)
##  (2,8)
##  (2,9)
##  (2,10)
##  (2,11)
##  (2,12)
##  (2,13)
##  (2,14)
##  (3,3)
##  (3,4)
##  (3,5)
##  (3,6)
##  (3,7)
##  (3,8)
##  (3,9)
##  (3,10)
##  (3,11)
##  (3,12)
##  (3,13)
##  (3,14)
##  (4,4)
##  (4,5)
##  (4,6)
##  (4,7)
##  (4,8)
##  (4,9)
##  (4,10)
##  (4,11)
##  (4,12)
##  (4,13)
##  (4,14)
##  (5,5)
##  (5,6)
##  (5,7)
##  (5,8)
##  (5,9)
##  (5,10)
##  (5,11)
##  (5,12)
##  (5,13)
##  (5,14)
##  (6,6)
##  (6,7)
##  (6,8)
##  (6,9)
##  (6,10)
##  (6,11)
##  (6,12)
##  (6,13)
##  (6,14)
##  (7,7)
##  (7,8)
##  (7,9)
##  (7,10)
##  (7,11)
##  (7,12)
##  (7,13)
##  (7,14)
##  (8,8)
##  (8,9)
##  (8,10)
##  (8,11)
##  (8,12)
##  (8,13)
##  (8,14)
##  (9,9)
##  (9,10)
##  (9,11)
##  (9,12)
##  (9,13)
##  (9,14)
##  (10,10)
##  (10,11)
##  (10,12)
##  (10,13)
##  (10,14)
##  (11,11)
##  (11,12)
##  (11,13)
##  (11,14)
##  (12,12)
##  (12,13)
##  (12,14)
##  (13,13)
##  (13,14)
##  (14,14)

\[\begin{eqnarray*} Null~model: y &=& \mu + \epsilon\\ Single~model: y &=& \mu + \beta_1q_1 + \epsilon \\ Additive~model: y &=& \mu + \beta_1q_1 + \beta_2q_2 + \epsilon\\ Full~model: y &=& \mu + \beta_1q_1 + \beta_2q_2 + \gamma(q_1 \times q_2) + \epsilon \\ \end{eqnarray*}\]
\[\begin{eqnarray*} LOD_1 (s) &=& l_1 (s) - l_0\\ LOD_a (s, t) &=& l_a (s, t) - l_0\\ LOD_f (s, t) &=& l_f (s, t) - l_0 \end{eqnarray*}\]
Maximum Likelihood \[\begin{eqnarray*} M_f (s, t) &=& LOD_{f.max} (s, t)\\ M_{fv1} (s, t) &=& LOD_{f.max} (s, t) - LOD_{1.max} (s)\\ M_i (s, t) &=& LOD_{f.max} (s, t) - LOD_{a.max} (s, t)\\ M_a (s, t) &=& LOD_{a.max} (s, t)\\ M_{av1} (s, t) &=& LOD_{a.max} (s, t) - LOD_{1.max} (s) \end{eqnarray*}\]

Table 1. LOD Score for the different models computed by scantwo.

table1<-data.frame(summary(out_mimulus))
knitr::kable(table1) 

chr1	chr2	pos1f	pos2f	lod.full	lod.fv1	lod.int	pos1a	pos2a	lod.add	lod.av1
1	1	0	80	6.372948	4.758527	3.7565362	0	50	2.6164112	1.0019910
1	2	95	55	8.041204	5.234967	3.8243210	75	45	4.2168828	1.4106456
1	3	75	95	8.294515	2.956925	1.6134148	55	95	6.6811004	1.3435103
1	4	65	110	5.055449	3.441029	1.7448630	55	210	3.3105862	1.6961660
1	5	130	20	4.022394	2.303933	0.7208028	75	10	3.3015917	1.5831307
1	6	130	5	5.539968	3.466080	1.9653507	55	35	3.5746171	1.5007291
1	7	50	160	6.743030	4.952041	3.0941551	55	160	3.6488751	1.8578856
1	8	55	90	4.660548	2.970651	1.0745116	55	90	3.5860369	1.8961390
1	9	200	40	4.948535	3.334115	2.9475765	55	35	2.0009583	0.3865381
1	10	95	135	4.283378	2.668958	1.5498400	55	135	2.7335384	1.1191182
1	11	110	165	6.387967	3.803961	2.5790576	55	160	3.8089097	1.2249035
1	12	120	25	5.682905	3.131572	1.5885686	60	60	4.0943367	1.5430035
1	13	55	115	4.868310	3.253890	2.7250139	55	95	2.1432963	0.5288761
1	14	115	115	11.615999	2.931903	1.2453314	60	115	10.3706672	1.6865720
2	2	25	145	7.733876	4.927639	3.8935532	45	130	3.8403231	1.0340859
2	3	50	95	8.287569	2.949979	0.4363755	50	95	7.8511939	2.5136038
2	4	105	215	8.572155	5.765918	4.0561284	50	210	4.5160266	1.7097894
2	5	120	85	7.540593	4.734356	2.8012993	50	90	4.7392938	1.9330566
2	6	165	200	6.637810	3.831573	1.8854373	50	15	4.7523727	1.9461356
2	7	130	5	5.714344	2.908106	0.9144877	50	165	4.7998559	1.9936187
2	8	125	90	6.512794	3.706557	2.4124695	50	90	4.1003242	1.2940870
2	9	45	70	6.348585	3.542348	3.1438642	50	40	3.2047207	0.3984835
2	10	40	100	6.070199	3.263962	2.3682145	50	135	3.7019848	0.8957477
2	11	45	155	8.040431	5.234194	2.3070262	50	160	5.7334046	2.9271674
2	12	130	60	7.245271	4.439033	1.8207699	50	60	5.4245007	2.6182635
2	13	0	60	4.833760	2.027522	1.5680322	50	125	3.2657274	0.4594902
2	14	40	115	11.557198	2.873103	0.2425833	130	115	11.3146145	2.6305194
3	3	30	125	7.601371	2.263781	1.3497598	80	105	6.2516108	0.9140207
3	4	90	100	8.250135	2.912545	1.4216918	95	100	6.8284436	1.4908534
3	5	90	95	8.737789	3.400199	1.8775325	95	10	6.8602569	1.5226667
3	6	95	140	7.919827	2.582237	1.4272034	95	10	6.4926235	1.1550334
3	7	95	160	7.197555	1.859965	0.5433560	95	160	6.6541993	1.3166092
3	8	95	0	7.663357	2.325767	0.9339274	95	90	6.7294298	1.3918397
3	9	95	200	7.776473	2.438883	1.8656412	95	0	5.9108318	0.5732417
3	10	90	90	8.141953	2.804363	1.5481610	95	100	6.5937922	1.2562020
3	11	90	110	8.286557	2.948967	0.7765364	95	160	7.5100205	2.1724303
3	12	95	65	9.425791	4.088201	0.9724884	95	60	8.4533025	3.1157123
3	13	90	60	6.806138	1.468548	1.0784343	95	125	5.7277037	0.3901135
3	14	80	110	14.294756	5.610661	2.1870317	85	115	12.1077247	3.4236295
4	4	210	240	5.902542	4.403230	3.5503288	210	240	2.3522132	0.8529010
4	5	180	85	7.599464	5.881003	4.2297546	210	15	3.3697093	1.6512483
4	6	105	10	4.754702	2.680814	1.1460697	210	35	3.6086317	1.5347438
4	7	75	255	4.644989	2.854000	1.5052802	210	160	3.1397089	1.3487194
4	8	160	20	4.270310	2.580412	1.2371491	210	90	3.0331608	1.3432628
4	9	100	185	5.239545	3.740233	3.2476863	210	160	1.9918590	0.4925467
4	10	210	170	5.894919	4.395607	3.0883345	210	135	2.8065847	1.3072725
4	11	100	20	6.032876	3.448870	2.0117977	210	165	4.0210787	1.4370725
4	12	215	70	7.192099	4.640766	2.9017055	210	70	4.2903939	1.7390607
4	13	95	40	4.401435	2.902123	2.3038380	210	125	2.0975969	0.5982847
4	14	210	115	11.817464	3.133369	1.5829835	210	115	10.2344802	1.5503851
5	5	15	90	4.110225	2.391764	0.4070998	15	90	3.7031250	1.9846640
5	6	65	35	5.219583	3.145695	1.6689851	90	35	3.5505980	1.4767100
5	7	90	250	5.616746	3.825756	2.1595630	10	160	3.4571825	1.6661930
5	8	5	90	6.622033	4.903572	2.7367469	10	90	3.8852862	2.1668251
5	9	150	150	5.796962	4.078500	3.7137957	10	140	2.0831658	0.3647048
5	10	95	160	3.838016	2.119555	0.8017632	10	135	3.0362527	1.3177917
5	11	85	100	5.464854	2.880848	1.1947614	10	160	4.2700930	1.6860868
5	12	90	70	5.916983	3.365650	1.3010801	90	60	4.6159033	2.0645701
5	13	90	105	3.755044	2.036583	1.6169280	10	0	2.1381160	0.4196549
5	14	90	115	10.791956	2.107861	0.9221243	90	115	9.8698317	1.1857366
6	6	15	30	5.512620	3.438732	2.4347159	35	235	3.0779041	1.0040162
6	7	35	165	6.415779	4.341891	2.2808855	35	160	4.1348934	2.0610054
6	8	35	180	5.808816	3.734928	2.1127675	35	90	3.6960484	1.6221604
6	9	85	35	3.850889	1.777001	1.2287430	15	0	2.6221464	0.5482584
6	10	5	140	6.546841	4.472953	3.1362517	35	135	3.4105894	1.3367014
6	11	35	100	5.884164	3.300158	1.3586980	10	160	4.5254659	1.9414597
6	12	0	60	5.659294	3.107960	1.0199801	35	60	4.6393136	2.0879804
6	13	30	110	4.057324	1.983436	1.3846947	35	125	2.6726298	0.5987418
6	14	90	115	11.211627	2.527532	1.1881122	35	115	10.0235147	1.3394195
7	7	95	165	5.028742	3.237753	1.9986793	165	170	3.0300631	1.2390736
7	8	250	90	6.246422	4.455433	2.6228682	160	90	3.6235539	1.8325644
7	9	165	90	5.256483	3.465494	3.0452589	165	115	2.2112246	0.4202351
7	10	160	135	5.377868	3.586879	2.3278078	165	100	3.0500602	1.2590707
7	11	160	160	7.378591	4.794585	2.9447981	165	160	4.4337927	1.8497865
7	12	210	70	6.687595	4.136262	2.7099330	160	60	3.9776620	1.4263288
7	13	245	35	4.129385	2.338395	1.7102354	165	125	2.4191492	0.6281597
7	14	85	115	12.355961	3.671866	2.1809299	160	115	10.1750310	1.4909359
8	8	90	120	4.721621	3.031723	1.4677660	90	105	3.2538545	1.5639565
8	9	65	150	4.141355	2.451457	2.0028699	90	0	2.1384849	0.4485870
8	10	95	40	3.743129	2.053231	0.8265317	90	135	2.9165969	1.2266989
8	11	95	200	7.679905	5.095899	3.8201808	90	160	3.8597246	1.2757184
8	12	95	60	7.216253	4.664920	3.1647851	90	60	4.0514677	1.5001345
8	13	95	120	6.450437	4.760539	3.9014381	90	125	2.5489986	0.8591006
8	14	90	115	10.939623	2.255528	0.7228954	90	115	10.2167277	1.5326325
9	9	5	40	3.433901	3.068338	1.9373418	5	35	1.4965593	1.1309965
9	10	35	135	4.749769	3.428203	3.0905171	40	135	1.6592514	0.3376861
9	11	235	90	5.306706	2.722700	2.2801645	35	160	3.0265415	0.4425353
9	12	185	60	5.103604	2.552271	2.2840969	40	60	2.8195072	0.2681740
9	13	145	60	4.734535	4.176722	3.7490690	140	125	0.9854659	0.4276528
9	14	65	115	10.291511	1.607416	1.1001934	140	115	9.1913177	0.5072226
10	10	110	135	4.813339	3.491774	1.2807118	120	130	3.5326275	2.2110622
10	11	160	205	5.733353	3.149347	2.0178333	135	160	3.7155194	1.1315132
10	12	0	90	5.012301	2.460968	1.3209884	135	60	3.6913129	1.1399796
10	13	105	125	4.299812	2.978247	2.2321502	135	125	2.0676622	0.7460969
10	14	135	115	13.065347	4.381252	2.8632751	135	115	10.2020720	1.5179769
11	11	140	160	6.758760	4.174754	2.6723989	100	195	4.0863613	1.5023551
11	12	165	60	6.148026	3.564020	1.3353766	160	60	4.8126498	2.2286436
11	13	260	120	5.398040	2.814034	2.3479120	160	125	3.0501281	0.4661219
11	14	30	115	12.620094	3.935999	2.5190260	100	115	10.1010681	1.4169730
12	12	20	50	4.783104	2.231771	1.1592046	0	70	3.6238993	1.0725661
12	13	50	115	5.361009	2.809676	2.1565007	60	125	3.2045089	0.6531756
12	14	100	115	12.801869	4.117774	1.9503750	60	115	10.8514938	2.1673987
13	13	85	95	3.917918	3.360105	1.8067348	95	105	2.1111832	1.5533701
13	14	125	115	9.910744	1.226649	0.5107573	125	115	9.3999869	0.7158918
14	14	115	255	10.692260	2.008165	0.6015512	60	115	10.0907090	1.4066138

plot(out_mimulus,col.scheme = "redblue")

plot(out_mimulus, lower="fv1",col.scheme = "redblue")

• Using permutations to assess significance level

out_mimulus_perm<-scantwo(mimulus_prob, method="hk", n.perm=1000, verbose=T,pheno.col=16)

## Doing permutation in batch mode ...

threshold <- summary(out_mimulus_perm)
summary(out_mimulus_perm,alpha=0.10)

##  (1000 permutations)
##     full  fv1  int add  av1  one
## 10% 8.82 6.74 5.71   6 3.23 3.32

• Permutation distribution plot: significant threshold in percentile 95 for full, fv1, add, av1, one

old.par <- par(mfrow=c(3, 2))
freq_perm_full<-hist(out_mimulus_perm$full, col="blue", main="full")
abline(freq_perm_full, v=threshold$full[1,1], col="red", lty=3, lwd=3)
freq_perm_fv1<-hist(out_mimulus_perm$fv1, col="yellow", main="fv1")
abline(freq_perm_fv1, v=threshold$fv1[1,1], col="red", lty=3, lwd=3)
freq_perm_add<-hist(out_mimulus_perm$add, col="pink", main="add")
abline(freq_perm_add, v=threshold$add[1,1], col="red", lty=3, lwd=3)
freq_perm_av1<-hist(out_mimulus_perm$av1, col="green", main="av1")
abline(freq_perm_av1, v=threshold$av1[1,1], col="red", lty=3, lwd=3)
freq_perm_one<-hist(out_mimulus_perm$one, col="orange", main="one")
abline(freq_perm_one, v=threshold$one[1,1], col="red", lty=3, lwd=3)
par(old.par)

• QTL interaction

summary_qtl1<- summary(out_mimulus,perms=out_mimulus_perm, alpha=c(0.1), pvalues = T)
summary_qtl1

##        pos1f pos2f lod.full pval lod.fv1  pval lod.int pval     pos1a
## c3:c14    80   110     14.3    0    5.61 0.616    2.19    1        85
##        pos2a lod.add pval lod.av1  pval
## c3:c14   115    12.1    0    3.42 0.055

summary(out_mimulus_perm,alpha=0.10)

##  (1000 permutations)
##     full  fv1  int add  av1  one
## 10% 8.82 6.74 5.71   6 3.23 3.32

Maximum Likelihood \[\begin{eqnarray*} M_f (s, t) &=& LOD_{f.max} (s, t)\\ M_{fv1} (s, t) &=& LOD_{f.max} (s, t) - LOD_{1.max} (s)\\ M_i (s, t) &=& LOD_{f.max} (s, t) - LOD_{a.max} (s, t)\\ M_a (s, t) &=& LOD_{a.max} (s, t)\\ M_{av1} (s, t) &=& LOD_{a.max} (s, t) - LOD_{1.max} (s) \end{eqnarray*}\]

Criteria for printing chromosome pair

• \(M_f (s, t) \ge T_f\) and \(M_{fv1}(s, t) \ge T{fv1}\\\)
• \(M_i(s, t) \ge T_i\\\)
• \(M_a(s, t) \ge T_a\) and \(M_{av1}(s, t) \ge T_{av1}\)

Table 2. LOD Score and p value for the interaction between C3 and C14 (envolving full model)

Linkage.Group1	Linkage.Group2	Position1	Position2	LOD.full	p.value.full	LOD.fv1	p.value.fv1	LOD.int	p.value.LOD.int
3	14	80	110	14.29476	0	5.610661	0.616	2.187032	1

Table 3. LOD Score and p-value for the interaction between linkage groups 3 and 14, considering additive models.

Linkage.Group1	Linkage.Group2	Position1	Position2	LOD.add	pvalue.add	LOD.av1	pvalue.av1
3	14	85	115	12.10772	0	3.42363	0.055

• QTL observation on the linkage map

qtl_mim <- makeqtl(mimulus_prob, chr=c(3,14), pos=c(85,112), what="prob")

plot(qtl_mim)

• Fitting the multiple-QTL model

out_fit1 <- fitqtl(mimulus_prob, pheno.col=16, qtl=qtl_mim, method="hk", get.ests=TRUE)
summary(out_fit1)

## 
##      fitqtl summary
## 
## Method: Haley-Knott regression 
## Model:  normal phenotype
## Number of observations : 276 
## 
## Full model result
## ----------------------------------  
## Model formula: y ~ Q1 + Q2 
## 
##        df        SS         MS      LOD     %var Pvalue(Chi2)    Pvalue(F)
## Model   4  577.2031 144.300770 12.63863 19.01292 6.917689e-12 1.041989e-11
## Error 271 2458.6439   9.072487                                            
## Total 275 3035.8469                                                       
## 
## 
## Drop one QTL at a time ANOVA table: 
## ----------------------------------  
##          df Type III SS   LOD   %var F value Pvalue(Chi2) Pvalue(F)    
## 3@85.0    2       150.8 3.567  4.966   8.309            0  0.000315 ***
## 14@112.3  2       364.1 8.277 11.994  20.068            0  7.46e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Estimated effects:
## -----------------
##               est      SE      t
## Intercept 15.7592  0.1850 85.203
## 3@85.0a    1.1395  0.2796  4.076
## 3@85.0d   -0.5221  0.4377 -1.193
## 14@112.3a  1.7384  0.2744  6.334
## 14@112.3d  0.1517  0.3675  0.413

• Calculate penalties

pen <- calc.penalties(out_mimulus_perm, alpha = 0.1)
summary(out_mimulus_perm)

## ww (1000 permutations)
##     full  fv1  int  add  av1  one
## 5%  9.12 7.08 6.10 6.44 3.44 3.71
## 10% 8.82 6.74 5.71 6.00 3.23 3.32

print(pen)

##     main    heavy    light 
## 3.320801 5.714057 3.422012

• Perform stepwise to confirm the model

out_stp <- stepwiseqtl(mimulus_prob,pheno.col=16, penalties = pen, max.qtl=4, method="hk", verbose=T, refine.locations = TRUE)

##  -Initial scan
## initial lod:  9.071802 
## ** new best ** (pLOD increased by 5.751)
##     no.qtl =  1   pLOD = 5.751001   formula: y ~ Q1 
##  -Step 1 
##  ---Scanning for additive qtl
##         plod = 5.997023 
##  ---Scanning for QTL interacting with Q1
##         plod = 5.257624 
##  ---Refining positions
##     no.qtl =  2   pLOD = 5.997023   formula: y ~ Q1 + Q2 
## ** new best ** (pLOD increased by 0.246)
##  -Step 2 
##  ---Scanning for additive qtl
##         plod = 5.535213 
##  ---Scanning for QTL interacting with Q1
##         plod = 3.91849 
##  ---Scanning for QTL interacting with Q2
##         plod = 4.571991 
##  ---Look for additional interactions
##         plod = 5.029571 
##  ---Refining positions
##     no.qtl =  3   pLOD = 5.535213   formula: y ~ Q1 + Q2 + Q3 
##  -Step 3 
##  ---Scanning for additive qtl
##         plod = 4.60837 
##  ---Scanning for QTL interacting with Q1
##         plod = 3.617351 
##  ---Scanning for QTL interacting with Q2
##         plod = 2.683523 
##  ---Scanning for QTL interacting with Q3
##         plod = 3.432879 
##  ---Look for additional interactions
##         plod = 4.246929 
##  ---Refining positions
##     no.qtl =  4   pLOD = 4.60837   formula: y ~ Q1 + Q2 + Q3 + Q4 
##  -Starting backward deletion
##  ---Dropping Q4 
##     no.qtl =  3   pLOD = 5.535213   formula: y ~ Q1 + Q2 + Q3 
##  ---Refining positions
##  ---Dropping Q3 
##     no.qtl =  2   pLOD = 5.997023   formula: y ~ Q1 + Q2 
##  ---Refining positions
##  ---Dropping Q2 
##     no.qtl =  1   pLOD = 5.751001   formula: y ~ Q1 
##  ---Refining positions
##  ---One last pass through refineqtl

out_stp

##   QTL object containing genotype probabilities. 
## 
##        name chr    pos n.gen
## Q1   3@85.0   3  85.00     3
## Q2 14@112.3  14 112.33     3
## 
##   Formula: y ~ Q1 + Q2 
## 
##   pLOD:  5.997

• Refine location

ref<-refineqtl(mimulus_prob,pheno.col=16,qtl_mim,chr=c(3,14),pos = c(80,110))

## pos: 85 112.3271 
## Iteration 1 
##  Q2 pos: 112.3271 -> 112.3271
##     LOD increase:  0 
##  Q1 pos: 85 -> 85
##     LOD increase:  0 
## all pos: 85 112.3271 -> 85 112.3271 
## LOD increase at this iteration:  0 
## overall pos: 85 112.3271 -> 85 112.3271 
## LOD increase overall:  0

ref

##   QTL object containing genotype probabilities. 
## 
##        name chr    pos n.gen
## Q1   3@85.0   3  85.00     3
## Q2 14@112.3  14 112.33     3

plotLodProfile(ref)

lodint(out_stp, chr = 3, qtl.index = 1)

##           chr pos      lod
## c3.loc20    3  20 1.895579
## c3.loc85    3  85 3.566824
## c3.loc140   3 140 1.884630

lodint(out_stp, chr = 14, qtl.index = 2)

##          chr      pos      lod
## BC498     14 108.3305 6.031331
## AAT374    14 112.3271 8.277387
## MgSTS388  14 117.0479 6.705145

out_fit_f <- fitqtl(mimulus_prob, pheno.col=16, qtl=ref, method="hk", get.ests=TRUE)
summary(out_fit_f)

## 
##      fitqtl summary
## 
## Method: Haley-Knott regression 
## Model:  normal phenotype
## Number of observations : 276 
## 
## Full model result
## ----------------------------------  
## Model formula: y ~ Q1 + Q2 
## 
##        df        SS         MS      LOD     %var Pvalue(Chi2)    Pvalue(F)
## Model   4  577.2031 144.300770 12.63863 19.01292 6.917689e-12 1.041989e-11
## Error 271 2458.6439   9.072487                                            
## Total 275 3035.8469                                                       
## 
## 
## Drop one QTL at a time ANOVA table: 
## ----------------------------------  
##          df Type III SS   LOD   %var F value Pvalue(Chi2) Pvalue(F)    
## 3@85.0    2       150.8 3.567  4.966   8.309            0  0.000315 ***
## 14@112.3  2       364.1 8.277 11.994  20.068            0  7.46e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Estimated effects:
## -----------------
##               est      SE      t
## Intercept 15.7592  0.1850 85.203
## 3@85.0a    1.1395  0.2796  4.076
## 3@85.0d   -0.5221  0.4377 -1.193
## 14@112.3a  1.7384  0.2744  6.334
## 14@112.3d  0.1517  0.3675  0.413

Table 4. Linkage Group, position, LOD score, QTLs effects (additive and dominant) and determination coefficient (%) for the Multiple Interval Mapping (CIM) for the phenotype corolla width, using the Haley-Knott regression.

	QTL	Linkage.Group	Position	LOD	Additive.effect	Dominance.effect	R2
3@85.0a	Q1	3	85.0000	3.5668	1.1395	-0.5221	4.9662
14@112.3a	Q2	14	112.3271	8.2774	1.7384	0.1517	11.9945

Comparison between the Composite Interval Mapping (CIM) and Multiple Interval Mapping (MIM)

A Composite Interval Mapping was performed in the same data set, as detailed in this page and in this vídeo

• Interval Mapping (IM) and Composite Interval Mapping (CIM)- entire genome

knitr::include_graphics('./mapchart/cimemim.png')

• Interval Mapping (IM) for chromosomes 3 and 14

knitr::include_graphics('./mapchart/3e14im.png')

• Composite Interval Mapping (IM) for chromosomes 3 and 14

knitr::include_graphics('./mapchart/cim.png')

• Multiple Interval Mapping (MIM) for chromosomes 3 and 14

plotLodProfile(ref)