从前做序列比对通常都是MEGA来做,然后偶然候会用genedoc来展示,但是常常碰到一个标题就是没办法得到矢量图。以是不停在想用R能不能展示多序列比对的效果,尤其是焦点关键domain的地区。
下面是paper中的一个效果,看着就很清晰。
library("phylotools")
library("ggplot2")
library("tidyr")
library("dplyr")
df <- read.fasta("data.fasta")
df
#用ggplot展示的话,我们就必要拆分成单个碱基的情势,有点雷同长矩阵酿成短矩阵的感觉。
new.df <- df %>% separate(seq.text,paste0("col",str_pad(1:28,2,side = "left",pad = "0")),'') %>% select(-col01) %>%mutate(seq.name=factor(seq.name,levels = rev(seq.name))) %>% pivot_longer(!seq.name)
下面就可以用ggplot画图了
p1<-ggplot(new.df,aes(x=name,y=seq.name))+
geom_tile(fill="white",alpha=0)+
geom_text(aes(label=value))+
theme_bw()+
theme(axis.text.x = element_blank(),
axis.ticks = element_blank(),
panel.grid = element_blank(),
panel.border = element_blank())
p1
下面添加两侧的数字:
df1<-data.frame(x=0,
y=11:1,
label=c(671,668,669,666,744,743,756,736,706,747,759)) #左侧
df2<-data.frame(x=28,
y=11:1,
label=c(695,692,693,690,770,665,777,761,727,768,780)) #右侧
p2 <- p1+
geom_text(data=df1,
aes(x=x,y=y,label=label),
inherit.aes = FALSE,hjust=1)+
geom_text(data=df2,
aes(x=x,y=y,label=label),
inherit.aes = FALSE,hjust=0)+
coord_equal(xlim = c(-1,29))
p2
下面就是添加颜色了:
我们把控制颜色的数据单独放在一个文件中,我只预备了部门。
df3 <- read.table("tmp.txt",sep="\t",header=T)
df3
p2 +
geom_tile(data=df3,aes(x=x,y=y,fill=group),
color="white",show.legend = FALSE)+
geom_text(data=df3,aes(x=x,y=y,label=label),
color="white")+
scale_fill_manual(values = c("#00adef","#ed1b24"))+
labs(x=NULL,y=NULL)
如许颠末控制,就可以得到文中的效果了。 |