使用 GCE 进行基因组大小评估

最近在尝试做基因组大小的评估。尝试了几款软件。今天介绍一下GCE的使用。该软件GCE(Genome Characteristics Estimation) 是华大基因用于基因组评估的软件。最早的版本（gce-1.0.0）发表于2012年, 其参考文献为：Estimation of genomic characteristics by analyzing k-mer frequency in de novo genome projects。时隔8年，终于更新啦！2020年更新版本为gce-1.0.2，软件下载地址 ftp://ftp.genomics.org.cn/pub/gce。

以前GCE 软件包中主要包含kmer_freq_hash 和 gce 两支程序。前者用于进行 kmer的频数统计，后者在前者的结果上进行基因组大小的准确估算。最近更新的版本（gce1.0.2）主要包含kmerfreq和gce两支程序。程序使用参数有所变动。

1.GCE 下载和安装

wget

ftp://ftp.genomics.org.cn/pub/gce/gce-1.0.2

tar -xzvf gce.tar.gz

cd gce-1.0.2

make

出现make: Nothing to be done for 'all'，上网查了一下，应该是已经编译好了，可以直接使用。.

温馨提示：记得添加环境变量

(1) gce

gce -h 可以看到其用法，如下图：

Usage:  gce(genomiccharactor estimator) [option]

Version: 1.0.2

Author:   BGIShenZhen

-f      depth frequency file with two columns: depthvalue and kmer species number #深度频率文件共两列：kmer重复次数和kmer种类数

-c      expected depth for unique kmer, which can be obtained by checking the data with human eyes # 唯一kmer的期望深度

-g      total kmer number, i.e. total number of kmerindividuals #全部kmer数量

-b     have bias(1) or not(0), default=0

-H     use hybrid mode(1) or not(0), default=0 #使用杂合模式（1），不使用（0），默认不使用。

-m    estimation mode: discrete mode(0) andcontinuous mode(1), default=0 #估算模型：离散型（0），连续型（1），默认离散型。

-M    max depth value, information for larger depthwill be ignored, default=1500 #最大深度值，默认1500，超过此数值的将被忽略，

-D     precision of expect value, default=1

-d     difference cut off, default=0.0001

-i      iterate cycle number cut off, default=10000

,-h  this help

Example:

(1) Before run gce, firstly get the total kmer number anddepth frequency file from the kmerfreq result file (example: AF.kmer.freq.stat)

     lessAF.kmer.freq.stat | grep "#Kmer indivdual number"

     less AF.kmer.freq.stat | perl-ne 'next if(/^#/ || /^\s/); print; ' | awk '{print $1"\t"$2}' >AF.kmer.freq.stat.2colum

(2) Run gce in homozygous mode, suitable for homozygousand near-homozygous genome (-g and -f must be set at the same time)

     gce-g 173854609857 -f AF.kmer.freq.stat.2colum >gce.table 2>gce.log

(3) Run gce in heterzygous mode, siutable forheterozgyous genome (-H and -c must be set at the same time)

     gce -g 173854609857 -fAF.kmer.freq.stat.2colum -c 75 -H 1 >gce2.table 2>gce2.log

 (2) kmerfreq 

kmerfreq  [options]

 Version 4.0

   -k  kmer size, recommand value13 to 19, default=17

#设置 kmer 的大小。推荐该值为 13~19，默认值为17

   -f  input file format: 1:fq|gz(one-line), 2: fa|gz(one-line), default=1 

#输入文件格式 1 fq|gz(one-line), 2: fa|gz(one-line),默认为1

   -p  output file prefix,default=reads_files.lib

#输出文件前缀，默认reads_files.lib

   -r  number of reads stored inbuffer memory, default=10000

   -t  thread number to use inparallel, default=10

#线程数 默认10

   -w  whether output kmer sequenceand frequency value, , 1:yes, 0:no, default=0

#是否输出kmer序列和频率值，1:yes, 0:no,默认不输出。

   -c  kmer frequency cutoff, equalor larger will be output, co-used with -w, default=5

#kmer频率终止值，等于或者大于默认值的会被输出。与-w一起使用，默认值为5。

   -m  whether output computermemory data, 1:yes, 0:no, default=0

#是否输出计算机内存数据，1:yes, 0:no,默认不输出。

   -q  kmer frequency cutoff, 0 forlower, 1 for equal and larger, co-used with -m, default=5

#kmer频率截止值，小于默认值为0，等于或者大于默认值为1，默认值为5.

   -h        get help information

#获取帮助信息

Example: kmerfreq reads_files.lib

        kmerfreq  -k 17 -t 10 -p Ecoli_K17reads_files.lib

        kmerfreq  -k 17 -t 10 -p Ecoli_K17-w 1 -c 5 reads_files.lib

        kmerfreq  -k 17 -t 10 -p Ecoli_K17-m 1 -q 5 reads_files.lib

GCE命令行：

kmerfreq  -k 17 -t 10 -p Ecoli_K17reads_files.lib

less AF.kmer.freq.stat | grep "#Kmer indivdual number"

less AF.kmer.freq.stat | perl-ne 'next if(/^#/ || /^\s/); print; ' | awk '{print $1"\t"$2}' >AF.kmer.freq.stat.2colum

gce -g 173854609857 -f AF.kmer.freq.stat.2colum >gce.table 2>gce.log   or

gce -g 173854609857 -f AF.kmer.freq.stat.2colum -c 75 -H 1 >gce2.table2>gce2.log

2.结果展示（以kmer=17为例）：

这里的基因组大小计算为： genome size=effective_kmer_individuals/coverage_depth=460468198.15143bp，即约460Mb。

参考：陈连福的生信博客：http://www.chenlianfu.com/?p=2335