RIdeogram 是用来展示 染色体组型 （idiogram）的一个R语言包，比如展示snp的密度分布；展示某类基因在染色体上的分布等等。
具体介绍可以参考
https://www.plob.org/article/16335.html
发表的论文是
RIdeogram: drawing SVG graphics to visualize and map genome-wide data on the idiograms
PeerJ 2020
附件中还提供了好几份数据集合代码，非常好的学习素材。
我在看这个包的帮助文档的时候也发现了展示2个或者3个基因组共线性分析的图。
需要准备的数据是两个数据框，
第一个是基因组的基本信息，包括
> head(karyotype_dual_comparison)
Chr Start End fill species size color
1 I 1 23037639 969696 Grape 12 252525
2 II 1 18779884 969696 Grape 12 252525
3 III 1 19341862 969696 Grape 12 252525
4 IV 1 23867706 969696 Grape 12 252525
5 V 1 25021643 969696 Grape 12 252525
6 VI 1 21508407 0ab276 Grape 12 252525
第二个是共线性分析的结果
> head(synteny_dual_comparison)
Species_1 Start_1 End_1 Species_2 Start_2 End_2 fill
1 1 12226377 12267836 2 5900307 5827251 cccccc
2 15 5635667 5667377 17 4459512 4393226 cccccc
3 9 7916366 7945659 3 8618518 8486865 cccccc
4 2 8214553 8242202 18 5964233 6027199 cccccc
5 13 2330522 2356593 14 6224069 6138821 cccccc
6 11 10861038 10886821 10 8099058 8011502 cccccc
了解了基本的数据输入格式，那么我们就可以将blast的比对结果转化成这种形式。
首先是使用blast比对
构建数据库
makeblastdb -in mt.fasta -dbtype nucl -out mt
比对
blastn -query cp.fasta -db mt -outfmt 6 > output.txt
构造RIdeogram的输入数据
df1.txt文件内容
Chr,Start,End,fill,species,size,color
I,1,131478,FF9D1E,chloroplast,12,252525
I,1,444567,FF9D1E,mitochondrion,12,252525
output6.txt文件内容
NC_044701 NC_044768 96.76 1603 43 7 54477 56074 375799 374201 0.0 2663
NC_044701 NC_044768 83.69 423 37 19 68192 68587 44003 44420 2e-101 370
NC_044701 NC_044768 85.28 326 38 5 66335 66654 360145 360466 1e-88 327
NC_044701 NC_044768 74.07 891 172 45 102938 103801 332888 333746 1e-83 311
NC_044701 NC_044768 90.54 148 11 3 36180 36324 64144 64291 4e-48 193
NC_044701 NC_044768 95.56 90 3 1 49 138 406657 406745 4e-33 143
NC_044701 NC_044768 96.51 86 2 1 110947 111032 420751 420667 2e-32 141
NC_044701 NC_044768 96.43 84 3 0 31738 31821 384266 384349 5e-32 139
NC_044701 NC_044768 94.94 79 4 0 53913 53991 2926 3004 2e-27 124
NC_044701 NC_044768 96.77 62 2 0 105594 105655 110375 110314 2e-21 104
NC_044701 NC_044768 87.34 79 6 3 88210 88284 240028 240106 2e-16 87.9
NC_044701 NC_044768 100.00 31 0 0 10365 10395 63403 63373 2e-07 58.4
NC_044701 NC_044768 96.97 33 0 1 66299 66330 360093 360125 2e-06 54.7
作图代码
df1<-read.csv("df1.txt",stringsAsFactors = F)
df1
df2<-read.csv("output6.txt",header=F,sep="\t",stringsAsFactors = F)
df3<-df2[,c(3,7,8,9,10)]
df3$fill<-ifelse(df3$V3>90,"0080cc",
ifelse(df3$V3<80,"0ab276","e64e60"))
df3$Species_1<-1
df3$Species_2<-1
head(df3)
df4<-df3%>%
select(Species_1,V7,V8,Species_2,V9,V10,fill)
colnames(df4)<-colnames(synteny_dual_comparison)
ideogram(karyotype =df1 ,
synteny=df4,output = "1.svg")
rsvg_pdf("1.svg",'3.pdf')
结果
image.png
但是这个图有一个缺点是不能体现出两条序列长度的差异，不知道能不能按长度的比列来显示，如何用代码实现暂时还不知道。
想到一个办法是出图后手动加上刻度线
今天先到这里了
欢迎大家关注我的公众号
小明的数据分析笔记本
网友评论