育种数据分析之放飞自我分享 http://blog.sciencenet.cn/u/yijiaobai 关注:生物统计,数量遗传,混合线性模型,生物信息,R,Perl,Python,GWAS,GS相关方法,文章及代码

博文

GBS hapmap 格式 转化为Plink格式:tassel vcftools方法汇总

已有 725 次阅读 2019-7-31 20:22 |个人分类:农学统计|系统分类:科研笔记

1. 参考资料

https://zhuanlan.zhihu.com/p/38590403

环境:linux系统

2. 安装软件

首先,要安装anaconda或者miniconda,然后使用conda install进行软件安装, 安装conda方法:
https://docs.anaconda.com/anaconda/install/linux/

2.1 安装tassel

git clone https://bitbucket.org/tasseladmin/tassel-5-standalone.git

将tassel文件中的文件, copy到home下的bin文件中, 不用设置路径了。

2.2 安装vcftools

https://anaconda.org/bioconda/vcftools

conda install -c bioconda vcftools

2.3 安装R语言

https://anaconda.org/r/r

conda install -c r r

3. 文件格式

3.1  hapmap格式:genotype.hmp.txt

行头:

rs#    alleles    chrom     pos     strand    assembly#    center    protLSID    assayLSID    panelLSID    QCcode    sample1 sample2 ...

内容:

rs#     alleles chrom   pos     strand  assembly#       center  protLSID        assayLSID       panelLSID       QCcode  Sample_YCX334/12Sample_ya>
1:1151  C       1       1151    +       NA      NA      NA      NA      NA      NA      N       N       N       C       C       C       N       C>
1:1203  T/C     1       1203    +       NA      NA      NA      NA      NA      NA      T       N       T       N       T       T       T       T>
1:1249  A/C     1       1249    +       NA      NA      NA      NA      NA      NA      A       N       A       N       A       A       A       A>
1:1266  G/A     1       1266    +       NA      NA      NA      NA      NA      NA      G       N       G       G       G       G       G       N>
1:1277  T/C     1       1277    +       NA      NA      NA      NA      NA      NA      T       T       T       T       T       T       T       N>
1:1325  C/T     1       1325    +       NA      NA      NA      NA      NA      NA      C       N       N       N       N       N       C       N>
1:1335  G/T     1       1335    +       NA      NA      NA      NA      NA      NA      G       G       G       G       G       G       G       G>
1:1362  G/A     1       1362    +       NA      NA      NA      NA      NA      NA      G       G       G       G       G       G       G       G>

3.2 plink格式:pedmap

plink格式是基因组选择中经常用到的文件类型, plink软件功能强大,运行速度快。

3.2.1 .map格式

格式说明链接: http://zzz.bwh.harvard.edu/plink/data.shtml#map

map格式的文件, 主要是图谱文件信息, 主要包括染色体名称, 所在的染色体和所在染色体的坐标.

1, map文件没有行头

2, map文件包括四列: 染色体, SNP名称, SNP位置,  碱基对坐标

  • 染色体编号为数字, 未知为0
  • SNP名称为字符或数字, 如果不重要, 可以从1编号, 注意要和bed文件SNP列一一对应
  • 染色体的摩尔未知(可选项, 可以用0)
  • SNP物理坐标

3, 如果只有SNP名称, 可以手动构建map文件, 第二列为SNP名称, 其它三列为0即可.

Example:

1 snp1 0 1
1 snp2 0 2
1 snp3 0 3
  • 这里有3个SNP, 分别名为snp1, snp3, snp3 (第二列)
  • 这三个SNP在第一个染色体上 (第一列)
  • 第三列为0
  • 第四列为SNP所在染色体的坐标
3.2.2 .ped格式

格式说明链接:http://zzz.bwh.harvard.edu/plink/data.shtml#ped

bed格式的文件, 主要包括SNP的信息, 包括个体ID, 系谱信息, 表型和SNP的分型信息.

1, 数据没有行头, 空格或者tab隔开的文件
2, 必须要有六列, 包括系谱信息, 表型信息

  • 第一列: Family ID # 如果没有, 可以用个体ID代替
  • 第二列: Individual ID # 个体ID编号
  • 第三列:  Paternal ID # 父本编号
  • 第四列:   Maternal ID # 母本编号
  • 第五列:   Sex (1=male; 2=female; other=unknown) # 性别, 如果未知, 用0表示
  • 第六列:   Phenotype # 表型数据, 如果未知, 用0表示
  • 第七列以后: 为SNP分型数据, 可以是AT CG或11 12, 或者A T C G或1 1 2 2

