macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

並列化に対応したアダプタートリミングツール AdapterRemoval 2

 

 Fossil material 由来などのごく短いDNA断片のハイスループットシーケンスでは、ライブラリーの調製中にインサートにライゲーションされたアダプター配列をシークエンシングする可能性がある[論文より ref.1]。このような汚染はよく知られた問題であり、下流分析にマイナスの影響を与える可能性がある[ref.2-5]。したがって、ワークフローの最初の部分には、通常、アダプタの汚染をフィルタリングまたは除去(トリミング)するステップが含まれている[ref.3]。(古代のクオリティの低いサンプルでも)ペアエンドリードののオーバーラップする領域を検出し、これらを品質認識的に崩壊(合併)させて鋳型分子全体を再構築することにより、精度をあげることができる[ref.6]。これは、短いインサートサイズが期待されるancient DNAシークエンシング、およびリード末端への死後のDNA修飾に起因する、すでに高い配列決定エラー率の減少が興味深い[ref.7 pubmed]場合に特に興味深い。

 AdapterRemoval v1 [ref.2]のオリジナルリリースでは、Needleman-Wunschアルゴリズムの修正版を使用して、アダプター配列と低品質塩基のトリミングのための使いやすいツールが提供されており、リード間のペアワイズアライメント(またはペア)および既知のアダプター配列を含まれる。シングルエンドリードの場合、このアライメントは、アダプター配列の5 '末端とrawシーケンシングリードの3'末端との間の最良のアラインメントを見出し、ペアエンドリードの場合には、mate1のリードと、mate2のリードの逆相補との間で最良のアラインメントを見出す。

 これらの凝集配列のアライメントは、各リードにおけるの3'末端の同定を可能にし、それにより外来性(アダプター)配列をトリミングすることを可能にする[ ref.2]。しかしながら、使用された改良Needleman-Wunschアルゴリズムはシーケンシング技術によって導入されたindelを考えないので、AdapterRemovalは主に、false indelsが少ないIllumina HTSプラットフォームから生成されたリード処理に適している。  AdapterRemovalはさらに、(オプションとして)重複するリードをマージするため、アダプターの位置合わせ手順の一部として重複フラグメントを使用する。これは、最高品質の塩基を選択し、重複位置で位置特異的スコアリングマトリックスとして扱うことで基本品質を再計算することによって達成される[ref.2]。

 しかし、AdapterRemoval v1は他の最新のツール[ref.8]と比較して、実行時間が比較的遅いという特徴がある。したがって、採用されたトリミング方法を変更することなく、スループットを向上させる目的で、AdapterRemoval v1の大幅な改訂を実施した。したがって、更新されたバージョンは、AdapterRemoval v1と同様の精度を示し、既存の分析パイプライン[ref.7]での使用に適切な互換品を提供する。

 AdapterRemoval v2は、AdapterRemoval v1で使用されるアライメントアルゴリズムを高速化するために、単一命令、複数データ(SIMD)命令(つまりコンシューマーグレードCPUで一般的にサポートされるSSE命令とSSE2命令)を使用することにより、複数のスレッドを使用してすべての操作のスループットをさらに向上させることができる。さらに、AdapterRemoval v2では、複数の異なる(ペアの)アダプタシーケンスを同時にトリミングし、リード(ペア)ごとに最適な一致を選択し、gzipおよびbzip2圧縮FASTQファイルを一貫して読み書きできる。 AdapterRemoval v2では、さらに、提供されたバーコードの単純な最大不一致数の比較を使用して、同時にdemultiplexingとアダプターのトリミングが可能になり、インサートの端を超えてシーケンスされるペアエンドリードも正しくトリミングされる。このようなリードは、ペアエンドの相方の逆相補バーコード配列に続いてアダプター配列によって終結され、両方とも除去されなければならない。最後に、AdapterRemoval v2は、インサートがオーバーラップするペアから推定アダプター配列を再構成することができる。これを使用して、実験エラーを検出し、文書化されていないデータセットを分析することができる。

 AdapterRemoval v2のパフォーマンスを評価するために、アダプターシーケンスを正しく調整する機能と、1つまたは複数のスレッドを使用する場合のスループットの点で、AdapterRemoval v1と現代のアダプタートリミングソフトウェアの選択肢を比較した。さらに、いくつかの他の重複シーケンスをマージするツールを使いAdapterRemovalの機能を比較した。後者については、比較をアダプター汚染の存在下でリードをマージするように設計されたプログラムに制限したが、アダプター汚染をほとんどまたはまったく含まないデータセットに適したいくつかの選択肢、例えばCOPE [9]、fastq-join [10]、FLASH [11]、およびXORRO(arXiv:1304.4620)に注意する。

 

 

 特徴(Githubより)

  • Trimming of adapters sequences from single-end and paired-end FASTQ reads.
  • Trimming of multiple, different adapters or adapter pairs.
  • Demultiplexing of single or double indexed reads, with or without trimming of adapter sequences.
  • Reconstruction of adapter sequences from paired-end reads, by the pairwise alignment of reads in the absence of a known adapter sequence.
  • Merging of overlapping read-pairs into higher-quality consensus sequences.
  • Multi-threading of all operations for increased throughput.
  • Reading and writing of gzip and bzip2 compressed files.
  • Reading and writing of interleaved FASTQ files.

 

