macでインフォマティクス

macでインフォマティクス

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

高速なロングリードのマッピング、エラー修正、アセンブリツール MECAT

2020 2/7 タイトル修正

 

 MECATは、1分子シークエンシング(SMRT)リードの超高速マッピング、エラー訂正、およびデノボアセンブリを行うツール。State of the artのアライナとエラー訂正ツールよりもはるかに効率的な、新しいアライメントとエラー訂正アルゴリズムを採用している。 MECATは、ラージゲノムの効率的なde novo アセンブリに使用できる。例えば、2.0GHz CPUを搭載した32スレッドコンピュータ環境下では、MECATは54xのSMRTヒトゲノムシーケンスデータを9.5日でアセンブリできる。これは現在のPBcR-Mhap pipelineの40倍速い。また、MECATを用いて、diploidのヒトゲノムの102x SMRTデータをわずか25日でアセンブリできる。後者のアセンブリは、54倍の一倍体SMRTデータから組み立てられた以前のゲノムの品質を大幅に改善するものである。 MECATの性能は、PBcR-Mhapパイプライン、FALCONおよびCanu(v1.3)と5つの実際のデータセットで比較した。 MECATによって作成されたコンティグの品質は、PBcR-MhapパイプラインおよびFALCONと同等以上だった。Githubの表に2.0GHzのCPUと512GBのRAMメモリを備えた32スレッドコンピュータでの上記ツールの比較がある(リンク)。

 

ランニングタイム

f:id:kazumaxneo:20180714182831p:plain

 Githubより転載。

 

"MECAT assembly"に関するツイート。

 

インストール

ubuntu16.04でソースからテストした。またdockerイメージもテストした。

依存

  • HDF5
  • dextract
#HDF5
wget https://support.hdfgroup.org/ftp/HDF5/releases/hdf5-1.8/hdf5-1.8.15-patch1/src/hdf5-1.8.15-patch1.tar.gz
tar xzvf hdf5-1.8.15-patch1.tar.gz
mkdir hdf5
cd hdf5-1.8.15-patch1
./configure --enable-cxx --prefix=/public/users/chenying/smrt_asm/hdf5
make
make install
cd ..

#dextract
git clone https://github.com/PacificBiosciences/DEXTRACTOR.git
cp MECAT/dextract_makefile DEXTRACTOR
cd DEXTRACTOR
export HDF5_INCLUDE=/public/users/chenying/smrt_asm/hdf5/include
export HDF5_LIB=/public/users/chenying/smrt_asm/hdf5/lib
make -f dextract_makefile
cd ..

export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/public/users/chenying/smrt_asm/hdf5/lib

本体 Github

git clone https://github.com/xiaochuanle/MECAT.git
cd MECAT
make -j 8
cd ..

#MECAT/Linux-amd64/binにパスを通す、dextractにもパスを通す

4つの代表的なコマンドがある。

  • mecat2pw, a fast and accurate pairwise mapping tool for SMRT reads

  • mecat2ref, a fast and accurate reference mapping tool for SMRT reads

  • mecat2cns, correct noisy reads based on their pairwise overlaps

  • mecat2canu, a modified and more efficient version of the Canu pipeline. Canu is a customized version of the Celera Assembler that designed for high-noise single-molecule sequencing

ここではdocker imageを使う。

https://hub.docker.com/r/robegan21/mecat/

docker pull robegan21/mecat

ホストからヘルプを表示。--rmでジョブ終了後コンテナを捨てる(ジョブを投げるたびに停止したコンテナが増殖してしまうので毎回消す)。

 > docker run --rm -it robegan21/mecat mecat2pw -h

$ docker run --rm -it robegan21/mecat mecat2pw -h

[parse_arguments, 142] unrecognised option 'h'

 

usage:

mecat2pw [-j task] [-d dataset] [-o output] [-w working dir] [-t threads] [-n candidates] [-g 0/1]

 

options:

-j <integer> job: 0 = seeding, 1 = align

default: 1

-d <string> reads file name

-o <string> output file name

-w <string> working folder name, will be created if not exist

-t <integer> number of cput threads

default: 1

-n <integer> number of candidates for gapped extension

Default: 100

-a <integer> minimum size of overlaps

Default: 2000 if x = 0, 500 if x = 1

-k <integer> minimum number of kmer match a matched block has

Default: 4 if x = 0, 2 if x = 1

-g <0/1> whether print gapped extension start point, 0 = no, 1 = yes

Default: 0

-x <0/x> sequencing technology: 0 = pacbio, 1 = nanopore

Default: 0

docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat mecat2ref -h

$ docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat mecat2ref -h

Error: argument to option 'h' is missing!

 

 

usage:

mecat2ref [-d reads] [-r reference] [-o output] [-w working dir] [-t threads]

 

options:

-d <string> reads file name

-r <string> reference file name

-o <string> output file name

-w <string> working folder name, will be created if not exist

-t <integer> number of cput threads

default: 1

-n <integer> number of of candidates for gap extension

default: 10

-b <integer> output the best b alignments

default: 10

-m <0/1/2> output format: 0 = ref, 1 = m4, 2 = sam

default: 0

-x <0/1> sequencing technology: 0 = pacbio, 1 = nanopore

default: 0

 

> docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat mecat2cns -h

$ docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat mecat2cns -h

input_type: 1

number of threads: 1

batch size: 100000

mapping ratio: 0.9

align size: 2000

cov: 6

min size: 5000

partition files: 10

tech: 0

USAGE:

mecat2cns [options] input reads output

 

OPTIONS:

-x <0/1> sequencing platform: 0 = PACBIO, 1 = NANOPORE

default: 0

-i <0/1> input type: 0 = candidte, 1 = m4

-t <Integer> number of threads (CPUs)

-p <Integer> batch size that the reads will be partitioned

-r <Real> minimum mapping ratio

-a <Integer> minimum overlap size

-c <Integer> minimum coverage under consideration

-l <Integer> minimum length of corrected sequence

-k <Integer> number of partition files when partitioning overlap results (if < 0, then it will be set to system limit value)

-h print usage info.

 

If 'x' is set to be '0' (pacbio), then the other options have the following default values: 

-i 1 -t 1 -p 100000 -r 0.9 -a 2000 -c 6 -l 5000 -k 10

 

If 'x' is set to be '1' (nanopore), then the other options have the following default values: 

-i 1 -t 1 -p 100000 -r 0.4 -a 400 -c 6 -l 2000 -k 10

docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat mecat2canu -h

$ docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat mecat2canu -h

 

