寻找网络中的hub节点

其实转录组走到现在我总觉得少了点什么东西，后来才想起来是cytospace寻找hub基因

这个其实很简单，网上的教程大把大把的，而且不需要写代码

参考：

RNA-seq入门实战（十）：PPI蛋白互作网络构建（下）——Cytoscape软件的使用 网络图分析工具Cytoscape及其插件的安装 Cytoscape教程|如何挖掘蛋白互作网络中的核心基因？（https://zhuanlan.zhihu.com/p/595214085） 寻找核心基因+子网络 （https://zhuanlan.zhihu.com/p/137653031）

但是随着我学习到的一些其他知识，让我加深了这方面的理解，所以还是决定写一下

所以这并不是一篇Cytoscape教程，而是一片探索性的文章，我也没有把推文标题写成寻找PPI中的hub基因，有需要的同学可以参考上面给出的其他内容

我这里一开始是用的是 初探mRNA、lncRNA联合分析之下游 的差异表达lncRNA共表达mRNA对应的基因

打算是进一步narrow down ：一是看看PPI网络和共表达网络有什么不同，二是找出来的PPI hub gene又是哪些

差异表达lncRNA共表达mRNA:

SWT1
AXIN2
HARS2
MMS19
HYLS1
PTRH1
MXRA7
PUF60
ASAH1
CLEC16A
NFIX
ZNF638
FIGNL1
BID
CNOT2
PTPN12
CD44
NCAM1
CANT1
DUSP14
REPS1
MFAP3L
EPB41
TCF4
IFNGR1
RNASEH2B

之所以说是进一步，是因为我们一般拿到共表达网络就可以了，不会再对其进一步做PPI找hub

在这篇推文中我们已经进行了

• 共表达网络
• （所有DEGs）PPI
• WGCNA

并做了比较

string

直接string检索可以发现有5个基因存在关联，特别是CD44处于中心位置，这和前面的lncRNA共表达mRNA的结果是一致的，可以发现CD44的连线很多

此外，还有EPB41 TCF4 AXIN2 NCAM1 ，其中EPB41 TCF4属于前面WGCNA、lncRNA共表达mRNA的overlap。而实际上前面WGCNA、lncRNA共表达mRNA的overlap的基因位置还比较特殊，属于右下角小局部网络中，这也进一步提示我们关注这个共表达网络中的小局部网络（不光WGCNA和lncRNA共表达mRNA的overlap集中在这里，而且所有共表达的mRNA的PPI hub中5个有两个位于这里，这同时反映了表达谱特征和真实生物网络）。

WGCNA、lncRNA共表达mRNA的overlap（黄框为overlap，图为共表达网络）：

有意思的是，AXIN2 NCAM1在string中的PPI关联是通过文本挖掘找到的，而并没有找到共表达信息，而在前面的lncRNA共表达mRNA网络中，这些基因（CD44 AXIN2 NCAM1）是通过了作者设置的与差异表达lncRNAs共表达（转录水平上）的阈值

The absolute values of correlation coefficients >0.95 and p < 0.001 were specified as screening criteria.

cytospace

若导到cytospace中

#node1 node2 node1_string_id node2_string_id neighborhood_on_chromosome gene_fusion phylogenetic_cooccurrence homology coexpression experimentally_determined_interaction database_annotated automated_textmining combined_score
AXIN2 CD44 9606.ENSP00000302625 9606.ENSP00000398632 0 0 0 0 0 0 0 0.546 0.546
CD44 EPB41 9606.ENSP00000398632 9606.ENSP00000345259 0 0 0 0 0.042 0 0 0.630 0.630
CD44 TCF4 9606.ENSP00000398632 9606.ENSP00000381382 0 0 0 0 0.104 0 0 0.527 0.558
CD44 NCAM1 9606.ENSP00000398632 9606.ENSP00000480132 0 0 0 0 0 0 0 0.712 0.712

太少了

拿参考文献提供的差异表达的mRNA做，顺便看看这样PPI cytospace 寻找hub基因的结果和作者基因水平的有哪里不一样

https://ftp.ncbi.nlm.nih.gov/geo/series/GSE209nnn/GSE209825/suppl/GSE209825_mRNA_ADS_vs_Normal.Differential_analysis_results.xls.gz

