Scissor实操（二）：逻辑回归分析流程解读
前言
在上一期推文：Scissor实操（一）：利用Cox回归模型识别特征细胞，Immugent介绍了如何使用Scissor对连续变量的表征(患者生存信息)进行关联分析，但我们知道在各种临床信息中还有一大类是分类变量。例如某个基因突变情况，不同阶段的肿瘤患者不同性别的肿瘤患者信息等等。那么在本篇推文中，Immugent将会介绍如何使用Scissor对分类变量的表型进行关联分析。
此外，由于Scissor还有最后一部分内容不值当单独写一篇推文（所以生信宝库还是很良心的，从不坑蒙粉丝，哈哈），Immugent也会在文末进行介绍。
代码展示
这次使用的也是来源于TCGA-LUAD的数据，只是本次演示的是如何使用Scissor对分类变量特征进行分析。
load(url(paste0(location, 'TCGA_LUAD_exp2.RData')))
load(url(paste0(location, 'TCGA_LUAD_TP53_mutation.RData')))
table(TP53_mutation)
## TP53_mutation
## 0 1
## 243 255
给定上述TP53突变= 1，野生型= 0的二元变量，下面在Scissor中使用family = 'binomial'来选择表型相关的细胞亚群:
phenotype <- TP53_mutation
tag <- c('wild-type', 'TP53 mutant')
infos4 <- Scissor(bulk_dataset, sc_dataset, phenotype, tag = tag, alpha = 0.5,
family = "binomial", Save_file = "Scissor_LUAD_TP53_mutation.RData")
infos5 <- Scissor(bulk_dataset, sc_dataset, phenotype, tag = tag, alpha = 0.2,
family = "binomial", Load_file = "Scissor_LUAD_TP53_mutation.RData")
总共有414个Scissor+细胞与TP53突变相关，318个Scissor-细胞与野生型相关，下面可以使用UMAP可视化这些选定的细胞:
Scissor_select <- rep(0, ncol(sc_dataset))
names(Scissor_select) <- colnames(sc_dataset)
Scissor_select[infos5$Scissor_pos] <- 1
Scissor_select[infos5$Scissor_neg] <- 2
sc_dataset <- AddMetaData(sc_dataset, metadata = Scissor_select, col.name = "scissor")
DimPlot(sc_dataset, reduction = 'umap', group.by = 'scissor', cols = c('grey','indianred1','royalblue'), pt.size = 1.2, order = c(2,1))
image.png
为了确定推断的表型-细胞关联是否可靠，可以使用reliability.test检验进行信度显著性检验。可靠性显著性检验的动机是:如果选择的单细胞和批量数据不适合表型与细胞的关联，相关性将是信息量较少，与表型标签的关联不佳。因此，相应的预测性能较差，与随机排列的标签没有显著区别。在本教程中，我们将测试上述应用程序中确定的关联，作为示例来展示如何运行reliability.test。
下面选择了205个存活较差或较好的细胞。保存在Scissor中的Rdata可以作为测试输入:
load('Scissor_LUAD_survival.RData')
numbers <- length(infos1$Scissor_pos) + length(infos1$Scissor_neg)
result1 <- reliability.test(X, Y, network, alpha = 0.05, family = "cox", cell_num = numbers, n = 10, nfold = 10)
names(result1)
## [1] "statistic" "p" "c_index_test_real" "c_index_test_back"
load('Scissor_LUAD_TP53_mutation.RData')
numbers <- length(infos5$Scissor_pos) + length(infos5$Scissor_neg)
result2 <- reliability.test(X, Y, network, alpha = 0.2, family = "binomial", cell_num = numbers, n = 10, nfold = 10)
最后，我们可以使用函数求值，为每一个选择的单元格获取一些支持信息。以具有存活信息的LUAD癌细胞为例，我们可通过运行以下代码对205个“Scissor+”细胞进行评估:
evaluate_summary <- evaluate.cell('Scissor_LUAD_survival.RData', infos1, FDR = 0.05, bootstrap_n = 100)
evaluate_summary[1:5,1:4]
all(evaluate_summary$`Mean correlation` & as.numeric(gsub('%','',evaluate_summary$`Correlation > 0`)) > 50)
## [1] TRUE
evaluate_summary[1:5,5:10]
image.png
总结
到目前为止，Scissor软件相关的三篇推文就已经完全更新完毕。但是按照生信宝库送返送到嘴里的一贯做法，未来可能还会更新一篇推文，通过已经发表的文章介绍如何将Scissor运用到实际的文章分析中，敬请期待。
总的来说，Scissor虽然有一些小bug，但是整体还是很稳健的，特别是对于解决单一表征为主的科学问题是非常合适的。未来可能还会出现比Scissor更好的分析软件，让我们拭目以待！
