今天看到了一个推送，研究的是分子分型，之前说要做这个东西，已经很长时间了，忙忙叨叨一直没有时间进行，终于找到时间进行分析了。
首先要加载这个包，这个包安装到加载真的是非常非常的慢。
1.1下载安装包
getOption('timeout')
options(timeout=6000000)
#这两步一定要运行，不然容易因为网速的问题安装不上
install.packages("Matrix", dependencies = TRUE)
library("devtools")
install_github("danro9685/CIMLR", ref = 'R')
library("CIMLR")
#安装包
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
if (!require("devtools"))
install.packages("devtools")
devtools::install_github("xlucpu/MOVICS")
1.2加载内置数据
安装好了R包之后，我们可以加载一下内部的数据，这里面做的是乳腺癌的
# load example data of breast cancer
load(system.file("extdata", "brca.tcga.RData", package = "MOVICS", mustWork = TRUE))
load(system.file("extdata", "brca.yau.RData", package = "MOVICS", mustWork = TRUE))
# 提取数据，了解数据特征
mo.data <- brca.tcga[1:4]
count <- brca.tcga$count
fpkm <- brca.tcga$fpkm
maf <- brca.tcga$maf
segment <- brca.tcga$segment
surv.info <- brca.tcga$clin.info
我加载完数据，其实我看了一下数据的结构： 关于brca.tcga
看起来是表达量的数据，有mRNA，miRNA，lnRNA等等
根据这些内容，我们可以继续提取count，fpkm，maf，segment，surv.info，其实这些在tcga中都是可以下载的。接下来我们按照流程拍一遍内置的数据看看可以得到什么东西。
1.3鉴定最佳的聚类数
# 鉴定最佳聚类数，耗时较多
optk.brca <- getClustNum(data= mo.data,
is.binary= c(F,F,F,T),
# note: the 4th data is somatic mutation which is a binary matrix
try.N.clust = 2:8,
# try cluster number from 2 to 8
fig.name= "CLUSTER NUMBER OF TCGA-BRCA")
这是第一个图，关于这个图的代码，我的理解是这样的：
data = mo.data: 这是用来进行聚类分析的多组学数据（如基因表达、DNA甲基化、突变等）。mo.data 是一个包含这些数据的变量。 is.binary = c(F, F, F, T): 这个参数用来指定每种组学数据是否是二进制格式（即只有0和1的值）。在这里，前面三组数据都不是二进制（F），而第四个数据（somatic mutation, 体细胞突变）是二进制数据（T）。 try.N.clust = 2:8: 表示尝试不同的聚类数量，从2到8。这个范围允许代码测试从2类到8类的数据分组，以确定最佳的聚类数量。 fig.name = "CLUSTER NUMBER OF TCGA-BRCA": 这是生成的图表名称，用于标识图的内容。在这里，图展示的是TCGA-BRCA（乳腺癌数据集）的最佳聚类数。
关于这个图的解释：图中展示了两个指标（Cluster Prediction Index和Gap-statistics），用于帮助选择合适的聚类数量。
Cluster Prediction Index（蓝色曲线，左轴）：这是用于评估聚类效果的一个指标。一般来说，这个值越高，聚类的效果越好。因此，我们希望找到聚类数下该值相对较高的点，以达到最佳聚类效果。 Gap-statistics（红色曲线，右轴）：Gap统计量是一种用来确定最佳聚类数的统计方法。它通过比较数据中的聚类情况和随机数据的聚类情况来判断聚类是否显著。通常，Gap统计量的值越高，表示聚类的显著性越高。 在图中，可以看到随着聚类数量的增加，这两个指标的变化趋势。我们通常会选择Cluster Prediction Index和Gap-statistics相对较高且稳定的聚类数。 从图中来看，聚类数量为3或5时，Cluster Prediction Index和Gap-statistics的值都相对较高和稳定，这可能表明这是合理的聚类数量。
iClusterBayes.res <- getiClusterBayes(data= mo.data,
N.clust= 5,
type= c("gaussian","gaussian","gaussian","binomial"),
n.burnin= 1800,
n.draw= 1200,
prior.gamma= c(0.5, 0.5, 0.5, 0.5),
sdev= 0.05,
thin= 3)
这段代码的解释是这样的：
data = mo.data: 输入的数据集，这里是 mo.data，即你的多组学数据。 N.clust = 5: 设定聚类的数量。这里设定为5，这个值应该是根据你之前的分析（比如上一个图中得到的最佳聚类数）来确定的。 type = c("gaussian","gaussian","gaussian","binomial"): 指定每个组学数据的分布类型。 "gaussian" 表示数据是连续的，服从正态分布（通常用于基因表达数据）。 "binomial" 表示二进制数据，适合用在体细胞突变等二值数据上。
在这里，你的多组学数据包含四种数据类型，前三种是连续的（正态分布），第四种是二值数据（如体细胞突变），所以这个参数的设置是合理的。
n.burnin = 1800: 表示在贝叶斯模型的马尔科夫链蒙特卡洛（MCMC）采样中，前1800次采样结果将被丢弃，作为“烧入期”（burn-in）。这个过程是为了让模型先稳定下来。 n.draw = 1200: 设定实际使用的采样次数为1200次。这些样本将用于最终的聚类分析。 prior.gamma = c(0.5, 0.5, 0.5, 0.5): 这是贝叶斯分析中的先验参数。这里设定了每种组学数据的gamma值为0.5，这表示我们在分析前认为每种数据对聚类的影响是均等的。 sdev = 0.05: 这个参数控制模型采样过程中的标准差，影响采样的步长。设置得较小有助于模型收敛，但也可能增加采样的时间。 thin = 3: 在采样过程中，每隔3个样本取一个数据点用于分析，这叫做“稀释”（thinning）。这样可以减少样本间的相关性，提高采样结果的独立性。
参数调整建议：
N.clust：这个值需要根据前面的聚类选择结果来设定，确保它是合理的。如果5是前面分析得到的最佳聚类数，那这个设置是合理的；如果另有最佳值，建议调整到对应的数量。 n.burnin 和 n.draw：这些参数的具体值会影响分析结果的准确性。一般来说，这两个参数的设置应该根据数据集的大小和复杂性来调整。如果数据量很大，可以适当增加 n.burnin 和 n.draw 的值，以保证采样稳定性；否则，当前的设置是合理的。 prior.gamma：如果你对不同数据类型对聚类的贡献有先验知识，可以调整这些值。例如，如果某个组学数据更重要，可以将其对应的 prior.gamma 值提高。 sdev 和 thin：如果模型收敛速度较慢，可以适当调低 sdev，或增大 thin 值。不过，当前的值一般是合理的默认设置，不建议随意更改。
图的解释
一致性矩阵：图中的蓝色方格代表样本间的相似性或一致性程度。蓝色越深，表示在多次聚类中，这些样本越常被分到同一类。理想情况下，每个类应该形成一个深蓝色的方块，与其他类有清晰的分界。这种模式说明聚类结果稳定且一致。 颜色编码的分组：图顶端和侧边的彩色条代表每个样本的类别（1到5类），对应于 N.clust=5 的设置。每种颜色表示一个不同的亚型。这些类别通过一致性矩阵确定，反映了数据中的自然分群结构。 层次聚类树：在图的顶部有一个树状结构（dendrogram），展示了样本的层次聚类关系。这个树图可以帮助更好地理解不同类别之间的关系。
通过上段的代码可以得到这个图，这个图也是反映了具体的聚类之间的关系
图的解释：
CDF 曲线：每条曲线代表不同的聚类数量（如2、3、4、5个聚类）。横轴是一致性指数（consensus index），表示样本间的一致性程度。纵轴是累积分布函数（CDF），反映了样本对该一致性指数的分布。 选择聚类数：通常情况下，CDF 曲线趋于平滑且接近1，表明该聚类数能较好地分开样本。通过比较不同聚类数的 CDF 曲线，可以判断哪一个聚类数的分群效果最佳。一般来说，曲线越平滑且曲线之间差异越小，表示聚类结果越稳定。
在这个图中，曲线颜色对应不同的聚类数（如红色代表2类，紫色代表5类）。可以看到，当聚类数增加时，曲线逐渐趋于平滑并向1靠近，这表明更大的聚类数可能提供更好的分群效果。在你的分析中，聚类数5可能是较优的选择，因为它的CDF曲线较平滑并接近1。
1.4取5个亚群进行分析
#为了和PAM50保持一致，取5个亚群（不同肿瘤，可以自己设定）
iClusterBayes.res <- getMOIC(data = mo.data,
N.clust = 5,
methodslist = "iClusterBayes",
# specify only ONE algorithm here
type = c("gaussian","gaussian","gaussian","binomial"),
# data type corresponding to the list
n.burnin = 1800,
n.draw = 1200,
prior.gamma = c(0.5, 0.5, 0.5, 0.5),
sdev = 0.05,
thin = 3)
这段代码通过 iClusterBayes 算法对多组学数据进行聚类，将所有数据整合到一起，分析它们之间的潜在分组关系。最终得到的 iClusterBayes.res 是包含5个聚类的结果，反映了数据中样本的结构。
参数调整建议
N.clust：这个值应根据之前的分析来设定。如果5类是最佳选择，可以保持不变。 n.burnin 和 n.draw：如果对结果的稳定性要求更高，可以适当增加这两个值。但目前的设置通常是合理的默认值。 prior.gamma：如果有先验知识，认为某种组学数据对聚类更重要，可以适当调整相应的 gamma 值。
prior.gamma 参数用来设置每种组学数据的先验权重，它反映了对每个组学数据在聚类中贡献的预期。默认值为 c(0.5, 0.5, 0.5, 0.5)，即假设每种数据对聚类的影响是均等的。 但是，如果你对各组学数据的重要性有先验知识，或者想引导模型更加关注某些特定的组学数据，可以通过调整这些数值来实现。建议的修改思路：
根据数据的质量和信噪比调整高质量的数据：如果某个组学数据的质量较高，信噪比好，且能较好地反映生物学特征，可以增加其权重，比如设为 0.7 或 0.8。低质量的数据：如果某个组学数据的噪声较大或测量误差较多，可以降低其权重，比如设为 0.3 或 0.2。 根据组学数据对生物学问题的相关性调整更重要的组学：如果某个组学数据与你的研究问题高度相关，例如基因表达数据在某些情况下比突变数据更能反映肿瘤的特征，你可以提升它的权重。次要的组学：如果某些组学数据对问题的相关性较低，可以降低其权重，使模型更多关注其他数据。 增量调整法初始设置：可以从均等权重（如当前的 c(0.5, 0.5, 0.5, 0.5)) 开始，然后尝试将某个数据类型的权重逐步增加或减少（比如将某个组学的数据权重调到 0.6 或 0.4）。观察结果：运行模型并观察聚类效果是否有显著提升。如果某个组学的数据对结果影响明显，可以逐步调整到更合适的值。 实验性调整 如果对各组学数据的重要性没有明确的先验知识，可以尝试多种权重组合，使用聚类结果的稳定性指标（如一致性矩阵或Gap统计量）来评估不同权重组合的效果，从而找到最佳的权重设置。根据你的数据特性和研究需求，prior.gamma 可以用来微调各组学数据的权重。合理的设置会让聚类结果更符合实际生物学意义。在不了解具体数据的情况下，建议先从均等权重开始，逐步调整。1.5继续完成聚类分析
moic.res.list <- getMOIC(data = mo.data,
methodslist = list("SNF", "PINSPlus",
"NEMO", "COCA",
"LRAcluster", "ConsensusClustering",
"IntNMF", "CIMLR", "MoCluster"), # 9种算法
N.clust = 5,
type = c("gaussian", "gaussian", "gaussian", "binomial"))
这段代码使用了getMOIC 函数对多组学数据进行集成聚类分析。
data = mo.data: 输入的多组学数据集，包含不同类型的数据，比如基因表达、DNA甲基化、突变等，存储在变量 mo.data 中。 methodslist = list("SNF", "PINSPlus", "NEMO", "COCA", "LRAcluster", "ConsensusClustering", "IntNMF", "CIMLR", "MoCluster"): 这是一个聚类算法列表，包含了9种常用的多组学聚类方法。通过组合多种算法，可以利用每种算法的优势，从而获得更可靠的聚类结果。 N.clust = 5: 设置聚类数为5。这个值可能是根据前面的聚类选择结果（比如一致性矩阵和一致性CDF分析）来确定的最佳聚类数。 type = c("gaussian", "gaussian", "gaussian", "binomial"): 指定每种组学数据的分布类型。 "gaussian" 表示数据是连续的，服从正态分布（适合用于基因表达数据）。 "binomial" 表示二进制数据，适用于体细胞突变等二值数据。
9种聚类的方法：
"SNF"：相似网络融合（Similarity Network Fusion），通过将多个组学数据融合为一个相似性网络来进行聚类。 "PINSPlus"：逐步识别非线性结构（Progressive Identification of Nonlinear Substructure），适用于处理稀疏和高维数据。 "NEMO"：集成的核方法，用于多组学数据的聚类分析。 "COCA"：基于共识的集群分析（Cluster Of Clusters Analysis），将多个聚类结果整合成一致的分型。 "LRAcluster"：低秩近似聚类（Low-Rank Approximation Clustering），适用于降维的分型方法。 "ConsensusClustering"：共识聚类方法，通过重复聚类得到一致性的结果。 "IntNMF"：集成非负矩阵分解（Integrated Non-negative Matrix Factorization），用于数据的矩阵分解聚类。 "CIMLR"：综合多层次高斯图模型（Comprehensive Integration of Multiple Layers of Reality）。 "MoCluster"：多组学数据的集成聚类方法。 通过集成多种算法，最终的聚类结果会更稳健。
这里指定的类型表示前三种数据是连续变量（如基因表达），第四种数据是二值变量（如突变数据）。这段代码执行多组学集成聚类分析，将不同类型的数据和多种聚类方法结合，以提高分型的稳定性和准确性。通过组合多种算法（如 SNF、ConsensusClustering 等），可以得到更稳健的分型结果，减少单一算法的偏差带来的不稳定性。主要目的是使用多种聚类算法对数据进行集成聚类分析，从而得到更可靠的5个亚型（分群）结果。输出的 moic.res.list 会包含所有算法的聚类结果，你可以进一步分析不同算法之间的相似性，并得出最终的聚类共识结果。
**图的解释
纵轴（k）：表示聚类的数量，从2到5（对应代码中的 N.clust = 5）。每一层表示在特定的聚类数下样本的分组情况。
横轴（samples）：表示每个样本。不同颜色的方块表示样本在不同聚类数下的类别分配。
颜色块：不同颜色表示不同的聚类类别。随着聚类数（k）的增加，每个样本可能会被分配到不同的类别，因此颜色也会发生变化。追踪图的主要作用是展示当聚类数从2增加到5时，每个样本在不同分群中的归属是否发生变化。颜色变化较少的样本说明它们的分群比较稳定，而颜色变化频繁的样本则说明它们在不同的k值下可能归属不同的类别。
如果图中大部分样本在 k=5 时形成了清晰的分群，这表明5个聚类可能是合理的分群选择。如果颜色变化过于频繁，说明样本分配不够稳定，可能需要重新评估合适的聚类数。
moic.res.list <- append(moic.res.list,
list("iClusterBayes" = iClusterBayes.res))
# 保存下结果
save(moic.res.list, file = "moic.res.list.rda")
然后我们把这个图保存一下，明天我们继续分析一下并且跑后面的图。拜拜。
