macでインフォマティクス

macでインフォマティクス

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

k-mer カウントして、配列も出力するツール jellyfish、BFCounter

K-merカウントを行うjellyfishと、k-merの全配列を書き出すBFCounterを紹介する。

 

2017.11追記

 

Jellyfish

公式サイト

JELLYFISH - Fast, Parallel k-mer Counting for DNA

Github

https://github.com/gmarcais/Jellyfish

 

ビルド

 ./configure 
make sudo
make install

またはbrewでも導入できる。

 

ラン

k-merカウントする。

jellyfish count -m 33 -s 100M -t 20 input.fasta/fastq -o mer_counts.jf
  • -m Length of mer
  • -s Initial hash size
  • -t Number of threads (1)

インプットファイルはfasta ,fastqとも対応している。

 mer_counts.jfができる。中身はbinaryである。このファイルを以後の解析に使う。

 

 

k-merの配列を全て出力する。

jellyfish dump mer_counts.jf > k33.fasta

中身を見てみる。

user-no-MacBook-Pro: user$ less -S k33.fa

>1

CCCAGGGTCTGGAACTATCGCCGGTGGAAAATATCGAAACGGAAGCGGGGCAAGGGGTTCAGGGTTGGTACCAAGGC

>15

AATGCCCTGGGGAGTAACTGAAAGGAAACGAAAAAATAATCTTCTTCCTATTCTATTCCTTCTGATTGAGTAACCAG

>1

CCTTGTTTGACAGTGCAGTCCGGGTTCCTAAAAGGTTGCAAGGGGATTTATTCATTCGCAGTACCCAAACGGACACC

>1

GGTGGAATAAATACCTGTGGTGAAGGTATTAACGCCAAAATATTGCACTGTGCCAAAATCGTTCAGAGTCCCCATCA

>1

AGAACGTCGTGGACGGCGTCCCCGTCGTGAGGGAACGGTGGAGGTCAAAGGGCGATCGGGAAGGGCGGCGGAACTAT

>16

TGGTGGACGGAAATCTTAAACATTGCCTAGAGAAAAGTTGGGAATTGACGGTGGTGGGGGCAGGGCAAGGAGTGGAA

>15

GGTTTCCTGGGCAAACTCGATGGAACCTGGACAATCGAGGAAATTTAAGCGCAAATTTTCATAATCAATCCCTGCCA

 

指定のk-merサイズの全塩基配列が出力されている。数値は登場回数である。重複させながら元の配列を分解しているので、元のfastqよりサイズははるかに大きくなる。

 

 

 

 

ヒストグラムデータを出力。

jellyfish histo mer_counts.jf > histogram.txt

出力の中身を見てみる。

user-no-MacBook-Pro: user$ less -S histogram.txt

1 3101726

2 8283

3 5424

4 12735

5 29386

6 61451

7 103516

8 164603

9 230629

10 291824

11 347466

Rに読み込んでグラフ化する。

 

純化されたバクテリアのシーケンスデータを分析する。

> R #Rに入る
input <- read.table("histogram.txt") #読み込み
plot(input[2:100,],type="l") #2列目が頻度情報。

f:id:kazumaxneo:20171120142312j:plain

 2つ目の山がゲノムのシングルコピーの領域に由来する33merのピークと思われる。ピーク座標を調べる。

input[0:20,]

13がピークと思われる(赤色)。

> input[0:20,]

   V1        V2

1   1 136238639

2   2   3362948

3   3    446016

4   4    142167

5   5    122950

6   6    159056

7   7    216472

8   8    287219

9   9    362505

10 10    430060

11 11    487906

12 12    520707

13 13    531598

14 14    531210

15 15    504837

16 16    472747

17 17    435696

18 18    390053

19 19    342325

20 20    298075

 

もっとちゃんと調べる。

 k15

f:id:kazumaxneo:20171120154345j:plain

k35

f:id:kazumaxneo:20171120154503j:plain

k55

f:id:kazumaxneo:20171120154620j:plain

k65

f:id:kazumaxneo:20171120154637j:plain

k77

f:id:kazumaxneo:20171120154650j:plain

k99

f:id:kazumaxneo:20171120155119j:plain

k127

f:id:kazumaxneo:20171120155356j:plain

 

 

kmergenuineでもbestなk-merを探す(参考にしたページ)。(kmergenuineはbrewで導入できる。)

kmergenie input.fq 

fitting model to histograms to estimate best k

table of predicted num. of genomic k-mers: histograms.dat

recommended coverage cut-off for best k: 2

best k: 65

kmergenuineでは65-merがベストと判定された。

 

 

65-merの入力テーブルを使い、ゲノムサイズを推定する。

sum(as.numeric(input[12:1130,1]*input[12:1130,2]))/65

リピートがあると、その分だけ少ない値が出るはずである。

 

 

k-mer coverageの貧弱なものを除きたいならkhmerを使う。khmerはk-merグラフの分布から低頻度k-mer配列を除去したり様々な機能をもつde novo assemblyをアシストするツールである。(pipでインストール可能

load-into-counting.py -T 20 -k 20 -x 3e6 output1 input.fastq
  • T スレッド数
  • k k-merサイズ
  • x 大まかなゲノムサイズ。
filter-abund.py -C 2 -T 20 output1 input.fastq

注 khmrは別途記事にします。

 

 

BFCounter

BFCounterは k-merを数えて、全部の配列を出力するツール。

 

Github

https://github.com/pmelsted/BFCounter

ビルドはソースコードをダウンロードしてmakeするだけである。パスを通しておく。

 

ラン

BFCounter count -n 1000 -k 127 -o single_out -t 12 input.fq
  • -n Estimated number of k-mers
  • -k Size of k-mers, at most 301
  • -t Number of threads to use (default 1)

出力はbinaryファイルである。ここからk-mer長の配列を抽出する。

BFCounter dump -i single_out -k 127 -o k-127_output
  • -i Filename for input file from count
  • -k Size of k-mers, at most 301
  • -o Filename for output

出力の中身を見てみる。

f:id:kazumaxneo:20170627225750j:plain

 このように全k-mer配列が出力されている。あとはスクリプトを使い自由にいじることができる。

このツールの機能はこれだけであるが、シンプルでなため使いやすい。

 

 

 

 引用

Guillaume Marcais and Carl Kingsford, A fast, lock-free approach for efficient parallel counting of occurrences of k-mers.

Bioinformatics (2011) 27(6): 764-770 (first published online January 7, 2011) doi:10.1093/bioinformatics/btr011)

 

Efficient counting of k-mers in DNA sequences using a bloom filter.

BMC Bioinformatics. 2011 Aug 10;12:333. doi: 10.1186/1471-2105-12-333.

Melsted P1, Pritchard JK.

 

K-mer analysis and genome size estimate