怎么使用R語言包RIdeogram展示blast雙序列比對結果

RIdeogram 是用來展示 染色體組型 （idiogram）的一個R語言包，比如展示snp的密度分布；展示某類基因在染色體上的分布等等。

發表的論文是

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')

結果

但是這個圖有一個缺點是不能體現出兩條序列長度的差異，不知道能不能按長度的比列來顯示，如何用代碼實現暫時還不知道。 想到一個辦法是出圖后手動加上刻度線。