如何解决如何在使用多变量对象或两个单独对象的组合 ggplot2 图之上添加额外统计信息
我有一个 ggplot2 图，它将两个单独的小提琴图绘制到一个图上，由本示例给出（感谢@jared_mamrot 提供）：
library(tidyverse)
data("Puromycin")
head(Puromycin)
dat1 <- Puromycin %>%
filter(state == "treated")
dat2 <- Puromycin %>%
filter(state == "untreated")
mycp <- ggplot() +
geom_violin(data = dat1,aes(x= state,y = conc,colour = "Puromycin (Treatment1)")) +
geom_violin(data = dat2,colour = "Puromycin (Treatment2)"))
mycp
我想添加箱线图或其他汇总统计信息，例如 http://www.sthda.com/english/wiki/ggplot2-violin-plot-quick-start-guide-r-software-and-data-visualization 和 https://www.maths.usyd.edu.au/u/UG/SM/STAT3022/r/current/Misc/data-visualization-2.1.pdf 中的那些，但尝试在这些地方建议的代码不会改变原始图。
mycp + geom_Boxplot()
感谢阅读，希望这是有道理的！
更新 ======================================== ====================================
所以上面的例子并没有完全反映我现在意识到的情况。本质上，我想将统计数据应用到使用两个单独对象作为其变量的组合 ggplot2 图上（此处为 TNBC_List1 和 ER_List1）这是一个示例（对不起，示例较长，我承认我在创建一个更简单的可重现示例时遇到了麻烦，而且我对编码非常陌生）：
# Libraries -------------------------------------------------------------
library(BiocManager)
library(GEOquery)
library(plyr)
library(dplyr)
library(Matrix)
library(devtools)
library(Seurat)
library(ggplot2)
library(cowplot)
library(SAVER)
library(Metap)
library(multtest)
# Loading Raw Data into RStudio ----------------------------------
filePaths = getGEOSuppFiles("GSE75688")
tarF <- list.files(path = "./GSE75688/",pattern = "*.tar",full.names = TRUE)
tarF
untar(tarF,exdir = "./GSE75688/")
gzipF <- list.files(path = "./GSE75688/",pattern = "*.gz",full.names = TRUE)
ldply(.data = gzipF,.fun = gunzip)
list.files(path = "./GSE75688/",full.names = TRUE)
list.files(path = "./GSE75688/",pattern = "\\.txt$",full.names = TRUE)
# full matrix ----------------------------------------------------------
fullmat <- read.table(file = './GSE75688//GSE75688_GEO_processed_Breast_Cancer_raw_TPM_matrix.txt',sep = '\t',header = FALSE,stringsAsFactors = FALSE)
fullmat <- data.frame(fullmat[,-1],row.names=fullmat[,1])
colnames(fullmat) <- as.character(fullmat[1,])
fullmat <- fullmat[-1,]
fullmat <- as.matrix(fullmat)
# BC01 ER+ matrix -----------------------------------------------------------
BC01mat <- grep(pattern =c("^BC01"),x = colnames(fullmat),value = TRUE)
BC01mat = fullmat[,grepl(c("^BC01"),colnames(fullmat))]
BC01mat = BC01mat[,!grepl("^BC01_Pooled",colnames(BC01mat))]
BC01mat = BC01mat[,!grepl("^BC01_Tumor",colnames(BC01mat))]
BC01pdat <- data.frame("samples" = colnames(BC01mat),"treatment" = "ER+")
# BC07 TNBC matrix -----------------------------------------------------------
BC07mat <- grep(pattern =c("^BC07"),value = TRUE)
BC07mat <- fullmat[,grepl(c("^BC07"),colnames(fullmat))]
BC07mat <- BC07mat[,!grepl("^BC07_Pooled",colnames(BC07mat))]
BC07mat <- BC07mat[,!grepl("^BC07_Tumor",!grepl("^BC07LN_Pooled",!grepl("^BC07LN",colnames(BC07mat))]
BC07pdat <- data.frame("samples" = colnames(BC07mat),"treatment" = "TNBC")
#merge samples together =========================================================================
joined <- cbind(BC01mat,BC07mat)
pdat_joined <- rbind(BC01pdat,BC07pdat)
#fdat ___________________________________________________________________________________
fdat <- grep(pattern =c("gene_name|gene_type"),value = TRUE)
fdat <- fullmat[,grepl(c("gene_name|gene_type"),colnames(fullmat))]
fdat <- as.data.frame(fdat,stringsAsFactors = FALSE)
fdat <- setNames(cbind(rownames(fdat),fdat,row.names = NULL),c("ensembl_id","gene_short_name","gene_type"))
rownames(pdat_joined) <- pdat_joined$samples
rownames(fdat) = make.names(fdat$gene_short_name,unique=TRUE)
rownames(joined) <- rownames(fdat)
# Create Seurat Object __________________________________________________________________
joined <- as.data.frame(joined)
sobj_pre <- CreateSeuratObject(counts = joined)
sobj_pre <-AddMetaData(sobj_pre,Metadata=pdat_joined)
head(sobj_pre@Meta.data)
#gene name input
sobj_pre[["RNA"]]@Meta.features<-fdat
head(sobj_pre[["RNA"]]@Meta.features)
#Downstream analysis -------------------------------------------------------
sobj <- sobj_pre
sobj <- FindVariableFeatures(object = sobj,mean.function = ExpMean,dispersion.function = LogVMR,nfeatures = 2000)
sobj <- ScaleData(object = sobj,features = rownames(sobj),block.size = 2000)
sobj <- RunPCA(sobj,npcs = 100,ndims.print = 1:10,nfeatures.print = 5)
sobj <- FindNeighbors(sobj,reduction = "pca",dims = 1:4,nn.eps = 0.5)
sobj <- FindClusters(sobj,resolution = 1,n.start = 10)
umap.method = 'umap-learn'
metric = 'correlation'
sobj <- RunUMAP(object = sobj,min.dist = 0.5,seed.use = 123)
p0 <- DimPlot(sobj,reduction = "umap",pt.size = 0.1,label=TRUE) + ggtitle(label = "Title")
p0
# ER+ score computation -------------------
ERlist <- list(c("CPB1","RP11-53O19.1","TFF1","MB","ANKRD30B","LINC00173","DSCAM-AS1","IGHG1","SERPINA5","ESR1","ILRP2","IGLC3","CA12","RP11-64B16.2","SLC7A2","AFF3","IGFBP4","GSTM3","ANKRD30A","GSTT1","GSTM1","AC026806.2","C19ORF33","STC2","HSPB8","RPL29P11","FBP1","AGR3","TCEAL1","CYP4B1","SYT1","COX6C","MT1E","SYTL2","THSD4","IFI6","K1AA1467","SLC39A6","ABCD3","SERPINA3","DEGS2","ERLIN2","HEBP1","BCL2","TCEAL3","PPT1","SLC7A8","RP11-96D1.10","H4C8","PI15","PLPP5","PLAAT4","galNT6","IL6ST","MYC","BST2","RP11-658F2.8","MRps30","MAPT","AMFR","TCEAL4","MED13L","ISG15","NDUFC2","TIMP3","RP13-39P12.3","PARD68"))
sobj <- AddModulescore(object = sobj,features = ERlist,name = "ER_List")
#TNBC computation -------------------
tnbclist <- list(c("FABP7","TSPAN8","CYP4Z1","HOXA10","CLDN1","TMSB15A","C10ORF10","TRPV6","HOXA9","ATP13A4","GLYATL2","RP11-48O20.4","DYRK3","MUCL1","ID4","FGFR2","SHOX2","Z83851.1","CD82","COL6A1","KRT23","GCHFR","PRICKLE1","GCNT2","KHDRBS3","SIPA1L2","LMO4","TFAP2B","SLC43A3","FURIN","ELF5","C1ORF116","ADD3","EFNA3","EFCAB4A","LTF","LRRC31","ARL4C","GPNMB","VIM","SDR16C5","RHOV","PXDC1","MALL","YAP1","A2ML1","RP1-257A7.5","RP11-353N4.6","ZBTB18","CTD-2314B22.3","galNT3","BCL11A","CXADR","SSFA2","ADM","GUCY1A3","GSTP1","ADCK3","SLC25A37","SFRP1","PRNP","DEGS1","RP11-110G21.2","AL589743.1","ATF3","SIVA1","TACSTD2","HEBP2"))
sobj <- AddModulescore(object = sobj,features = tnbclist,name = "TNBC_List")
#ggplot2 issue ----------------------------------------------------------------------------
sobj[["ClusterName"]] <- Idents(object = sobj)
sobjlists <- FetchData(object = sobj,vars = c("ER_List1","TNBC_List1","ClusterName"))
library(reshape2)
melt(sobjlists,id.vars = c("ER_List1","ClusterName"))
p <- ggplot() + geom_violin(data = sobjlists,aes(x= ClusterName,y = ER_List1,fill = ER_List1,colour = "ER+ Signature"))+ geom_violin(data = sobjlists,y = TNBC_List1,fill = TNBC_List1,colour="TNBC Signature"))
扩展 ======================================== ================================
如果你想这样做但有两个对象（例如sobjlists1和sobjlists2）而不是我的例子显示的（两个变量但一个对象），rbind这两个然后按照@StupidWolf 说的做
library(reshape2)
sobjlists1= melt(sobjlists1,id.vars = "treatment")
sobjlists2= melt(sobjlists2,id.vars = "treatment")
combosobjlists <- rbind(sobjlists1,sobjlists2)
然后使用 combosobjlists 继续他们的代码：
ggplot(combosobjlists,y = value)) +
geom_violin(aes(fill=variable)) +
geom_Boxplot(aes(col=variable),width = 0.2,position=position_dodge(0.9))
希望这个帖子有帮助！
解决方法
尽量只包含显示问题的最少代码。就像在您的示例中一样，无需从整个 seurat 处理开始。您可以只提供带有 dput() 的 data.frame ，我们可以看到 ggplot2 的问题，请参阅 this post。
创建一些示例数据：
library(Seurat)
library(ggplot2)
genes = c(unlist(c(ERlist,tnbclist)))
mat = matrix(rnbinom(500*length(genes),mu=500,size=1),ncol=500)
rownames(mat) = genes
colnames(mat) = paste0("cell",1:500)
sobj = CreateSeuratObject(mat)
sobj = NormalizeData(sobj)
添加一些虚构的集群：
sobj$ClusterName = factor(sample(0:1,ncol(sobj),replace=TRUE))
添加您的模块分数：
sobj = AddModuleScore(object = sobj,features = tnbclist,name = "TNBC_List",ctrl=5)
sobj = AddModuleScore(object = sobj,features = ERlist,name = "ER_List",ctrl=5)
我们得到了数据，您需要做的是正确旋转它。用 ggplot2 绘制两次会导致各种问题：
sobjlists = FetchData(object = sobj,vars = c("ER_List1","TNBC_List1","ClusterName"))
head(sobjlists)
ER_List1 TNBC_List1 ClusterName
cell1 -0.05391108 -0.008736057 1
cell2 0.07074816 -0.039064126 1
cell3 0.08688374 -0.066967324 1
cell4 -0.12503649 0.120665057 0
cell5 0.05356685 -0.072293651 0
cell6 -0.20053804 0.178977042 1
应该是这样的：
library(reshape2)
sobjlists = melt(sobjlists,id.vars = "ClusterName")
ClusterName variable value
1 1 ER_List1 -0.05391108
2 1 ER_List1 0.07074816
3 1 ER_List1 0.08688374
4 0 ER_List1 -0.12503649
5 0 ER_List1 0.05356685
6 1 ER_List1 -0.20053804
现在我们绘制：
ggplot(sobjlists,aes(x= ClusterName,y = value)) +
geom_violin(aes(fill=variable)) +
geom_boxplot(aes(col=variable),width = 0.2,position=position_dodge(0.9))
为了能够在不指定的情况下使用图中的数据（如 geom_boxplot() ），您需要将数据放入 ggplot() 函数调用中。那么下面的函数就可以继承它们了。
您也不需要每种颜色的额外小提琴图
library(tidyverse)
data("Puromycin")
head(Puromycin)
mycp <- ggplot(Puromycin,aes(x= state,y = conc,colour=state))+geom_violin()
mycp + geom_boxplot(width=0.1,color= "black") +
scale_color_discrete(
labels= c("Puromycin (Treatment1)","Puromycin (Treatment2)")
)
结果：
版权声明：本文内容由互联网用户自发贡献，该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请发送邮件至 dio@foxmail.com 举报，一经查实，本站将立刻删除。
