References
http://zzz.bwh.harvard.edu/plink/tutorial.shtml
https://www.cog-genomics.org/plink/1.9/
软件的输入:
.map文件:SNP(marker)的名字和位置信息
image.png空格分割
- 四列分别是染色体(1-22,X,Y,XY,MT,0表示没有定位到染色体上的marker/SNP),marker/SNP的ID,摩尔根遗传距离(如果这列没有的话使用--map3标记来排除,e.g., plink --file mydata --map3),碱基的位置
- 一般来说,常规的关联分析不需要遗传距离,可以设置为0
- PLINK的输出中X,Y,XY,MT会变成23,24,25,26
- X表示X染色体中不与Y重叠的部分(所以这区域的SNP在计数的时候,纯和的情况下女性中计2个,男性中只计1个;杂合的情况下,女性中两个allele各记一个,男性中不计数,因为这种情况实际上是不存在的,我也不知道实际的genotyping是不是会出现这种数据和为什么会出现)。
- Y表示Y染色体中不与X重叠的部分(所以计数规则是女性全部当作缺失值不计数,男性中杂合也当作缺失值,纯和只计一个)
- XY表示XY重叠的区域,可以称作伪常染色体区域,在计数的时候男女没有差别,所以跟常染色体一样全部正常计数就可以
- MT表示线粒体,由于是单倍体,而且没有男女之分,所以杂合情况下不计,纯和情况下的allele只计1个
.ped文件: 表型+基因型文件
image.png空格分隔
- 前六列是必须的,分别是家庭ID、个人ID、父亲ID、母亲ID、性别(1男2女)、表型(分类和定量变量都可以)。
表型如果是分类变量,默认的编码是-9、0代表缺失,1是control、2是case,如果是0/1编码control/case,用--1指定(e.g., plink --file mydata --1),缺失值的编码可以用--missing-phenotype指定(e.g., plink --file mydata --missing-phenotype 99) - 第七列以及往后是基因型,成对出现(e.g.,genotype1-allele1,genotype1-allele2,genetype2-allele1,genotype2-allele2 ...),可以是ATGC或者1234或者其它任何字符,0是默认的缺失值,可以用--missing-genotype指定(e.g., plink --file mydata --missing-genotype N),输出的时候可以指定缺失值的格式--output-missing-phenotype NA,--output-missing-genotype NA
- 输入文件中如果前六列有缺失,也无妨,输入的时候通过--no-fid、--no-parentes、--no-sex、--no-pheno分别指定缺失的列就可以
.phe文件:表型文件(表型可以是分类变量也可以是连续变量,可以放到一个文件也可以分开放)
image.pngimage.png
一般操作流程:
- 把PED文件转换为压缩格式的文件(会生成三个文件.bed .bim .fam,输入需要.ped和.map)
plink --file [prefix of .ped/.map] --out outbfile
# or
plink --map [a.map] --ped [b.ped] --out outbfile
后面的所有操作都用压缩处理后的文件作为输入,用--bfile指定
- 各种统计量
- 缺失率
plink --bfile [prefix of .bed/bim/fam] --missing --out outfile
# or
plink --bed [x.bed] --bim <x.bim> --fam <x.fam> --missing --out outfile
输出中.imiss是个人的基因型缺失数目 (N_MISS)和缺失率 (F_MISS, i.e., genotyping rate of individual),.lmiss是对于每个基因型来说,缺失的人数和比例
- allele频率
# 输出minor allele frequency (MAF)
plink --bfile [prefix] --freq --out outfile
# 如果只对某一个SNP感兴趣(e.g., rs1891905)
plink --bfile [prefix] --freq --snp rs1891905 --outfile
- LD
待补充
- 简单的关联分析(貌似大家一般不用PLINK来做,一般都转到rvtest, RAREMETALWORKER, seqMeta等去做)
plink --bfile [prefix] --assoc --out outfile
格式转换:
plink --bfile [bfile prefix] --recode vcf --out [vcf file prefix]
VCF生成后可以转到rvtest或者raremetalworker去做association test
网友评论