はじめての環境DNA
研究室から現場まで、実践で使える知識を!
データの解析方法

SRA-toolkitを使ってシーケンスファイルを取得する

自分が扱ったことのない分析機器で分析されたデータを解析してみたいなと思うことありませんか?

やっぱり、所属機関にある分析機器にも制限がありますし、ほしいと思ってもなかなか手に入れられないのが普通だと思います。

技術者としてデータ解析は実際に動かしてみたことがあるかどうか、調べてみたことがあるかが一番大切だと思っています。とりあえず環境DNA分析の過程で必要となりそうなものに関して取り組めるようにしておきたいですね。

私が主に扱うのはハイスループットシーケンサー(HTS)から出力されるクオリティスコア付きの塩基配列データです。そして、HTSのデータは論文掲載時にNCBIやDDBJなどの公共データバンクに登録する必要があります。

公共データバンクはInternational Nucleotide Sequence Database Collaboration (INSDC) が運用しており、DDBJ (DNA Data Bank of Japan)と欧州のEBI、米国のNCBIによって構成されています。

ISNDCについて -DDBJ

DDBJ (DNA Data Bank of Japan)は、遺伝研が運営する塩基配列データベース管理拠点で、日本からの登録の99%以上はDDBJを通じて行われています。

上記記事で紹介しているMiFishプライマーの開発論文には、以下のようにDRR030411–030428とあります。

Data accessibility

Custom Ruby scripts used in in silico evaluation of interspecific variation are available from http://dx.doi.org/10.5061/dryad.54v2q. Raw reads from the MiSeq sequencing are available from the DDBJ Sequence Read Archive (DRR030411–030428). The bioinformatic pipeline from data pre-processing through taxonomic assignment (including Perl scripts) is available from http://dx.doi.org/10.5061/dryad.n245j

https://royalsocietypublishing.org/doi/10.1098/rsos.150088

これは、DDBJから発行されるSubmission(DRA), Experiment(DRX), Run(DRR), Analysis(DRZ)のアクセッション番号の一つです。日本語サイトだと以下の記事が分かりやすいです。

次世代シーケンスデータベース(SRA)の見かた

DDBJ Sequence Read Archive Handbook

DRAの中にDRXやDRR, DRZの情報が含まれているので、情報を遡っていくとMiFishプライマーの開発論文に関する分析情報はDRA003066というアクセッションナンバーにまとめられているということが分かります。

DRAの番号が分かれば、含まれているFastq.gzファイルや分析の情報などを知ることができます。

NCBI SRAからFastq.gzファイルをダウンロード

Fastq.gzをDDBJからダウンロードしようとした際、SRAファイルがなかったり、ダウンロードがうまくいかない場合があります。

DDBJのデータはNCBIにもミラーリングされているため、基本的にNCBI SRAからデータを取得するようにしておくと良いと思います。

先ほどの論文のDRA Number(DRA003066)をNCBI SRAで検索して確認してみましょう。

DRR003066の情報は18個に分かれており、1つ目のデータを見てみると以下のような感じになっています。

黄色下線部の部分を見るとPAIREDになっているので、このデータはPaired-endシーケンスデータになります。

以降は配列データを扱うprefetchを利用してSRAファイルを取得後、fasterq-dump(fastq-dumpでマルチスレッド処理することで高速化)でFastqファイルへと変換する方法と、直接Fastqファイルをダウンロードする方法を解説していきたいと思います。

動作環境

  • windows10
  • WSL
  • Ubuntu LTS 20.04
  • fastq-dump 2.9.3-1
  • prefetch 2.9.3-1

SRA Toolkit

SRA Toolkitは公共データバンク(INSDC)のデータを扱うためのツールとライブラリパッケージです。windows, max os, linuxそれぞれに対応しています。

SRA-Toolkitのインストール

sra-toolkitはaptコマンドで入手することができますが、Ubuntu LTS 20.04ではまだサポートされていないようです。

Ubuntu LTS 18.04以前の方は、

updateコマンドを実行して、パッケージリポジトリを更新し、最新のパッケージ情報を取得してから、

でインストールすることができます。

Ubuntu LTS 18.04以降のバージョンの方は、古いファイルをインストールすることで使用できます。

updateコマンドを実行して、パッケージリポジトリを更新し、最新のパッケージ情報を取得しておきます。

もう少し楽な方法が見つかったので更新しました。以下は実行しないでください↓


ダウンロードはLinuxコマンドのwgetを利用します。


SRA-Toolkitをインストールする先のフォルダを作成して、パスを通します。永続的にパスを反映させるのであればホームディレクトリにある.bashrcにパスを記載しておくと楽なのでそこまで解説します。

ここにコンパイル済のパッケージをwgetコマンドでダウンロードします。

※ https://ftp ~ 以降は実行時のsra-toolkitのバージョンに合わせてください。
  link : https://github.com/ncbi/sra-tools/wiki/01.-Downloading-SRA-Toolkit

ダウンロードできたら圧縮ファイルをtarコマンドで解凍します。

2022年1月現在ではsratoolkit.2.11.3-ubuntu64フォルダが作成されます。この中のbinフォルダに必要な実行コマンドが入っているので、そこにパスを通します。