インストール

mac os10.12でテストした。

本体 Github

https://github.com/MikkelSchubert/adapterremoval

git clone https://github.com/MikkelSchubert/adapterremoval.git 
cd adapterremoval
make
sudo make install

$ sudo make install

Password:

Building AdapterRemoval with gzip support: yes

Building AdapterRemoval with bzip2 support: yes

Building AdapterRemoval with pthreads support: yes

Linking executable 'build/AdapterRemoval'

clang: warning: argument unused during compilation: '-ansi' [-Wunused-command-line-argument]

Constructing man-page 'build/AdapterRemoval.1' from 'AdapterRemoval.pod'

Installing AdapterRemoval ..

  .. binary into /usr/local/bin/

  .. man-page into /usr/local/share/man/man1/

  .. README into /usr/local/share/adapterremoval/

  .. examples into /usr/local/share/adapterremoval/examples/

 

AdapterRemoval

$ AdapterRemoval

AdapterRemoval ver. 2.1.7

 

This program searches for and removes remnant adapter sequences from

your read data.  The program can analyze both single end and paired end

data.  For detailed explanation of the parameters, please refer to the

man page.  For comments, suggestions  and feedback please contact Stinus

Lindgreen (stinus@binf.ku.dk) and Mikkel Schubert (MikkelSch@gmail.com).

 

If you use the program, please cite the paper:

    Schubert, Lindgreen, and Orlando (2016). AdapterRemoval v2: rapid

    adapter trimming, identification, and read merging.

    BMC Research Notes, 12;9(1):88.

 

    http://bmcresnotes.biomedcentral.com/articles/10.1186/s13104-016-1900-2

 

 

Arguments:                           Description:

  --help                             Display this message.

  --version                          Print the version string.

 

  --file1 FILE                       Input file containing mate 1 reads or single-ended reads [REQUIRED].

  --file2 FILE                       Input file containing mate 2 reads [OPTIONAL].

 

FASTQ OPTIONS:

  --qualitybase BASE                 Quality base used to encode Phred scores in input; either 33, 64, or solexa

                                       [current: 33].

  --qualitybase-output BASE          Quality base used to encode Phred scores in output; either 33, 64, or solexa. By

                                       default, reads will be written in the same format as the that specified using

                                       --qualitybase.

  --qualitymax BASE                  Specifies the maximum Phred score expected in input files, and used when writing

                                       output. ASCII encoded values are limited to the characters '!' (ASCII = 33) to

                                       '~' (ASCII = 126), meaning that possible scores are 0 - 93 with offset 33, and

                                       0 - 62 for offset 64 and Solexa scores [default: 41].

  --mate-separator CHAR              Character separating the mate number (1 or 2) from the read name in FASTQ

                                       records [default: '/'].

  --interleaved                      This option enables both the --interleaved-input option and the

                                       --interleaved-output option [current: off].

  --interleaved-input                The (single) input file provided contains both the mate 1 and mate 2 reads, one

                                       pair after the other, with one mate 1 reads followed by one mate 2 read. This

                                       option is implied by the --interleaved option [current: off].

  --interleaved-output               If set, trimmed paired-end reads are written to a single file containing mate 1

                                       and mate 2 reads, one pair after the other. This option is implied by the

                                       --interleaved option [current: off].

 

