直接下载编译好的linux可执行文件就可以了,代码如下
wgethttps://github.com/gmarcais/Jellyfish/releases/download/v1.1.12/jellyfish-linuxmv jellyfish-linux jellyfishchmod +x jellyfish
jellyfish有很多个子命令,下面介绍常用的子命令
1. 计算kmer分布
count子命令用来计算kmer分布,用法如下
jellyfishcount-m31-t10-s1Gtest.fq
-m参数指定kmer的长度,-t指定并行的线程数,-s指定内存中hash的大小,这个参数可以根据基因组的大小适当调整,比如人类基因组3G,这里就设置成3G;test.fq是输入的序列文件。
需要注意的是,jellysifh 只能接受fastq或者fasta文件作为输入,而且不能是压缩文件。
默认情况下会生成名为mer_counts.jf的文件,该文件是一个二进制文件,可以通过其他命令来查看该文件中的内容。
2. 生成kmer 统计表
dump子命令可以将上述步骤生成的二进制文件转换成纯文本的文件,用法如下
jellyfishdumpmer_counts.jf>kmer_count.fasta
kmer_count.fasta文件中每一条序列就是一个kmer,序列标识符是该kmer出现的次数,示意如下
>1150
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
>20
GCTACCATGATAGCCAAGGAAATCCCACAAA
也支持输出成表格的形式,只需要添加-c和-t两个参数,用法如下
jellyfishdump-c-tmer_counts.jf>kmer_count.tsv
输出内容如下
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA 1150
GCTACCATGATAGCCAAGGAAATCCCACAAA 20
TATTTGCGGAGCTTGTTAGATAACAATCAAA 18
第一列为kmer,第二列为该kmer频数。
3. 查询某个kmer出现的次数
query子命令可以快速查询某个kmer的频数,用法如下
jellyfishquerymer_counts.jfGCTACCATGATAGCCAAGGAAATCCCACAAA
输出结果:
GCTACCATGATAGCCAAGGAAATCCCACAAA 20
4. 统计kmer基本信息
stats子命令会给出kmer的基本统计信息,用法如下
jellyfishstatsmer_counts.jf
会在命令行输出如下统计结果
Unique: 130512636
Distinct: 260032104
Total: 3502352808
Max_count: 43395
unique 代表只出现了1次的kmer个数;Distinct 代表kmer的个数,Total 代表所有kmer频数总和,Max count 代表kmer的最大频数。
5. 统计kmer频数分布
histo子命令可以给出kmer的频数分布,用法如下
jellyfishhistomer_counts.jf>kmer_hist.tsv
输出内容如下
1 130512636
2 24879103
3 12766220
4 8746042
第一列代表kmer的频数,第二列代表出现该频数的kmer的总数。利用这个数据,可以画出kmer频数分布曲线,对应的R语言代码如下
x <- read.table(input, header = F, sep =" ", stringsAsFactors = F)pdf("kmer_hist.pdf")plot(x,type="l")dev.off()
最终的效果图如下
网友评论