晶诚所至 生命所能

Engage to Life Energy

 
原创文章|一个个性化oncoplot的实现
发布日期:2020-06-12浏览:

 

 

 


问题的提出

 


肿瘤研究的相关文章中经常出现oncoplot,这种表现形式几乎成为肿瘤研究文章中的必备元素,见诸于上至TCGAThe Cancer Genome Atlas Programm)发表于Nautre、Science、Cell的综合性研究 (Comprehensive profiling)文章,下至零零散散的研究者发表的oncology领域的文章等等。

oncoplot的样子看上去像是heatmap,或者更简单点,类似于个Excel表格,有行(row)、列(column)和单元格(cell)。这种table或者grid的布局形式经常见到,比如网页设计领域里BootstrapGrid System,用于网页元素的排版布局,又比如Rggplot2作图系统对图形元素的布局。利用这种平面的、二维的、行列的形式可以把布局做的很简单,也可以异常复杂。类比于Excel表格,你可以想象成通过填充单元格设置单元格格式合并单元格 设置行列宽度高度等等操作制作一个复杂的、美观的表格。说回oncoplot和heatmap, 它们的单元格里填充颜色(当然可以可以在单元格里填充任何图像要素,比如文字、几何图形等,但最简单的形式就是填充颜色)。 heatmap中经常用颜色深浅代表某个连续型变量数值大小,oncoplot则用颜色代表某个离散型变量的类型,最常见的是表征某些样本的某些oncogene的突变类型的形式。
有专门的、现成的工具可以用来做oncoplot图,比如R的maftoolsComplexHeatmap包等,但其函数参数繁杂,要实现灵活性着实不易。何不用ggplot2来实现这些图形呢?ggplot2用户众多,其基于的grid包又有足够的灵活性(实际上ComplexHeatmap,pheatmap等亦是基于grid来写的),ggplot2加上grid再加上两者的中间层gtable(gtable亦是基于grid),实现这类图,基本上迎刃而解,无往而不利。下面我们先提出一个个性化的问题,然后用ggplot2的方案来解决。最终给出代码。

?


问题如下


