ブラックボックス解析を卒業するRNA-Seq実践ガイド——AI時代に求められる「考えて進める解析者」へ

RNA-Seq解析を学び始めると、多くの人が最初にぶつかるのが、Rや各種パッケージの扱いです。
コードを書くことやエラー対応に時間を取られ、本来向き合うべきデータの解釈や、生物学的な意味づけまで手が回らないことは少なくありません。

一方で、近年はGUIベースの自動解析ツールも充実してきました。手軽に結果を得られる反面、前処理や正規化、バッチエフェクトの影響を十分に吟味しないまま解析が進んでしまうこともあります。RNA-Seq解析で本当に重要なのは、計算を終えることではなく、各ステップでデータがどう変化したかを自分の目で確認し、妥当な判断を重ねることです。

Subio Platform は、そのための“見ながら考える解析環境”を目指して設計されています。
視覚化に重点を置き、必要十分な統計機能を組み合わせることで、解析の意図を持って進められるようにしています。

本チュートリアルでは、RNA-Seq解析の基本フローをなぞりながら、
「どこを見て、どう判断するか」という解析の思考プロセスを、Subio Platform上で実践的に学びます。
「実行してなんとなく結果が出せる」状態から、「データを理解し、解析を説明できるようになる」状態へ。

【まず結果を出したい方へ】

このページは、自分の手でRNA-Seq解析を理解しながら進めたい方向けの完全ガイドです。
一方で、次のような状況の方には、より短いルートもあります。

すぐに自分のデータを解析したい
 → データ解析サービス:解析済み結果だけでなく、自由に探索できるSSAファイルとして納品します。

自分のデータで操作を最短で覚えたい
 → オンライントレーニング:実データを使って必要な解析フローだけを集中的に学べます。録画して手元の操作マニュアルにもできます。

▶まずはデモデータで基礎から理解したい
 → このままチュートリアルへ進んでください。

1-a.FASTQファイルのインポート

要点:解析の第一歩は、測定した実験データの取り込みです。生データ(FASTQ)から処理する場合(1-a)と、すでに遺伝子ごとの発現量 (gene counts, FPKM, TPMなど) に変換された場合(1-b)に分けて説明します。

オミクスデータ解析を学ぶ、効率的な方法とは?」でも触れたように、FASTQファイルの処理は、多くの初心者が最初につまずきやすく、時間を取られる工程です。
しかし、この工程は研究の核心ではなく、手順に従って進める実行プロセスに近いものです。

コマンドラインを学び、エラーを乗り越えてパイプラインを完了させることには、確かに達成感があります。ですが、これはあくまで準備段階に過ぎません。研究の核心である「データの理解」や「生物学的な洞察」はこの先にあります。

本来最も時間をかけるべきは、その先の解析と解釈です。前段の工程に過度に時間と労力を費やしてしまうと、全体として非効率になりがちです。

もしこの工程で手が止まっているのであれば、この部分は他の人やサービス、あるいはAIに任せてしまうのも一つの合理的な選択です。

一方で、「FASTQファイルを入れてボタンを押すだけ」で完結するツールや、その操作だけを説明するガイドは手軽ではありますが、解析の理解を深める機会を得にくい側面もあります。

すでにGene Counts(あるいは、必ずしも最適とは言えないRPKM/FPKM、TPM、CPMなど)に変換されたデータをお持ちの場合は、無理にこのステップから始める必要はありません。1-bから進めることを推奨します。

本セクション(1-a)では、この前処理をSubio Platform上でどのように完結させるかを具体的に解説します。ただし、これはあくまで準備段階であり、本質的な解析はこの先にあることを意識しながら進めていきましょう。

Subio PlatformにRNA-SeqのFASTQファイルをインポートするのは、  fastp HISAT2 StringTie という3つのツールからなるパイプラインを実行させ、発現量を計算することになります。ただし、ワークステーションも、UNIXのコマンドラインによる操作スキルも必要ありません。普通のWindowsまたはMacのコンピューターで動かせます。コマンドラインを使える方にとっても、Subio Platformを使う方が操作がずっとラクです。

