1.介绍

CancerSubtypes是一个用于癌症亚型分析的包，包括从数据集处理到结果验证的各种功能。在CancerSubtypes软件包，我们提供了从原始数据到结果可视化分析癌症亚型的统一框架。主要功能包括基因组数据预处理、癌症亚型识别、结果验证、可视化和比较。CancerSubtypes为基因组数据预处理提供了常用的数据归一化方法。同时，有四种特征选择方法来筛选基因组数据集中的关键特征。常见的癌症亚型识别方法集成在这个包中，例如共识聚类(CC)[来自R包ConsensusClusterPlus]，非负矩阵分解(CNMF)[来自R包NMF]，综合集群(iCluster)[来自R包iCluster]，相似网络融合[来自R包SNFtool]，联合SNF和CC (SNF.CC)和加权相似网络融合。我们以统一的输入和输出数据格式实现这些癌症亚型识别方法。分析癌症亚型的过程可以很容易地在标准工作流程中进行。CancerSubtypes提供了最有用的特征选择方法和亚型验证方法，帮助用户专注于他们的癌症基因组数据，不同方法的结果可以很容易地进行可视化的比较和评估。

2数据处理

对于基本数据处理，CancerSubtypes提供了数据分布检查、归一化和特征选择的方法。研究中有四种特征选择方法(方差- var、中值绝对偏差- mad、COX模型、主成分分析- pca)CancerSubtypes包中。所有数据处理方法都具有相同的输入和输出数据格式。

2.1基本处理

2.1.1通过检查数据分布对原始数据进行分析。

准备一个TCGA基因表达数据集进行分析。library(CancerSubtypes) library("RTCGA.mRNA") rm(list = ls()) data(BRCA.mRNA) mRNA=t(as.matrix(BRCA.mRNA[，-1]) colnames(mRNA)=BRCA. mRNA)mRNA[，1] ###观察数据集的均值、方差和中位数绝对偏差分布，有助于用户获得数据的分布特征，例如评估数据集是否符合正态分布。data.checkDistribution (mRNA)

2.1.2缺失值特征的数据赋值

原始基因组数据集总是包含缺失的观测数据，特别是在微阵列基因表达数据中。在很少的样本中删除所有缺失观测的特征是不明智的，因为有用的信息将被丢弃。常用的方法是为缺失的观测值赋值。CancerSubtypes集成了基因组数据集的三种常用imputation方法。

index=which(is.na(mRNA)) res1=data.imputation(mRNA,fun="median") res2=data.imputation(mRNA,fun="mean") res3=data.imputation(mRNA,fun="microarray")

2.1.3数据归一化。

result1 = data.normalization编写此表达式(mRNA,类型=“feature_Median log2 = FALSE) result2 = data.normalization (mRNA、类型=“feature_zscore log2 = FALSE)

2.2特征选择

2.2.1基于最大方差的特征选择。

###前1000个方差最大的特征将被选择。data1=FSbyVar(mRNA, cut.type="topk"，value=1000) ###选择(方差>0.5)的特征。data2 = FSbyVar (mRNA, cut.type =“截止”,值= 0.5)

2.2.2基于最大变异中值绝对偏差(MAD)的特征选择。

(mRNA, cut.type="topk"，value=1000) data2=FSbyMAD(mRNA, cut.type="cutoff"，value=0.5)

2.2.3基于主成分分析的特征降维与提取。

mRNA1=data.imputation(mRNA,fun="microarray") data1=FSbyPCA(mRNA1, PC_percent=0.9,scale = TRUE)

2.2.4基于Cox回归模型的特征选择。

data(GeneExp,time,status) data1=FSbyCox(GeneExp,time,status,cutoff=0.05)

3.肿瘤亚型鉴定的聚类方法。

3.1肿瘤亚型鉴定的共识聚类

共识聚类(CC, 2003)作为一种无监督的亚型发现方法，是许多基因组研究中常用的和有价值的方法，并有许多成功的应用。

输入数据集为单基因表达矩阵。data(GeneExp) result=ExecuteCC(clusterNum=3,d=GeneExp,maxK=10,clusterAlg="hc"，distance="pearson"，title="GBM") ###输入数据集是多基因组数据作为列表数据(GeneExp) data(miRNAExp) GBM=list(GeneExp=GeneExp,miRNAExp=miRNAExp) result=ExecuteCC(clusterNum=3,d=GBM,maxK=10,clusterAlg="hc"，distance="pearson"，title="GBM")

