欢迎关注微信公众号"生信小王子"！
QTL分析是进行基因精细定位和克隆的基础，今天小编教大家使用R包" qtl "进行QTL分析。
在开始分析前，我们需要准备两个输入文件：基因型和表型文件。
基因型文件：
表型文件：
基因型和表型文件均保存为逗号分隔的csv文件。
准备好两个输入文件后，我们就可以开始分析啦！
## 安装R包
install.packages("qtl")
## 加载R包
library("qtl")
## 导入基因型和表型数据
sug <- read.cross("csvs", ".", "gen.csv", "phe.csv")
## 查看输入文件相关信息
summary(sug)
此外，还有一些函数可以统计对应的信息。
## 样本数
nind(sug)
## 染色体数
nchr(sug)
## 标记数
totmar(sug)
## 每个染色体上的标记数
nmar(sug)
## 表型数
nphe(sug)
除了文字信息外，我们还可以用图来展示这些信息。
plot(sug)
这三张图分别展示了缺失的基因型数据，遗传图谱和表型数据分布。
也可以单独展示这三张图。
## 展示缺失基因型数据（黑色为缺失的基因型）
plotMissing(sug)
## 绘制遗传图谱
plotMap(sug)
## 绘制表型分布直方图
plotPheno(sug, pheno.col=1)
## 计算基因型概率
sug <- calc.genoprob(sug, step=1)
## 使用默认方法进行single-QTL全基因组扫描
out.em <- scanone(sug)
## 查看扫描结果
summary(out.em)
## 挑选LOD>3的结果
summary(out.em, threshold=3)
## 展示结果
plot(out.em)
我们可以看到7号和15号染色体上各有一个显著的峰。
## 使用Haley-Knott回归方法进行全基因组扫描
out.hk <- scanone(sug, method="hk")
## 使用Multiple imputation法进行全基因组扫描
sug <- sim.geno(sug, step=1, n.draws=64)
out.imp <- scanone(sug, method="imp")
## 比较三种方法结果的差异
plot(out.em, out.hk, out.imp, col=c("blue", "red", "green"))
我们可以看到，三种方法的结果并没有明显差异。
## 进行1000次Permutation test
operm <- scanone(sug, method="hk", n.perm=1000)
## 获得显著性阈值
summary(operm, alpha=c(0.05, 0.2))
## 从扫描结果中挑选显著的位点
summary(out.hk, perms=operm, alpha=0.2, pvalues=TRUE)
接下来，我们需要估计QTL区间。因为我们通过LOD值过滤后的QTL位点位于7号和15号染色体上，所以我们首先对7号染色体上的QTL区间的进行估计。
## 获得7号染色体1.5倍LOD区间和95%贝叶斯区间
lodint(out.hk, chr=7, drop=1.5)
bayesint(out.hk, chr=7, prob=0.95)
第一行和第三行是区间的范围，第二行是预测QTL的位置。
## 获得区间两侧最近的标记
lodint(out.hk, chr=7, expandtomarkers=TRUE)
bayesint(out.hk, chr=7, expandtomarkers=TRUE)
## 获得离QTL最近的标记
max(out.hk)
mar <- find.marker(sug, chr=7, pos=47.7)
## 统计不同基因型个体的表型
plotPXG(sug, marker=mar)
红色的点表示基因型进行过填补的个体。
## 统计不同基因型个体的效应
effectplot(sug, mname1=mar)
Tips：我们常说的LOD值=log10 (L1/L0) ，L1指该位点含QTL的概率，L0指该位点不含QTL的概率。LOD值为3表示该位点含QLT的概率是不含QTL概率的1000倍。
参考资料：
https://rqtl.org/tutorials/rqtltour2.pdf
网友评论