探秘全基因组关联分析第二期：数据质控

 

上一期我们了分享了GWAS分析需要的数据格式，以及不同格式之间的转换。现在我们已经准备好了表型数据和基因数据，是不是就想马上进行关联分析了？心急吃不了热豆腐，为了提高关联分析结果的准确性，需要对数据进行质控，去掉不合格的样本和变异数据。

 

 

1、SNP及个体缺失过滤 

 

人工采集的数据，可能存在位点基因型和个体基因数据缺失（表型缺失的直接去掉），这些缺失数据影响关联分析的准确性，需要将缺失率控制在一定标准以下。建议首先以宽松的阈值（0.2；> 20%）过滤SNP和个体，从而过滤掉缺失程度很高的SNP和个体；再使用更严格的阈值过滤（(0.02；> 2%）。

 

# SNP缺失过滤

$plink --noweb --bfile $project.raw.mark --geno 0.2 --allow-no-sex --make-bed --out ${project}.filter.mds1

# 个体缺失过滤

$plink --noweb --bfile ${project}.filter.mds1 --mind 0.2 --allow-no-sex --make-bed --out ${project}.filter.mds2

 

注意：以上步骤更换更严格的参数再过滤一遍。

 

 

2、性别和亲缘关系检测（可选） 

 

性别检测基于X染色体近交系（纯合子性）估计，一般女性受试者的F值 < 0.2，男性受试者的F值 > 0.8，不满足这些要求的被标记为“PROBLEM”。

 

# 性别检测

$plink --noweb --bfile ${project}.raw.mark --check-sex

# 输出结果保存在plink.sexcheck文件中，提取性别异常个体

$grep "PROBLEM" plink.sexcheck | awk '{print $1,$2}' >sex_removelist.txt

# 删除性别异常个体（不建议删除，除非明确该样本数据有污染）

$plink --noweb --bfile ${project}.raw.mark --remove sex_removelist.txt --make-bed --out ${project}.raw.mark2

 

亲缘关系检测基于遗传信息，判断样本亲缘关系的指标分为状态同源（identical by state，IBS）和血缘同源（Identity By Descent，IBD），通常IBD无法直接观察，但IBS可以通过两个个体基因型算出（如下图），再根据IBS以及等位基因频率的分布推断IBD。

 

 

# 亲缘关系检测

$plink --noweb --bfile ${project}.raw.mark --genome

# 输出文件保存在plink.genome文件中，提取亲缘关系异常的样本

sed 's/^\s\+//' plink.genome | sed 's/\s\+/\t/g' | awk -v dst=0.85 'NR>2 {if($12 > dst) {print $1,$2; print $3,$4}}' | sort | uniq >genome_removelist.txt

# 删除亲缘关系异常个体（不建议删除）

$plink --noweb --bfile ${project}.raw.mark --remove genome_removelist.txt --make-bed --out ${project}.raw.mark2

 

 

3、哈温平衡过滤 

 

哈迪-温伯格（Hardy-Weinberg）法则是群体遗传中最重要的原理，提出在一个不发生突变、迁移和选择的无限大的随机交配的群体中（理想状态下），基因频率和基因型频率将逐代保持不变。一对等位基因的3种基因型分布比例符合以下规律：

 

(p + q)^2 = 1 等价于 p^2 + 2pq + q^2 = 1

注：p和q分别表示两个等位基因频率，且p + q = 1。

$plink --noweb --bfile ${project}.raw.mark --hwe 1e-10 --hwe-all --make-bed --out ${project}.filter.haw

 

 

4、最小等位基因频率过滤 

 

最小等位基因频率（MAF）通常是指在给定人群中的不常见的等位基因发生频率。

MAF如果非常小，比如低于0.02，那么意味着大部分位点都是相同的基因型，这些位点贡献的信息非常少，增加假阳性；更有甚者MAF为0，即所有位点只有一种基因型，这些位点没有贡献信息，放在计算中增加计算量，没有意义，所以要根据MAF进行过滤。

 

# 最小等位基因频率过滤（这里MAF阈值设为0.05）

$plink --noweb --bfile ${project}.raw.mark --maf 0.05 --allow-no-sex --make-bed --out ${project}.filter.maf

 

 

5、群体分层 

 

群体分层（Population stratification）：是最常见的差异来源，指的是case/control组的样本来自于不同的祖先群体，其分型结果自然是有差异的。

不同群体SNP频率不一样，导致后面做关联分析的时候可能出现假阳性位点（不一定是显著信号位点与该表型有关，可能是与群体SNP频率差异有关），因此我们需要在关联分析前对群体分层校正。

 

# 主成分分析

$plink --noweb --bfile ${project}.raw.mark --pca 10 --out pca

# 提取离群样本

根据主成分分析结果，绘图展示，确定离群样本，写入pca_removelist.txt文件

# 删除离群个体（可选）

$plink --noweb --bfile ${project}.raw.mark --remove pca_removelist.txt --make-bed --out ${project}.filter.pc

 

 

6、杂合性过滤 

 

杂合性是指某一个位点上含有一对及其以上的不同的等位基因。包括同系合性和同种5408性。群体遗传多态性的均匀度的度量常采用杂合度作为参数。杂合性是在同源染色体上的一个或多个位点上有不同等位基因存在的状态，是种群的基本属性之一。

 

# 连锁过滤（LD），得到不连锁的SNP

$plink --noweb --bfile ${project}.raw.mark --indep-pairwise 50 5 0.2 --out indepSNP

# 提取不连锁的SNP进行杂合性分析

$plink --noweb --bfile ${project}.raw.mark --extract indepSNP.prune.in --het --out hetSNP

# 提取杂合度较高的个体

sed 's/^\s\+//' hetSNP.het | sed 's/\s\+/\t/g' | awk -v f=0.35 'NR>1 {if(($5-$3)/$5 > f) {print $1,$2}}' >hetSNP_removelist.txt

# 删除杂合度高的个体（可选）

$plink --noweb --bfile ${project}.raw.mark --remove hetSNP_removelist.txt --make-bed --out ${project}.filter.het

 

 

以上就是本期分享的内容，下一期我们将讲解GWAS关联分析。