ネットワーク解析wgcna

WGCNA blockwiseModules の主要パラメータ:deepSplit・minModuleSize・mergeCutHeight の意味と決め方

WGCNAでモジュールを検出する blockwiseModules() は、相関 → TOM → 階層クラスタ → モジュール統合までを一気にやってくれる便利な関数です。ただ便利な反面、出てくるモジュールの数やサイズは引数の設定で大きく変わります。「モジュールが多すぎる」「逆に少なすぎる」「灰色(未割当)の遺伝子だらけ」——こうした悩みは、たいてい数個のパラメータで解決します。

この記事では blockwiseModules() の主要パラメータについて、それぞれが何を左右し、どちらに動かすと何が起きるかを整理し、最後に症状別のチューニング早見表にまとめます。

なお本記事は WGCNA解析シリーズの番外編(逆引きリファレンス) です。part 1 / part 2 の続きというより、パラメータで詰まったときにいつでも引ける早見集として使ってください。

WGCNA自体の流れがまだ曖昧な方は、先にこちらの概念編から読むとスムーズです。

関連記事WGCNA解析 part 1:RNA-seqデータを使って、疾患に絡む遺伝子ネットワークを可視化する

まず全体像:blockwiseModules の最小呼び出し

細かい話の前に、主要な引数だけを並べた呼び出し例を見ておきます。この記事では、ここに出てくる引数を順番に解説していきます。

net <- blockwiseModules(
  datExpr,                 # サンプル × 遺伝子 の発現行列
  power          = 7,      # softPower(β)。pickSoftThreshold で決める
  networkType    = "signed",
  TOMType        = "signed",
  corType        = "bicor",
  maxBlockSize   = 20000,  # 全遺伝子を1ブロックに収めるよう大きめに
  deepSplit      = 2,      # モジュールの細かさ(0〜4)
  minModuleSize  = 30,     # モジュールの最小サイズ
  mergeCutHeight = 0.25,   # 似たモジュールを統合する閾値
  numericLabels  = TRUE,
  saveTOMs       = TRUE,
  randomSeed     = 1234    # 再現性のため固定
)

大まかに言うと、

  • power / networkType / TOMType / corType が「ネットワークの作り方(土台)」

  • maxBlockSize が「計算のしかた(ブロック分割)」

  • deepSplit / minModuleSize / mergeCutHeight が「モジュールの切り出し方(数とサイズ)」

を決めます。実際にコードを動かしながら確認したい方は、全Rコード付きの完全版にこの呼び出しの前後がすべて載っています(前半無料)。

関連記事【前半無料|全Rコード付き】WGCNA完全パイプライン🔒 一部有料

ネットワークの土台:power と相関の取り方

まず「ネットワークそのものをどう作るか」を決めます。ここは下流のすべてに効いてくるので、最初に固めておきます。

power(softPower / β)

相関に掛ける重み=softPower(β) です。WGCNAで最初に決める、そして最も影響の大きいパラメータです。

  • pickSoftThreshold() を実行し、スケールフリー性(R²)が 0.8 以上になる最小のβを選ぶのが基本

  • signed ネットワークでは unsigned より高めのβが必要になりやすい(目安:unsigned 6前後/signed 12前後)

βの決め方そのものは part 1 で図つきで解説しています。

関連記事WGCNA解析 part 1(softPowerの決め方)

βが低すぎるとノイズまで拾ってモジュールがぼやけ、高すぎると連結が切れて灰色(未割当)が増えます。後述の「灰色が多すぎ」問題の第一容疑者はだいたい power です。

networkType / TOMType / corType

ネットワークの「符号」と「相関の種類」を決めます。

  • networkTypeunsigned(相関の絶対値)/ signedsigned hybrid。負の相関を別物として扱いたい(生物学的に素直)なら signed 推奨TOMType も基本これに揃えます

  • corTypepearson(既定)/ bicorbicor(biweight midcorrelation)は外れ値に強いので、サンプルにばらつきがある実データではこちらが無難。ただし bicor を使うときは maxPOutliers = 0.05 あたりを併用するのが定石

signed か unsigned かでβの最適値も変わるので、この2つは power を決める前に固定しておきましょう。

計算のしかた:maxBlockSize とブロック分割の罠

地味ですが、最も事故が多いパラメータです。

blockwiseModules は、遺伝子数が maxBlockSize を超えると、遺伝子を複数のブロックに自動分割してから別々にネットワークを組みます(メモリ節約のため)。問題は、このとき本来つながっているはずの遺伝子が別ブロックに割り振られ、同じモジュールに入らなくなることがある点です。

  • 既定値は maxBlockSize = 5000多くの解析では遺伝子数がこれを超えるため、知らずに使うと意図せずブロック分割されている

  • できれば 全遺伝子を1ブロックに収めるよう、遺伝子数より大きい値を指定するのが理想

  • ただしメモリと相談。ざっくりの目安:

搭載メモリ1ブロックに入れられる遺伝子数の目安
8 GB〜13,000
16 GB〜20,000
32 GB〜30,000 以上

