RでWGCNA解析～TOMの計算と妥当性の検討～

さて、前回の記事からだいぶ時間が空いてしまいましたが次のステップです。前回はSoft Thresholdを決定したので、その値に基づいて冪乗し、TOMを計算しその結果が妥当かどうかをＲ上で検討します。データはwgcna　tutorialにあるものを用いるのと、前回の記事の続きからコードを書いていくので注意して下さい。TOMについては以前書いた記事を参照してください。

#まずは相関係数を重み付きで計算。前回の記事で10に設定したのでそれで行う。
softPower = 10;
adjacency = adjacency(Data2, power = softPower);

#データによってはここの計算に結構時間がかかる。気長に待ちましょう。
TOM = TOMsimilarity(adjacency);
dissTOM = 1-TOM

さて、この結果をwgcna解析に使っていくわけですが、Soft Thresholdの選択はそもそも恣意的であるので一旦検討を挟む必要があります。具体的にはスケールフリー性を持つかどうかをR上で検証していくことになります。まずはネットワーク中のノードの次数の分布をヒストグラムで確認してみます。

k=as.vector(apply(adjacency, 2, sum, na.rm=T))
hist(k)

出力は以下のようになるはずです。kが増えるほど指数関数的に減少しているので、スケールフリー性を持つと言えそうです。

次はWGCNAパッケージのscaleFreePlot関数で次数kがパレート分布に従っているかを確認します。

scaleFreePlot(k, main="Check scale free topology\n")

出力は以下の通り。R^2値が0.94と高い値を示しているので、結果は信用できそうです。

次回はdissTOMに基づいたクラスタリングとモジュール検出を行っていきます。

前回の記事↓

次回の記事↓