PRINSEQ++(plusplus)でRNA-seqのRAWデータのQC(クォリティコントロール)をする

これまでの記事ではfastq-dumpを用いてNCBI上からRNA-seqのデータをダウンロードするなどの「データの入手」に焦点を当ててきた。

参考: SRA Toolkitの使い方 ~fastq-dumpでSRAファイルをダウンロード~

今回の記事ではそれら、もしくは自分がRNA-seqに解析委託を出して返ってきたデータに対してまず行う「クォリティチェック」のあとの段階で行う「クォリティコントロール(QC)」について解説を行う。クォリティチェックについては後々記事にする。

クォリティコントロール(QC)はなぜするのか？

様々なサイトでこれについては解説がなされているが要約すると

• ポリAテール(末端)を除去
• 信頼度の低いリードの末端を除去
• ゴミリード(短すぎるリード)を除去

などが挙げられる。

そして今回はそれを行っていく。

PRINSEQについて

PRINSEQはこういったQCのために開発されたプログラムで、多種多様な機能を持ち合わせるPerlで作られた……そう、Perlで作られたプログラムだ。

そのとおり、結構解析には時間がかかる。

また、私が行った解析においてはめちゃくちゃなメモリ使用と、バグによって結局解析がうまくいかなかった(メモリ64GBでは足りない？)。

使い方はこちらのサイトに詳しいので、もしオリジナルのPRINSEQを使う場合はこちらを参照したほうがいいだろう。

こうした遅い、高メモリ要求、マルチスレッド非対応といった従来のPRINSEQを踏まえて、新たに今年(2019年)になってPRINSEQ++(plusplus)が発表された。C++で書かれているから”++”。

Cantu, V. A., Sadural, J., & Edwards, R. (2019). PRINSEQ++, a multi-threaded tool for fast and efficient quality control and preprocessing of sequencing datasets. PeerJ Preprints, 7, e27553v1. リンク

このPRINSEQ++は従来のPRINSEQと比べて5倍ほど高速になっているという。実際使ってみたところ、かなり早くエラーなく解析を終えることができた非常に優秀なソフトウエアである。

PRINSEQ++のインストール

さて、では早速このPRINSEQ++をインストールしていこう。想定している環境はUbuntuである。

GitHubから念の為確認してほしいが、現時点(2021/05/25)での最新のバージョン(v1.2)のダウンロードは以下の2通り。

GitHubからClone

2021/05/25追記: こっちのほうがかんたん。

$ git clone https://github.com/Adrian-Cantu/PRINSEQ-plus-plus.git
$ cd PRINSEQ-plus-plus
$ ./autogen.sh
$ ./configure
$ make
$ make test 
$ sudo make install

従来の方法

$ wget https://github.com/Adrian-Cantu/PRINSEQ-plus-plus/releases/download/v1.2/prinseq-plus-plus-1.2.tar.gz

ダウンロードし終わったら以下の要領で解凍、コンパイル、インストール。

$ tar -xvf prinseq-plus-plus-1.2.tar.gz
$ cd prinseq-plus-plus-1.2
$ ./configure
$ make
$ make test
$ sudo make install

全部パスも通った状態になるのですぐに使えるようになる。

$ prinseq++ -h

で、

        PRINSEQ++ 1.2

PRINSEQ++ is a C++ implementation of the prinseq-lite.pl program. It can be used
to filter, reformat or trim genomic and metagenomic sequence data. It is 5X faster
than prinseq-lite.pl and uses less RAM thanks to the use of multi-threading
and the cboost libraries. It can read and write compressed (gzip) files, drastically
reducing the use of hard drive.


Option:
    -h | -help
        Print the help page; ignore other arguments.

    -v | -version
        Print version; ignore other arguments.

    -threads
        Nuber of threads to use. Note that if more than one thread is used, output


    -noiupac
        Filter sequence with characters other than A, C, G, T or N.

    -derep
        Filter duplicated sequences. This only remove exact duplicates.

    -lc_entropy=[float]
        Filter sequences with entropy lower than [float]. [float] should be in
        the 0-1 interval. (Default=0.5)

    -lc_dust=[float]
        Filter sequences with dust_score lower than [float]. [float] should be in
        the 0-1 interval. (Default=0.5)

    ***** TRIM OPTIONS *****

    -trim_left <integer>
        Trim <integer> bases from the left (5'->3').

    -trim_right <integer>


    -trim_qual_left=[float]
        Trim recursively from the 3'-end chunks of length -trim_qual_step if the
        mean quality of the first -trim_qual_window bases is less than [float].
        (Default=20)

    -trim_qual_right=[float]
        Trim recursively from the 5'-end chunks of length -trim_qual_step if the
        mean quality of the last -trim_qual_window bases is less than [float].
        (Default=20)

    -trim_qual_window [int]
        Size of the window used by trim_qual_left and trim_qual_right (Default=5)

    -trim_qual_step [int]
        Step size used by trim_qual_left and trim_qual_right (Default=2)

    -trim_qual_type <string>
        Type of quality score calculation to use. Allowed options are min,
        mean, max and sum. [default= min]

と出てくる。

これを見てやればいいだけなのだが、今回は実際に前回落としたようなペアエンドのリードのクオリティコントロールを行ってみる。

使用例

$ prinseq++ -fastq SRR8137461_1.fastq -fastq2 SRR8137461_2.fastq -out_name trim_SRR8137461 -trim_left 5 -trim_tail_right 5 -trim_qual_right 30 -ns_max_n 20 -min_len 30 -threads 24

モロバレのアクセッション番号で研究内容をばらしていくスタイル

何らかの形で手に入れたペアエンドのfastqファイルをダウンロードし、用意されている状態から始まる。

参考: SRA Toolkitの使い方 ~fastq-dumpでSRAファイルをダウンロード~

まずインプットオプション“-fastq”、”-fastq2″でインプットfastqファイルを指定する。ペアエンドの場合は両方指定してやる。

FASTQ形式以外にもFASTA形式も使うことが可能だが、「出力はデフォルトでFASTQ」「クォリティ値は一律31 (A)として扱う」となっている。

“-out_name”で出力の際のファイル名を指定。FASTQ形式で出てくる。

“-trim_left” でleft方向 (5′->3′) からいくつの塩基を削除するかしていする。

“-trim_tail_left” で3’末端に付加されたポリA/T鎖を削除する。

“-trim_qual_right” では、指定したクォリティ値と比較して5’末端から再帰的に-trim_qual_windowの範囲でクォリティが低い部分を削っていく。-trim_qual_windowのデフォルトは5である。

“-ns_max_n” はそのリードにおいて指定した値以上Nが出てきたリードは切り捨てるというもの。オリジナルのPRINSEQでは割合だったが、こちらは絶対値。ちょっと不便な気がする。

“-min_len” で指定した値未満の長さのリードを切り捨てる。短すぎる配列は意味をなさないので切ってしまう。

“-threads” で使用するスレッド数を指定してやる。

といった例を紹介した。

Ryzen5 3400Gでもそう時間はかからなかったので、かなり早いほうだろう。

関連記事

RNA-seqデータをSTARでゲノムデータにマッピングする

参考文献リスト

• Schmieder R and Edwards R: Quality control and preprocessing of metagenomic datasets. Bioinformatics 2011, 27:863-864. [PMID: 21278185]
• Cantu, V. A., Sadural, J., & Edwards, R. (2019). PRINSEQ++, a multi-threaded tool for fast and efficient quality control and preprocessing of sequencing datasets. PeerJ Preprints, 7, e27553v1. リンク