甲基化数据分析QTL定位重测序下游分析

plink 进行PCA分析

2020-09-23  本文已影响0人  斩毛毛

当我们进行群体分析时,获得vcf文件后,可以根据变异位点对这些样本进行PCA分析,现简单介绍

1、软件安装

conda install plink

2、简单操作

本次使用train.vcf.gz 作为例子

plink --vcf F_M_trans.recode.vcf.gz --recode --out testacc --const-fid --allow-extra-chr
  
 # --vcf vcf 或者vcf.gz
 # --recode 输出格式
 # --out 输入前缀
 # --const-fid  添加群体信息
 # --allow-extra-chr 允许非标准染色体编号

得到3个以.map, .nosex, .ped结尾的文件。

plink --allow-extra-chr --file testacc --noweb --make-bed --out testacc
 
# --file .ped + .map 文件前缀
# --make-bed 建立一个新的二进制文件

得到2个以.bim, .bed结尾的文件。


plink --allow-extra-chr --threads 20 -bfile testacc --pca 20 --out testacc

 
# --threads 线程数
# --pca 主成分

得到2个以.eigenval, .eigenvec结尾的文件;其中.eigenval代表每个PCA所占的比重,另外一个记录特征向量,用于坐标轴的绘制

pop     name    pca1    pca2
B201    B201-1  -0.621307       -0.0873675
B201    B201-2  -0.374507       -0.00563367
B201    B201-3  -0.422489       -0.00238888
B202    B202-1  0.228835        -0.172511
B202    B202-2  0.21213 -0.142105
B202    B202-3  0.247804        -0.202722
B203    B203-1  0.216946        -0.141285
B203    B203-3  0.204785        -0.0981921
B204    B204-1  0.136368        0.930181
B204    B204-2  0.169868        -0.0802494
(ggplot2)
pca <- read.table("pca_info.txt",sep = "\t",header = T)
ggplot(pca, aes(x=pca1,y=pca2)) +
        geom_point(aes(color=pop, shape=pop),size=1.5)+
        labs(x="PC1",y="PC2")+theme_bw()+theme(legend.title = element_blank())

结果


** 若想分析部分样本,则可以使用--remove参数,后接一个文件,其格式为: 第一列:群体编号, 第二列:样本名称,在这个例子中

echo '0\tSP23' > remove.txt
plink --remove remove.txt --allow-extra-chr -bfile testacc --pca 20 --out testacc_dele 
上一篇下一篇

猜你喜欢

热点阅读