3, 上面六列, 必须要有, 如果没有相关数据, 用0表示.

Example:

1 1 0 0 1  0  G G  2 2  C C
1 2 0 0 2  0  A A  0 0  A C
1 3 1 2 1  2  0 0  1 2  A C
2 1 0 0 1  0  A A  2 2  0 0
2 2 0 0 2  2  A A  2 2  0 0
2 3 1 2 1  2  A A  2 2  A A
  • 数据包括两个家系 (第一列)
  • 每个家系有三个个体 (第二列)
  • 第三列父本编号
  • 第四列母本编号
  • 第五列性别
  • 第六列表型值
  • 第七列, 第八列为一个基因型
  • 第九列, 第十列为第二个基因型
  • 第十一列, 第十二列为第三个基因型

4. 测试数据

将下面文件保存为:hmp.txt
注意, 知乎中是.,但数据应该是#,下面这个代码是正确的。

rs#    alleles    chrom    pos    strand    assembly#    center    protLSID    assayLSID    panelLSID    QCcode    Box302.YX.2017.06.002    Box302.YX.2017.06.003    Box302.YX.2017.06.004    Box302.YX.2017.06.005    Box302.YX.2017.06.11    Box302.YX.2017.06.013    Box302.YX.2017.06.018    Box302.YX.2017.06.019    Box302.YX.2017.06.020    Box302.YX.2017.06.21
10000235    A/C    6    80388997    +    NA    NA    NA    NA    NA    NA    AC    AA    AC    AA    AA    AA    AA    AA    AC    AC
10000345    G/A    8    17865885    +    NA    NA    NA    NA    NA    NA    GG    GG    AA    GG    GG    GG    GG    AA    AA    GG
10004575    G/T    7    32792755    +    NA    NA    NA    NA    NA    NA    GG    GG    GG    GG    GG    GT    GG    GG    GG    GG
10006974    T/C    1    14418789    +    NA    NA    NA    NA    NA    NA    CT    CC    TT    CT    CT    CT    CT    TT    TT    CC
10006986    A/G    2    76416246    +    NA    NA    NA    NA    NA    NA    AA    AA    AA    AA    AA    AG    AG    GG    GG    AA
10007074    G/C    14    105043500    +    NA    NA    NA    NA    NA    NA    GG    GC    GC    GC    GC    GC    GC    GC    GG    GC
10007097    C/A    15    118863849    +    NA    NA    NA    NA    NA    NA    CC    CC    CC    CC    CC    CC    AC    CC    CC    AC
10007113    G/A    9    127852320    +    NA    NA    NA    NA    NA    NA    GG    GG    GG    AG    GG    GG    AG    AA    AA    GG
10007117    C/T    2    150034851    +    NA    NA    NA    NA    NA    NA    CC    CT    CC    CT    CT    CT    TT    CC    CT    CT
10007153    A/C    6    145256960    +    NA    NA    NA    NA    NA    NA    AA    AC    AA    AC    AA    AA    AA    AC    AA    AA
10007668    G/A    10    19391507    +    NA    NA    NA    NA    NA    NA    GG    GG    GG    AG    AG    AA    GG    GG    GG    GG
12784072    G/A    17    55229968    +    NA    NA    NA    NA    NA    NA    AG    AG    AG    GG    AG    GG    GG    AG    AA    GG
12784632    T/C    6    28132948    +    NA    NA    NA    NA    NA    NA    TT    TT    TT    TT    TT    TT    TT    CT    TT    TT
12784633    T/C    6    28132948    +    NA    NA    NA    NA    NA    NA    TT    TT    TT    TT    TT    TT    TT    CT    TT    TT
12784634    T/C    6    28132948    +    NA    NA    NA    NA    NA    NA    TT    TT    TT    TT    TT    TT    TT    CT    TT    TT
12784635    T/C    6    28132948    +    NA    NA    NA    NA    NA    NA    TT    TT    TT    TT    TT    TT    TT    CT    TT    TT
12784636    G/A    1    160773437    +    NA    NA    NA    NA    NA    NA    AG    GG    GG    GG    AG    AG    GG    AG    AG    GG
12784637    G/A    1    160773437    +    NA    NA    NA    NA    NA    NA    AG    GG    GG    GG    AG    AG    GG    AG    AG    GG
12784638    G/A    1    160773437    +    NA    NA    NA    NA    NA    NA    AG    GG    GG    GG    AG    AG    GG    AG    AG    GG
12784639    G/A    1    160773437    +    NA    NA    NA    NA    NA    NA    AG    GG    GG    GG    AG    AG    GG    AG    AG    GG

