GWAS分析-说人话(20)-单倍体关联分析
前言
来到这里,已经渐渐不是人话了...
估计这辈子我也没有想到,我把我最厌恶的东西写出来,
居然是最多人看的......
分析方法很多,这一小节吐槽官网的这一行代码,大神请点赞后直接略过~
plink --bfile mydata --hap-window 3 --hap-assoc
1.背景
我一般不交代“科学”背景,但是还得“装”一下。
什么是“Haplotype”?
就是你有一堆SNPs,"近朱者赤,近墨者黑",就是要“在一起”遗传的。
所以Haplotype分析就是要找到这一堆SNPs。
和GWAS分析中常提到的“连锁”有什么关系?
“连锁”是虚的!(概率相关的概念)
单体型是实的! (物理存在的附近)。
人话:
“连锁”法找女朋友全靠算命,星座风水;
“单倍型”法找女朋友,直接就是你邻座、邻居,青梅足马啊~。
2.今天的主角:单倍型分析(Haplotype testing)
道理很简单,我们都知道单纯靠一个SNP决定遗传是不靠谱的~
(别问我为什么知道,问就是100个位点,你做101个都未必有阳性实验结果,在其他人中也验证不出来)
本质上,就像看警匪片,要揪,就把整个犯罪团伙揪出来!
(非人话地说:单倍型分析在定位疾病和性状有关的基因方面具有更好的功效)
3.今天的吐槽大会:
plink --bfile mydata --hap-window 3 --hap-assoc
前面的:plink --bfile mydata 这里看过之前的都没有问题了,除非你没有看~
后面的这个 --hap-assoc 也就是单倍型关联分析嘛~
至于中间的这个“ --hap-window 3” 是什么鬼?
官网写的是这样的:
to form all 3-SNP haplotypes across the entire dataset (respecting chromosome boundaries, however). In this case the windows will be automatically named WIN1, WIN2, etc.
然后又是这样的:
This command can take a comma-delimited list of values, e.g.
--hap-window 1,2,3
to perform all single SNP tests (1-SNP haplotypes) as well as sliding windows of all 2-SNP and 3-SNP haplotypes.
每个单词都懂,但是就是看不懂啊~~~~~
放心,这锅不该英语老师来背~
要读懂这里,首先要了解一下:“sliding window”
因为冲动的后果就是这样:
ERROR: Problem with specification of haplotype sliding window
sliding:滑行的;window:窗 (简单、总体、多数情况来说,就是你规定的基因范围)
sliding window相当于就是:
我设定SNP个数一个一个玩!~~~
即:1就是一个一个玩,2就是两个两个玩,3就是三个三个玩!!!!
吐槽:
知道为什么程序员写的“说明”都“看不懂”了吧!以来就三个!~~~那么猛!
具体如下:(注意,不是随机组合抽取,是按固定个数一直往下走,“slide”~,直到走不下去)
6个SNPs,一个一个玩(可以玩6次,6个windows):
--hap-window 1的情况既结果,看SNP一列就是1个1个SNP来
6个SNPs,三个三个玩(可以玩4次,4个windows):
--hap-window 1的情况既结果,看SNP一列就是3个3个SNP来6个SNPs,五个五个玩(只能玩2次,2个windows):
--hap-window 1的情况既结果,看SNP一列就是5个5个SNP来
例子够了吧?还不点赞?
其实,吐槽完了,plink还是挺“银杏”的~
只要你在 --hap-window 中用逗号这样分开“1,2,3”( --hap-window 1,2,3),就会给出1,2,3各种情况时的结果,如下:
1个,2个,3个,各种情况的玩~所以,
plink --bfile mydata --hap-window 3 --hap-assoc,
你学废了吗?
后记:
目前文献,把基本关联做完了以后,就常规做一个LD检验,之后就继续做单体型与疾病关联分析,然而这个后续的单体型与疾病关联分析,并不以前面的LD所分辨出来的BLOCK为基础,而是他做了多少个SNP,就用多少个SNP一起做单体型与疾病关联分析。
实现这里的“用多少个SNP一起做单体型与疾病关联分析”,就是在这一小节之前使用extract提取你要的SNP文件(具体参考:),如6个SNPs(注意这个snps.txt,就是一个一列,没有列名的,想要提取snp名字的txt文件!不要想复杂了 ):
plink --noweb --bfile data --extract snps.txt --recode --make-bed --out 6snps
然后这个6snp文件就是你的mydata文件。
不写了,反正你们也不会点赞的~
参考文献:
1.提升SNP和单倍体认识逼格:http://doc.aporc.org/attach/Course001/Bioinformatics-3.pdf
2.大神讨论:https://www.dxy.cn/bbs/newweb/pc/post/20879971