Scissor：联合表型数据,Bulk-seq和scRNA(2)

前面一个帖子讲了scissor的原理以及paper中的一些应用实例。几天我们来测试这个工具。

========安装========

devtools::install_github('sunduanchen/Scissor')

devtools::install_github("jinworks/scAB")

注：因为我们还要用到scAB工具中的例子，所以顺便安装一下。

library(Scissor)

library(Seurat)

library(preprocessCore)

library(scAB)

=======加载数据======

data("data_survival")

dim(sc_dataset)

sc_dataset

head(bulk_dataset[,1:10])

head(phenotype)

=======测试数据=======

sc_dataset <- run_seurat(sc_dataset,verbose = FALSE)

sc_dataset

UMAP_celltype <- DimPlot(sc_dataset, reduction ="umap",group.by="celltype",label = T)

UMAP_celltype

bulk2 <- as.matrix(bulk_dataset)  #我也不清楚我在linux下面报错了，所以强制转换了一下

infos1 <- Scissor(bulk2, sc_dataset, phenotype, alpha = 0.05,family = "cox", Save_file = './Scissor_survival.RData')    #phenotype也没让人选，哎。。。

运行大概几分钟，在本地文件夹下即可生成结果Rdata，检查结果可见共识别429个Scissor+细胞（即与表型呈正相关的细胞）和1001个Scissor阴性细胞（即与表型呈负相关的细胞）。

========可视化结果===========

Scissor_select <- rep("Background", ncol(sc_dataset))

names(Scissor_select) <- colnames(sc_dataset)

Scissor_select[infos1$Scissor_pos] <- "Scissor+"

Scissor_select[infos1$Scissor_neg] <- "Scissor-"

sc_dataset <- AddMetaData(sc_dataset, metadata = Scissor_select, col.name = "scissor")

UMAP_scissor <- DimPlot(sc_dataset, reduction = 'umap',

        group.by = 'scissor',

        cols = c('grey','royalblue','indianred1'),

        pt.size = 0.001, order = c("Scissor+","Scissor-"))

patchwork::wrap_plots(plots = list(UMAP_celltype,UMAP_scissor), ncol = 2)

下面，我们统计一下scissor细胞在各个cluster的分布情况：

table(sc_dataset$scissor,sc_dataset$celltype)

从结果可以发现，Scissor+ 细胞主要为肿瘤细胞，Scissor阴性细胞主要为T细胞和CAFs：

用下面的图展示会更直观一点：

library(gplots)

balloonplot(table(sc_dataset$scissor,sc_dataset$celltype))

但是，问题来了，这个出来的结果，到底如何与生物学问题挂钩和解释呢？我看phenotype里面有time和status两个表型，是与这2个表型都密切相关吗？