的印R包实现了一种检测生物通路选择的方法。总体思路是使用启发式方法(模拟退火)在生物通路中搜索呈现不寻常特征的基因子网络。
总体思想很简单:我们考虑一个具有预先定义分数的基因列表(例如，像Fst这样的分化测量)，我们希望找到在零假设下表现出高于预期分数的基因网络。
为此，我们将使用转化为基因网络的生物途径数据库，并在这些图中搜索高分子网络。
关于算法的详细信息可以在Gouy等人(2017)．
请引用这篇论文印对于您的项目:
印输入一个基因得分的数据帧。第一列必须对应基因ID(例如Entrez)，第二列是基因得分(每个基因一个值)。
另一个输入是生物通路(基因网络)的列表graphNEL格式。我们建议使用包装石墨获取路径数据:
library(graphite) # pathwayDatabases() #查看路径和可用的物种#获取路径列表:paths <- graphite:: paths ("hsapiens"， "kegg") #将前3条路径转换为图形:kegg_human <- lapply(paths[1:3]， graphite::pathwayGraph) head(kegg_human)注意，基因得分数据帧和路径列表(例如entrez)之间的基因标识符必须相同。石墨提供转换基因标识符的函数。
提供了来自Daub等人(2013,MBE)以及人类KEGG途径的示例数据集:
库(图章)数据(daub13)头(分数)#基因分数基因评分## 1 1 0.9200665 ## 2 10 1.5974385 ## 3 100 1.6885589 ## 4 1000 3.3314333 ## 5 10000 1.6668512 ## 6 10001 1.3529425我们首先必须在提供的生物通路中搜索高分子网络，使用模拟退火:
#在前3条KEGG路径上运行模拟退火:HSS <- searchSubnet(kegg_human, scores)对于每个路径，该函数返回所找到的得分最高的子网络、它的得分以及模拟退火运行的其他信息。
然后，为了测试高分子网络的显著性，我们生成一个高分的零分布:
#生成经验null分布null <- nullDist(kegg_human, scores, n = 1000)注意零对象是一个简单的零高分向量(这里是1000)。因此，如果您想计算更精确的p值，您可以在之后运行其他迭代，并将输出与前一个向量连接起来。
这个分布最后被用来计算p值和更新印对象:
HSS <- testSubnet(HSS, null)当p值计算完成后，您可以生成一个汇总表(每个路径一行):
#结果:生成汇总表选项卡<- summary(HSS) head(选项卡)##路径网。子网大小。子网大小。得分## 1急性髓系白血病53 8 3.26023952470392 ## 2粘附结66 8 3.69037387285983 ## 3脂肪细胞因子信号通路62 7 3.79223744028381 ## p.val子网。基因## 1 0.033 673 2885 3815 3845 5291 5605 6654 8503 ## 2 0.017 1496 1499 1956 2534 5787 5797 51176 83439 ## 3 0.014 32 2538 4852 51422 53632 92579 126129#你可以像下面这样写汇总表:Table (tab， # file = "signet_output. "tsv"， # sep = "\t"， # quote = FALSE， # row.names = FALSE)请注意，搜索高分子网络和生成空分布可能需要几个小时。然而，这些步骤很容易在集群上并行化，因为不同的迭代是相互独立的。
Cytoscape (www.cytoscape.org)是一个专门用于网络可视化的外部软件。印允许生成一个XGMML文件在Cytoscape中加载(文件>导入>网络>文件…)。由于该函数，可以将该文件写入工作目录writeXGMML．
如果函数的输入是单个的印对象时，将表示整个路径，属于得分最高的子网络(HSS)的节点将用红色突出显示。
writeXGMML(HSS[[1]]， filename = "cytoscape_input.xgmml")如果提供了路径列表(signetList)，则所有p值低于给定阈值(默认值:0.01)的子网络都将被合并和表示。注意，在这种情况下，只保留属于HSS的节点进行表示。
writeXGMML(HSS, filename = "cytoscape_input. txt ")Xgmml "， threshold = 0.01)然后可以在Cytoscape中对表示进行精细定制。