上图(摘选自https://www.nature.com/articles/nature13385)意义明了,内容显而易见。但有一个问题,比如某个样本的TP53基因如果既发生了Missense突变又发生了Nonsense突变,那么该单元格里应当包含两个信息了,图形里怎么表示呢?这种情况并不罕见。

 

 

 

问题的解决





为了说明问题,我们从COSMIC(https://cancer.sanger.ac.uk/cosmic) 下载CosmicMutantExportCensus.tsv表格,这是一个综合性的数据,包含了若干万样本的oncogene的变异情况, 我们对数据做了个筛选,保留一个很小的数据集,筛选后的数据可以从github得到(见代码部分),最终保留内容如下(仅显示头部几行):


gene
sample
mutation
EGFR
1004872
Insertion - In frame
EGFR
P-0000670-T01-IM3
Insertion - In frame
RBM10
P-0002758-T01-IM3
Deletion - Frameshift
STK11
P-0006311-T02-IM5
Deletion - In frame
RB1
P-0001171-T01-IM3
Substitution - Missense
ATM
P-0011041-T01-IM5
Deletion - Frameshift /// Substitution - Missense

由上表可见,P-0011041-T01-IM5的ATM基因发生了两种突变。

 


第一种方法


一个简单的考虑,我们把“某个个体的某个gene发生了2个或2个以上的突变”归结为 “multipe hit”即可,作图如下,我们用红色突出显示该类别:‍

第一种方法的改进

上图中multipe hit”单独作为一类展示,但不太够突出,利用gtable稍微做些改动,让这个类型更突出一些。图形如下:


第二种方法


第一种方法的弊端在于,你不能通过“multipe hit”这一个类别知道某个样本的某个基因到底是发生了几个类型的突变,各自是什么类型的突变。我们考虑第二种方法,具体显示突变类型,简单来讲就是把一个单元格分隔,图形如下:


通过上图,可以看出P-0011041-T01-IM5的ATM基因发生了两种突变,分别为Deletion - Frameshift 和Substitution - Missense,某样本在TP53基因中甚至发生了三种突变。‍

 


第三种方法


上述图片都是静态的(statitc)图片,完全可以利用动态的(interactive)的图片来表示更为丰富的内容。比如鼠标hover到某个单元格,弹出详细信息等。Javascript最适合来做这种工作,比如大家熟悉的D3.js,R里可以利用htmlwidgetsl类的包(比如plotly)来实现,当然这是另一个话题了,在此不再赘述。

 

 

比较





‍点击图片查看大图


 

代码




 

代码和数据文件见https://github.com/lilvxi/heatmapCellSplit,代码附录于此,仅仅利用data.table处理数据,利用ggplot2+gtable+grid来作图:
	
library(data.table) library(RColorBrewer) library(grid) library(gridExtra) library(ggplot2) library(gtable) load('todo.RData') ###prepare the data (color, coorinate etc) using data.table manipulation; todo[,.N,sample][order(-N)][,sample]->sams todo[,.N,gene][order(-N)][,gene]->genes melt(dcast(todo,gene~sample,value.var='mutation'),id.var='gene')->todo names(todo)<-c('gene','sample','mutation') todo[is.na(todo)]<-'None' todo[,mut:=mutation] todo[grep("///",mut),mut:='multipe hit'] mut.Cate<-c("Substitution - Missense","Substitution - Nonsense", "Insertion - In frame","Insertion - Frameshift", "Deletion - In frame","Deletion - Frameshift", "Nonstop extension","Whole gene deletion") setNames(brewer.pal(n=10,'Paired')[0-5:6],mut.Cate)->mut.Cate c(mut.Cate,'None'='gray','multipe hit'='red')->mut.Cate todo[,mut:=ordered(mut,levels=names(mut.Cate))] inter.x<-setNames(ppoints(length(sams),1/2),sams) inter.y<-setNames(rev(ppoints(length(genes),1/2)),genes) todo[,x:=inter.x[as.character(sample)]] todo[,y:=inter.y[gene]] todo[,w:=1/length(sams)] todo[,h:=1/length(genes)] ################################the ggplot graph tu<-ggplot(todo,aes(x=x,y=y,height=h,width=w,fill=mut)) +    geom_tile(colour='white',show.legend=F)+    scale_x_continuous(expand=expansion(),breaks=inter.x,labels=names(inter.x))+    scale_y_continuous(expand=expansion(),breaks=inter.y,labels=names(inter.y),position = "right")+    scale_fill_manual(values=mut.Cate) +    theme_minimal() %+replace% theme(axis.text.x=element_text(angle=90,hjust=1,vjust=0.5),legend.position='top') +    labs(x='',y='') gai<-ggplotGrob(tu) #######################################ok1 zz<-todo[mut=='multipe hit'] apply(zz[,.(x,y,w,h)],1,function(xxx) { segmentsGrob(x0=c(0,0),x1=c(1,1),y0=c(0,1),y1=c(1,0),gp=gpar(col='#2B1B17'),vp=viewport(xxx[1],xxx[2],xxx[3],xxx[4])) })->rs append(rs,list(rectGrob(zz$x,zz$y,zz$w,zz$h,gp=gpar(col='#2B1B17'))))->rs addOn<-gTree(children=Reduce(gList,rs)) gtable_add_grob(gai,addOn,t=7,l=5)->ok1 ##################################ok2 apply(zz[,.(x,y,w,h,mutation)],1,function(xxx) { strsplit(xxx[5]," /// ")[[1]]->jus grid.rect(x=ppoints(length(jus),0.5),width=1/length(jus),vp=viewport(xxx[1],xxx[2],xxx[3],xxx[4]),gp=gpar(fill=mut.Cate[jus],col='#2B1B17')) })->rs addOn<-gTree(children=Reduce(gList,rs)) gtable_add_grob(gai,addOn,t=7,l=5)->ok2 ###################################combination grid.arrange(gai,ok1,ok2)

 

上一条:客户文章 | 复旦大学任国栋课题组揭示HYL1拮抗外切体调控miRNA加工生成新机制
下一条:客户文章|东亚人群精神分裂症的遗传结构基础
返回
网站地图 | 法律声明 | 联系我们

地址:上海市松江区中心路1158号5幢5楼

电话:400-9200-612  传真:+86 21 6090 1207/1208-8154

晶能生物技术(上海)有限公司 Copyright 2012 Genergy Inc. 沪ICP备10017363号

友情链接: