macでインフォマティクス

macでインフォマティクス

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

HiFiロングリードアセンブリのためのリピートを認識したポリッシングツール NextPolish2

Preprint Pacbio polish 2023 Trio binning assembly misassembly

 

　PacBio社が開発した高忠実度（HiFi）ロングリードシーケンス技術により、ゲノムアセンブリの塩基レベルの精度は大幅に向上したが、これらのアセンブリには、特にHiFiロングリードのエラーが発生しやすい領域内に、塩基レベルのエラーが残っている。しかし、既存のゲノム研磨ツールは、HiFiロングリードからアセンブルされたゲノムのエラーを修正する際に、通常、オーバーコレクションやハプロタイプスイッチエラーが発生する。NextPolish2は、HiFiロングリードからアセンブルされた高精度ゲノムに残存する塩基エラーを、過剰補正やハプロタイプスイッチエラーを発生させずに修正することができる。NextPolish2は、テロメア間（T2T）ゲノムの精度をさらに向上させるために大きな意義があると考えられる。NextPolish2は、https://github.com/Nextomics/NextPolish2で自由に利用できる。

ここでは、主にHiFiロングリードを用いて構築されたT2Tゲノムのエラー訂正のために、バージョンアップしたゲノム研磨ツールNextPolish2を紹介する。NextPolish2は、CHM13のヒトT2Tゲノムアセンブリの研磨に採用された最新の研磨パイプライン（Racon + Merfin [10]、以下RM）と比較して、「高精度」ドラフトアセンブリにおいて、繰り返しエレメントの多い領域でも過剰な修正を起こさずに塩基エラーを修正できる。また、内蔵のフェージングモジュールにより、エラー塩基を修正するだけでなく、元のハプロタイプの一貫性を維持することができる。実際、著者らの評価では、ヘテロ接合領域でのスイッチエラーをわずかに減少させることができた。

 

インストール

ubuntu18LTSでビルドした（rustc 1.69.0）。

  • NextPolish2 is written in rust
#Rustがない場合導入する
curl --proto '=https' --tlsv1.2 -sSf https://sh.rustup.rs | sh

#本体
git clone --recursive git@github.com:Nextomics/NextPolish2.git
cd NextPolish2 && cargo build --release
cd targt/release/

> ./nextPolish2 

Repeat-aware polishing genomes assembled using HiFi long reads

 

Usage: nextPolish2 [OPTIONS] <HiFi.map.bam> <genome.fa[.gz]> <short.read.yak>...

 

Arguments:

  <HiFi.map.bam>       HiFi-to-ref mapping file in sorted BAM format.

  <genome.fa[.gz]>     genome assembly file in [GZIP] FASTA format.

  <short.read.yak>...  one or more k-mer dataset in yak format.

 

Options:

  -o, --out <FILE>               output file. [default: stdout]

  -u, --uppercase                output in uppercase sequences.

  -k, --min_kmer_count <INT>     filter kmers in k-mer dataset with count <= INT. [default: 1]

  -t, --thread <INT>             number of threads. [default: 1]

  -i, --iter_count <INT>         number of iterations to attempt phasing. [default: 2]

  -m, --model <STR>              phasing model. [default: ref] [possible values: ref, len]

  -l, --min_read_len <INT>       filter reads with length <= INT. [default: 1000]

  -s, --use_supplementary        use supplementary alignments.

  -S, --use_secondary            use secondary alignments, consider setting `min_map_qual` to -1 when using this option.

  -a, --min_map_len <INT.FLOAT>  filter alignments with alignment length <= min(INT, FLOAT * read_length). [default: 500.5]

  -q, --min_map_qual <INT>       filter alignments with mapping quality <= INT. [default: 1]

  -c, --max_clip_len <INT>       filter alignments with unaligned length >= INT. [default: 100]

  -r, --use_all_reads            use all unfiltered reads, reads with different haplotypes from the reference assembly are discarded by default.

  -h, --help                     Print help (see more with '--help')

  -V, --version                  Print version

 

 

 

実行方法

１、HiFiリードをゲノムにマッピングする。Winnowmap（紹介）を使う場合、まずmeryl（紹介）を使ってリファレンスゲノムのk-merカウントを行う必要がある。ここではminimap2を使う。

#winnowmap
meryl count k=15 output merylDB asm.fa.gz
meryl print greater-than distinct=0.9998 merylDB > repetitive_k15.txt
winnowmap -t 5 -W repetitive_k15.txt -ax map-pb asm.fa.gz hifi.fasta.gz|samtools sort -o hifi_sort.bam -
samtools index -@ 8 hifi_sort.bam

#minimap2
minimap2 -ax map-hifi -t 20 asm.fa.gz hifi.fasta.gz |samtools sort -@ 8 -o hifi_sort.bam -
samtools index -@ 8 hifi_sort.bam

 

２、yak（紹介）でショートリードのk-merカウント。シングルトンのみカウントするなら-b 37を外す（-b 37はヒトゲノムで推奨されている）

#21-mer
yak count -o k21.yak -k 21 -b 37 <(zcat sr.R*.fastq.gz) <(zcat sr.R*.fastq.gz)

#31-mer
yak count -o k31.yak -k 31 -b 37 <(zcat sr.R*.fastq.gz) <(zcat sr.R*.fastq.gz) 
  • -k    k-mer size [31]
  • -p    prefix length [10]
  • -b    set Bloom filter size to 2**INT bits; 0 to disable [0]
  • -H   use INT hash functions for Bloom filter [4]
  • -t     number of worker threads [4]
  • -o    dump the count hash table to FILE []
  • -K    chunk size [100m]

 

３、nextpolish2でポリッシングする。

nextPolish2 -t 20 hifi_sort.bam asm.fa.gz k21.yak k31.yak > asm.np2.fa
  • -t    number of threads. [default: 1]
  • -i     number of iterations to attempt phasing. [default: 2]

asm.np2.faが出力される。

 

レポジトリより

  • NextPolish2が修正できるのは、HiFiリードがマッピングされた領域のみ。HiFiリードがマッピングされていない領域（通常、エラー率が高い）については、マッピングパラメータの調整で改善できる可能性がある。
  • NextPolish2の性能はショートリードの品質に大きく依存する。
  • NextPolish2が修正できるのは一部の構造ミスアセンブリのみ。
  • RaconやNextPolish1などのロングノイズリード用に設計されたゲノムポリッシングツールは、修正した内容よりも多くのエラーを導入することが判明したため（論文表1）、HiFiロングリードを用いてアセンブルしたゲノムのエラー修正には推奨できない（論文より）。

 


Trio binningで得られたハプロタイプゲノムアセンブリなら、異なるハプロタイプを読んだリードをポリッシング前に除外することで精度を上げることができるとされる。

手順はこちらで説明されています。

https://github.com/Nextomics/NextPolish2/blob/main/doc/benchmark3.md

 

引用

NextPolish2:a repeat-aware polishing tool for genomes assembled using HiFi long reads
Jiang Hu, Zhuo Wang, Fan Liang, Shanlin Liu,  Kai Ye, De-Peng Wang

bioRxiv, Posted April 27, 2023 

 

関連