いろいろ倉庫

KNIME、EXCEL、R、Pythonなどの備忘録

【R】ネットワーク解析を学んでみたい③

・お題:生体内の様々な分子は相互作用し、絶妙なバランスを保っている。生体内の分子の役割を理解するのに、ネットワーク解析が有向らしい。学んでみたい。

 

・私は素人なので、参考サイトをなぞるような形で学んでいきたい。今回は以下のサイトを参考に学ばせていただく。九州工大の竹本先生のHP。日本語でここまで詳しく解説されているHPは他に見たことがないので、ぜひこちらのHPをご覧いただきたい。

sites.google.com

・竹本先生の貴著「生物ネットワーク解析」も名著だと思うので、気になる方はぜひチェックしてみてほしい。

www.coronasha.co.jp

 

・前回、グラフの特徴を求めて遊んでみた。今回はネットワークの処理から。

・まず、連結という言葉に関して、とあるノードから別のノードに到達可能かどうか、ということを表す言葉らしい。有効グラフの場合、矢印を辿ってどちらからでも到達可能なら強連結、どちらか一方から到達可能なら弱連結と呼ぶみたい。連結されているノードの最も大きい塊を、「連結成分」と呼ぶらしい。

・適当にグラフを作成して試してみる。

A B 1
B C 1
C A 1
D E 1
D F 1
F G 1
E G 1
G H 1
F H 1
H I 1

をtextファイルとして保存。

>library(igraph)
>library(RColorBrewer)

> d <- read.table("20230325_igraph.txt")#読み込み
> g <- graph.data.frame(d[1:2],directed=F)#グラフ化
> plot(g)#図示

・大きく二つの塊(連結成分)がある。連結成分の様子を見てみる。

> cls <- clusters(g,"weak")

>cls

$membership
A B C D F E G H I 
1 1 1 2 2 2 2 2 2 

$csize
[1] 3 6

$no
[1] 2

・csizeが最大の連結成分以外の連結成分に所属しているノードを消すみたいな形で。もっとも大きい塊(最大連結成分)を求めることができるらしい。

> g_comp <- delete.vertices(g,
+                           subset(V(g),
+                                  cls$membership!=which(cls$csize==max(cls$csize))1))
> g_comp

IGRAPH 41f22c7 UN-- 6 7 -- 
+ attr: name (v/c)
+ edges from 41f22c7 (vertex names):
[1] D--E D--F F--G E--G G--H F--H H--I

 

・この後、元サイトではネットワークのランダム化に関して、igraphでの適用例を解説されている。ネットワークの構造指標が有意であるかどうかを検証する際に必要とのことだけれど、私の理解が追い付いていないので、今回はスキップする。ちょっとずつ勉強を進めていきたい。。

 

・次に、ネットワークのグループ分け(コミュニティ抽出)をやってみる。生物系のネットー枠にいてグループ分けは、ネットワークの生物学的意義を見出すのに重要らしい。データセットとして、元サイトで提供されているデータセットをダウンロードし、解凍してeco_EM+TCA.txtをカレントディレクトリに保存した。

> d_eco <- read.table("eco_EM+TCA.txt")

> g_eco <- simplify(graph.data.frame(d_eco,directed=F),
+                   remove.multiple=T,
+                   remove.loops=T)
> layout_eco <- layout.fruchterman.reingold(g_eco)
> plot(g_eco, layout=layout_eco)

・糖の代謝クエン酸回路に関するネットワークらしい。

・ここからどうやってグループを設定していくかというと、Q値(モジュラリティ)というものを指標に塊に分けるのが一般的みたい。コミュニティというのは教師なし学習クラスタリングに似ており、何等か指標を最適化する分け方を見つけようという操作になる。直感的には、クラスタリングはギュッと集まっているところがクラスタっぽいし、ネットワークでいうとノード同士がめっちゃエッジで繋がっていればコミュニティっぽいということになる。これを定量化した指標がQ値であり、Q値を最大にするグループの分け方を何とかして求めることになる。

・これにはいろいろな手法が提案されており、計算リソースや性能などから、必要に報じて手法を選択することになる(目的に応じて手法の選択が必要)。今回は、貪欲法、固有ベクトル法、焼きなまし法を試してみたい。元サイトにはさらに多くの手法が紹介されているので、気になる方は元サイトをご確認いただきたい。

#貪欲法。高速。

> data_fg <- fastgreedy.community(g_eco)

> data_fg$membership
 [1] 4 1 1 5 2 2 1 5 2 1 2 3 2 2 1 1 4 3 5 4 5 5 2 2 3 3 4

・貪欲法は階層型にどんどんコミュニティを形成して最適化していくアルゴリズムなので、樹形図を描くことができる。

> dendPlot(data_fg)#plot_dendrogram(data_fg)でも同じ樹形図を描けるみたい。

> V(g_eco)$color <- data_fg$membership#色を設定。
> plot(g_eco, layout = layout_eco)#図示

 

#固有ベクトル

> data_ev <- leading.eigenvector.community(g_eco, options = list(maxiter=1000000, ncv=5))

> data_ev$membership
 [1] 1 3 3 3 2 2 5 3 2 5 5 4 2 2 5 5 1 4 3 1 3 3 2 2 4 4 1

> V(g_eco)$color <- data_ev$membership
> plot(g_eco, layout = layout_eco)

 

#焼きなまし法。検出精度が高い。

> data_sp <- spinglass.community(g_eco)

> data_sp$membership
 [1] 4 6 6 5 2 2 6 5 2 1 1 3 2 2 1 1 4 3 5 4 5 5 2 2 3 3 4

> V(g_eco)$color <- data_sp$membership
> plot(g_eco, layout = layout_eco)

・元サイトではこのあとlinkcommライブラリを用いて重なりのあるコミュニティ検出をやる。今回はスキップする。興味のある方は元サイトをご確認いただきたい。

 

つづく。