非负矩阵分解用于癌症亚型鉴定

非负矩阵分解(CNMF, 2004)作为一种有效的降维方法，被用于区分高维基因组数据的分子模式，为类发现提供了一种强有力的方法。我们应用NMF包执行癌症基因组数据集的非负矩阵分解。因此，该方法允许用户输入并行处理的核心cpu数量。

输入数据集为单基因表达矩阵。data(GeneExp) result=ExecuteCNMF(GeneExp,clusterNum=3,nrun=30) ###输入数据集是多基因组数据作为列表数据(GeneExp) data(miRNAExp) GBM=list(GeneExp=GeneExp,miRNAExp=miRNAExp) result=ExecuteCNMF(GBM,clusterNum=3,nrun=30)

3.3肿瘤亚型的整合聚类鉴定

整合聚类(iCluster, 2009)使用联合潜变量模型对多种类型的组学数据进行整合聚类。

data(GeneExp) data(miRNAExp) data1=FSbyVar(GeneExp, cut.type="topk"，value=1000) data2=FSbyVar(miRNAExp, cut.type="topk"，value=300) GBM=list(GeneExp=data1,miRNAExp=data2) result=ExecuteiCluster(datasets=GBM, k=3, lambda=list(0.44,0.33,0.28))

3.4相似网络融合用于癌症亚型识别

相似网络融合(SNF, 2014)是一种融合相似网络的计算方法，用于聚合多组数据。

data(GeneExp) data(miRNAExp) GBM=list(GeneExp=GeneExp,miRNAExp=miRNAExp) result=ExecuteSNF(GBM, clusterNum=3, K=20, alpha=0.5, t=20)

3.5 SNF与CC的集成方法用于癌症亚型鉴定

我们提出将SNF和CC结合起来，形成一种新的癌症亚型识别方法。

data(GeneExp) data(miRNAExp) data(time) data(status) data1=FSbyCox(GeneExp,time,status,cutoff=0.05) data2=FSbyCox(miRNAExp,time,status,cutoff=0.05) GBM=list(GeneExp=data1,miRNAExp=data2) result=ExecuteSNF。CC(GBM, clusterNum=3, K=20, alpha=0.5, t=20,maxK = 10, pItem = 0.8,reps=500, title =" GBM"， plot =" png"， finalLinkage ="average")

3.6加权相似网络融合用于癌症亚型识别

WSNF是一种借助基因调控网络信息进行癌症亚型识别的方法。它利用miRNA-TF-mRNA调控网络来考虑特征的重要性。

data(GeneExp) data(miRNAExp) data(Ranking) ####检索基因的特征排名gene_Name=rownames(GeneExp) index1=match(gene_Name,Ranking$mRNA_TF_miRNA.v21_SYMBOL) gene_ranking=data.frame(gene_Name,Ranking[index1，]，stringsAsFactors=FALSE) index2=which(is.na(gene_ranking$ranking_default)) gene_ranking$ranking_default[index2]=min(gene_ranking$ranking_default,na。rm =TRUE) ####检索there特征排序基因miRNA_ID=rownames(miRNAExp) index3=match(miRNA_ID, ranking$ mRNA_TF_miRNA_ID) miRNA_ranking=data.frame(miRNA_ID, ranking [index3，]，stringsAsFactors=FALSE) index4=which(is.na(miRNA_ranking$ranking_default)) miRNA_ranking$ranking_default[index4]=min(miRNA_ranking$ranking_default,na。rm =TRUE) ###集群ranking1=list(gene_ranking$ranking_default,miRNA_ranking$ranking_default) GBM=list(GeneExp,miRNAExp) result=ExecuteWSNF(datasets=GBM, feature_ranking=ranking1, beta = 0.8, clusterNum=3, K = 20,alpha = 0.5, t= 20, plot =TRUE)

4已确定的癌症亚型的结果验证、解释和可视化。

计算方法鉴定的癌症亚型应符合生物学意义，并揭示其独特的分子模式。

4.1轮廓宽度

