内容

1简介

的GeneExpressionSignaturePackage利用基因表达谱来测量不同生物状态之间的相似性。给出了两种相似度度量算法:文中提到的GSEA算法(Iorio et al. 2010)的PGSEA算法PGSEA包中。基于基因表达特征的测量方法的进一步描述可以在Lamb中找到(Lamb et al. 2006),胡(Hu and Agarwal 2009)和人工(Iorio et al. 2010)．

本手册是一个简单的介绍结构，功能和使用GeneExpressionSignature包中。它展示了如何通过一系列计算步骤确定生物相似性，以及如何将这些信息用于进一步的聚类分析。

当前版本GeneExpressionSignature只能用于来自同一平台的数据，例如HG-U133A平台。

2开始

一个完整的分析程序接受一组表示不同生物状态的基因表达谱作为输入，并生成相似矩阵作为输出。它可以分为三个步骤:1)数据排序，2)秩归并，3)相似度测量。

首先，我们通过在R会话中输入以下命令来加载包:

库(GeneExpressionSignature)

2.1数据排名

基因表达谱在分析前应进行适当的预处理，包括背景校正、归一化和总结。下面的程序使用基因表达水平的等级，而不是精确的值。首先对微阵列探针集标识符按照不同的表达值(计数或比值)进行排序，得到排序的基因列表。需要注意的是，数据预处理没有标准的方法，有一个函数getRLs ()其中采用C-MAP中的方法进行数据预处理，仅供参考。我们可以通过调用来获得排序列表矩阵getRLs ()．

你的实验数据可以用来分析，或者用户可以用R包从GEO数据库下载基因表达谱GEOquery．用户可以在GEOquery欲知详情。

作为一个例子，我们从GEO数据库下载数据包GEOquery．然后将处理表达值组合成处理矩阵和对照表达值。

#如果您有网络访问#GSM118720 <- getGEO('GSM118720') # GSM118721 <- getGEO('GSM118721') If (require(GEOquery)){#处理基因表达配置文件GSM118720 <- getGEO(filename = system. zip)。文件(“extdata / GSM118720。GSM118721 <- getGEO(filename=system. exe);文件(“extdata / GSM118721。control <- as.matrix(as.numeric(Table(GSM118721)[， 2])) treatment <- as.matrix(as.numeric(Table(GSM118720)[， 2])) ranked_list <- getRLs(control, treatment)}

2．2合并排序

通过秩归并，将多个排序列表合并为一个排序列表，称为原型排序列表(prototype ranking list, PRL)，表示某种生物状态。该程序主要在相似性测量之前进行，适用于将多个排位表分配给一个单一生物状态，具有不同细胞类型或实验条件的特定情况。

但应考虑两种不同的情况:1)所有具有相同生物状态的排名表都同等重要;2)每个排名列表都有自己的排名权重。这个包提供了两种常用的算法:一种利用Iorio提出的Kruskal算法(Iorio et al. 2010)对于前者和另一种情况，平均排名技术是一种简单但相当有用的方法。函数RankMering ()用于根据其表型数据将排序列表聚合为一个或多个prl。我们要做的就是构造aExpressionSet对象作为输入，排序列表作为分析数据，相应的生物状态作为表型数据。

为方便起见，将数据存储为ExpressionSet类eset对象作为输入数据，带有排序列表(通过调用getRLs ())作为分析数据，相应的生物学状态作为表型数据。作为一个例子，我们从加载培养的exampleSet数据开始，它是C-MAP的一个子集(Lamb et al. 2006)作为样本数据，这是一个大型参考目录的基因表达数据从培养的人类细胞扰动许多化学物质和遗传试剂。子数据集由50对基因表达谱组成，涉及22283个基因。这些剖面分别从处理过15种化合物的细胞中获得，这些化合物的值已经转换为等级顺序。

