macでインフォマティクス

macでインフォマティクス

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

CCSリードの精度を推定する Yak

Githubより

Yakは当初、2つの特定のユースケースのために開発された。1) CCSリードとアセンブリコンティグの塩基精度をロバストに推定すること、2) CCSリードの系統的なエラー率を調査することである。ショートリードのk-merスペクトラムと配列を比較したり、スペクトラムを比較したりすることで目標を達成する。リファレンスゲノムgrand truthデータは必要ない。

ベース精度の推定がトリッキーであることは注目に値する。精度がQ50に近づくと、ショートリード中のサンプルされていないk-merと誤ったk-merの両方がナイーブな推定器に干渉する可能性がある。Yakはこの問題に対処するために経験的モデルを導入した。その推定値は、カバレッジとショートリードの質の影響をあまり受けない。

 

インストール

macos10.14でビルドした。

Github

git clone https://github.com/lh3/yak
cd yak && make

#conda、ここでは高速なmambaを使う
mamba install -c bioconda yak -y

yak

$ yak

Usage: yak <command> <argument>

Command:

  count     count k-mers

  qv        evaluate quality values

  triobin   trio binning

  trioeval  evaluate phasing accuracy with trio

  inspect   k-mer hash tables

  version   print version number

yak count

$ yak count

Usage: yak count [options] <in.fa> [in.fa]

Options:

  -k INT     k-mer size [31]

  -p INT     prefix length [10]

  -b INT     set Bloom filter size to 2**INT bits; 0 to disable [0]

  -H INT     use INT hash functions for Bloom filter [4]

  -t INT     number of worker threads [4]

  -o FILE    dump the count hash table to FILE []

  -K INT     chunk size [100m]

Note: -b37 is recommended for human reads

yak qv

$ yak qv

Usage: yak qv [options] <kmer.hash> <seq.fa>

Options:

  -l NUM      min sequence length [0]

  -f FLOAT    min k-mer fraction [0.5]

  -e FLOAT    false positive rate [4e-05]

  -p          print QV for each sequence

  -t INT      number of threads [4]

  -K NUM      batch size [1g]

yak triobin -h

$ yak triobin -h

Usage: yak triobin [options] <pat.yak> <mat.yak> <seq.fa>

Options:

  -c INT     min occurrence [2]

  -d INT     mid occurrence [5]

  -t INT     number of threads [8]

yak trioeval -h

$ yak trioeval -h

Usage: yak trioeval [options] <pat.yak> <mat.yak> <seq.fa>

Options:

  -c INT     min occurrence [2]

  -d INT     mid occurrence [5]

  -n INT     min streak [2]

  -t INT     number of threads [8]

  -e         print error positions (out of order)

yak inspect -h

$ yak inspect -h

Usage: yak inspect [options] <in1.yak> [in2.yak]

Options:

  -m INT    max count (effective with in2.yak) [20]

Notes: when in2.yak is present, inspect evaluates the k-mer QV of in1.yak and

  the k-mer sensitivity of in2.yak.

yak version

$ yak version

0.1-r58-dirty

 

 

実行方法

1、k-merハッシュテーブルの構築

#assembly
yak count -K 1.5g -t 16 -o assembly.yak assembly.fa.gz

#single fastq(githubではCCS readsを例にしている,シングルトンk-merは除外)
yak count -b 37 -t 16 -o reads.yak reads.fq.gz

#paired-end fastq (シングルトンk-merは除外)
yak count -b37 -t32 -o sr.yak <(zcat sr*.fq.gz) <(zcat sr*.fq.gz)
  • -K   chunk size [100m]
  • -k   k-mer size [31]
  • -t    number of worker threads [4]
  • -b   set Bloom filter size to 2**INT bits; 0 to disable [0]
  • -o   dump the count hash table to FILE []

 

2、精度推定

yak qv -t32 -p -K3.2g -l100k sr.yak asm.fa.gz > asm-sr.qv.txt

 

3、k-merヒストグラムやk-merの保存

#k-mer histgram
yak inspect sr.yak > sr.hist

#k-mer配列プリント
yak inspect -p sr.yak > sr.kmers

 

引用

https://github.com/lh3/yak

 

関連