遺伝子数が多くてどうしても1ブロックに収まらない場合は分割もやむなしですが、「なぜか関連遺伝子が別モジュールに分かれる」時はまず maxBlockSize を疑う と覚えておくと事故が減ります。

モジュールの切り出し:数とサイズを決める

ここからが、検出されるモジュールの数とサイズを直接いじる「3兄弟」と微調整です。まず1つだけ動かして結果を見るのが鉄則です。

deepSplit:モジュールの「細かさ」

deepSplit は動的樹形図カットの感度で、0〜4 の整数を取ります。

  • 値を上げるほど敏感になり、モジュールは多く・小さくなる

  • 値を下げるほど大まかになり、モジュールは少なく・大きくなる

  • 迷ったら既定の 2 から始め、結果を見て ±1 する

minModuleSize:モジュールの最小サイズ

1つのモジュールに最低何遺伝子必要か、の下限です(既定 30)。

  • 下げると小さなモジュールが乱立する(ノイズ的なモジュールも増える)

  • 上げると小さい塊は統合・除外され、モジュールは大きく・少なくなる

  • 遺伝子数が少ない解析では 20〜30、多い解析では 30〜50 あたりが扱いやすい

mergeCutHeight:似たモジュールの統合

検出後、固有ベクトル(ME)が似ているモジュール同士を統合する閾値です(既定 0.15)。

内部的には「ME相関 > 1 − mergeCutHeight のモジュールを統合」します。つまり既定の 0.15 なら、相関 0.85 以上のモジュールがマージ されます。

  • 値を上げる(例:0.25)と統合が進み、モジュールは減って大きくなる

  • 値を下げると統合が緩み、似たモジュールが分かれたまま残る

「似たようなモジュールが多くて解釈しづらい」ときは、まずここを上げるのが手っ取り早いです。

その他の微調整:reassignThreshold / pamStage / randomSeed

普段はデフォルトのままで構いませんが、知っておくと便利な引数です。

  • reassignThreshold:遺伝子をより相関の高いモジュールへ再割り当てする際のp値閾値(既定 1e-6)

  • pamStage / pamRespectsDendro:PAM段階。TRUE(既定)だと未割当の遺伝子をできるだけモジュールに割り振る=灰色が減る

  • numericLabels = TRUE:モジュールラベルを色名でなく番号で受け取る(扱いが楽)

  • randomSeed:ブロック分割時の前処理に乱数を使うため、再現性のために必ず固定 する

【早見表】ケース別チューニング

ここまでを「困ったケース → どの引数をどっちに動かすか」でまとめます。まず1つだけ動かして結果を見るのが鉄則です。

ケース対処(動かす引数)
モジュールが多すぎるdeepSplit ↓ / minModuleSize ↑ / mergeCutHeight
モジュールが少なすぎるdeepSplit ↑ / minModuleSize ↓ / mergeCutHeight
灰色(未割当)が多すぎるpower を見直す(高すぎ?)/ minModuleSize ↓ / pamStage = TRUE を確認
似たモジュールが分かれすぎmergeCutHeight
関連遺伝子が別モジュールに割れるmaxBlockSize を大きく(ブロック分割の解消)
結果が毎回変わるrandomSeed を固定

kME(signedKME)でモジュールの質を確認する

パラメータを調整したら、モジュールが「中身として」きちんとまとまっているかを確認します。ここで使うのが kME(module membership=signedKME) です。

kME は各遺伝子とモジュール固有ベクトル(ME)との相関で、その遺伝子がどれだけそのモジュールらしいかを表します。signedKME() で計算でき、kMEが高い遺伝子ほどモジュールの中心(ハブ候補)です。良いパラメータ設定なら、各モジュールに高kMEの遺伝子がしっかり存在します。

kME・kWithin を使ったハブ遺伝子の絞り込みは part 2 で具体的に解説しています。

関連記事WGCNA解析 part 2:疾患関連モジュールの機能を読み解き、ターゲット候補遺伝子を絞り込む

よくあるハマりどころ

  • maxBlockSize を既定のまま大規模データに使う → 知らぬ間にブロック分割。まず疑うべき定番の事故

  • randomSeed を固定し忘れる → 再現できず、レビューで詰む

  • power をいい加減に決める → 下流のすべてがブレる。pickSoftThreshold を省略しない

  • 一度に複数の引数を動かす → 何が効いたか分からなくなる。1つずつ

まとめ

blockwiseModules のモジュール数・サイズは、

  • 土台power / networkType / TOMType / corType

  • 計算maxBlockSize(ブロック分割の罠)

  • 切り出しdeepSplit / minModuleSize / mergeCutHeight

で決まります。困ったら早見表に戻り、1つずつ動かして kME で質を確認する——この繰り返しで、納得のいくモジュールに近づけます。

ポイント

パラメータを踏まえて、実際に最初から回したい方へ — GEOデータ取得 → 正規化 → モジュール検出 → ハブ遺伝子の絞り込みまで、通しで動かせる全Rコードを実装記事にまとめています。前半無料なので、まず動かして試せます(後半 ¥500)。

関連記事WGCNA を R で実装:GEOデータ取得〜ハブ遺伝子まで全Rコード付き(前半無料)🔒 一部有料