EMMAX是一个统计测试软件。根据样本结构进行大规模的人类或者模式生物的体的关联映射。可以通过EMMA算法获得计算效率,EMMAX还利用了每个基因座仅解释一小部分复杂性状的这一事实,避免使用混合模型进行关联映射和重复的方差分量估计过程而导致计算时间的显着增加。
安装包下载网址
https://genome.sph.umich.edu/wiki/EMMAX
Tar解压:
tar xvf emmax-beta-07Mar2010.tar
这个安装超级简单,解压后直接加绝对路径引用就OK啦,当然嫌绝对路径麻烦的小伙伴可以加个环境变量,按个人习惯来就好
之前我在格式那篇教程里曾经介绍过.tped、.tfam、 *.pheno格式,这三个格式都是作为EMMAX的输入格式使用的
所以在介绍EMMAX用法前需要给大家介绍一下怎么把filter和beagle处理后的vcf文件转化为EMMAX的输入文件
方法一:
当然是自己写脚本处理啦,知道文件格式,相当于自己处理数据,而且我总觉得自己写的靠谱一点,我也倾向于自己写
方法二:
踩在前辈的肩膀上,也有现成的转化软件—-plink—-
安装包下载网址
http://www.cog-genomics.org/plink/1.9/index
Unzip解压
unzip -d plink_linux_x86_64.zip
PS(plink_linux_x86_64(不知道为什么直接用unzip解压到目录下产生的文件不能运行所以加一个-d参数解压到目标文件夹下)。) -version查看plink是否安装成功
/data/home/mjchen/app/plink/plink_linux_x86_64/plink -version
PLINK v1.90b5.4 64-bit (10 Apr 2018)
废话不多,直接上代码,我不怎么用plink转,所以只实验过一次,其他小伙伴们要是有什么其它更方便的转化方法请多指教
1.vcf转ped/map:
/data/home/mjchen/app/plink/plink_linux_x86_64/plink --vcf E19.vcf --allow-extra-chr --recode --out E19
2.准备输入的基因型文件
使用PLINK软件将基因型文件(bed或ped格式)转换为tped / tfam格式
plink --bfile [bed_prefix] (or --file [ped_prefix]) --recode12 --output-missing-genotype 0 --transpose --out [tped_prefix]
example
/data/home/mjchen/app/plink/plink_linux_x86_64/plink --file F2-Mo17vsB73_beagle --recode12 --output-missing-genotype 0 --transpose --out F2-Mo17vsB73_beagle
此处需要注意的是—file输入的文件不要包含.ped的扩展名,不然命令运行会报错,另—out也一样只要输入需要的文件名即可输出文件会直接带上.tped的扩展名
—bfile {prefix}:输入bed文件
—file:输入ped或者map文件
—recode:生成文本文件集。类型有很多:<01 |="" 12=""> <23 |="" a="" |="" a-transpose="" |="" ad="" |="" beagle="" |="" beagle-nomap="" |="" bimbam="" |="" bimbam-1chr="" |="" compound-genotypes="" |="" fastphase="" |="" fastphase-1chr="" |="" hv="" |="" hv-1chr="" |="" lgen="" |="" lgen-ref="" |="" list="" |="" oxford="" |="" rlist="" |="" structure="" |="" transpose="" |="" vcf="" |="" vcf-fid="" |="" vcf-iid="">
‘12’修饰符导致A1(通常是次要)等位基因编码为’1’,A2等位基因编码为’2’,而’01’则映射A1→0和A2→1。 (PLINK迫使你将’01’与 - {output-} missing-genotype结合起来)23>01>
准备输入的表型文件 按照.tfam文件的相同顺序重新设置表型文件的格式。 表型文件在每一行有三个条目,FAMID,INDID和表型值。 缺失的表型值应表示为“NA”。 表型文件的示例行。 (制表符或空格分隔)
创建基于标记的亲属关系矩阵 使用emmax-kin创建亲属关系矩阵(IBS或BN,BN是首选)。 确保.tped和.tfam文件存在相同的前缀
IBS matrix % emmax-kin-intel64 -v -s -d 10 [tped_prefix] (will generate [tped_prefix].aIBS.kinf) BN (Balding-Nichols) matrix % emmax-kin-intel64 -v -d 10 [tped_prefix] (will generate [tped_prefix].aBN.kinf)
好了,输入文件都准备好了,可以运行EMMAX了
% emmax -v -d 10 -t [tped_prefix] -p [pheno_file] -k [kin_file] -o [out_prefix]
运行结束后将生成以下文件:
[out_prefix].reml : REML 输出文件有六列
Log-likelihood with variance component
Log-likelihood without variance component
\delta = \sigma_e^2 / \sigma_g^2 (Ratio between variance parameters)
\sigma_g^2 (genetic variance parameter)
sigma_e^2 (residual variance parameter)
The pseudo-heritability estimates . (Explained variance by the kinship matrix)
[out_prefix].ps : 每行包括以下内容
SNP ID
Beta (1 is effect allele)
SE(beta)
p-value.
最后呢,大家可以根据p-value做曼哈顿图进行定位了。
联系客服