Plink提取tagSNPs
Plink提取tagSNPs
1)识别haplotype block
在提取tagSNPs之前需要先识别haplotype block,这一步的代码很多教程里都给了出来:
plink --file mydata --blocks
这里的mydata是由vcf文件经过plink转换来的,
如果是ped/map格式的文件(mydata.ped/mydata.map),就用**–file**
如果是bed/fam格式的文件(mydata.bed/mydata.fam),就用**–bfile**
当我自己运行的时候出现了小问题(好像是因为表型的问题),文件不能正确的输出,根据软件给出的的提示加入参数:no-pheno-req
修改后的代码为:
plink --file mydata --blocks no-pheno-req
这一步输出两个文件,分别是plink.blocks和plink.blocks.det
plink.blocks文件(上图)
plink.blocks.det文件(上图)
2)提取tagSNPs
网上的教程和官方文档也给出了相应的代码:
plink --bfile mydata --show-tags mysnps.txt
这里的mysnps.txt是从plink.blocks文件中提取出来的snp编号(一行一个):
但是运行结束后只得到了一个plink.tag文件,以为是哪里出错了,于是查阅官网,官方的解释:
plink.tags
that lists all the SNPs in the dataset that tag the SNPs in mysnps.txt (including the SNPs in the original file). A
message is also written to the LOG file that indicates how many new SNPs were added
If the option
–list-all
is also added, then an additional file is generated that gives some more details for each target SNP (i.e. each SNP listed in mysnps.txt, in the above example) regarding how many and which tags were set for it. The file is named
plink.tags.list
就是说如果用上述代码,最后只会输出一个plink.tag文件,如果再加上–list-all,则会再输出一个plink.tags.list文件,修改后代码为:
plink --file mydata --list-all --show-tags mysnps.txt
plink.tag文件和mysnp.txt文件类似,plink.tag,list文件结果如下: