R语言技巧

R包:OmicCircos

2020-01-03  本文已影响0人  小潤澤

前言:

在R的世界里面,画圈图的工具除了RCircos和Zuguang Gu的circlize以外,还有个OmicCircos也可以画圈图
当然,这个我是学习大神@ nnlrl 的推送:https://www.jianshu.com/p/99e9b7ab5731
在这里先注明出处。
这里我只谈谈我的学习记录

关于这个包:

这个R包收录在了Bioconductor,所以只用输入下列命令就可安装:

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("OmicCircos")

import R包自带数据

#载入R包
library(OmicCircos)
options(stringsAsFactors = FALSE)
#载入测试数据
data("TCGA.BC.fus")
data("TCGA.BC.cnv.2k.60")
data("TCGA.BC.gene.exp.2k.60")
data("TCGA.BC.sample60")
data("TCGA.BC_Her2_cnv_exp")
image.png

任何绘图,第一重要的永远是数据结构和数据类型,下面,我就一一展示下这些数据

1. TCGA.BC_Her2_cnv_exp

image.png

这份数据包括了染色体号,基因位置,基因名和统计推断中的相关特征值(t-value,p-value,FDR等)

2. TCGA.BC.cnv.2k.60

image.png

这一个表格包括了染色体号,基因位置信息,基因名和不同组织拷贝数变异的情况(可能叙述的不是很准确,我没弄清楚这个value代表什么)

3. TCGA.BC.sample60

image.png

这个表格显示了不同的sample的不同肿瘤亚型

4. TCGA.BC.gene.exp.2k.60

image.png

这一个表格包括了染色体号,基因位置信息,基因名和不同组织基因表达的情况

5.TCGA.BC.fus

image.png

该表介绍了融合基因的相关情况

筛选目的信息:

#cnv的p值进行log10转化
pvalue <- -1 * log10(TCGA.BC_Her2_cnv_exp[,5])
pvalue <- cbind(TCGA.BC_Her2_cnv_exp[,c(1:3)], pvalue)

#筛选Her2亚型
Her2.i <- which(TCGA.BC.sample60[,2] == "Her2")
Her2.n <- TCGA.BC.sample60[Her2.i,1]
Her2.j <- which(colnames(TCGA.BC.cnv.2k.60) %in% Her2.n)

#筛选Her2亚型患者cnv信息
cnv    <- TCGA.BC.cnv.2k.60[,c(1:3,Her2.j)] 
cnv.m  <- cnv[,c(4:ncol(cnv))]
cnv.m[cnv.m >  2] <- 2
cnv.m[cnv.m < -2] <- -2
cnv <- cbind(cnv[,1:3], cnv.m)

#筛选Her2亚型患者的基因表达信息
Her2.j   <- which(colnames(TCGA.BC.gene.exp.2k.60) %in% Her2.n)
gene.exp <- TCGA.BC.gene.exp.2k.60[,c(1:3,Her2.j)]

筛选完成的表长这个样子:

  1. cnv:
    此时仅有Her2 患者的cnv信息


    image.png

2.gene.exp
此时仅含Her2 患者的基因表达信息


image.png

绘图:

colors <- rainbow(10, alpha=0.5);

#画图,首先建立一个空白图层
pdf("OmicCircos4vignette10.pdf", 8,8);
par(mar=c(2, 2, 2, 2));
plot(c(1,800), c(1,800), type="n", axes=FALSE, xlab="", ylab="", main="");


#选择需要放大的坐标
zoom <- c(1, 22, 939245.5, 154143883, 0, 180);

#开始画图
circos(R=400, cir="hg18", W=4,   type="chr", print.chr.lab=TRUE, scale=TRUE, zoom=zoom);#最外层染色体信息
circos(R=300, cir="hg18", W=100, mapping=gene.exp, col.v=4,  type="heatmap2", cluster=TRUE, col.bar=TRUE, lwd=0.01, zoom=zoom);#次外层基因表达信息
circos(R=220, cir="hg18", W=80,  mapping=cnv,      col.v=4,   type="ml3", B=FALSE, lwd=1, cutoff=0, zoom=zoom);#第三层cnv信息
circos(R=140, cir="hg18", W=80,  mapping=pvalue,   col.v=4,    type="l",   B=TRUE, lwd=1, col=colors[1], zoom=zoom);#第四层p值曲线
circos(R=130, cir="hg18", W=10,  mapping=TCGA.BC.fus, type="link", lwd=2, zoom=zoom);#基因共表达连线


