DEComplexDisease旨在寻找复杂疾病的差异表达基因(DEGs),这些疾病的特征是表达谱的异质性。与现有的DEG分析工具不同,它不假设复杂疾病的患者具有共同的DEG。通过应用双聚类算法,DEComplexDisease发现了许多患者共享的deg。通过这种方式,DEComplexDisease以一种新的语法描述了复杂疾病的deg,例如,一个由200个基因组成的基因列表在30%的研究复杂疾病患者中存在差异表达。
使用DEComplexDisease分析结果,用户可以根据共享签名找到受相同机制影响的患者。这可以通过对特征患者或基因的双聚类分析和富集分析的断点建模来实现,例如携带相同突变的患者。DEComplexDisease还为几乎所有的输出提供可视化工具。
如果(!requireNamespace("BiocManager", quiet =TRUE)) install.packages("BiocManager")::install("DEComplexDisease")
DEComplexDisease需要强大的计算能力来完成整个分析,原因如下:(1)它是为大数据设计的,通常有数百名患者;(2)实现了多步双聚类算法。因此,强烈建议使用多核并行计算。
在DEComplexDisease
包装,我们使用mclapply
的工具平行
包做多进程计算。中设置“cores = n”DEComplexDisease
功能要利用多核硬件,其中“n”为核心数。
在分析之前,最好开始一个新的会话或运行
rm(列表= ls ())
清除当前会话。否则,整个会话将被复制到每个进程,这可能会消耗巨大的ROM内存。一旦内存用完,DEComplexDisease
可能报告错误的结果或其他问题。
在本手册中,所有示例都假设用户通过设置拥有一台4核计算机核
= 4。
DEComplexDisease
接受两个强制输入:表达式矩阵(exp)和示例注释向量(cl)。表达矩阵可以是RNA-seq计数数据或归一化微阵列表达数据。样本注释载体(cl)是表示样本疾病状态的载体。强制要求,cl
的列的长度和顺序是一样的经验值
.cl
只有两个可能的值:1和0。1为患者对应样本,0为对照组或正常样本。可选的输入是临床注释数据,可用于Plot()函数中的可视化。临床注释应该有row.names
的列名匹配经验值
.
库(“DEComplexDisease”)
##警告:软件包“DEComplexDisease”已弃用,将从## Bioconductor 3.17版中删除
#load RNA-seq计数矩阵数据(exp) #load样本注释向量数据(cl) #load样本ER状态注释数据(ann.er) exp[1:5,1:5]
# # TCGA.D8.A27G。01 TCGA.A2.A04N。01a tcga.a8.a096.01a tcga.bh.a0bz。01A ## POLRMTP1 59 5 80 91 ## C6orf52 93 59 59 25 40 ## NDUFS6 3356 2160 3529 6269 ## KCTD17 761 438 988 972 ## CEP126 510 317 1303 497 ## TCGA.BH.A0E2.01A ## POLRMTP1 58 ## C6orf52 196 ## NDUFS6 2684 ## KCTD17 437 ## CEP126 199
头(cl, 4)
# # TCGA.D8.A27G。01 TCGA.A2.A04N。01a tcga.a8.a096.01a tcga.bh.a0bz。01a ## 1 1 1 1 1
DEComplexDisease适用bi.deg
函数来估计描述对照样本的基因表达谱的分布参数。正常样品采集依据cl
.bi.deg
有三种方法估计对照样本的参考表达剖面:“edger”,“deseq2”或“normalized”。“edger”或“deseq2”用于RNA-seq计数数据,通过实现由刨边机或DESeq2.估计每个基因的离散度(disp)和平均值(mu)。\(x_{i,j}\)的计数数由
P =pnbinom(x, size=1/disp, mu=mu, lower。尾巴= F)
“normalized”用于规范化的RNA-seq或微阵列数据。估计每个基因的平均值(mu)和标准偏差(sd)。z分数和p值由
z = (xμ)/ sd p = pnorm (z, lower.tail = F)
使用用户定义的p值截断值,bi.deg
赋值1或-1表示上调或下调基因,其中1为上调基因,-1为下调基因。其他基因被赋值为0。在这个例子中,
Deg =bi.deg(exp, cl, method="edger", cutoff=0.05, cores=4)
##使用经典模式。
bi.deg
返回一个度
对象,它是一个由1,0和-1组成的矩阵。用户可以使用Plot()查看分析结果。
情节(度,安=安。呃,show.genes = row.names(度)[1:5])
##所有突变类型:向上,向下。
假设' alter_fun '是可向量化的。如果它没有生成正确的##情节,请在“oncopprint()”中设置“alter_fun_is_vectorized = FALSE”。
在这里,用户可以使用show.genes
显示所选基因。否则,max.n
将显示富集程度最高的基因。在本例中,我们还使用了a安
显示患者ER状态的标注。
来自bi.deg
有噪声,例如,当一个样本再次测试时,与疾病无关的deg参考。他的分析背后的想法是,与疾病相关的deg应该在其他患者身上观察到。换句话说,当有足够的样本时,患者可以用于高信度DEGs的交叉验证。deg.specific
对二进制DEG矩阵进行双聚类分析,以找到高置信疾病相关的DEG。
deg.specific
用于寻找交叉验证的deg。用户只需通过设置选择合适的最小基因和最小患者数min.genes
而且min.patients
.请注意min.patients
包括种子病人本身。
Res.deg = degree .specific(deg, min.genes=30, min.patients=5, cores=4)
如果用户只想看到某些选定患者的结果,则使用test.patients
减少计算时间。
res. degree .test= degree .specific(deg, test.patients=colnames(deg)[1:10], min.genes=50, min.patients=8, cores=4)
deg.specific
返回一个deg.specific
对象,它是一个由交叉验证的deg和支持邻居组成的列表。如果test.patients
, adeg.specific.test
对象返回。用户可以使用Plot()来可视化deg.specific
而且deg.specific.test
对象。在图2中,我们显示了所有使用的乳腺癌患者中验证过的DEGs。
情节(res.deg安=安。呃,马克斯。n = 5)
##所有突变类型:下,上。
假设' alter_fun '是可向量化的。如果它没有生成正确的##情节,请在“oncopprint()”中设置“alter_fun_is_vectorized = FALSE”。
为res.deg.test
,还可以绘制:
情节(res.deg.test安=安。呃,max.n = 5)
##所有突变类型:向上,向下。
假设' alter_fun '是可向量化的。如果它没有生成正确的##情节,请在“oncopprint()”中设置“alter_fun_is_vectorized = FALSE”。
在这项工作中,DEG模块是指许多患者共享的DEG列表。DEG模块应该是患者共享的签名相似的因果机制。seed.module
而且cluster.module
都是为查找此类模块而设计的。
就像deg.specific
,seed.module
的输出进行双聚类分析bi.deg
把每个病人都当做种子。不同的是,该函数具有更复杂的设置和步骤来预测许多患者共有的DEGs模块。
种子。Mod1 =种子。模块(deg, res.deg=res.deg, min.genes=50, min.patients=20, overlap=0.85, cores=4)
或
种子。Mod2 = seed。模块(deg, test.patients=colnames(deg)[1:10], min.genes=50, min.patients=20, overlap=0.85, cores=4)
无论参数设置是什么,seed.module
首先尝试找到一个所有患者共享的模块,其中最终的患者数量可能小于min.patients
而基因数量可能小于min.genes
.如果这样的模块存在,它将被命名为M0
的模块基因M0
将从模块基因的进一步发现中移除。也就是说,模块基因的M0
不会包含在其他模块中。
双聚类分析由一个DEG种子开始,由一个患者的DEG组成。如果res.deg
时,仅将交叉验证过的deg作为种子,并将交叉验证过的deg初始化种子。否则,所有的患者都将被使用,所有的deg都将被用作种子。将患者的deg逐渐去除,检查是否有残留的种子min.patients
保持相似度时不小于重叠
.seed.module
将记录双聚类分析中基因-患者数的轨迹,存储在曲线
每个病人的输出键。如果test.patients
是设定了,只有病人在列吗test.patients
用作双聚类分析的种子。
seed.module
将记录三种场景下的双聚类结果:
max.genes
记录种子被观察时的患者和基因信息min.patient
.max.patients
存储患者和基因信息时min.genes
,这也是双聚类分析的终点。模型
当基因-患者编号时,存储基因/患者信息曲线
满足'model定义的一些建模标准。方法’,这可能表明分子机制的包合/排除。在当前版本中,'model。方法有四个可能的值:聚类”、“max。“平方”、“最小斜率”和“最小相似度”,分别表示四种不同的建模方法:
slope.clustering
斜率变化最大,可能是分子机制的包合/排除;max.square
是具有最大产物的基因患者数;min.slope
基因-患者数曲线斜率最小;min.similarity
相似度分数最小的点吗seed.module
将会返回一颗‘种子’。模块的对象。这是一张单子。它有一个前缀为“decd”的键:
它可能有一个键“M0”:
其他键是患者id,这是用患者的DEG种子预测的模块。每一个都有几个键:
曲线
;“种子。模块的对象可以用Plot()可视化。
情节(种子。mod1安=安。Er, type="model", max.n=5)
所有突变类型:>= 0.85。
假设' alter_fun '是可向量化的。如果它没有生成正确的##情节,请在“oncopprint()”中设置“alter_fun_is_vectorized = FALSE”。
对于大数据,通常有很多患者样本,而DEComplexDisease可能会预测太多的患者种子模块。cluster.module
通过基于患者和基因特征相似性的改进k-means对患者种子模块进行聚类。将同一聚类内的患者按照连接程度进行排序,找出具有代表性的患者。如果vote.seed
设为false,则将有代表性患者的双聚类分析结果作为模块的最终结果。否则,将通过投票的方法生成一个通用的种子,并使用新种子进行双聚类分析来预测最终结果。所有模块将标记为“M1”,“M2”,“M3”…
集群。Mod1 <- cluster.module(seed. module)mod1核心= 4)
或
集群。Mod2 <- cluster.module(seed. module)mod1,投票。种子= TRUE,核= 4)
cluster.module
返回一个cluster.module
对象。
排序(名称(cluster.mod1),减少= TRUE)
## [1] "decd。输入“decd。clustering" "M4" "M3" ## [5] "M2" "M1" "M0"
名称(cluster.mod1[[“decd.input”]])
##[1] "重叠" "度" "测试。病人“min.genes”##[5]“min.patients””模型。方法”“vote.seed”
名称(cluster.mod1[[“decd.clustering”]])
##[1]“组”“代表”
名称(cluster.mod1[[“M1”]])
##[1]“max。基因”“马克斯。病人”“基因。删除“病人。添加“##[5]”曲线“种子”“模型”
一个cluster.module
有一个前缀为“decd”的键:
“其他键”前缀为“M”,表示模块。每个模块有几个键:
min.patients
在双聚类分析中观察到;min.genes
在双聚类分析中达到”;曲线
;的deg.module
可以通过Plot()进行可视化。
情节(集群。mod1安=安。Er, type="model", max.n=5)
所有突变类型:>= 0.85。
假设' alter_fun '是可向量化的。如果它没有生成正确的##情节,请在“oncopprint()”中设置“alter_fun_is_vectorized = FALSE”。
DEG模块可能与基因或患者有部分重叠。使用“模块。重叠” to check the gene and patient overlap among modules. This function is useful to check the relationship of modules and to choose the proper modules during the exploratory discovery steps.
module.overlap(集群。mod1 max.n = 5)
module.compare
用于比较来自不同研究的模块,例如不同疾病或同一疾病的不同数据。在下面的例子中,我们计算了ER+和ER-乳腺癌患者的模块。对各模块进行了比较module.compare
.
Res.mod1 <- seed。模块(deg[,1:26], min.genes=50, min.patients=10, overlap=0.85, cores=4) res.mod1 <- cluster.module(res.mod1) res.mod2 <- seed. modules (,1:26)Module (deg[,27:52], min.genes=50, min.patients=10, overlap=0.85, cores=4) res.mod2 <- cluster.module(res.mod2)
Module.compare (res.mod1, res.mod2, max. mod2)n1 = 5, max.n2 = 5)
请注意,虽然这个函数是设计来进行模块比较的,但它不能在条件之间找到组成不同的模块。
在模块发现步骤中,对二进制DEG矩阵采用了双聚类算法。种子患者的deg逐渐去除到一些min.genes
,这可能会导致患者数量的增加min.patients
.的输出deg.module
工具,记录基因-患者数变化轨迹为曲线
对于每个模块。
名称(cluster.mod1 [[" M1 "]][[“曲线”]])
“不。基因”“不。病人”“分数”
头(cluster.mod1[[“M1”]][[“曲线”]][[“no.gene”]])
## [1] 84 82 81 75 72 71
头(cluster.mod1[[“M1”]][[“曲线”]][[“no.patient”]])
## [1] 10 11 12 13 14 15
module.curve
能以简单的方式显示患者基因数量。
module.curve(集群。mod1,“M1”)
在这条曲线中,module.curve
突出“max”的要点。基因”、“模型”、“max.patients”。
使用默认设置,deg.module
可能无法找到的最佳断点模型
一些模块。用户可以修改模型
所有或部分模块的结果。module.modelling
提供此类工具。它提供了两种方式来修改建模结果keep.gene.num
或者改变model.method
.
用户可以更改模型
结果通过手动设置keep.gene.num
的基因号模型
结果是选择。keep.gene.num
可以是整数值或矢量。如果是整数,则所有模块的值都相同keep.gene.num
这将改变所有模块的建模结果。如果它是一个向量,它的元素应该以模块名作为它们的名称。否则,将只使用第一个元素,并设置所有模块。当keep.gene.num
是向量,它的长度不必与模块相同。可以只更改部分模块。左边的模块将使用默认设置。
在当前版本中,model.method
有四个可能的值:" slope "。聚类”、“max。“平方”、“最小斜率”和“最小相似度”,分别表示四种不同的建模方法:
slope.clustering
具有最大的斜率变化,这可能表明分子机制的包含/排除。max.square
是具有最大产物的基因患者数;min.slope
为基因-患者数曲线中斜率最小的点;min.similarity
是基于相似度得分,选择相似度得分最小的点。x=c(50,40) names(x)<-c("M1","M3") new.cluster.mod1=module.modeling(cluster. modeling)Mod1, keep.gene.num = x, model.method='slope。clustering', cores=4) #here, only "M1" and "M3" are modified new.cluster.mod1=module.modeling(cluster.mod1, keep.gene.num = 50) # here, all the modules are modified module.curve(new.cluster.mod1, "M1")
在复杂疾病中,许多患者具有相似的临床特征或携带相同的基因突变。这些特征患者应该受到共同机制的影响。model.screen
用于寻找可能与特征患者或基因相关的模块。
module.screen(集群。mod1 feature.patients =样本(colnames(度),10))
#搜索模块模块。Mod1, feature.patients=sample(colnames(deg),10), method="fisher.test")