アメリエフの技術ブログ

Amelieff Staff Blog

Reactome Pathwayの最上位カテゴリを自動的に検索する方法

バイオインフォマティクス

こんにちはtsuyuhです。
寒さ・乾燥・花粉の三重苦に打ちひしがれております。

今回はReactome Pathwayに関する話題です。

Reactome Pathwayといえば、
「RNA-seq解析によって発現変動遺伝子群を同定し、それらのエンリッチメント解析をおこない、有意にエンリッチしたReactome Pathwayの一覧を取得する」ということがよくやられていると思います。

実際、弊社の受託解析サービスでも、バルクのRNA-seq解析やシングルセルRNA-seq解析でReactome Pathwayのエンリッチメント解析をよくおこなっています。

ところが、取得したReactome Pathwayの一覧を見ただけでは、カテゴリの種類が多岐にわたりすぎて傾向がいまいちよくわからない、なんてことがよくあります。

Reactome Pathway には上位概念・下位概念が存在し、それらの関係が階層的に定義されています。
そのため、それぞれのPathwayの最上位にあたるカテゴリ (Top Level Pathway) を調べればより傾向を掴みやすくはなりますが、リストアップされた大量のPathwayを地道に検索していくのは嫌すぎます。

そこで、そんな作業を自動的におこなえる便利な関数をご紹介します。

f:id:tsuyuh:20211203142356p:plain:w150
pathway

※ Reactome Pathwayのエンリッチメント解析については過去のブログをご参照ください。
　 ReactomePAを使ってみた - アメリエフの技術ブログ

概要

Rパッケージ "ReactomeContentService4R" の getPathway() 関数を用いればOK。

詳細

作業環境

  • R v4.1.1
  • ReactomeContentService4R v1.2.0

ReactomeContentService4R のインストール

Bioconductorからインストール可能。
※ R4.1以上が必要です。
　 (私はRをインストールし直す羽目になりました。)

# パッケージのインストール
> if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
> BiocManager::install("ReactomeContentService4R")

# パッケージのロード
> library(ReactomeContentService4R)

getPathways() 関数の基本的な使い方

ReactomeContentService4Rパッケージに含まれるgetPathways()関数は、パスウェイIDを入力し、top.level = TRUE オプションをつけて実行することで、入力したパスウェイの最上位のカテゴリ名 (Pathway IDなど) を返してくれます。

> getPathways("パスウェイID", top.level = TRUE)


例えば、R-MMU-69656 (Cyclin A:Cdk2-associated events at S phase entry) を入力すると、

> getPathways("R-MMU-69656", top.level = TRUE)

     dbId displayName          stId       name     stIdVersion      oldStId     schemaClass  speciesName releaseDate hasDiagram hasEHLD diagramWidth diagramHeight inferred inDisease
1 9859918  Cell Cycle R-MMU-1640170 Cell Cycle R-MMU-1640170.1 REACT_288106 TopLevelPathway Mus musculus  2021-09-15       TRUE   FALSE          651           136     TRUE     FALSE

のように出力され、最上位のカテゴリが R-MMU-1640170 (Cell Cycle) であることがわかりました。

getPathways() 関数を繰り返し実行する

getPathways() は複数の引数を扱えないようなので、for文で回しつつ出力をデータフレームに格納していく方法をとります。

例えば、こういうデータがあったとします (ReactomePAパッケージによるエンリッチメント解析結果をイメージしています)。

> head(res[,2:3])
           ID                     Description
1    R-MMU-69656    Cyclin A:Cdk2-associated events at S phase entry
2    R-MMU-75105    Fatty acyl-CoA biosynthesis
3    R-MMU-8953750  Transcriptional Regulation by E2F6
・・・

これを1行ずつgetPathways() にかけていきます。
ついでに、色々試していて不便だった下記ポイントを無理やり解消します。

  • 入力したパスウェイIDも出力結果上で参照できるようにしたい
  • 何故か "oldStId" 列が出力されないパスウェイがあるので、rbind() が上手くいかない

結論

> res_top <- data.frame()
> x <- res[, 2]
> for (i in x) {
    res_top_each <- cbind(i, getPathways(i, top.level = TRUE))
    res_top_each <- res_top_each[, colnames(res_top_each) != "oldStId"]
    res_top <- rbind(res_top, res_top_each)
}
# 出力結果
> res_top
              i    dbId                     displayName          stId                            name     stIdVersion     schemaClass  speciesName releaseDate hasDiagram hasEHLD diagramWidth diagramHeight inferred inDisease
1   R-MMU-69656 9859918                      Cell Cycle R-MMU-1640170                      Cell Cycle R-MMU-1640170.1 TopLevelPathway Mus musculus  2021-09-15       TRUE   FALSE          651           136     TRUE     FALSE
2   R-MMU-75105 9860010                      Metabolism R-MMU-1430728                      Metabolism R-MMU-1430728.1 TopLevelPathway Mus musculus  2021-09-15       TRUE   FALSE          558           586     TRUE     FALSE
3 R-MMU-8953750 9859976 Gene expression (Transcription)   R-MMU-74160 Gene expression (Transcription)   R-MMU-74160.1 TopLevelPathway Mus musculus  2021-09-15       TRUE   FALSE         1046           355     TRUE     FALSE
・・・

これですべてのPathwayについて最上位のカテゴリを検索することができました。

※ なお、解析には1行につき数秒かかります。

まとめ

getPathways() 関数を使えば、エンリッチしたReactome Pathway の傾向がわかりやすくなる！