R语言R语言可视化

箱线图走一个

2019-02-11  本文已影响5人  小洁忘了怎么分身
# 从内置数据集的表达矩阵中找TP53基因的表达量
rm(list=ls())
options(stringsAsFactors = F)
if(!suppressMessages(require("hgu95av2.db")))BiocManager::install("hgu95av2.db")
if(!suppressMessages(require("CLL")))BiocManager::install("CLL")
if(!suppressMessages(require(tidyverse)))install.packages("tidyverse")
if(!suppressMessages(require(ggpubr)))install.packages("ggpubr")
suppressMessages(library(hgu95av2.db))
suppressMessages(library(CLL))
suppressMessages(library(ggpubr))

data(sCLLex)
# sCLLex
exprSet <- exprs(sCLLex) #探针的表达量
exprSet[1:4,1:4]
#>           CLL11.CEL CLL12.CEL CLL13.CEL CLL14.CEL
#> 1000_at    5.743132  6.219412  5.523328  5.340477
#> 1001_at    2.285143  2.291229  2.287986  2.295313
#> 1002_f_at  3.309294  3.318466  3.354423  3.327130
#> 1003_s_at  1.085264  1.117288  1.084010  1.103217
pd <- pData(sCLLex) #sampleID与disease的对应关系
head(pd)
#>           SampleID  Disease
#> CLL11.CEL    CLL11 progres.
#> CLL12.CEL    CLL12   stable
#> CLL13.CEL    CLL13 progres.
#> CLL14.CEL    CLL14 progres.
#> CLL15.CEL    CLL15 progres.
#> CLL16.CEL    CLL16 progres.

p2s <- toTable(hgu95av2SYMBOL) #探针与symbol的对应关系
head(p2s)
#>    probe_id  symbol
#> 1   1000_at   MAPK3
#> 2   1001_at    TIE1
#> 3 1002_f_at CYP2C19
#> 4 1003_s_at   CXCR5
#> 5   1004_at   CXCR5
#> 6   1005_at   DUSP1

p3 <- filter(p2s,symbol=='TP53')

# boxplot [find TP53 has 3 probe IDs]
probe_tp53 <- p3$probe_id
for(i in 1:3){
boxplot(exprSet[probe_tp53[i],] ~ pd$Disease)
}
#用ggpubr作图
#http://www.sthda.com/english/articles/24-ggpubr-publication-ready-p#lots/
expd <- rownames_to_column(as.data.frame(exprSet))
expd[1:4,1:4]
#>     rowname CLL11.CEL CLL12.CEL CLL13.CEL
#> 1   1000_at  5.743132  6.219412  5.523328
#> 2   1001_at  2.285143  2.291229  2.287986
#> 3 1002_f_at  3.309294  3.318466  3.354423
#> 4 1003_s_at  1.085264  1.117288  1.084010
expd2 <- gather(expd,
                   key = 'sample',
                   value = 'exp',-1)
pd <- rownames_to_column(pd)
head(pd)
#>     rowname SampleID  Disease
#> 1 CLL11.CEL    CLL11 progres.
#> 2 CLL12.CEL    CLL12   stable
#> 3 CLL13.CEL    CLL13 progres.
#> 4 CLL14.CEL    CLL14 progres.
#> 5 CLL15.CEL    CLL15 progres.
#> 6 CLL16.CEL    CLL16 progres.
expd3 <- inner_join(expd2,pd,by=c('sample'='rowname'))
i=1 ###可换1,2
for(i in 1:3){
p <- ggboxplot(filter(expd3,rowname==probe_tp53[i]), 
               x = 'Disease',
               y = 'exp',
               color = "Disease", palette =c("#00AFBB", "#E7B800", "#FC4E07"),
               add = "jitter", shape = "Disease")
print(p)
}


上一篇 下一篇

猜你喜欢

热点阅读