5. 处理代码

5.1 排序

run_pipeline.pl -SortGenotypeFilePlugin -inputFile hmp.txt -outputFile test.sort.hmp.txt -fileType Hapmap

生成一个test.sort.hmp.txt

5.2 生成 vcf4.0 格式的文件

run_pipeline.pl -fork1 -h test.sort.hmp.txt -export -exportType VCF

生成一个test.vcf文件

5.3 使用vcftools生成plink文件

vcftools --vcf test.vcf --plink --out  tassel.test.vcf2plink

生成tassel.test.vcf2plink文件

5.4 使用plink将vcf文件, 变为bed文件

plink --file tassel.test.vcf2plink --make-bed --out tassel.test.vcf2plink

在这里插入图片描述

5.5 使用plink将bed文件转化为map和ped文件

plink --bfile tassel.test.vcf2plink --recode --out result

结果生成:result.pedresult.map文件:

(base) [dengfei@ny01 hmp_2_plink]$ cat result.ped 
Box302.YX.2017.06.002 Box302.YX.2017.06.002 0 0 0 -9 C T A G A G A G A G A A C C T T T T T T T T C A A A G G G G G G G G G G C C A G
Box302.YX.2017.06.003 Box302.YX.2017.06.003 0 0 0 -9 C C G G G G G G G G A A T C T T T T T T T T A A C A G G G G G G G G C G C C A G
Box302.YX.2017.06.004 Box302.YX.2017.06.004 0 0 0 -9 T T G G G G G G G G A A C C T T T T T T T T C A A A G G A A G G G G C G C C A G
Box302.YX.2017.06.005 Box302.YX.2017.06.005 0 0 0 -9 C T G G G G G G G G A A T C T T T T T T T T A A C A G G G G A G A G C G C C G G
Box302.YX.2017.06.11 Box302.YX.2017.06.11 0 0 0 -9 C T A G A G A G A G A A T C T T T T T T T T A A A A G G G G G G A G C G C C A G
Box302.YX.2017.06.013 Box302.YX.2017.06.013 0 0 0 -9 C T A G A G A G A G G A T C T T T T T T T T A A A A T G G G G G A A C G C C G G
Box302.YX.2017.06.018 Box302.YX.2017.06.018 0 0 0 -9 C T G G G G G G G G G A T T T T T T T T T T A A A A G G G G A G G G C G A C G G
Box302.YX.2017.06.019 Box302.YX.2017.06.019 0 0 0 -9 T T A G A G A G A G G G C C C T C T C T C T A A C A G G A A A A G G C G C C A G
Box302.YX.2017.06.020 Box302.YX.2017.06.020 0 0 0 -9 T T A G A G A G A G G G T C T T T T T T T T C A A A G G A A A A G G G G C C A A
Box302.YX.2017.06.21 Box302.YX.2017.06.21 0 0 0 -9 C C G G G G G G G G A A T C T T T T T T T T C A A A G G G G G G G G C G A C G G
(base) [dengfei@ny01 hmp_2_plink]$ cat result.map 
1    10006974    0    14418789
1    12784636    0    160773437
1    12784637    0    160773437
1    12784638    0    160773437
1    12784639    0    160773437
2    10006986    0    76416246
2    10007117    0    150034851
6    12784632    0    28132948
6    12784633    0    28132948
6    12784634    0    28132948
6    12784635    0    28132948
6    10000235    0    80388997
6    10007153    0    145256960
7    10004575    0    32792755
8    10000345    0    17865885
9    10007113    0    127852320
10    10007668    0    19391507

公众号:育种数据分析之放飞自我



http://blog.sciencenet.cn/blog-2577109-1191913.html

上一篇:GGE Biplot 双标图如何看?
下一篇:因为一个序言爱上一本书《无穷的开始--世界进步的本源》

0

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备14006957 )

GMT+8, 2019-10-21 21:25

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部