パイプラインを実行する前に、これらのツールをインストールする必要があります。もし難しいようでしたら、FASTQ処理の問題解決サービスをご注文ください。

FASTQファイルは、.gz圧縮形式のままにしておいてください。同じフォルダー内に、途中結果のファイルが大量に生成されますので、ディスクに十分な空きスペースがあること(FASTQファイルの合計サイズの5倍以上)をご確認ください。もし、実行途中にディスクが足りなくなった場合は、途中で止まります。FASTQファイルは、外付けディスクに置いてあっても大丈夫です。

Paired-end サンプルのFASTQファイルをインポートする場合は、ファイル名の付け方のルールに従ってください。また、何らかのトラブルで実行できないときは、トラブルシューティングをお試しください。

RNA-Seq FASTQ ファイルをインポートする。

目次に戻る

1-b.遺伝子レベルの発現量データのインポート

すでに遺伝子レベルの発現量データ(gene counts, FPKM, TPMなど)に変換されたデータをお持ちでしたら、FASTQファイルを処理するために頑張る必要はありません。そのテキストファイルをSubio Platformにインポートしてください。 

ここでは、例としてGSE49110  の 発現量のテキストファイルをダウンロード してインポートしてみましょう。または、SSAファイルをダウンロードして次のステップに飛ばすこともできます。

まずはダウンロードしたファイルをエクセルで編集しますが、ここにあるとおりちょっとしたコツがあります。いらない列と行を削除して、IDである遺伝子名と発現量のテーブルに整形します。

この編集したテキストファイルを Subio Platform にインポートします。Import Samples ウィザードを開始したら、最初のページで “Multiple Samples in One File” オプションを使い、次のページで “Create A New Platform” を使って下さい。

インポートが完了したら、SOFT formatted family file からサンプルの属性情報を取り込んでおきます。これでより多くの情報を見ることができ、また、キーワードでフィルターをかけられるようになります。

それでは、この8個のSampleから成るSeriesを作って、データの視覚化および解析をしましょう。

RNA-Seq データ解析チュートリアル (01) - RNA-SeqのCountsデータのインポート

上のムービーは、日本語字幕を表示できます。

【最短で結果を出したい方】
データ解析サービスを見る

【自分で操作できるようになりたい方】
オンライントレーニングを見る

目次に戻る

2.Series の作成と設定

要点:統計解析において「実験デザイン(どのサンプルがどの群か)」を正確に定義しなければ解析を始められません。属性情報を整理し、比較の「軸」を明確にすることで、迷いのない解析環境を整えます。

Seriesがロードされたら、画面左のSeriesパネルにMeasurement ListやDataSetなどのオブジェクトが表示されます。Analysis Browserの上段ではデフォルトでScatter Plotが描画されます。一方、下段ではSetup Seriesタブが開きます。
Setup Seriesタブでは、基本的に画面の左から順にボタンを押していけばいいですので、最初のボタンEdit Parametersを押します。通常は、Sample Informationから情報をインポートするのが早いです。このデータセットは、4つの状態(コントロールと3種類のsiRNA処理)があり、それぞれ2回の反復があります。

Setup DataSetタブに移動し、DataSetの編集と作成を行います。先ほど設定したParameterの優先順位を設定することでSample Groupを定義・整理するのがこのタブの役割です。このチュートリアルでは2つのDataSetを設定します。さらに、Sample Info.タブに情報を記録し、関連する添付ファイルを保存します。

RNA-Seq データ解析チュートリアル (02) - Series の作成と設定

上のムービーは、日本語字幕を表示できます。

目次に戻る

3.正規化と前処理

要点:正規化は、総リード数などの違いによって生じる系統誤差を補正する不可欠な工程ですが、盲信は禁物です。(Case Study 379参照)プリセットに頼り切らず、ヒストグラムで補正前・補正後の分布を視覚的に確認し、データの歪みが正しく解消されているかを見極める目を養います。

正規化および前処理は、これをどうするかによってその後の解析結果が変わってしまう大事なステップです。しかも、実際のデータは教科書が想定しているようなものとは限りません。目の前のデータの特徴を正しく理解し、それにあった処理をしなければ、その後の解析で誤った結論を導き出してしまうので十分に注意して下さい。