data(exampleSet) show(exampleSet) #> ExpressionSet (storageMode: lockedEnvironment) #> assayData: 22283 features, 50 samples #> element names: exprs #> protocolData: none #> phenoData #> sampleNames: 1 2…50 (50 total) #> varLabels: state #> varMetadata: labelDescription #> featureData: none #> experimentData: use 'experimentData(object)' #>注释:exprs(exampleSet)[c(1:10)， c(1:3)] #> 123 #> 1 11264 14408 13919 #> 2 12746 12365 3080 #> 3 8267 5630 13060 #> 4 2193 16694 16084 #> 5 9556 6044 8294 #> 6 279 5120 4826 #> 7 15381 10225 10883 #> 8 9452 10777 13359 #> 9 6149 6213 6800 #> 10 4943 12760 3444 levels(as(phenoData(exampleSet)， "data.frame")[，]1]) # >[1]“alsterpaulone”“氮胞苷”“喜树碱”“白菊素”# >[5]“柔红霉素”“阿霉素”“ellipticine”“乙丙酸”# >[9]“非瑟酮”“harmine”“木犀草素”“mitto蒽醌”#>[13]“孤雌素”“staurosporine”“硫代链霉素”

Rank归并过程将从50对表达谱中生成15个PRL的合并集，每个PRL分别对应15种化合物中的一种。

MergingSet <- rankmerged (exampleSet， "Spearman"，加权= TRUE) show(MergingSet) #> ExpressionSet (storageMode: lockedEnvironment) #> assayData: 22283 features, 15 samples #> element names: exprs #> protocolData: none #> phenoData #> sampleNames: alsterpaullone azacitidine…thiostrepton (15 total) #> varLabels: state #> varMetadata: labelDescription #> featureData: none #> experimentData: use 'experimentData(object)' #>注释:

2．3相似度测量

通过秩归并得到一个状态的单一组合PRL。这些prl被用来通过评分函数来测量不同生物状态下基因特征的相似性ScoreGSEA ()而且ScorePGSEA ()．并不是所有的基因都参与相似性测量，只有被称为基因特征的基因子集，其组合表达模式是生物状态的唯一特征。在相似性评分过程中，作为基因标记的基因通常是由先验知识预先定义的。人工(Iorio et al. 2010)提出了一种以上调基因最多和下调基因最多作为基因签名的“最优签名”方法。基因特征的大小需要考虑，这是在相似性度量中除了prl之外的另一个参数。在大多数情况下，基因签名的默认大小是250全基因组表达谱。

假设N是prl的数量(也与生物状态的数量相同)，则N × N通过相似度测量生成距离矩阵。对于mergingSet，我们将得到一个15 x 15矩阵对应于这些化合物之间的相似距离。

ds <- ScoreGSEA(MergingSet, 250， "avg") ds[1:5, 1:5] #> alsterpaullone azacitiine喜树碱白苷daunorubicin #> alsterpaullone 0.0000000 0.6176992 0.4669311 0.6896005 0.6895031 0.8515960 0.6413233 0.5372661 #>喜树碱0.4669311 0.6125031 0.0000000 0.7897938 0.5372661 #> daunorubicin 0.5288110 0.6413233 0.5372661 0.7443612 0.0000000 #> daunorubicin 0.5288110 0.6413233 0.5372661 0.7443612

2.4签名的距离

如上所述，四种算法作为函数实现getRLs ()，RankMerging ()，ScoreGSEA (),ScorePGSEA ()其中，一种算法用于数据预处理，一种叫做Iorio算法用于秩归并，另外两种算法叫做GSEA和PGSEA用于相似度度量。此外,函数SignatureDistance ()是提供作为一个单一的入口和方便的访问点，排名合并和相似度测量，这贯穿
包括排名合并和评分，建议在大多数情况下使用。由于没有标准的数据预处理方法和基因表达数据类型的不确定性，没有将数据排序集成到该函数中。此外，还没有有效的方法来整合来自不同平台的数据。函数getRLs ()其中采用C-MAP中的方法进行数据预处理，仅供参考。