usage: canu [-correct | -trim | -assemble] \

            [-s <assembly-specifications-file>] \

             -p <assembly-prefix> \

             -d <assembly-directory> \

             genomeSize=<number>[g|m|k] \

             errorRate=0.X \

            [other-options] \

            [-pacbio-raw | -pacbio-corrected | -nanopore-raw | -nanopore-corrected] *fastq

 

  By default, all three stages (correct, trim, assemble) are computed.

  To compute only a single stage, use:

    -correct  - generate corrected reads

    -trim     - generate trimmed reads

    -assemble - generate an assembly

 

  The assembly is computed in the (created) -d <assembly-directory>, with most

  files named using the -p <assembly-prefix>.

 

  The genome size is your best guess of the genome size of what is being assembled.

  It is used mostly to compute coverage in reads.  Fractional values are allowed: '4.7m'

  is the same as '4700k' and '4700000'

 

  The errorRate is not used correctly (we're working on it).  Set it to 0.06 and

  use the various utg*ErrorRate options.

 

  A full list of options can be printed with '-options'.  All options

  can be supplied in an optional sepc file.

 

  Reads can be either FASTA or FASTQ format, uncompressed, or compressed

  with gz, bz2 or xz.  Reads are specified by the technology they were

  generated with:

    -pacbio-raw         <files>

    -pacbio-corrected   <files>

    -nanopore-raw       <files>

    -nanopore-corrected <files>

 

Complete documentation at http://canu.readthedocs.org/en/latest/

 

ERROR:  Invalid command line option '-h'.  Did you forget quotes around options with spaces?

ERROR:  Assembly name prefix not supplied with -p.

ERROR:  Directory not supplied with -d.

ERROR:  Required parameter 'genomeSize' is not set

 

 

 

 

ラン

テスト1、Pacbio

1、 テストシーケンスデータをダウンロード。

cd /home/kazu/MECAT/
wget http://gembox.cbcb.umd.edu/mhap/raw/ecoli_filtered.fastq.gz
gzip -dv ecoli_filtered.fastq.gz

ecoli_filtered.fastqができる。 

 

dockerを使い、ホスト側からジョブを投げる。見にくいのでMECATのコマンド部分手前でエスケープして改行し、さらにMECATコマンド部は太字にする。

2、mecat2pw: a fast and accurate pairwise mapping tool for SMRT reads

2、mappingによるオーバーラップ検出。ecoli_filtered.fastqを指定する。

#現在のパス /home/kazu/MECAT/
sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
mecat2pw -j 0 -d ecoli_filtered.fastq -o ecoli_filtered.fastq.pm.can -w wrk_dir -t 16

ecoli_filtered.fastq.pm.canが出力される。

 

3、mecat2cns: correct noisy reads based on their pairwise overlaps

エラー訂正。ecoli_filtered.fastq.pm.canとecoli_filtered.fastqを指定する。

sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
mecat2cns -i 0 -t 16 ecoli_filtered.fastq.pm.can ecoli_filtered.fastq corrected_ecoli_filtered.fasta

 corrected_ecoli_filtered.fastaが出力される。

 

4、extract the longest 25X corrected reads。推定ゲノムサイズとカバレッジを指定する(リンク)。

#extract_sequences usage 
extract_sequences inputReads outputReads-prefix genomeSize coverage

sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
extract_sequences corrected_ecoli_filtered.fasta corrected_ecoli_25x 4800000 25

corrected_ecoli_25x.fastaが 出力される。

 

5、mecat2canu: a modified and more efficient version of the Canu pipeline. Canu is a customized version of the Celera Assembler that designed for high-noise single-molecule sequencing

アセンブリ

sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
mecat2canu -trim-assemble -p ecoli -d ecoli genomeSize=4800000 ErrorRate=0.02 maxMemory=40 maxThreads=16 useGrid=0 Overlapper=mecat2asmpw -pacbio-corrected corrected_ecoli_25x.fasta

 

 

スト2、Nanopore

1、ダウンロード。

wget http://nanopore.s3.climb.ac.uk/MAP006-PCR-1_2D_pass.fasta

2、mappingによるオーバーラップ検出。

#現在のパス /home/kazu/MECAT/
sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
mecat2pw -j 0 -d MAP006-PCR-1_2D_pass.fasta -o candidatex.txt -w wrk_dir -t 16 -x 1

3、エラー訂正。

sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
mecat2cns -i 0 -t 16 -x 1 candidates.txt MAP006-PCR-1_2D_pass.fasta corrected_ecoli.fasta

4、extract the longest 25X corrected reads

sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
extract_sequences corrected_ecoli.fasta corrected_ecoli_25x.fasta 4800000 25

5、アセンブリ

sudo docker run --rm -it -v /home/kazu/MECAT/:/root robegan21/mecat \
mecat2canu -trim-assemble -p ecoli -d ecoli genomeSize=4800000 \ ErrorRate=0.06 maxMemory=40 maxThreads=16 useGrid=0 \
Overlapper=mecat2asmpw -nanopore-corrected corrected_ecoli_25x.fasta

 

 

マッピングにはmecat2refを使う。

 

 

引用

MECAT: fast mapping, error correction, and de novo assembly for single-molecule sequencing reads

Xiao CL, Chen Y, Xie SQ, Chen KN, Wang Y, Han Y, Luo F, Xie Z.

Nat Methods. 2017 Nov;14(11):1072-1074