このチュートリアルでは、プリセットの「RNA-Seq (Counts)」というシナリオから出発して、データに合わせて調整します。その過程で、何を知るためにどこを見るか(操作するか)、そしてどのように判断するかの一例を学んでください。繰り返しますが、ここで紹介するやり方がどんなデータにも通用するわけではありません。実際のデータ解析では、そのデータに合わせてその都度調整が必要です。自信を持てない場合は、どうぞデータ解析サービスをご利用ください。

こちらを併せてご覧ください。

RNA-Seq データ解析チュートリアル (03) - 正規化と前処理

上のムービーは、日本語字幕を表示できます。

【最短で結果を出したい方】
データ解析サービスを見る

【自分で操作できるようになりたい方】
オンライントレーニングを見る

目次に戻る

4.フィルタリング(Quality Control)

要点:「ゴミを入れれば、ゴミが出る(Garbage In, Garbage Out)」。低発現で信頼性の低い遺伝子や、変化のないノイズを適切に除外することで、真に生物学的な意味を持つシグナルだけを抽出するテクニックを学びます。

この前の正規化の説明で見たとおり、測定値があったとしてもすべての測定値が信頼できるわけではありません。解析する前に、解析に値する測定値を抽出し、解析不能の値は解析の邪魔にならないように操作する必要があります。ここで使うのが、Basic Plug-inに含まれるFilterツールです。プラグインをお持ちでない方は、5日間の無料お試しをご利用ください

フィルターを使う際の基本的な考え方は、「信頼できる遺伝子を抽出する」ではなく、「解析に値しない遺伝子を抽出して除く」です。この二つの違いは解析経験がないとまず分からないでしょう。前者は、例えばコントロールでは発現していなかったのに、何らかの処理をすると発現してくる遺伝子を逃してしまうことになります。この違いを、フィルターツールを実際に使ってみて確かめてください。

もう一つの基本的な考え方は、2段階で行うということです。まず、値が低すぎる遺伝子を除きます。次に、発現変動しない遺伝子を除きます。

たとえば、RDESeq2を使った解析手順を独学していると、次のようなサンプルコードがあります。

dds <- dds[rowSums(counts(dds)) >= 10,]

こういうとき、「どんなデータでも同じ閾値を適用していいのかな?」と迷いながらも、よくわからないので同じいいやとコピペしてはいませんか?適切な閾値はどうやって決めたらいいのでしょうか?答えは、「データを見ないと分からない」です。だから100サンプル未満の小さなデータセットであれば、このチュートリアルのように「データを見ながら進める」のが正解だし、100サンプル以上の大きなデータセットであれば、ランダムに100サンプルくらいピックアップして行う仮解析で閾値を見定めたうえで、全体の解析パイプラインを動かすのが正解です。データを見もしないでサンプルコードをそのまま走らせると、とんでもない解析結果が出てくるのは仕方のないことです。

RNA-Seq データ解析チュートリアル (04) - フィルタリング (Quality Control)

上のムービーは、日本語字幕を表示できます。

【5日間無料トライアルはこちら】
ノイズを除去するフィルタリングを試す

目次に戻る

5.PCAとクラスタリング

要点:数万個の遺伝子の動きは人間の理解を超えますので、理解可能なものにする必要があります。サンプルのバラツキ(実験の質)や群間の類似性を直感的に把握し、「このデータで次に何を調べるべきか」という全体像を俯瞰する視点を手に入れます。

前段までで準備が整いましたので、ここからいよいよ、さまざまな解析手法を駆使してデータから何らかの意味を抽出する段階に入ります。すべての遺伝子を対象とするのでなく、フィルターを使ってQuality Controlを通過した遺伝子群だけを対象とすることにご注意ください。