轮廓宽度用于衡量一个样本与其识别的子类型相比与其他子类型匹配的相似程度，高值表示样本匹配良好。每条水平线代表剪影图中的一个样本。线的长度是样本的轮廓宽度。

data(GeneExp) data(miRNAExp) GBM=list(GeneExp=GeneExp,miRNAExp=miRNAExp) result=ExecuteSNF(GBM, clusterNum=3, K=20, alpha=0.5, t=20,plot = FALSE) ###相似性样本矩阵sil=廓形相似矩阵(result$group, result$distanceMatrix) plot(sil)

注:如果输入矩阵是样本间的不相似矩阵，请使用轮廓()在集群包中计算轮廓宽度，否则会产生错误的结果。

sil1=silhouette(result$group, result$distanceMatrix) plot(sil1) ##错误的结果

所有样本的轮廓宽度均为负。

4.2生存分析

生存分析用于判断亚型之间的不同生存模式。

data(GeneExp) data(miRNAExp) data(time) data(status) data1=FSbyCox(GeneExp,time,status,cutoff=0.05) data2=FSbyCox(miRNAExp,time,status,cutoff=0.05) GBM=list(GeneExp=data1,miRNAExp=data2) ####ExecuteSNF result1=ExecuteSNF(GBM, clusterNum=3, K=20, alpha=0.5, t=20,plot = FALSE) group1=result1$group distanceMatrix1=result1$distanceMatrix p_value=survAnalysis(mainTitle="GBM1"，time,status,group1, distanceMatrix1,similarity=TRUE)
## ## ***************************************************** ## GBM1集群= 3叫:# # survdiff(公式= Surv(时间、状态)~集团)# # # # N观察预期(执着)^ 2 / E(执着)^ 2 / V # #集团46 44 = 1 41.3 0.173 0.311 # #集团40 39 27.3 4.990 7.387 = 2 # #集团14 13 27.4 7.529 11.661 = 3 # # # # 2自由度Chisq = 14.2, p = 8 e-04

这是一个由三部分组成的组合图:生存曲线、样本相似矩阵热图和已识别癌症亚型的廓形宽度图。所有图中的样本都已根据确定的癌症亚型进行了重组。这种图提供了容易评估的可见结果。

# # # # 2. executesnf。CC result2 = ExecuteSNF。CC(GBM, clusterNum=3, K=20, alpha=0.5, t=20, maxK =5, pItem = 0.8,reps=500, title =" GBM2"， plot =" png"， finalLinkage ="average")
group2=result2$group distanceMatrix2=result2$distanceMatrix p_value=survAnalysis(mainTitle="GBM2"，时间，状态，group2, distanceMatrix2，相似性=TRUE)
## ## ***************************************************** ## GBM2集群= 3叫:# # survdiff(公式= Surv(时间、状态)~集团)# # # # N观察预期(执着)^ 2 / E(执着)^ 2 / V 40 39 # #集团= 1 27.3 4.990 7.387 # #集团44 41.3 0.173 0.311 = 2 46 # #集团14 13 27.4 7.529 11.661 = 3 # # # # 2自由度Chisq = 14.2, p = 8 e-04

4.3聚类的统计学意义

聚类的统计显著性是一种纯统计方法来检验子类型之间的显著性差异数据分布。不同表达是测试每个子类型与一个参考组(总是一组正常样本)之间的表达差异。

data(GeneExp) data(miRNAExp) data(time) data(status) GBM=list(GeneExp=GeneExp,miRNAExp=miRNAExp) result=ExecuteSNF(GBM, clusterNum=3, K=20, alpha=0.5, t=20,plot = FALSE) group=result$group sigclust=sigclustTest(miRNAExp,group, nsim=1000, nrep=1, icovest=1) sigclust
子类型1子类型2子类型3 ##子类型1 1.000 0.016 0.00 ##子类型2 0.016 1.000 0.23子类型3 0.000 0.230 1.00

SigClust汇总图显示了聚类指数(CI)的分布。蓝色的点表示模拟的ci，用随机的垂直抖动绘制，以便更好地可视化。实线和虚线对应于估计的非参数密度，高斯密度与模拟ci相吻合。p值表示两个子类型之间的显著性水平。更多信息请参考[5]sigclust。

4.4差异表达分析

差异表达分析是检测每个亚型与参考组(总是一组正常样本)之间的表达差异。这里我们应用limma包进行各亚型与正常样本之间的不同表达分析。