## 局部放大操作
the.col1=rainbow(10, alpha=0.5)[1];
highlight <- c(140, 400, 11, 282412.5, 11, 133770314.5, the.col1, the.col1);
circos(R=110, cir="hg18", W=40, mapping=highlight, type="hl", lwd=2, zoom=zoom);
the.col2=rainbow(10, alpha=0.5)[6];
highlight <- c(140, 400, 17, 739525, 17, 78385909, the.col2, the.col2);
circos(R=110, cir="hg18", W=40, mapping=highlight, type="hl", lwd=2, zoom=zoom);

## highlight link
highlight.link1 <- c(400, 400, 140, 376.8544, 384.0021, 450, 540.5);
circos(cir="hg18", mapping=highlight.link1, type="highlight.link", col=the.col1, lwd=1);
highlight.link2 <- c(400, 400, 140, 419.1154, 423.3032, 543, 627);
circos(cir="hg18", mapping=highlight.link2, type="highlight.link", col=the.col2, lwd=1);

## zoom in chromosome 11
zoom <- c(11, 11, 282412.5, 133770314.5, 180, 270);
circos(R=400, cir="hg18", W=4,   type="chr", print.chr.lab=TRUE, scale=TRUE, zoom=zoom);
circos(R=300, cir="hg18", W=100, mapping=gene.exp, col.v=4,  type="heatmap2", cluster=TRUE, lwd=0.01, zoom=zoom);
circos(R=220, cir="hg18", W=80,  mapping=cnv,      col.v=4,   type="ml3", B=FALSE, lwd=1, cutoff=0, zoom=zoom);
circos(R=140, cir="hg18", W=80,  mapping=pvalue,   col.v=4,    type="l",   B=TRUE, lwd=1, col=colors[1], zoom=zoom);

## zoom in chromosome 17
zoom <- c(17, 17, 739525, 78385909, 274, 356);
circos(R=400, cir="hg18", W=4,   type="chr", print.chr.lab=TRUE, scale=TRUE, zoom=zoom);
circos(R=300, cir="hg18", W=100, mapping=gene.exp, col.v=4,  type="heatmap2", cluster=TRUE, lwd=0.01, zoom=zoom);
circos(R=220, cir="hg18", W=80,  mapping=cnv,      col.v=4,   type="ml3", B=FALSE, lwd=1, cutoff=0, zoom=zoom);
circos(R=140, cir="hg18", W=80,  mapping=pvalue,   col.v=4,    type="l",   B=TRUE, lwd=1, col=colors[1], zoom=zoom);
image.png

几个重要的参数:

“chr”:染色体或片段的绘图
“chr2”:染色体或部分基因组片段图
“heatmap”:热图
“heatmap2”:基因组坐标的heatmaps
“hightlight.link”:用于缩放的链接线
“hl”:突出显示
“label”:基因标签或文本注释
“label2”:具有相同圆周坐标的基因标签或文本注释
“link.pg”:基于贝塞尔曲线的链接多边形
“link”:基于贝塞尔曲线的链接线
“link2”:具有较小染色体内弧的链接线
“l”:折线
“ls”:阶梯状线图
“lh”:水平线
“ml”:多折线
“ml2”:多水平线
“ml3”:多阶梯状线图
“box”:箱线图
“hist”:多个样本的直方图
“ms”:多个样本的点图(类似箱线图)
“h”:柱状图
“s”:点图
“b”:条形图
“quant75”:75%分位数线
“quant90”:90%分位数线
“ss”:与数值成比例的点
“sv”:与方差成比例的点
“s.sd”:与标准偏差成比例的点
“ci95”:95%置信区间线
“b2”:条形图(双向)
“b3”:相同高度的条形图
“s2”:固定半径的点
“arc”:半径可变的弧(表示基因片段)
“arc2”:具有固定半径的弧

参考文献:https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3921174/

上一篇下一篇

猜你喜欢

热点阅读