签名距离(exampleSet，签名长度= 250,MergingDistance = "Spearman"，评分方法= "GSEA"，评分距离= "avg"，加权= TRUE) #> alsterpaullone偶胞嘧啶喜树碱白菊柔红霉素#> alsterpaullone 0.0000000 0.6176992 0.4669311 0.6896005 0.8515960 0.8597938 0.56413233 #>喜树碱0.4669311 0.6125031 0.0000000 0.7897938 0.5372661 #>柔红霉素0.56896005 0.8515960 0.7897938 0.0000000 0.7443612 #>多柔红霉素0.4449537 0.223770 0.5590938 0.4805674 #>椭圆素0.61471760.6958627 0.6060621 0.8399921 0.6013959 0.9625380 0.9150898 0.8846840 0.9174653 #>非瑟汀0.6191321 0.7401457 0.7258576 0.9056164 0.7204894 #>哈明0.7381854 0.9011707 0.8082408 0.6264392 0.7486181 #>木犀草素0.6601723 0.8357249 0.65869045 0.4627429 0.6882700 #>米托蒽醌0.5351687 0.6326825 0.6586904 0.8367069 0.5450045 #>单品内酯0.9183664 0.8581793 0.8943531 0.8865647 0.8487120 #>星孢素0.6984201 0.6982204 0.7120952 0.9037265 0.7625792 #>硫链蛋白0.9258501 0.8624173 0.8523135 0.9235488 0.8449146 #>阿霉素椭圆药乙烯酸非塞丁harmine #> alsterpaullone 0.4449537 0.6147176 0.9546259 0.6191321 0.7401457 0.9011707 #>喜树碱0.5590938 0.6060621 0.9150898 0.7258576 0.8082408 #>白菊0.8152383 0.8399921 0.8846840 0.9056164 0.6264392 #>阿霉素0.0000000 0.5737439 0.95905890.6130763 0.8146797 #>椭圆素0.5737439 0.0000000 0.9130729 0.8065379 0.9558315 0.9745611 #>非甾体素0.6130763 0.8065379 0.9558315 0.0000000 0.9077328 #>木犀草素0.7326149 0.7415050 0.84137980 0.8872799 0.5954765 #>米托蒽醌0.4266442 0.6073227 0.9880406 0.7107602 0.8455741 #>单胞内酯0.9058101 0.8260994 0.6586242 0.9955009 0.9361217 #> staurosporine 0.67291080.7785568 0.9677592 0.7667645 0.9525701 #>柔红霉素0.8801426 0.8887309 0.7367122 0.9806517 0.9818760 #> alsterpaullone 0.6601723 0.5351687 0.9183664 0.6984201 0.9258501 #>氮胞苷0.8357249 0.6326825 0.8581793 0.6982204 0.8624173 #>喜树碱0.6559045 0.6586904 0.8943531 0.7120952 0.8523135 #>白菊0.4627429 0.8367069 0.8865647 0.9037265 0.9235488 #>柔红霉素0.6882700 0.5450045 0.8487120 0.7625792 0.8449146 #>阿霉素0.7326149 0.4266442 0.9058101 0.6729108 0.8801426 #>椭圆药0.7415050 0.6073227 0.8260994 0.7785568 0.8887309 #>乙烯酸0.8413584 0.9880406 0.6586242 0.9677592 0.7367122 #>非甾体素0.8872799 0.7107602 0.9955009 0.9667645 0.9806517 #>木犀草素0.5954765 0.8455741 0.9361217 0.9525701 0.9842047 #>米托蒽醌0.7964905 0.0000000 0.9181195 0.7469812 0.8984673 #>单胞内酯0.8139365 0.9181195 0.0000000 0.98642200.6180173 #> staurosporine 0.9017486 0.7469812 0.9864220 0.0000000 0.9521839 #> thiostrepton 0.8642047 0.8984673 0.6180173 0.9521839 0.0000000

3.实现细节

3．1自适应加权秩归并

Iorio的秩归并算法采用了Kruskal算法(Cormen, Leiserson, and Rivest 1990)将同一生物状态对应的排序列表合并。首先必须计算这些排序表之间的距离，使用“Spearman”算法或“Kendall tau”算法计算两个排序表之间的距离度量。需要注意的是，用Kendall tau距离进行秩归并比较耗时，所以我们建议选择“Spearman”距离。接下来，使用“Borda合并”算法合并两个或多个具有相同生物状态的排名列表。

根据“Kruskal”算法方法(Cormen, Leiserson, and Rivest 1990)时，该排序归并算法首先寻找Spearman’s Footrule距离最小的两个排序表，然后使用Borda归并法将它们归并，得到一个新的排序表。最后，新列表替换两个未合并的列表。这个过程直到只剩下一个列表时才会终止。

