使用Plink对SNP数据进行质量控制

前

在做和基因相关的分析时，拿到基因型数据，首先要进行质量控制。

在使用Plink进行质量控制时，一般包含以下几步：

1. individual and SNP missingness,（筛选个体和SNP）
2. inconsistencies in assigned and genetic sex of subjects (see sex discrepancy),（性别是否出错）
3. minor allele frequency (MAF),（次等位基因频率筛选）
4. deviations from Hardy–Weinberg equilibrium (HWE),（哈代平衡检测）
5. heterozygosity rate,（杂合率）
6. relatedness,（相关性）
7. ethnic outliers (see population stratification).（人口分层）

当然一般也要根据数据进行步骤的选择，因为我们处理数据是为了后续的相关性分析，所以这里没有使用Plink做相关分析，同时因为人口分层暂时没有特别了解，所以只记录使用Plink进行上述前五步的工作。

在Plink中做这几步的质量控制，可以直接使用命令+参数完成筛选，也可以使用Plink生成检测结果来手动筛选。当然，一般对数据要求严格的实验都是先生成结果，经过查看以后确定参数再使用Plink进行筛选。

中

具体来讲，根据我的习惯一般先进行个体和SNP的筛选，次等位基因频率、哈代平衡的检测，之后做性别的检测，剔除性别出错的样本，再做一遍上述的筛选。

对个体和SNP的检测和筛选

根据个体和SNP的缺失率进行挑选。

–missing 生成个体和SNP缺失的检测结果。

个体缺失位点结果文件为plink.imiss：
第一列为家系ID，第二列为个体ID，第三列是否表型缺失，第四列缺失的SNP个数，第五列总SNP个数，第六列缺失率。

单个SNP缺失的结果文件为plink.lmiss：
第一列为染色体，第二列为SNP名称，第三列为缺失个数，第四列为总个数，第五列为缺失率

（这里有一点更正，之前SNP和个体写反了，感谢评论区提醒。 2020/08/06）
–geno 筛选SNP
–mind 筛选个体

–mind 0.05：如果一个SNP5%的个体中都是缺失的，那么就删掉该个体。
–geno 0.05：如果一个个体有5%的SNP都是缺失的，那么就删掉该SNP。

哈代平衡

根据基因型频率进行筛选。对于哈代平衡可以自行查找资料理解。

–hardy 计算所有位点的哈温检测结果

结果文件为plink.hwe：

第一列为染色体，第二列为SNP的ID，第三列为TEST类型，第四列A1为 minor 位点，第五列A2为 major 位点，第六列为基因型分布，对应A1A1, A1A2, A2A2的个数，第六列为观测杂合度频率，第七列为期望杂合度频率，第八列为哈温平衡的卡方检验P-value值。

–hwe 直接过滤

次等位基因频率

根据基因频率进行挑选。

这里谈一点个人对次等位基因的理解，某一位点可能为AA，AT，TT，甚至出现更大的变异导致CC，CG，所以才有次等位基因频率，而不是最小，很多人将minor allele frequency（MAF）解释为最小等位基因频率，感觉是不够准确的。

某一位点可能有A，T，C，G四种基因情况，四种情况的频率和应该为1，第二常见的基因型为次等位基因，次等位基因出现的基因频率就称为次等位基因频率。

而使用MAF进行过滤，是因为一般在人群中肯定是少数患病，那么导致少数人患病的就可能是这个次等位基因，而当MAF小于一定值时，意味着大部分位点都是相同的基因型，这些SNP位点贡献的信息很少（和所研究关系很小），会增加假阳性。极端情况MAF为0，该位点只有一种基因型，那么致病基因肯定和这个位点无关，故而过滤掉。

–freq 计算每个SNP位点的次基因频率

结果文件为plink.frq：
第一列为染色体号，第二列为SNP ID，第三列为 minor 位点，第四列为 major 位点，第五列为次等位基因频率，第六列为等位基因观察数。

–maf 直接过滤

性别检测

–check-sex 可以生成性别检测的结果，然后根据结果构造removesample.txt（名称自由），再使用 --remove removesamp.txt 将性别出错的样本剔除。
removesample 文本第一列为FID，第二列为IID。

结果文件为plink.sexcheck：
第一列为家系ID，第二列为个体ID，第三列为原始数据中所给性别，第四列为根据SNP推断的性别，第五列是否正常，第六列为F值。

女性受试者的F值必须小于0.2，男性受试者的F值必须大于0.8。这个F值是基于X染色体近交（纯合子）估计。不符合这些要求的受试者第五列为“PROBLEM”，符合则为“OK”。

后

将上述质量控制的完整流程总结成可用命令如下：

plink --bfile inputname --geno 0.05 --mind 0.1 --maf 0.05  --hwe 0.000001 --make-bed --out outputname1

plink --bfile outputname1 --check-sex 

plink --bfile outputname1 --remove removesamp.txt --make-bed --out outputname2

plink --bfile outputname2 --geno 0.05 --mind 0.1 --maf 0.05  --hwe 0.000001 --make-bed --out outputname3

上述inputname、outputname均为文件名称。

上图是某次质量控制Plink日志文件的说明。

其实除了上述质量控制，根据实验有时也需要做亲缘关系的排查，结果如下。
还有–het，杂合率的检测，因为本人实验并不需要这些，有时间再写。

上述均根据Plink1.90完成。Plink基本操作可参考PLink常用命令总结，更多有关Plink命令的使用可以在参考中去查找。

参考：
http://www.cog-genomics.org/plink2/