Sangshin | 基因组组装行动（三）：Kmer 评估基因组

写在前面

• 以下内容均来自我在菲沙基因（Frasergen）暑期生信培训班上记录的课堂笔记

1.Kmer定义

• 杂合Kmer
  

  杂合Kmer

• 杂合Kmer峰
  

  杂合Kmer峰

• 重复Kmer峰
  

  重复Kmer峰

• 杂合变化模拟Kmer图
  

  杂合变化模拟Kmer图

2.基因组大小预估方式

• 也就是基因组大小=（基因切割的Kmer数目）/（主峰深度）
• 杂合和重复会影响统计，而一般的基于Kmer预估基因组的软件会对此做处理,让结果更贴近真实。

3.Kmer分析实战

• 软件：GCE

3.1 下载安装gce软件

wget ftp://ftp.genomics.org.cn/pub/gce/gce-1.0.2.tar.gz
tar -zxvf gce-1.0.2.tar.gz

3.2 使用gce软件中的kmerfreq脚本切割kmer并统计频率深度表格

gce-1.0.2/kmerfreq -k 17 -t 10 -p freq list_of_clean
#-k 17：切割kmer的长度
#list_of_clean是质控后的文件名
#cat list_of_clean 
#/local_data1/pop_clean_1P.fastq.gz
#/local_data1/pop_clean_2P.fastq.gz

less freq.kmer.freq.stat|perl -ne 'next if(/^#/ || /^\s/); print; ' | awk '{print $1"\t"$2}' > freq.stat.2colum
#total kmer number, i.e. total number of kmer individuals
c=`less freq.kmer.freq.stat| grep "#Kmer indivdual number"|awk '{print $4}'`
echo $c

3.3 将上一行的$c也就是总的kmer种类数目和freq.stat.2colum也就是频率深度表格输入，得到基因组大小，杂合，重复信息 ，-H和-c的选择很重要

#纯合模式
gce -g $c -f freq.stat.2colum 2>gce.log

#杂合模式
gce -H 1 -g $c -c 60 -f freq.stat.2colum 2>heterozgyousgce.log

-g：kmer总数, 从kmerfreq分析结果获取
-f freq.stat.2colum：Kmer频率分布，从kmerfreq分析结果获取
-H 1：是否启动杂合模式(1是杂合模式，推算出杂合率, 0是非杂合模式没有杂合度)
-c 60： Kmer主峰深度，由gce自己选，或者根据情况自己选择(峰的选择很重要)

• GCE输出结果说明

3.4 Kmer作图

c=`awk '$1==60' freq.stat.2colum|awk '{print $2}'`
echo $c
#选取合理的深度范围
head -n 500 freq.stat.2colum > freq.stat.2colum.500
#作图
Rscript distribution.r freq.stat.2colum.500 ./ $c
convert kmer_distribution.svg kmer_distribution.png
sz kmer_distribution.png

• R作图脚本（上面的distribution.r）

library(ggplot2)
#1. data # 读入 深度-Kmer种类数频率 表格
args <- commandArgs()
file=args[6] 
a<-read.table(file,sep="\t")
#2. output
setwd(args[7])
#3. ylim 峰值大小，就是Kmer的种类数峰值大小，作为y的max值
peak=args[8]
peak<-as.numeric(peak)
#4. plot 作图，
svg("kmer_distribution.svg", width=10) 
ggplot(a,aes(x=V1,y=V2),col="red")+geom_line(color="green")+geom_point(color="red")+xlim(0,200)+ylim(0,peak)+xlab("
Depth of Kmer Species")+ylab("Frequency of Kmer Species")+theme_bw()+theme(axis.title=element_text(size=20))
dev.off()

• 结果

总结

• Survey分析内容回顾
  

  Survey分析内容回顾

• Tips
  

  Tips