library("RTCGA.mRNA") #require(TCGAbiolinks) rm(list = ls()) data(BRCA.mRNA) mRNA=t(as.matrix(BRCA.mRNA[，-1])) colnames(mRNA)=BRCA. mRNA。mRNA[，1] mRNA1=data.imputation(mRNA,fun="microarray") mRNA1=FSbyMAD(mRNA1, cut.type="topk"，value=5000)
###拆分正常和肿瘤样本index=which(as.numeric(substr(colnames(mRNA1)，14,15))>9) mRNA_normal=mRNA1[，index] mRNA_tumor=mRNA1[，-index] ###删除重复样本index1=which(as.numeric(substr(colnames(mRNA_tumor)，14,15))>1) mRNA_tumor=mRNA_tumor[，-index1] #####识别癌症亚型result=ExecuteCC(clusterNum=5,d=mRNA_tumor,maxK=5,clusterAlg="hc"，distance="pearson"，title="BRCA") group=result$groupres = DiffExp.limma (Tumor_Data = mRNA_tumor Normal_Data = mRNA_normal组=组,topk = NULL, RNAseq = FALSE)
1亚型头部不同表达基因(res[[1]])
## ID logFC AveExpr t P.Value adj.P.Val ## 2439 COL10A1 4.712350 5.0043839 38.57960 3.234769e-151 1.617385e-147 ## 3725 MMP11 3.528733 2.6563911 34.61589 3.206291e-134 8.015157e -131 ## 3797 CA4 -4.357770 -1.0831256 -27.32458 7.230942e-101 1.205157e-97 ## 3235 CAV1 -3.006418 -1.5612746 -26.90504 7.051163e-99 8.813954e-96 ## 4586 CXCL3 -2.499337 -2.9294699 - 26.675861e -96 ## 4133 RYR3 -2.575299 -0.4590268 -25.94906 2.103670e-91 ## B ## 2439 334.8796 ## 3725295.9708 ## 3797 219.5287 ## 3235 214.9674 ## 4586 207.7722 ## 4133 204.5236
2亚型头部不同表达基因(res[[2]])
## ID logFC AveExpr t P.Value adj.P.Val B ## 637 ITM2A -6.173926 0.8088468 -12.356546 1.091109e-18 5.455546e-15 31.11131 ## 4942 CKAP2 4.229372 -2.1539566 9.959475 1.093416e-14 2.267182e-11 22.64778 # 4596 RCBTB2 -4.492770 -0.2484435 -9.904508 1.360309e-14 2.267182e-11 22.44505 ## 3297 GPR19 4.272492 -2.2195806 9.791601 2.132269e-14 2.665337e-11 22.02758 ## 4829 CASP4 -3.767096 0.7556694 -9.613507 4.342293e-11 21.36631 ## 2065 MCM4 5.079992 -3.2200565 9.518389 6.355510e-145.296258 e-11 21.01180

5的结论

的CancerSubtypesR包提供了一套癌症亚型分析工具，并将分析嵌入到标准化的工作流程中。它为在全基因组范围内分析癌症亚型提供了有力的方法。

6.参考文献

[1] Monti, Stefano等，“共识聚类:一种基于重采样的基因表达微阵列数据的类发现和可视化方法。”机器学习52.1-2(2003):91-118。

[2] Brunet, Jean-Philippe等，“利用矩阵分解发现元基因和分子模式。”国家科学院学报101.12(2004):4164-4169。

[3]沈荣来，亚当·b·奥尔申，马克·拉达尼。“使用联合潜在变量模型对多种基因组数据类型进行综合聚类，并应用于乳腺癌和肺癌亚型分析。”生物信息学25.22(2009):2906-2912。

[4] Wang, Bo等，“在基因组规模上聚合数据类型的相似网络融合。”自然方法11.3(2014):333-337。

[5]刘玉峰等。“高维、低?c样本量数据。”美国统计协会杂志(2012)。

[6]卢梭，彼得J.“剪影:对聚类分析的解释和验证的图形辅助。”计算与应用数学杂志20(1987):53-65。

[7]徐涛，乐天东，刘林，等。基因表达与肿瘤分型的关系[J]。公共科学学报，2016,11(4):e0152792。