上一部分介绍了数据导入部分，有需要的可以移步：
理论上来说，在数据导入之后，第一项工作就应该是识别和纠正数据中的错误，确保没有太明显的问题。
所以这一块就基于 R/qtl 中的自带文件介绍一下几种常见的问题。
1. 表型数据
检测具有异常表型的个体，这种异常多是由于数据输入过程中存在错误，需要仔细甄查。以示例数据集 ch3a 为例：
library(qtlbook)
data(ch3a)
该数据集包括 5 个表型，首先可以绘制每个表型的直方图分布（未展示）。
par(mfrow=c(3,2))
for(i in 1:5)
plot.pheno(ch3a, pheno.col=i)
但这个图中不是很容易看到异常表型个体，相较而言散点图可能更明显一些。
pairs(jitter( as.matrix(ch3a$pheno) ), cex=0.6, las=1)
在散点图中，可以明显看到 phe4 中有一个个体表型为0，较为异常，需要检查一下是数据缺失或输入错误造成的。
如果为 0 的异常表型是输入错误导致的，可以将其赋值为缺失，如下命令：
ch3a$pheno[ch3a$pheno == 0] <- NA
2. 偏分离检测
遗传标记如果展示出明显的偏分离，则可能是基因分型过程出现问题造成的。
以示例数据 ch3b 为例，geno.table 可以检测每个标记的基因型频率，最后一列是其是否符合孟德尔比例的卡方检验的 p 值。
data(ch3b)
gt <- geno.table(ch3b)
gt[ gt$P.value < 1e-7, ]
上述命令可以调出显著偏分离的遗传标记（p < 1e-7）。
有一些文献是会过滤掉偏分离的标记的，不过大家的操作没有很一致，这块我后面再多看些文献。师兄曾有尝试，去掉和不去掉差别不大，目前我是倾向于过滤掉这些异常位点的。
3. 个体基因型比较
有时候个体可能混淆，R/qtl 还提供了判断是否存有着非常相似的基因型的个体。在一些群体内，如果出现基因型高度甚至完全一致的个体，可以是错误导致的。
同样以示例数据及 ch3a 为例，构建两两之间相同遗传标记所占比例的直方图：
data(ch3a)
cg <- comparegeno(ch3a) # 比较个体之间遗传标记的相似度
hist(cg, breaks=200, xlab="Proportion of identical genotypes") # 绘制直方图
rug(cg) # 在直方图下绘制轴须图
从轴须图中可以明显看到右端存在基因型异常相似的个体，可以调出来看看是哪些个体，做进一步的检查：
which(cg > 0.9, arr.ind=TRUE)
可以看出，138 和 5 之间，以及 12 和 55 之间，异常相似，值得怀疑。
4. 标记顺序
标记在染色体上的顺序如果发生错误，会直接影响 QTL 的分析结果。另外，即使是基于高质量物理图谱得到的结果，也可能发生标记位置的错误。这一部分是重点。
4.1 成对标记之间的重组率
一般来说，不同染色体上的标记应该是不相关的，同一条染色体上的标记相对来说连接会更紧密。举个例子，如果 1 号染色体上的某个标记和 3 号染色体上的某个标记连锁，那可能是哪里出错了。
函数 est.rf 可用于估计所有标记两两之间的重组率。
data(ch3c)
ch3c <- est.rf(ch3c)
plotRF(ch3c, alternate.chrid=TRUE)
上面的命令能够得到上图所示的不同标记之间的重组率关系。可以发现，1、7、12、13、15 号染色体上存在重组率异常的标记。
plot.rf(ch3c, chr=c(1, 7, 12, 13, 15))
细看的话，可以看清楚具体哪些标记可能存在问题。
还有方法进一步确定这些异常的标记，就是估计标记之间的遗传距离并绘图。
nm <- est.map(ch3c, error.prob=0.001)
plot(nm)
从上面的两个图中，可以清楚地看到哪些标记可能出现了问题。存在一种可能，1 号染色体上的那个异常标记属于 7 号染色体；13 号染色体上的异常标记属于 12 号染色体。R/qtl 包提供了改变标记位置的选项，不过也可以直接去掉这些标记。
4.2 染色体内标记的顺序检查
和基因组上标记的重组率类似，单条染色体上的标记顺序也需要推定（个人推测，这些还是主要基于以前没有物理图谱时的分析）。
这里推定的原理是最大似然和最少的 obligate crossovers，在一个窗口内列举所有标记的排列组合以及响应的 crossover，推定最少 crossovers 的为最佳顺序（简单来说，就是在上面的图种，把染色体内的方块也捋顺了。但现在如果有高质量物理图谱的情况下，这一步是否十分必要，待进一步阅读和思考。故未详细介绍。）
补充：还是有一定必要的，如果相较参考基因组存在较大的结构变异（如：倒位），这一块的分析还是有必要看看的。
4.3 估计遗传图谱
确定好标记的顺序后，可以进行遗传图谱的估计，并将前后的遗传图谱进行比较。
nm <- est.map(ch3c, error.prob=0.001, verbose=FALSE)
plot.map(ch3c, nm)
而后，完成遗传图谱的替换即可。
ch3c <- replace.map(ch3c, nm)
5. 可能的基因分型错误
简单来说，就是 crossover 出现的距离很近，可能是测序错误造成的。但标记太稀疏的话就无法判断了。
data(hyper)
newmap <- est.map(hyper, error.prob=0.01)
hyper <- replace.map(hyper, newmap)
hyper <- calc.errorlod(hyper)
top <- top.errorlod(hyper, cutoff=5)
top
plot.geno(hyper, 16, top$id[top$chr==16], cutoff=5)
图中标注 × 的位置未可能的 crossover，红圈内可能是基因分型错误。
少量错误对结果影响不大，如果类似的标记很多，就需要重新审视一下数据的情况了。
6. crossovers 数目统计
函数 countXO 可以统计每个个体上观察到的 crossovers 数目。如果观察到异常多或者异常少的个体，需要进一步检查确定没有搞错。
nxo <- countXO(hyper)
plot(nxo, ylab="No. crossovers")
7. 缺失基因型信息统计
这一块就是看一下哪些个体内缺失标记的数目异常高，以及在基因组中的哪些位置。
（对我来说，先是知道这回事就可以了，后面有需要再进一步学习）
这一部分看起来好像没什么用，但做科研嘛，严谨还是必要的。后面就要开始介绍 QTL 定位分析的计算过程了。
要是觉得有用可以登录一下账号，点个赞，以表支持！
网友评论