[GWAS]杂合度质控_个体杂合度-优快云博客

本文链接：https://blog.youkuaiyun.com/qq_63511424/article/details/131757741

文章介绍了在自然群体中，过高或过低的杂合度可能是样品污染或近亲繁殖的信号。通过Plink工具，可以进行杂合度过滤，如删除偏离平均值±3SD的个体。在非自然群体如自交系和F1杂交种中则无需此步骤。文章详细展示了如何使用Plink命令行工具进行数据处理，包括生成不同报告并利用R语言进行杂合度计算和可视化，以及如何识别并移除超出标准差范围的个体。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

在一群自然群体种，基因型个体的杂合度过高或者过低都是不正常的，我们需要根据杂合度来进行过滤。这里的偏差可能表明样品收到污染，或者近亲繁殖。
我们采取的做法是删除样品杂合率平均值中偏离±3 SD（sd就是标准差）的个体

非自然群体中，比如自交系，杂交种F1，这些群体不需要过滤杂合度

参数过滤合手动过滤：
plink有个特点，所有的过滤标准，都可以生成过滤之前的文件，然后可以手动过滤，也可以用参数进行过滤。

比如：–missing生成结果，也可以用–geno用–mind过滤
比如：–hardy生成结果，可以使用–hwe过滤
比如：–freq生成结果，可以用–maf过滤，当时杂合度 --het，没有用来过滤的函数，所以只能通过编程来实现。

首先我们提取ID，然后用–remove来实现过滤
->计算杂合度：

plink --bfile HapMap_3_r3_9 --het --out R_check

我们来分析结果文件，首先解读一下各个列的信息：

FID 👉Family ID #家系ID
D 👉Within-family ID #个体ID
HOM) Observed number of homozygotes #实际纯合个数
HOM) Expected number of homozygotes #期望纯合个数
NM) Number of non-missing autosomal genotypes #总个数
Method-of-moments F coefficient estimate

杂合度=（N-O）/N
我们将杂合度可视化

# 先要设置工作路径

het = read.table("R_check.het", head=TRUE)
het$HET_RATE = (het$"N.NM." - het$"O.HOM.")/het$"N.NM."
hist(het$HET_RATE, xlab="Heterozygosity Rate", ylab="Frequency", main= "Heterozygosity Rate")

pdf("heterozygosity.pdf")
het$HET_RATE = (het$"N.NM." - het$"O.HOM.")/het$"N.NM."
hist(het$HET_RATE, xlab="Heterozygosity Rate", ylab="Frequency", main= "Heterozygosity Rate")
dev.off()

接着计算杂合度三倍标准差以外的个体
首先我们来查看哪些个体在3倍的标准差之外

# 先要设置工作路径

het = read.table("R_check.het", head=TRUE)
het$HET_RATE = (het$"N.NM." - het$"O.HOM.")/het$"N.NM."
het_fail = subset(het, (het$HET_RATE < mean(het$HET_RATE)-3*sd(het$HET_RATE)) | (het$HET_RATE > mean(het$HET_RATE)+3*sd(het$HET_RATE)));
het_fail$HET_DST = (het_fail$HET_RATE-mean(het$HET_RATE))/sd(het$HET_RATE);

write.table(het_fail, "fail-het-qc.txt", row.names=FALSE)

我们从结果可以看出，这两个个体杂合度在3倍标准差之外
去掉这两个个体之后，先对数据进行清洗，去掉引号，然后提取家系合个体ID

我们使用remove去掉这两个个体

plink --bfile HapMap_3_r3_9 --remove het_fail_ind.txt --make-bed --out HapMap_3_r3_10