OUTPUT FILES:

  --basename BASENAME                Default prefix for all output files for which no filename was explicitly set

                                       [current: your_output].

  --settings FILE                    Output file containing information on the parameters used in the run as well as

                                       overall statistics on the reads after trimming [default: BASENAME.settings]

  --output1 FILE                     Output file containing trimmed mate1 reads [default: BASENAME.pair1.truncated

                                       (PE), BASENAME.truncated (SE), or BASENAME.paired.truncated (interleaved PE)]

  --output2 FILE                     Output file containing trimmed mate 2 reads [default: BASENAME.pair2.truncated

                                       (only used in PE mode, but not if --interleaved-output is enabled)]

  --singleton FILE                   Output file to which containing paired reads for which the mate has been

                                       discarded [default: BASENAME.singleton.truncated]

  --outputcollapsed FILE             If --collapsed is set, contains overlapping mate-pairs which have been merged

                                       into a single read (PE mode) or reads for which the adapter was identified by

                                       a minimum overlap, indicating that the entire template molecule is present.

                                       This does not include which have subsequently been trimmed due to low-quality

                                       or ambiguous nucleotides [default: BASENAME.collapsed]

  --outputcollapsedtruncated FILE    Collapsed reads (see --outputcollapsed) which were trimmed due the presence of

                                       low-quality or ambiguous nucleotides [default: BASENAME.collapsed.truncated]

  --discarded FILE                   Contains reads discarded due to the --minlength, --maxlength or --maxns options

                                       [default: BASENAME.discarded]

 

OUTPUT COMPRESSION:

  --gzip                             Enable gzip compression [current: off]

  --gzip-level LEVEL                 Compression level, 0 - 9 [current: 6]

  --bzip2                            Enable bzip2 compression [current: off]

  --bzip2-level LEVEL                Compression level, 0 - 9 [current: 9]

 

TRIMMING SETTINGS:

  --adapter1 SEQUENCE                Adapter sequence expected to be found in mate 1 reads [current:

                                       AGATCGGAAGAGCACACGTCTGAACTCCAGTCACNNNNNNATCTCGTATGCCGTCTTCTGCTTG].

  --adapter2 SEQUENCE                Adapter sequence expected to be found in mate 2 reads [current:

                                       AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTATCATT].

  --adapter-list FILENAME            Read table of white-space separated adapters pairs, used as if the first column

                                       was supplied to --adapter1, and the second column was supplied to --adapter2;

                                       only the first adapter in each pair is required SE trimming mode [current:

                                       <not set>].

 

  --mm MISMATCH_RATE                 Max error-rate when aligning reads and/or adapters. If > 1, the max error-rate

                                       is set to 1 / MISMATCH_RATE; if < 0, the defaults are used, otherwise the

                                       user-supplied value is used directly. [defaults: 1/3 for trimming; 1/10 when

                                       identifing adapters].

  --maxns MAX                        Reads containing more ambiguous bases (N) than this number after trimming are

                                       discarded [current: 1000].

  --shift N                          Consider alignments where up to N nucleotides are missing from the 5' termini

                                       [current: 2].

 

  --trimns                           If set, trim ambiguous bases (N) at 5'/3' termini [current: off]

  --trimqualities                    If set, trim bases at 5'/3' termini with quality scores <= to --minquality value

                                       [current: off]

  --minquality PHRED                 Inclusive minimum; see --trimqualities for details [current: 2]

  --minlength LENGTH                 Reads shorter than this length are discarded following trimming [current: 15].

  --maxlength LENGTH                 Reads longer than this length are discarded following trimming [current:

                                       4294967295].

  --collapse                         When set, paired ended read alignments of --minalignmentlength or more bases are

                                       combined into a single consensus sequence, representing the complete insert,

                                       and written to either basename.collapsed or basename.collapsed.truncated (if

                                       trimmed due to low-quality bases following collapse); for single-ended reads,

                                       putative complete inserts are identified as having at least

                                       --minalignmentlength bases overlap with the adapter sequence, and are written

                                       to the the same files [current: off].

  --minalignmentlength LENGTH        If --collapse is set, paired reads must overlap at least this number of bases to

                                       be collapsed, and single-ended reads must overlap at least this number of

                                       bases with the adapter to be considered complete template molecules [current:

                                       11].

  --minadapteroverlap LENGTH         In single-end mode, reads are only trimmed if the overlap between read and the

                                       adapter is at least X bases long, not counting ambiguous nucleotides (N); this

                                       is independant of the --minalignmentlength when using --collapse, allowing a

                                       conservative selection of putative complete inserts while ensuring that all

                                       possible adapter contamination is trimmed [current: 0].

 

