LD衰减图

LD衰减距离指的是,当平均LD系数衰减到一定大小(最大值的一半/0.5以下)的时候,对应的物理距离。通常用LD衰减距离来描述LD衰减速度。LD衰减速度越快,即衰减距离越小,说明该群体遗传多样性越高;LD衰减速度越慢,通常驯化程度越高,选择强度越大,导致遗传多样性下降。
LD系数衰退速度会受到不同因素的影响而有所不同。常见的因素包括:

1)物种类型LD存在的本质是两个位点的连锁遗传导致的相关性。但这种相关性理论上会随着世代的增加、重组次数的增加而不断下降。所以,那些繁殖力强、时代间隔短的物种(例如,昆虫),其LD衰减的速度是非常快的。例如在家蚕和野蚕群体中,LD系数下降到最大值的1/2仅仅需要46bp和7bp的距离

2)群体类型相同物种的不同群体,由于其遗传背景不同,LD衰减速度也存在很大的差异。驯化选择,会导致群体遗传多样性下降,位点间的相关性(连锁程度)加强。所以,通常驯化程度越高,选择强度越大的群体,LD衰减速度是最慢的。例如,栽培稻比野生稻通常更大的LD衰减距离。类似的,自然选择、遗传漂变导致的群体遗传多样性下降,也会减慢LD衰减的速度。

3)在染色体的位置染色体不同区域的LD衰减距离而是不同的。通常着丝粒区更难重组,所以LD衰减更慢。而基因组上那些受选择的区域相比普通的区域,LD衰减速度也是更慢的。

下载软件:PopLDdecay

git clone https://github.com/BGI-shenzhen/PopLDdecay.git
chmod 755 configure
./configure

计算每个亚群的LD

~/PopLDdecay/bin/PopLDdecay -InVCF Filter.snp.vcf -SubPop 1-population.txt -OutStat p1.stat -MaxDist 1000
~/PopLDdecay/bin/PopLDdecay -InVCF Filter.snp.vcf -SubPop 2-population.txt -OutStat p2.stat -MaxDist 1000
~/PopLDdecay/bin/PopLDdecay -InVCF Filter.snp.vcf -SubPop 3-population.txt -OutStat p3.stat -MaxDist 1000
~/PopLDdecay/bin/PopLDdecay -InVCF Filter.snp.vcf -OutStat all.stat -MaxDist 2000
#-SubPop 群体ID信息,一行一个ID
#-OutStat 输出文件
#-MaxDist 两个SNP间的最大距离,默认是300kb,最大不超过5mb

创建draw.list

#上一步生成的文件的绝对路径 以及去掉后缀的文件名
/home/ug0797/data/test-fst/test2/p1.stat.gz p1
/home/ug0797/data/test-fst/test2/p2.stat.gz p2
/home/ug0797/data/test-fst/test2/p3.stat.gz p3

画图

~/PopLDdecay/bin/Plot_MultiPop.pl -inList draw.list -output draw -bin1 5000 -bin2 50000
# -bin1 -bin2设置的大一点,曲线会更平滑!
image.png

发现不止1000kb

用ggplot画

用poplddecay画图时会返回有以下内容的文件,前两列是画图所需的的数据


image.png

其中列名#Dist、Mean_r^2在R语言中无法识别,所以首先要改列名

修改列名(draw_p*就是poplddecay画图后生成的文件)

colnames(draw_p1)[1]<-'Dist'
colnames(draw_p1)[2]<-'Mean_r.2'
colnames(draw_p2)[1]<-'Dist'
colnames(draw_p2)[2]<-'Mean_r.2'
colnames(draw_p3)[1]<-'Dist'
colnames(draw_p3)[2]<-'Mean_r.2'

删除最后一行数据(只画到1000kb)

draw_p1 <- draw_p1[-22,]
draw_p2 <- draw_p2[-22,]
draw_p3 <- draw_p3[-22,]

添加分组信息

Group=("G1")
data1=data.frame(draw_p1,Group)
View(data1)
Group=("G2")
data2=data.frame(draw_p2,Group)
View(data2)
Group=("G3")
data3=data.frame(draw_p3,Group)

将三组数据合并

library(dplyr)
total_data<-dplyr::bind_rows(data1,data2,data3)

ggplot画图

p <- ggplot(data=total_data,aes(x=Dist/1000,y=Mean_r.2,colour=Group))+geom_line()+labs(y=expression(r^{2}),title="LD decay")+
xlab("Distance")+theme_bw()+xlim(0,1000)+ylim(0,0.6)#xlim,ylim设置轴范围
p<- p+theme(plot.title = element_text(hjust = 0.5))#使主标题居中
p

image.png

基础知识参考文章:
https://links.jianshu.com/go?to=http%3A%2F%2Fwww.omicshare.com%2Fforum%2Fthread-878-1-1.html

©著作权归作者所有,转载或内容合作请联系作者
平台声明:文章内容(如有图片或视频亦包括在内)由作者上传并发布,文章内容仅代表作者本人观点,简书系信息发布平台,仅提供信息存储服务。

推荐阅读更多精彩内容