为了方便，用户可以通过函数直接获取每个状态的PRLRankMerging ()，该算法使用“sprarman”、“BordaMerging”和“Kruskal”算法，将获得的具有相同生物状态的排名列表进行聚合。例如，我们将50个样本数据合并为15个样本。

MergingSet <- rankmerged (exampleSet， "Spearman"，加权= TRUE) show(MergingSet) #> ExpressionSet (storageMode: lockedEnvironment) #> assayData: 22283 features, 15 samples #> element names: exprs #> protocolData: none #> phenoData #> sampleNames: alsterpaullone azacitidine…thiostrepton (15 total) #> varLabels: state #> varMetadata: labelDescription #> featureData: none #> experimentData: use 'experimentData(object)' #>注释:

3.1.1等加权秩归并

解决这个问题的一个简单但相当有用的方法是平均排名技术。当我们假设每个排名列表的重要性是相同的权重时，这个技术是一个两步的过程。第一步是计算每个排名列表的平均排名，然后第二步是构建他们的最终排名。

3.1.2相似度测量

将具有相同生物状态的排序列表合并为单个PRL后，采用基因集富集分析(Gene Set Enrichment Analysis, GSEA)和参数基因集富集分析(Parametric Gene Set Enrichment Analysis, PGSEA)来衡量这些PRL之间的相似性。

GSEA算法(Subramanian et al. 2005)是一种非参数的、基于秩的相似性度量方法，用于确定先验定义的基因集是否在两个生物状态之间显示出统计上显著的、和谐的差异，而PGSEA算法是一种基于参数统计分析模型的改进的基因集富集分析方法。这两个函数都给出了相应的p值ScoreGSEA ()从蒙特卡罗程序中计算经验p值(North, Curtis, and Sham 2002)．

ds <- ScoreGSEA(MergingSet,250，"avg") ds[1:5,1:5] #> alsterpaullone azacitiine喜树碱白苷daunorubicin #> alsterpaullone 0.0000000 0.6176992 0.4669311 0.6896005 0.6895031 0.8515960 0.6413233 0.5372661 #>喜树碱0.4669311 0.6125031 0.0000000 0.7897938 0.5372661 #> daunorubicin 0.5288110 0.6413233 0.5372661 0.7443612 0.0000000 #> daunorubicin 0.5288110 0.6413233 0.5372661 0.7443612
ds <- ScorePGSEA(MergingSet,250，"avg") ds[1:5,1:5] #> alsterpaullone aztothecin daunorubicin #> alsterpaullone 0.0000000 0.5477505 0.4136914 0.6340643 0.4182223 #> azacitidine 0.5477505 0.0000000 0.5646674 0.8402056 0.7438953 0.4305080 #> daunorubicin 0.4182223 0.5478599 0.4305080 0.7067854 0.0000000 #>

4进一步分析

来说明如何使用GeneExpressionSignature在基因表达特征分析中，亲和传播聚类可以利用基因特征的相似性对这些生物状态进行分组。亲和传播聚类算法通过最大化一个称为净相似度的目标函数来迭代搜索最优聚类。这里，我们用in函数apcluster包装将15种生物状态分为3组。在这一步，R包apcluster也应该安装在你的电脑上。

如果(要求(apcluster)){库(apcluster) clusterResult < - apcluster (1 - ds)显示(clusterResult)} # >加载所需的包:apcluster # > # >附加包:“apcluster”# >以下对象是蒙面的包:统计数据:# > # # # > > >的热图APResult对象# > # >样本的数量= 15 # > = 122 #的迭代次数> = 0.2561047 # >输入偏好和相似之处= 6.432731 # >首选项之和= 0.768314 # = 7.201045 # > >净相似集群的数量= 3 # > # >范本:星团:b>星团1，样例阿霉素:> alsterpaulone氮胞苷camptothecin柔红霉素多柔比星> ellipticine非瑟酮mitto蒽醌staurosporine >星团2，样例木甲素:>白菊harmine木甲素#>星团3，样例木甲素:#>乙丙酸孤雌素硫代链素