まずは発現プロファイルの全体像を俯瞰することで、データの概観を把握し、解析のポイントを明らかにするのがよいでしょう。この目的で役に立つのが主成分分析(PCA)です。PCAの結果を見るときの要点は3つです。まず、点と点の距離が近ければ発現プロファイルが似ている、そして離れていれば発現プロファイルも大きく異なることを表します。次に、原点からの方向です。同じ方向であればある遺伝子群が似たような動きをしていることを表し、原点から遠いほど変動幅が大きいことを表します。また、原点を挟んで点があるということは、発現変動が逆向きであることを示しています。そして、最後に主成分と寄与率です。縦方向と横方向は、あるいはそれぞれの主成分は異なる遺伝子群の動きを代表していると言えます。寄与率が大きいほど、そのような動きをする遺伝子の数が多く、概観を捉えていると言っていいでしょう。しかし、寄与率が大きいことが必ずしも生物学的に重要というわけではありません。生物学では、ごく少数の遺伝子群の動きが将来の方向性を決定づけることはよくあります。つまり、寄与率の小さな主成分が決定的な動きを捉えているかもしれません。

以上を踏まえて、このデータのPCAの結果を見てみましょう。繰り返しの2点間の距離は、siRNA間の距離に比べてずっと小さく、ばらつきの小さい良質の実験データであることが分かります。原点のコントロールと比べて、各siRNA処理群は同一方向(右下)への動きと、個別の動き(右と下)があるようです。

次に、階層型クラスタリングを適用してみます。ヒートマップを概観すると、すべてのsiRNAで共通で発現上昇(赤)する遺伝子は多くあるのと対照的に、発現が下降(青)する遺伝子はsiRNAごとに異なっていることが分かります。解析者がすべきことは、発現差のある遺伝子リストや図を作ることではなく、このような違いや特徴に気づき、細胞内で何が起きているかを洞察することです。

また、PCAと階層型クラスタリングは、基本的には同様の解析結果を違う見方をしているにすぎません。二つの結果を見比べながら、データの理解を深めましょう。

RNA-Seq データ解析チュートリアル (05) - PCAとクラスタリング

上のムービーは、日本語字幕を表示できます。

【5日間無料トライアルはこちら】
PCAやクラスタリングでデータの全体像を可視化してみる

目次に戻る

6.発現差のある遺伝子の抽出

要点:単なるP値のカットオフだけでなく、ベン図などを活用して多角的に遺伝子を絞り込みます。異なる条件間で「共通して動くもの」と「固有に動くもの」を明確にし、次の考察のターゲットを決定します。

このケースでは、前章のクラスタリングで見た通り、発現が上昇する遺伝子群は3種のsiRNAで共通なのに対し、発現が下降する遺伝子群の共通項はあまりに少ないという偏りがありそうです。Basic Plug-inのツールを使って発現差のある遺伝子を抽出したら、Venn Diagramツールで確認してみましょう。

RNA-Seq データ解析チュートリアル (06) - 発現差のある遺伝子の抽出

上のムービーは、日本語字幕を表示できます。

【5日間無料トライアルはこちら】
発現差解析(DEG)機能を試す

目次に戻る

7.遺伝子アノテーションと、エンリッチメント解析

要点:「数字のリスト」を「生物学的な言葉」に翻訳します。統計的に有意な変動が、細胞内のどのような機能(GO)や代謝経路(Pathway)に集約されているかを突き止め、発見のヒントを得ます。

RNA-Seqのデータの場合、countまたはFPKMのテーブルにGene Symbol(ENSGやENSTのID、Entrez Gene IDの場合もある)しか付いていないことも多いです。そこで、遺伝子アノテーションのテーブルをデータベースサイトから取得して埋める必要があります。ここではNCBI FTP Siteサイトを使いますが、ENSGやENSTのIDがついている場合は Ensembl BioMartから取得してください。それ以外のIDの場合は、そのIDと遺伝子アノテーションを管理しているデータベースサイトを使ってください。

遺伝子アノテーションのインポートが終わると、Subio Platformでの検索機能や、Advanced Plug-inに含まれるEnrichment Analysisツールが活用できるようになります。

用語が紛らわしいのですが、Gene Ontology (GO) 解析、パスウェイ解析、ネットワーク解析などのキーワードで尋ねてこられる場合、その多くがエンリッチメント解析のことをおっしゃっています。Gene Set Enrichment Analysis(GSEA)、Ingenuity Pathway Analysis (IPA)、David Functional Annotation、Meatscapeなどはエンリッチメント解析を行うツールです。

