recos軟件繪制不同類型TE在染色體的分布熱圖

1、前期準(zhǔn)備

第一步, 獲取repeat注釋結(jié)果,使用RepeatMasker對基因組進(jìn)行repeat分析,可以得到以out結(jié)尾的文件,命令行如下:

#genome.fa:為進(jìn)行repeat分析的基因組文件
#denovo.lib:從頭repeat分析得到的repeat庫
RepeatMasker genome.fa  -lib denovo.lib -s  -nolow  -norna  -gff  -engine  ncbi  -parallel  8  -no_is  -dir  ./ 
# 輸出文件會在-dir指定的路徑下生成 genome.fa.out

第二步,獲取關(guān)心TE類型在基因組的分布

#通過模式匹配獲取位置分布,下面命令行以LTR/Gypsy為例進(jìn)行展示,
#   后面我們還依次得到LTR/Copia、LINE、DNA的分布
sed 's/^\s*//' genome.fa.out|awk 'BEGIN{OFS="\t"}{if($11~/LTR\/Gypsy/){print $5,$6,$7,$11}}' | sort -k1,1 -k2,2n >LTR_Gypsy.bed
bedtools merge  -i  LTR_Gypsy.bed >LTR_Gypsy.merged.bed

第三步,生成基因組的配置文件

#計算染色體長度
samtools faidx genome.fa -o genome.fai
#刪除冗余列,只保留序列名稱和長度列,并且只保留染色體水平的序列
grep "^Chr" genome.fai|cut -f1,2 >genome.length.txt
#以100Kb為窗口,生成基因組的bed文件
bedtools makewindows -g genome.length.txt -w 100000 >genome.bed
#生成繪圖用到的基因組信息文件
cat genome.length.txt|while read chr len;do index=$[index + 1]; echo -e "$chr\t$len\t12\t1\t1\t0\tblack\tblack\tblack\t$index";done >genome.list 

第四步,計算不同類型TE在基因組以100Kb為窗口中的覆蓋度

#還是以LTR/Gypsy為例(該文件來自第二步),計算窗口中的覆蓋度,
#其中g(shù)enome.bed來自第三步,其他類型操作類似依次得到LTR_Copia.BED、LINE.BED、DNA.BED
bedtools coverage -a genome.bed -b LTR_Gypsy.merged.bed| awk 'BEGIN{OFS="\t"}{print $(1),$(2),$(3),$(7)}' >LTR_Gypsy.BED

2、生成配置文件

首先從整體上看下配置文件(example.ini)的內(nèi)容,后面我們詳細(xì)說明ref_in*中color_file 、lab、lab_color 、lab_size 、min、max參數(shù)。

[canvas]                                                                                                                                       
width  = 1000
height = 1500
direction = horizontal
axis_ratio = 0.05
name_ratio = 0.1
margin = 10,10,5,20
inner_ratio = 0.1,0.85,0.05,0,0

[axis]
canvas_position = bottom
ticks_minor = 2Mb
ticks_major = 10Mb
ticks_minor_len = 5
ticks_major_len = 10
axis_line = 0.4
axis_color = rgb(0,0,0)
axis_label = 0.85
axis_label_size = 12
axis_label_color = rgb(0,255,255)
axis_width = 1
axis_opacity = 1
label_unit = Mb

[chromosome]
canvas_position = left
chromosome_list = min_ref.list
chroms = 1,2,3,4,5,6,7,8,9,10
name_position = topleft
round = 0

[ref_in1]
file = LTR_Gypsy.BED
type = heatmap
pos0 = 0
pos1 = 0.2
color_file = color.txt
lab = LTR/Gypsy
lab_color = black
lab_size = 10
min = 0
max = 1

[ref_in2]
file = LTR_Copia.BED
type = heatmap
pos0 = 0.2
pos1 = 0.4
color_file = color.txt
lab = LTR/Copia
lab_color = black
lab_size = 10
min = 0
max = 1

[ref_in3]
file = DNA.BED
type = heatmap
pos0 = 0.4
pos1 = 0.6
color_file = color.txt
lab = DNA
lab_color = black
lab_size = 10
min = 0
max = 1

[ref_in4]
file = LINE.BED
type = heatmap
pos0 = 0.6
pos1 = 0.8
color_file = color.txt
lab = LINE
lab_color = black
lab_size = 10
min = 0
max = 1

如果type為heatmap,我們可以使用內(nèi)置生成顏色的方案,這就需要設(shè)置low_color和high_color ,也可以自己生成漸變色列表。如果使用自定義的顏色方案,就需要配置color_file 參數(shù),里面是漸變色碼。顏色和值的對應(yīng)關(guān)系是線性的,文件開頭的顏色值對應(yīng)值較小的顏色,后面的顏色值對應(yīng)值較大時的顏色。如果使用自己的顏色方案,可以使用下面R代碼生成漸變色列表文件,只需要修改里面的顏色值即可:

#!/usr/bin/Rscript                                                                                                                             
code <- colorRampPalette(c("blue", "orange", "red"),space = "rgb")(100)
write.table(code, file="color.txt", quote=FALSE, col.names=FALSE, row.names=FALSE)
library("RColorBrewer")
display.brewer.pal(11,"PiYG")
data<-brewer.pal(11,"PiYG")
write.table(rev(data), file="color2.txt", quote=FALSE, col.names=FALSE, row.names=FALSE)

這里的lab參數(shù)用于設(shè)置標(biāo)簽名稱,來標(biāo)注每一部分的TE類型;lab_color設(shè)置標(biāo)簽顏色;lab_size設(shè)定標(biāo)簽的字體大小。
對于一個窗口,如果沒有覆蓋,那么其覆蓋度就是0,而如果完全被覆蓋,即這個窗口都是某一類型的TE,那么其覆蓋度就是1,所以我們設(shè)置min=0、max=1。

3. 運(yùn)行軟件生成圖片

recos下載安裝完成后,按照上面描述生成相關(guān)配置文件,就可以運(yùn)行軟件了:

#其中example.ini為配置文件名稱,-o后面的sample表示輸出文件名稱前綴,下例輸出文件為sample.svg
perl recos_wrapper.pl -c example.ini -o sample
最終生成的不同TE類型的分布熱圖
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時請結(jié)合常識與多方信息審慎甄別。
平臺聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點(diǎn),簡書系信息發(fā)布平臺,僅提供信息存儲服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容