性别错误 (Sex discrepancy)：有时数据录入会有错误，我们应基于X染色体的杂合性，确定是否有性别录入错误的个体。PLINK默认对于男性，X染色体的纯合性（homozygosity）估计值应高于0.8，女性应低于0.2，但实际情况更为复杂，本文简要对性别检查 与 阈值选择进行讲解。
本文使用数据来源：https://onlinelibrary.wiley.com/doi/full/10.1002/mpr.1608
PLINK中提供了check-sex选项，通过计算X染色体的纯合性来判断样本性别是否有错误：—check-sex [女性阈值] [男性阈值] 默认的参数为 0.2 0.8，但该标准过于严格，并不适用与实际情况，下面会讲解。首先介绍基本操作：
#首先进行LD-pruning
plink \
--bfile HapMap_3_r3_1 \
--indep-pairwise 50 5 0.2 \
--out HapMap_3_r3_1
#估计X染色体的纯合性
plink \
--bfile HapMap_3_r3_1 \
--extract HapMap_3_r3_1.prune.in \
--check-sex \
--out HapMap_3_r3_1
注意：
- 要提前确认X染色体上PAR区域已经被分离 （可以使用—split-x 选项来进行分离）
- 默认的阈值是，男性>0.8， 女性<0.2
log文件如下：
PLINK v1.90b6.21 64-bit (19 Oct 2020) www.cog-genomics.org/plink/1.9/
(C) 2005-2020 Shaun Purcell, Christopher Chang GNU General Public License v3
Logging to HapMap_3_r3_1.log.
Options in effect:
--bfile HapMap_3_r3_1
--check-sex
--extract HapMap_3_r3_1.prune.in
--out HapMap_3_r3_1
191875 MB RAM detected; reserving 95937 MB for main workspace.
Allocated 3037 MB successfully, after larger attempt(s) failed.
1457897 variants loaded from .bim file.
165 people (80 males, 85 females) loaded from .fam.
112 phenotype values loaded from .fam.
--extract: 166240 variants remaining.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 112 founders and 53 nonfounders present.
Calculating allele frequencies... done.
Warning: 35 het. haploid genotypes present (see HapMap_3_r3_1.hh ); many
commands treat these as missing.
Total genotyping rate is 0.996691.
166240 variants and 165 people pass filters and QC.
Among remaining phenotypes, 56 are cases and 56 are controls. (53 phenotypes
are missing.)
--check-sex: 3993 Xchr and 0 Ychr variant(s) scanned, 1 problem detected.
Report written to HapMap_3_r3_1.sexcheck .
–check-sex: 3993 Xchr and 0 Ychr variant(s) scanned, 1 problem detected.从日志文件中可以看出我们检测出了1个不一致的样本。上述代码运行结果如下：
head HapMap_3_r3_1.sexcheck
FID IID PEDSEX SNPSEX STATUS F
1328 NA06989 2 2 OK 0.02151
1377 NA11891 1 1 OK 1
1349 NA11843 1 1 OK 1
1330 NA12341 2 2 OK 0.009813
1444 NA12739 1 1 OK 1
1344 NA10850 2 2 OK 0.0533
1328 NA06984 1 1 OK 0.9989
1463 NA12877 1 1 OK 1
1418 NA12275 2 2 OK -0.0747
前两列为FID与IID，三四列为PED文件中的性别，与通过X染色体估计的性别，第五列表示是否一致，不一致时会表示为PROBLEM，第六列是估计的F系数。可以提取出不一致的样本：
cat HapMap_3_r3_1.sexcheck | grep PROBLEM
1349 NA10854 2 1 PROBLEM 0.9933
#将有问题的样本FID与IID输出到sex_error_to_remove.id文件中
cat HapMap_3_r3_1.sexcheck | grep PROBLEM | awk '{print $1,$2}'> sex_error_to_remove.id
后续PLINK分析可以使用 —remove sex_error_to_remove.id 来剔除性别错误的样本.以上结果的F系数分布如图所示：
实际情况：上面的例子是一种很理想的分布，可以看到样本被清晰地分开，但实际情况要复杂很多，例如在千人基因组项目Phase1的样本中，一些女性样本的F系数也大于0.8，即使进行LD-pruning后，女性样本最大的F系数也超过了0.66。所以一般来说，我们会进行两次—check-sex，第一次用于观察分布，在看到一个右侧紧凑的峰（男性），与左侧平缓的峰后（女性），根据其分界线来决定判定阈值，（真实的数据的分布大多如下图所示）
这里常用的阈值为 -check-sex 0.7 0.8
#如果使用真实数据，通过目视决定阈值后，更改参数，重新判定性别
plink \
--bfile HapMap_3_r3_1 \
--extract HapMap_3_r3_1.prune.in \
--check-sex 0.7 0.8 \
--out HapMap_3_r3_1
参考：
GWASLab - GWAS实验室 