利用细胞图对聚类结果进行可视化。在网络中，节点表示不同的化合物(用不同化合物处理的细胞状态)，边缘表示这两种化合物之间的相似距离低于阈值，这里为0.68。不同的颜色表示不同的基团，如化合物的分类。我们注意到，最大的组编号为9个节点，另外两个组由每个组3个节点组成。

会话信息

sessionInfo() #> R正在开发中(不稳定)(2022-10-25 r83175) #>平台:x86_64-pc-linux-gnu(64位)#>运行在:Ubuntu 22.04.1 LTS #> #>矩阵产品:默认#> BLAS: /home/biocbuild/bbs-3.17-bioc/R/lib/libRblas。so #> LAPACK: /usr/lib/x86_64-linux-gnu/ LAPACK /liblapack.so.3.10.0 #> #> locale: #> [1] LC_CTYPE=en_US。UTF-8 LC_NUMERIC= c# > [3] LC_TIME=en_GB LC_COLLATE= c# > [5] LC_MONETARY=en_US。utf - 8 LC_MESSAGES = en_US。UTF-8 #> [7] LC_PAPER=en_US。UTF-8 LC_NAME= c# > [9] LC_ADDRESS=C LC_TELEPHONE= c# > [11] LC_MEASUREMENT=en_US。UTF-8 LC_IDENTIFICATION= c# > #>附加基础包:#> [1]stats graphics grDevices utils datasets methods base #> #>其他附加包:#> [1]apcluster_1.4.10 GEOquery_2.67.0 #> [3] Biobase_2.59.0 BiocGenerics_0.45.0 #> [5] GeneExpressionSignature_1.45.0 BiocStyle_2.27.0 # b> #>通过命名空间加载(且未附加):#> [1] sass_0.4.2 utf8_1.2.2 generics_0.1.3 #> [31] assertthat_0.2.1 #> [7] stringi_1.7.8 hms_1.1.2 digest_0.6.30 #> [10] magrittr_2.0.3 grid_4.3.0 evaluate_0.17 #> [13] bookdown_0.29 fastmap_1.1.0 jsonlite_1.8.3 #> [19] Matrix_1.5-1 DBI_1.1.3 limma_3.55.0 #> [19] BiocManager_1.30.19 purrr_0.3.5 fansi_1.0.3 #> [22] jquerylib_0.1.4 cli_3.4.1 rlang_1.0.6 #> [25] ellipsis_0.3.2 cachem_1.0.6 #> [31] assertthat_0.2.1 #> [31] assertthat_0.10 #> [31] assertthat_0.2.1 #vctrs_0.5.0 R6_2.5.1 #> [34] lifecycle_1.0.3 string_1 .4.1 pkgconfig_2.0.3 #> [37] pillar_1.8.1 bslib_0.4.0 data.table_1.14.4 #> [40] glue_1.6.2 Rcpp_1.0.9 highr_0.9 #> [43] xfun_0.34 tibble_3.1.8 tidyselect_1.2.0 #> [46] knitr_1.40 htmltools_0.5.3 rmarkdown_2.17 #> [49] readr_2.1.3 compiler_4.3.0

参考文献

科门，T. T.， C. E.雷瑟森，R. L.里维斯特，1990。“算法导论”。

胡，G.和P.阿加瓦尔，2009。“基于基因组表达谱的人类疾病-药物网络”《公共科学图书馆•综合》4 (8): e6536。

Iorio, F.， R. Bosotti, E. Scacheri, V. Belcastro, P. Mithbaokar, R. Ferriero, L. Murino等，2010。从转录反应发现药物作用模式和药物重定位美国国家科学院院刊107(33): 14621。

兰姆，J. E. D.克劳福德，D.佩克，J. W.莫德尔，I. C.布拉特，M. J.罗贝尔，J.勒纳等。2006。连通性地图:使用基因表达签名连接小分子、基因和疾病。科学313(5795): 1929。

诺斯，BV, D.柯蒂斯和PC .沙姆，2002。蒙特卡罗程序中经验P值计算的注记美国人类遗传学杂志71(2): 439。

萨勃拉曼尼亚，A.， P. Tamayo, V. K. Mootha, S. Mukherjee, B. L. Ebert, M. A. Gillette, A. Paulovich等，2005。基因集富集分析:解释全基因组表达谱的一种基于知识的方法美国国家科学院院刊102(43): 15545。