エンリッチメント解析については、別に詳しいチュートリアルがありますので、併せてご覧ください。

RNA-Seq データ解析チュートリアル (07) - 遺伝子アノテーションと、エンリッチメント解析

上のムービーは、日本語字幕を表示できます。

【5日間無料トライアルはこちら】
遺伝子リストからGO・Pathwayを見つけてみる

目次に戻る

8.ゲノム上の位置特異的に発現制御されている遺伝子と、モチーフ配列

要点:物理的な位置関係から、新たな制御メカニズムの可能性を探ります。特定のゲノム領域に依存した変動や、転写因子が結合する「モチーフ配列」との関連を調べることで、発現変動の背景にある「要因」を推察する手がかりを得ます。

発現差のある遺伝子を抽出したり、クラスタリングで特定の発現パターンを持つ遺伝子群を特定したら、それらがゲノム上に偏って存在していないかを確認してみるのも一つの手です。

もし特定の領域にある遺伝子群が一斉に発現上昇または下降していれば、その領域のエピジェネティックな状態変化や、染色体の構造変化が発現変動に関与している可能性が考えられます。特に偏りが見られなければ、転写因子による制御の可能性がより高いと推測できるかもしれません。

また、モチーフ配列の検索機能を活用することで、転写開始点(TSS)近傍に特定の配列を持つ遺伝子を特定し、それらの発現パターンを確認することも可能です。こうした位置情報を活用した解析を組み合わせることで、データに対する理解の解像度をさらに高めることが期待できます。ちなみに、モチーフ配列の検索は、 IUPAC nucleotide code に対応しています。

このケースでは、ゲノム上の位置に偏った発現制御はないように見えます。siRNAによって抑制したERR alphaの結合モチーフ配列はWikipediaで見つけることができたので、これを転写開始位置近傍に持つ遺伝子を抽出することができました。これと、3種すべてのsiRNAで発現抑制された遺伝子リストとの重複を調べると、一つの遺伝子に辿り着きました。しかし、上述のとおり、3種すべてで抑制されている必要はないかもしれないですし、発現差解析の閾値がきつすぎだった可能性もあります。これらの条件を緩めることで、ERR alphaが直接制御する遺伝子の候補はもう少し広がるでしょう。

RNA-Seq データ解析チュートリアル (08) - ゲノム上の位置特異的に発現制御されている遺伝子と、モチーフ配列

上のムービーは、日本語字幕を表示できます。

【5日間無料トライアルはこちら】
コンセンサス配列を特定し、TSS近傍に持つ遺伝子を探してみる

目次に戻る

9.おわりに:ツールを使いこなし、「決断」できる解析者へ

チュートリアルはあくまで「手法(Methods)」の紹介に過ぎません。実際の解析現場において、どの閾値が適切か、バッチエフェクトをどう処理すべきかを判断するのは、ツールではなく解析者自身です。客観的な数字を超えた、「主観を伴う意思決定」の重要性を再確認しましょう。

本チュートリアルでは、RNA-Seqデータ解析の標準的な流れを解説してきました。しかし、実際のデータ解析が一本道で進むことはありません。「正解」はデータの特性や研究の目的によって常に変化します。

手法を学び終えたとき、次に来るのが最も本質的なステップです。それは、「主観に基づいた適切な判断を下せる解析者」へと進化することです。

現実の測定データは完璧ではありません。だからこそ、不足している部分は解析者の「主観」で補う必要があります。ただし、その主観はデータと研究目的に基づいた、妥当なものでなければなりません。

それこそが、AI時代においても価値を失わない「生き残れる解析者」の姿ではないでしょうか。

Subio Platformは、あなたが「ツールに使われる作業者」から、「自らの主観で解析を支配し、意思決定を行える研究者」へと進化するための解析基盤です。
ソフトウェアの提供に加え、解析サービスとオンライントレーニングによってそのプロセスを一貫してサポートしています。

次は、あなた自身のデータで、この違いをぜひ体感してみてください。

▶ 今すぐ自分のデータを解析したい方はこちら(データ解析サービス)
▶ 自分のデータで操作を習得したい方はこちら(オンライントレーニング)

目次に戻る