一時的でいい場合は、

永続的にパスを通しておきたい場合は、nanoコマンドを使用し、.bashrcファイルを編集します。

編集モードになったら一番下まで下矢印でいって、下記を記述してください。

編集内容を保存して出る場合は、Ctrl + xyです。

完了したらfasterq-dumpとprefetchのバージョンを確認して、インストールできたか確認してみましょう。

2022/1月現在でfasterq-dump: 2.11.3, prefetch: 2.11.2と返ってきました。

prefetch

開発者githubはこちら。オプションは以下の通り。

prefetchはSRAデータをダウンロードしてくるコマンドです。一度.sraファイルをダウンロードしておくとfastqへの変換に失敗してしまったりした際にやり直しが楽です。

DRR030411–030428 のRunデータがまとめられているので、それらを一括でダウンロードします。

まず、DRR030411- 030428までをリストにしたテキストファイルを用意します(以下のファイルを使っていただいても結構です)。

ダウンロードしたいRunデータのリスト(SRR_Acc_list.txt)が用意できたら、下記コマンドを実行して、.sraファイルをダウンロードします (10分くらい 計約4GB) 。

./の部分は例えばResultsフォルダなどを作ってもらって./Resultsにするとそのフォルダにファイルが保存されます。

また、listファイルを作成したりするのが面倒で、DRRの番号が連番なら以下のようなコマンドで同様の結果を得ることもできます。

次はfastqerq-dumpを使ってダウンロードした.sraファイルからFastqファイルを生成してみましょう。

fasterq-dump

開発者githubはこちら。オプションは以下の通り。

基本的に環境DNAメタバーコーデイングのシーケンスデータはpaired-endで読まれていることが多いです。Illuminaシーケンサーでのpaired-endシーケンスは解読したい150bp ~ 500bp程度の領域の配列を150bp ~ 300bpずつ解読できるようなキットで両側から配列を読んで、オーバーラップしている部分を繋ぎ合わせて配列決定するような方法です。

single-endシーケンスと比べて以下のように長い配列長で読むことができるのが特徴です。

.sraファイルからFastqファイルを取得

ダウンロードした.sraファイルをfastqファイルに変換します。

1行のコマンドで変換は終了です。

オプションを2つ使用しています。-pは進捗状況を表示、-eは複数のコアを使用しての並列処理オプションです。*2022/2/23 追加

_1.fastqと_2.fastq はイルミナのHTSから出力された状態だとR1, R2にあたるファイルで、fowardプライマーとReverseプライマーで読まれた配列になります。

  • DRR0304XX.sra_1.fastq
  • DRR0304XX.sra_2.fastq

.fastqファイルはgzipで圧縮して名前を_1 →_R1、_2 → _R2に変換しておきます。

これで、それぞれの.sraファイルのfastqファイルを取得することができました。

paired-end sequenceをfasterq-dumpでダウンロード

.sraファイルからfastqファイルを作成するのに使用したfasterq-dumpは直接SRAから.fastqファイルをダウンロードすることもできます。その際、注意点としてpaired-endでシーケンスされた.sraデータは--split-filesオプションを使用しないと、つながった状態のものがダウンロードされるようです。

そのため、分割してR1とR2ファイルをダウンロードするには--split-filesオプションを使用します。DRA003066の中にあるFastqファイルはDRR030411–030428と複数ファイルに分かれて登録されているので、forを使ってそれぞれのファイルをダウンロードします。

これでDRR0030XXのXXの部分を11から28まで変化させて自動的にファイルをダウンロードすると以下のようなファイル名のFastqファイルが出力されます。

  • DRR0304XX_1.fastq
  • DRR0304XX_2.fastq

fastqファイルがダウンロードできればprefetchの章でも解説したように、圧縮しつつファイル名を変更しておきます。

まとめスクリプト

prefetchを使用して、.sraファイルをダウンロードしたのちに_R1.fastq.gz, _R2.fastq.gzファイルを取得するパターン。

.sraファイルはsraフォルダを、fastqファイルはfastqフォルダを、.fastq.gzはfastqgzフォルダを作成してそこに保存するようにしました。

また、以下はアクセッションナンバーのリストが記載されたtxtファイルを用意していない前提のスクリプトです。

以下は、fastqerq-dumpを使ってSRAから直接fastqファイルをダウンロードするまとめスクリプトです。

オプション : PMiFish piplineで解析してみる

PMiFish piplineについては以下の記事を参考に進めてもらえればと思います。

MiFishプライマーの開発論文では、MiFish-UとMiFish-Eの2種類のプライマーセットを使用しているので、PMiFish piplineのDatabaseフォルダにあるPrimer.txtは以下のように記載します。

複数のプライマーを使用しているので、setting.txtのDivideはNo→Yesに変更しておきます。

あとはそのままの設定でもOKです。パイプラインを動かして結果を出力してみましょう。

これ以降の解析についてはまた別記事で解説したいと思います。最後まで読んでくださりありがとうございます。

参考

ABOUT ME
しばた
分析・データ解析なんでもします。環境DNA使って研究したい方、気軽に相談ください。