1.质控和LD过滤
1.1 质控:删除信息缺失严重的个体和位点
plink2 --vcf zy_bsz-bi.vcf --mind 0.10 --maf 0.05 --geno 0.05 --hwe 0.0001 \
--out zy_bsz_fi --recode vcf-iid \
--allow-extra-chr --set-missing-var-ids @:# \
--keep-allele-order
--mind 0.10 去除基因型丢失率大于10%的个体样本：
--geno 0.05 大于95%的个体都具有的变异位点才保留，其他去除：
--maf 0.01 次等位基因频率，频率较低的第二等位基因的频率（防止假阳性）。
--hwe 0.0001 保留符合Hardy—Weinbery 的变异位点
1.2 LD过滤:降低连锁程度
plink2 --vcf zy_bsz_fi.vcf \
--indep-pairwise 50 10 0.2 \ #滑窗口50 bp、步长10bp 删除LD大于0.2的SNP对中一个,可换
--out tmp.ld --allow-extra-chr \
--set-missing-var-ids @:#
plink2 --vcf zy_bsz_fi.vcf \
--extract tmp.ld.prune.in \
--out all.LDfilter --recode vcf-iid \
--keep-allele-order --allow-extra-chr \
--set-missing-var-ids @:#
2.用plink 进行pca分析
plink2 --vcf LD_filter_56_0.2.vcf --pca 10 --double-id --allow-extra-chr -out 56_02pca
###--pca 10表示降维为10维
### 个体少时，强制执行 --bad-freqs
plink2 --vcf zy_bsz.vcf --pca 3 --double-id --allow-extra-chr -out bsz_zy_pca --bad-freqs
2.1. 生成两个文件
(1) XX.eigenval 特征值（比如10维会生成10个数字：A1、A2、A3......A10)。
第一主成分的解释变量 = A1 / (A1+A2+A3+.....A10)
(2) XX.eigenvec 用于画图的文件
XX.eigenvec 需手动加上第一行，name group PC1 PC2 PC3......（可在notepad+里将短空格转换为/t（即tab键），然后将其copy到csv文件里。）
4.R作图
####PCA
setwd("E:/yangling02/East_Asia/PCA/56_pca_02/")
df=read.csv("56_02_pca.csv",header=T)
library(ggplot2)
library(scatterplot3d)
color <- c("orange","palevioletred3","pink","black","blue","grey","purple")
df$group=factor(df$group)
colors <- color[as.numeric(df$group)]
p=ggplot(data=df,mapping=aes(PC1,PC3,color=colors))
p+geom_point(color=colors)#+geom_text(label=paste(df$name),color=colors)
p=df
pp=scatterplot3d(p[,3:5], color=colors, angle=65, pch=16,type="h",axis = TRUE,label.tick.marks = TRUE,xlab="PC1(22.2%)",ylab="PC2(15.4%)",zlab="PC3(11.1%)")
legend("right", legend=levels(p$group), col=c("orange","palevioletred3","pink","black","blue","grey","purple"), pch=16, xpd=T, horiz=F)#
text(pp$xyz.convert(p[,3:5]),labels=p$name,cex=0.8,col = colors,adj=c(1,-1,font=2))#给散点加文字，adj让文字在点上
3 GCTA法