The absolute logarithmic fold change (LFC) ≥ 1 (|log2 fold change|≥ 1) and p-value <0.05 were adopted as the screening criteria

   55  awk '$(NF-1) < 0.05 && $(NF-2) < -1 {print $3}' GSE209825_mRNA_ADS_vs_Normal.Differential_analysis_results.xls | sed '1d' > DEG_list.txt
   56  awk '$(NF-1) < 0.05 && $(NF-2) > 1 {print $3}' GSE209825_mRNA_ADS_vs_Normal.Differential_analysis_results.xls | sed '1d' >> DEG_list.txt

$ wc -l  DEG_list_uniq.txt
2492 DEG_list_uniq.txt

超过2000，直接导入cytospace

使用MCODE寻找子网络（模块）：

原文献：

In order to construct a visual network map, interaction relationships for the list of differential genes were extracted from the STRING protein interaction database (https://www.string-db.org/). The network data files were imported directly into Cytoscape software for visual editing. The size of a node in a PPI map is proportional to the degree of the node. Among them, PPI nodes with relatively high connectivity include UBA52, AKT1, SUPT20H, RPL19, EGF, and MYC

relatively high connectivity 指的是什么算法呢？

看起来有点像最简单的degree

反正作者应该是直接用对整个网络找的hub，没在这之前划分子模块，所以我们后面还是对整个网络进行hub的寻找

cytoHubba: identifying hub objects and sub-networks from complex interactome （https://bmcsystbiol.biomedcentral.com/articles/10.1186/1752-0509-8-S4-S11）

cytoHubba提供了12种算法，虽然我在其一开始发表的文章中只找到11种，这里还是归纳一下：

Local-based Methods

1. Degree: 度中心性算法根据节点的连接数来评估其重要性。节点的度越高，表示其拥有更多的连接，因此在网络中更加中心。

1. Maximum Neighborhood Component (MNC): 最大邻域组件算法计算节点邻域中最大子图的大小。该算法寻找具有最大子图的节点，这些节点在网络中具有重要的连接。
2. Density of Maximum Neighborhood Component (DMNC): 最大邻域组件密度算法是针对 MNC 的改进，它计算最大邻域组件中边的数量与节点数量的比值（密度）。该算法用于衡量节点在邻域中的连接密度。
3. Maximal Clique Centrality (MCC): 最大团中心性算法通过计算节点所属的最大团数量来评估节点的中心性。最大团是指一个子图，其中的每个节点都与其他所有节点直接相连。

Global-based methods

1. Closeness: 接近中心性计算从一个节点到网络中所有其他节点的最短路径长度的倒数和。接近中心性较高的节点更接近所有其他节点。
2. EcCentricity: 离心度计算从一个节点到网络中所有其他节点的最短路径长度。离心度较低的节点更为中心。
3. Radiality: 这段是AI给的，我看原文公式看看不懂

在网络分析中，径向度是一种用于衡量节点在网络中的中心性的指标。它关注的是节点与其可达邻居之间的距离关系。具体来说，径向度将高中心性赋予那些在其可达邻域中离其他节点的距离较短的节点。 为了计算节点的径向度，首先需要确定节点的可达邻域。可达邻域包括节点本身以及所有与该节点距离不超过其直径的节点。直径是指可达邻域中节点之间的最大距离。 然后，对于每个节点，计算它与可达邻域内所有其他节点之间的平均距离。距离越短，节点的径向度越高。因此，径向度高的节点可以被认为是网络中的中心节点，它们在网络的结构和连接中起着重要的作用。 通过计算节点的径向度，我们可以获取关于网络中节点的重要性和中心性的信息，这对于理解网络结构和节点的功能扮演了重要的角色。

1. Bottleneck: 这个我也看不太懂，感觉和后面的基于通过最短路径的stress（数量）、betweenness（频率）都很像

可以将瓶颈定义为具有高介数中心性的蛋白质（即，具有许多“最短路径”的网络节点，类似于高速公路地图上的主要桥梁和隧道）。 参考：The Importance of Bottlenecks in Protein Networks - NCBI （https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1853125/）

该方法基于网络在信息流传递过程中的瓶颈原理，寻找在信息流中具有关键作用的基因。Bottleneck方法关注的是基因在网络中连接不同模块之间的重要性，例如，一个基因连接了多个模块，其在信息传递中起到了关键的桥梁作用。因此，Bottleneck方法寻找的是在网络中起到关键枢纽作用的基因。

1. Stress: 压力中心性衡量节点对网络内信息流的影响。它考虑通过特定节点的最短路径的数量。具有高压力中心性的节点是指许多最短路径通过的节点，表明它们在网络中的重要性。
2. Betweenness: 中介中心性基于网络中基因之间的最短路径数量，通过计算一个基因在所有最短路径中被作为中间节点的频率来评估其重要性。Betweenness方法关注的是基因在不同基因之间传递信息的能力。一个具有高Betweenness值的基因意味着它在网络中连接了许多最短路径，因此在信息传递中具有重要性。
3. Edge Percolated Component (EPC): EPC算法通过迭代地移除基于边权重weight和边中介数Betweenness的最不重要的边来识别显著

上面带了引用没给参考的部分都是AI帮我解释但我不确定的，需要甄别着看，也希望懂的老师能在评论区说说，或者笔记投稿，我在检索上面这11个方法的时候其实无论中外互联网感觉相关资料都很少

挑选一些方法寻找top10的hub nodes：

其实可以发现全局算法的结果较为一致，最简单的局部算法degree也差不多，但是MCC我们发现找到的都是核糖体大小亚基蛋白

随着调整不同算法出图，我发现作者的PPI hub中始终有一些没有出现在我的top10中

检查：

> PPI_hub<-c("MYC","EGF","AKT1","SUPT20H","UBA52","RPL19")
> PPI_hub %in%  deg_list_uniq
[1]  TRUE FALSE  TRUE FALSE  TRUE  TRUE

可以发现EGF SUPT20H在转录水平的表达差异并不显著，而在作者的DNA水平定量上显著，所以我们构建的PPI网络中没有这两个节点。其实，这也在侧面说明了在真实生物世界中网络调控处于中心位置的蛋白并不一定差异表达非常显著，可能仅仅是略微的差别就有可能引起一系列变化从而影响表型。

img

img

在检索学习如何鉴定网络中心节点的时候，我检索到了一个整理的很好的网页，可以对上面11种方法补充着看，这里也贴出来

网络中心性

网络中心性 (https://wiki.swarma.org/index.php/%E7%BD%91%E7%BB%9C%E4%B8%AD%E5%BF%83%E6%80%A7#.E6.B8.97.E6.BB.A4.E4.B8.AD.E5.BF.83.E6.80.A7_Percolation_centrality.EF.BC.88PC.EF.BC.89)

在图论和网络分析中，中心性 centrality 指标用于识别图中最重要的顶点。其应用包括在社交网络中识别出最有影响力的个人，在因特网或城市网络中识别出最为关键的基础设施节点，以及识别疾病的超级传播者。中心性的概念最初是在社交网络分析中发展起来的，许多用于衡量中心性的术语都反映出了它们的社会学起源。中心性不应与节点影响度相混淆，后者意在量化网络中每个节点的影响。

其中包括了一些其它中心性，如调和中心性、特征向量中心性、网页排名中心性等

其中，网页排名中心性 PageRank centrality

PageRank算法是搜索引擎Google所采用的评估网页重要性的算法，GeneRank算法在PageRank的基础上通过生物网络的拓扑结构和节点的初始值对节点分数进行迭代，可以成功地对网络中的基因进行优先排序

original GeneRank：

Modified GeneRank：

我们这里不做公式的深入解析，只是借此引入节点的初始值

虽然我们在前面发现

EGF SUPT20H在转录水平的表达差异并不显著，而在作者的DNA水平定量上显著，所以我们构建的PPI网络中没有这两个节点。其实，这也在侧面说明了在真实生物世界中网络调控处于中心位置的蛋白并不一定差异表达非常显著，可能仅仅是略微的差别就有可能引起一系列变化从而影响表型。

但是初始分数（log2FC的绝对值或者-lgPvalue）其实还是能在一定程度上说明节点的在真实的生物世界网络中的重要性

私以为其实可以两个都做做，先不带初始值看看hub节点是哪些，这时的hub节点都是因为在网络中占据了重要地位。然后带上初始值，就可以看看哪些节点是因为带上了初始值而变成hub

思路打开

另外这些无初始值的算法并不一定需要用cytospace不可，R包 centiserve已经整合了

所以不一定非要是PPI，只要是网络结构就可以做

思路打开

现在来看看有初始值的PPI如何找到hub节点

string（初始值）

string最开始是没有这个功能的，后面给加上了

但是这个功能和我一开始想的并不一样

我以为的是可以加上初始值来确定hub基因了

其实只是一个类似GSEA的功能，关于基因集分析可以参考我前面的这篇系统介绍的推文基因集分析的前世今生（附进行通路富集分析的9个tips）

测试数据结果：

看来string和cytospace还是分割的，功能联系没那么紧密啊

我手上有21个基因及其初始值的列表输进去看看

“具有值/等级的蛋白质”功能 基于大基因列表（覆盖整个基因组或至少其大部分）中的值的分布进行统计富集检测。 考虑到我们进行的测试数量，几乎不可能在如此少量的值上产生具有统计意义的信号。然而，您仍然可以继续并将您的数据集提交给基于常规基因集的分析，我们将您的蛋白质集与整个蛋白质组进行比较

外面的环就代表了初始分数

相当于没用初始分数，只是让你看看在对应的网络中哪些节点处于中心位置并且初始分高（倒也是一种思路）

但我还想看看初始分数最后对hub基因的影响

想要找hub基因还是需要导出来去cytospace，但是这个过程就丢失掉了nodes初始分数，只是edges发挥作用

我们刚刚介绍的GeneRank也是需要输入初始值的

我们可以自己手搓一下上面那个GeneRank算法，并不难

鉴于R的运行速度，还是推荐用相对来说快一点的python

如果熟练可以自己用基本数据结结构比如列表、字典完成整个任务，此外我发现NetworkX这个第三方库用起来还可以，如果用python画PPI可以考虑一下这个

https://networkx.org/

总之，在完成算法迭代前，我的初始分数top10为

CCR1, IFNAR2, MAPT, OAS1, HLA-C, ABO, NSF, PDE4A, CAT, FUBP1

完成算法迭代后top10为

max_dict = {}
for index, value in enumerate(ls_b):
    max_dict[value] = index
  
sorted_dict = sorted(max_dict.items(), reverse=True)[:10]
top_10_indexes = [index for _, index in sorted_dict]
print(top_10_indexes)
[13627, 3767, 1609, 2927, 6738, 482, 452, 1953, 11160, 3862]
selected_elements = [nodes[index] for index in top_10_indexes]
print(selected_elements)
['ABO', 'CCR1', 'IFNAR2', 'MAPT', 'OAS1', 'HLA-C', 'JAK1', 'NSF', 'CYP4B1', 'PDE4A']

可以发现，网络的拓扑结构特征被融入到了基因中，特别是JAK1和CYP4B1这些一开始并不在top10的节点

这样的top节点作为hub基因会不会更有意义呢？

那是那些句话，这里就刚好放到最后作为小结：

EGF SUPT20H在转录水平的表达差异并不显著，而在作者的DNA水平定量上显著，所以我们构建的PPI网络中没有这两个节点。其实，这也在侧面说明了在真实生物世界中网络调控处于中心位置的蛋白并不一定差异表达非常显著，可能仅仅是略微的差别就有可能引起一系列变化从而影响表型。 但是初始分数（log2FC的绝对值或者-lgPvalue）其实还是能在一定程度上说明节点的在真实的生物世界网络中的重要性 私以为其实可以两个都做做，先不带初始值看看hub节点是哪些，这时的hub节点都是因为在网络中占据了重要地位。然后带上初始值，就可以看看哪些节点是因为带上了初始值而变成hub

无论是WGCNA、PPI+cytospace的各种方法、调控序列共表达网络还是其它考虑了初始权值的方法，其实都可以根据自己的实际情况来，思路打开没有坏处，narrow down将会是高通量数据的主旋律，就算是学到后面用机器学习（回归模型、随机森林、支持向量机等等）也是同样的思路