macでインフォマティクス

macでインフォマティクス

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

異なるk-merの割合を計算してエラー率推定やゲノムサイズ推定を行う KmerStream

2019 5/27インストール方法追記

2019 7/29タイトル修正

 

KmerStreamは異なるk-merの数を計算する方法論。シーケンス業者のクオリティに依存せず純粋にk-merの頻度からエラー率を見積もることができるため、上手く使えばシーケンスの品質管理などに使用することができる。サンプリングを行うためメモリ使用量が少ないのも特徴となっている。

 

公式HP 

https://bioweb.pasteur.fr/packages/pack@KmerStream@1.1


 

インストール

Github


gitからソースをダウンロードしてmakeする。 brewでもインストールできる。

git clone https://github.com/pmelsted/KmerStream.git
cd
KmerStream/
make

#homebrew
brew
install Kmerstream

> ./KmerStream

$ ./KmerStream

Need to specify files for input

KmerStream 1.1

 

Estimates occurrences of k-mers in fastq or fasta files and saves results

 

Usage: KmerStream [options] ... FASTQ files

 

-k, --kmer-size=INT      Size of k-mers, either a single value or comma separated list

-q, --quality-cutoff=INT Comma separated list, keep k-mers with bases above quality threshold in PHRED (default 0)

-o, --output=STRING      Filename for output

-e, --error-rate=FLOAT   Error rate guaranteed (default value 0.01)

-t, --threads=INT        SNumber of threads to use (default value 1)

-s, --seed=INT           Seed value for the randomness (default value 0, use time based randomness)

-b, --bam                Input is in BAM format (default false)

    --binary             Output is written in binary format (default false)

    --tsv                Output is written in TSV format (default false)

    --verbose            Print lots of messages during run

    --online             Prints out estimates every 100K reads

    --q64                set if PHRED+64 scores are used (@...h) default used PHRED+33

 

 

実行例 

付属のスクリプトを使いゲノムサイズを推定。 phread quality score≥30以上のk-merだけ使い計算した。

KmerStream -k 33,55,77,91,111,127 -t 18 -q 30 --tsv -o output.tsv
KmerStreamEstimate.py output.tsv #出力されたtsvファイルを指定

出力

Vendor:  continuum

Product: anaconda

Message: trial mode expires in 30 days

Q k F0 f1 F1 F0-f1 G e lambda

33 33 4281284 387290 76765368 3893994.0 3893799.43151 0.00505018692713 19.7147719985

33 55 4468082 546551 64067064 3921531.0 3921306.1254 0.00853804385099 16.3381949665

33 77 4564542 651000 53357825 3913542.0 3913387.55198 0.0122080429515 13.6346897135

33 91 4604518 690472 47100340 3914046.0 3914179.73801 0.0146620835559 12.0332593679

33 111 4630474 720842 38730707 3909632.0 3911960.09364 0.0185627390272 9.90058847045

33 127 4587192 717597 32455542 3869595.0 3879373.50588 0.0218474819721 8.36618127922

F0、F1、の他に、k=33,55,77,91,111,127merrから推定されたゲノムサイズがプリントされた。ここではバクテリアのシーケンスデータ(fastq)を使用した。リピート遺伝子を100程度持つ3.95Mのゲノムになる。Gのカラムが推定されたゲノムサイズだが、だいたいあってる。

 

F0、F1, F2について(論文より)

F0 , F1 and F2 . F0 is the number of distinct k -mers in the set of reads and F1 is simply the number of k -mers in the reads, counted with repetition. The second moment, F2 , puts a higher weight on the number of k -mers that have high abundance values.

 

引用

KmerStream: streaming algorithms for k-mer abundance estimation.

Melsted P, Halldórsson BV.

Bioinformatics. 2014 Dec 15;30(24):3541-7.

 

参考資料

http://tech.plaid.co.jp/streaming-algorithm/