DEMULTIPLEXING:

  --barcode-list FILENAME            List of barcodes or barcode pairs for single or double-indexed demultiplexing.

                                       Note that both indexes should be specified for both single-end and paired-end

                                       trimming, if double-indexed multiplexing was used, in order to ensure that the

                                       demultiplexed reads can be trimmed correctly [current: <not set>].

  --barcode-mm N                     Maximum number of mismatches allowed when counting mismatches in both the mate 1

                                       and the mate 2 barcode for paired reads.

  --barcode-mm-r1 N                  Maximum number of mismatches allowed for the mate 1 barcode; if not set, this

                                       value is equal to the '--barcode-mm' value; cannot be higher than the

                                       '--barcode-mm value'.

  --barcode-mm-r2 N                  Maximum number of mismatches allowed for the mate 2 barcode; if not set, this

                                       value is equal to the '--barcode-mm' value; cannot be higher than the

                                       '--barcode-mm value'.

 

MISC:

  --identify-adapters                Attempt to identify the adapter pair of PE reads, by searching for overlapping

                                       reads [current: off].

  --seed SEED                        Sets the RNG seed used when choosing between bases with equal Phred scores when

                                       collapsing. Note that runs are not deterministic if more than one thread is

                                       used. If not specified, a seed is generated using the current time.

  --threads THREADS                  Maximum number of threads [current: 1]

 

——

 

condaでも導入できます(Github)。

 

 

ラン

/usr/local/share/adapterremoval/examples/にあるテストシーケンスデータをランする。リンクを貼っておく。

ln -s /usr/local/share/adapterremoval/examples/* .

 

 

1、アダプター配列を探す。

AdapterRemoval --identify-adapters --file1 reads_1.fq --file2 reads_2.fq 
  • --identify-adapters    Attempt to identify the adapter pair of PE reads, by searching for overlapping  reads [current: off].
  • --file1    Input file containing mate 1 reads or single-ended reads [REQUIRED].  
  • --file2    Input file containing mate 2 reads [OPTIONAL].

 

2、シングルエンドのアダプタートリミング。(アダプターは自動で検出される。low qualityとNのトリミングも実行される)

AdapterRemoval --file1 reads_1.fq --basename output_single --trimns --trimqualities --gzip --threads 8
  • --basename    Default prefix for all output files for which no filename was explicitly set
  • --trimns    If set, trim ambiguous bases (N) at 5'/3' termini [current: off]
  • --gzip     Enable gzip compression [current: off]
  • --threads   Maximum number of threads [current: 1]

3つのファイルが出力される。アダプタートリミングされたリードは~truncated~というファイルに保存される。settings~にはシーケンスデータのstatisticsや実行条件がまとめられている。

 

3、ペアエンドのアダプタートリミング。

AdapterRemoval --file1 reads_1.fq --file2 reads_2.fq --basename output_paired --gzip --trimqualities --collapse --threads 8
  • --trimqualities    If set, trim bases at 5'/3' termini with quality scores <= to --minquality value   [current: off]

 7つのファイルが出力される。ペアエンドの両方のリードのアダプタートリミングがうまくいったリードはtruncated~に保存される。ペアエンドの片方がdiscardされ順番が壊れたリードはcollapse~に保存される。ペアエンドがマージされた場合はoutput_paired.collapsed~に保存され、--trimns や --trimqualitiesオプションの影響で トリムされてマージされたリードはoutput_paired.collapsed.truncatedに保存される。

f:id:kazumaxneo:20180516150936j:plain

 

4、ペアエンドのinterleave fastqのアダプタートリミング。

AdapterRemoval --interleaved --file1 interleaved.fq --basename output_interleaved

 

注:デフォルトではPhred+33 としてトリミングを実行する。昔のPhred+64のfastqなら"--qualitybase 64"をつけて実行する。"--qualitybase-output 33"もつけて実行するとPhred+33に修正して出力してくれる。

 

この他、Githubではv2.1で実装されたDemultiplexing機能を使い、Demultiplexingとアダプタートリミングを同時実行する方法も紹介されています。さらに、double indexのDemultiplexingや、アダプターが複数ある特殊なシーケンスデータのアダプター除去も可能です。詳細はGithubで確認してください。

引用

AdapterRemoval v2: rapid adapter trimming, identification, and read merging

Mikkel Schubert, Stinus Lindgreen, Ludovic Orlando

BMC Res Notes. 2016; 9: 88.

 

AdapterRemoval: easy cleaning of next-generation sequencing reads.

Lindgreen S

BMC Res Notes. 2012 Jul 2;5:337.