基因組survey——K-mer頻譜

Kmer是從測序數(shù)據(jù)中滑窗提取出的長度為k的寡聚核苷酸序列淡喜,可以評估基因組大小秕磷、雜合度、重復序列比例等炼团。
在測序reads均勻分布的前提下澎嚣,有以下公式:
基因組長度 = 總堿基數(shù) / 平均測序深度 = 總kmer數(shù) / 平均kmer深度
根據(jù)該公式,可以使用總kmer數(shù)和平均kmer深度估算基因組長度瘟芝。K值使用滿足以下公式計算的最大奇數(shù):4 ^ K / genome > 200

做kmer有幾種軟件:SOAPdenovo2的Kmerfreq模塊易桃,jellyfish和KAT(K-mer Analysis Toolkit)

jellyfish:簡書1簡書2GitHub
產(chǎn)生的kmer頻數(shù)分布數(shù)據(jù)需要用R包GenomeScopefindGSE來進行統(tǒng)計估計锌俱。(http://qb.cshl.edu/genomescope/)

這里我還是用一種后來出現(xiàn)但整合了jellyfish和其他分析的軟件KAT( is accomplished through an integrated and modified version of Jellyfish2's counting method)
https://doi.org/10.1093/bioinformatics/btw663
Github: https://github.com/TGAC/KAT
Document:https://kat.readthedocs.io/en/latest/using.html

安裝

  1. 從bioconda里安裝最新版本的kat:
bioconda install kat
#這里我用以下這種簡單快捷的方法:
conda install -c bioconda kat
conda activate
kat -h
  1. 從GitHub安裝:
    需要先安裝很多依賴包晤郑,不然可能很多功能用不了:

-GCC V4.8+
-make
-autoconf V2.53+
-automake V1.11+
-libtool V2.4.2+
-pthreads (probably already installed)
-zlib
-Python V3.5+ with the tabulate, scipy, numpy and matplotlib packages and C API installed. This is optional but highly recommended, without python KAT functionality is limited: no plots, no distribution analysis, and no documentation.
-Sphinx-doc V1.3+ (Optional: only required for building the documentation.

然后是一系列的安裝編譯過程:

git clone https://github.com/TGAC/KAT.git
cd KAT
./build_boost.sh
./autogen.sh
./configure
make

二、運行

  1. 給定一個k值(kmer = 17贸宏,小基因組一般用17造寝,大基因組才用默認的27),產(chǎn)生不同kmers數(shù)量的直方圖
    usage:kat hist [options] (<input>)+
    options: -o 輸出文件 [kat.hist]锚赤;-t [1]線程數(shù)匹舞; -l [1] histogram的最低值;-h [10000] histogram的最高值线脚;-m [27] kmer的長度
nohup kat hist -o KAT_result/kat.hist -t 16 <(gzip -d -c NGS_10G_?.fq.gz) 2> kat_hist.log & #不支持壓縮格式赐稽,需先解壓
nohup kat hist -m 17 -o KAT_result/kat.hist  -t 16 NGS_10G_1.clean.fq NGS_10G_2.clean.fq &  #kat1.hist.png
nohup kat hist -m 17 -o KAT_result/kat.hist  -h 400 -t 16 NGS_10G_1.clean.fq NGS_10G_2.clean.fq &  #kat2.hist.png叫榕,縮短x軸看得更清楚測序深度
kat1.hist.png

kat2.hist.png
  1. Kmer的GC覆蓋圖
    usage: kat gcp (<input>)+
    options: -o 輸出文件 [kat.hist];-t [1]線程數(shù); -x [1] 當創(chuàng)建污染矩陣時姊舵,用于gc數(shù)據(jù)的bins的數(shù)量; -y [1000] 當創(chuàng)建污染矩陣時晰绎,用于coverage數(shù)據(jù)的bins的數(shù)量; -m [27] kmer的長度
nohup kat gcp -m 17 -o KAT_result/kat.gcp -t 12 NGS_10G_1.clean.fq NGS_10G_2.clean.fq &
kat.gcp.mx.png

SOAPdenovo2的Kmerfreq模塊

ls NGS_10G_?.clean.fq > NGS10.clean.lst
KmerFreq_AR -k 17 -t 12 -m 1 -p KmerFreq.10G -q 33 NGS10.clean.lst > NGS.10G.kmer.count 2>NGS.10G.kmerfreq.cerr

freq.cz和.freq.cz.len文件是用于Corrector_AR對reads進行校正分析

Corrector_AR -k 17 -t 12 -l 3 -Q 33 KmerFreq.10G.freq.cz KmerFreq.10G.freq.cz.len NGS10.clean.lst >corr.cout 2>corr.cerr ##沒有校正的reads就不需要再跑下一步
ls NGS_10G_1.clean.fq.cor.pair_1.fq.gz NGS_10G_2.clean.fq.cor.pair_2.fq.gz > kmer.input.fil.data
KmerFreq_AR -k 17 -t 12 -m 1 -p KmerFreq.10G.cor -q 33 kmer.input.fil.data > NGS.10G.kmer.cor.count 2>NGS.10G.kmerfreq.cor.cerr

jellyfish和genomescope

#/usr/bin/jellyfish    jellyfish 1.1.10
nohup /usr/bin/jellyfish count -m 17 -s 1G -t 12 -C -o jellyfish.out <(gzip -d -c NGS_10G_?.fq.gz) jellyfish.out.log & 
nohup /usr/bin/jellyfish count -m 17 -s 1G -t 12 -C -o jellyfish.out NGS_10G_1.fq NGS_10G_2.fq 2> jellyfish.out.log & 
/usr/bin/jellyfish histo -t 12 jellyfish.out_0 > jellyfish.histo 
~/Software/genomescope/genomescope.R jellyfish.histo 17 149 genomescope.result
plot.png
?著作權歸作者所有,轉載或內(nèi)容合作請聯(lián)系作者
  • 序言:七十年代末,一起剝皮案震驚了整個濱河市括丁,隨后出現(xiàn)的幾起案子荞下,更是在濱河造成了極大的恐慌,老刑警劉巖史飞,帶你破解...
    沈念sama閱讀 216,324評論 6 498
  • 序言:濱河連續(xù)發(fā)生了三起死亡事件尖昏,死亡現(xiàn)場離奇詭異,居然都是意外死亡构资,警方通過查閱死者的電腦和手機内地,發(fā)現(xiàn)死者居然都...
    沈念sama閱讀 92,356評論 3 392
  • 文/潘曉璐 我一進店門落追,熙熙樓的掌柜王于貴愁眉苦臉地迎上來孽亲,“玉大人对妄,你說我怎么就攤上這事〖旱ィ” “怎么了唉窃?”我有些...
    開封第一講書人閱讀 162,328評論 0 353
  • 文/不壞的土叔 我叫張陵,是天一觀的道長纹笼。 經(jīng)常有香客問我纹份,道長,這世上最難降的妖魔是什么允乐? 我笑而不...
    開封第一講書人閱讀 58,147評論 1 292
  • 正文 為了忘掉前任矮嫉,我火速辦了婚禮,結果婚禮上牍疏,老公的妹妹穿的比我還像新娘蠢笋。我一直安慰自己,他們只是感情好鳞陨,可當我...
    茶點故事閱讀 67,160評論 6 388
  • 文/花漫 我一把揭開白布昨寞。 她就那樣靜靜地躺著,像睡著了一般厦滤。 火紅的嫁衣襯著肌膚如雪援岩。 梳的紋絲不亂的頭發(fā)上,一...
    開封第一講書人閱讀 51,115評論 1 296
  • 那天掏导,我揣著相機與錄音享怀,去河邊找鬼。 笑死趟咆,一個胖子當著我的面吹牛添瓷,可吹牛的內(nèi)容都是我干的梅屉。 我是一名探鬼主播,決...
    沈念sama閱讀 40,025評論 3 417
  • 文/蒼蘭香墨 我猛地睜開眼鳞贷,長吁一口氣:“原來是場噩夢啊……” “哼坯汤!你這毒婦竟也來了?” 一聲冷哼從身側響起搀愧,我...
    開封第一講書人閱讀 38,867評論 0 274
  • 序言:老撾萬榮一對情侶失蹤惰聂,失蹤者是張志新(化名)和其女友劉穎,沒想到半個月后咱筛,有當?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體搓幌,經(jīng)...
    沈念sama閱讀 45,307評論 1 310
  • 正文 獨居荒郊野嶺守林人離奇死亡,尸身上長有42處帶血的膿包…… 初始之章·張勛 以下內(nèi)容為張勛視角 年9月15日...
    茶點故事閱讀 37,528評論 2 332
  • 正文 我和宋清朗相戀三年眷蚓,在試婚紗的時候發(fā)現(xiàn)自己被綠了鼻种。 大學時的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片。...
    茶點故事閱讀 39,688評論 1 348
  • 序言:一個原本活蹦亂跳的男人離奇死亡沙热,死狀恐怖,靈堂內(nèi)的尸體忽然破棺而出罢缸,到底是詐尸還是另有隱情篙贸,我是刑警寧澤,帶...
    沈念sama閱讀 35,409評論 5 343
  • 正文 年R本政府宣布枫疆,位于F島的核電站爵川,受9級特大地震影響,放射性物質(zhì)發(fā)生泄漏息楔。R本人自食惡果不足惜寝贡,卻給世界環(huán)境...
    茶點故事閱讀 41,001評論 3 325
  • 文/蒙蒙 一、第九天 我趴在偏房一處隱蔽的房頂上張望值依。 院中可真熱鬧圃泡,春花似錦、人聲如沸愿险。這莊子的主人今日做“春日...
    開封第一講書人閱讀 31,657評論 0 22
  • 文/蒼蘭香墨 我抬頭看了看天上的太陽辆亏。三九已至风秤,卻和暖如春,著一層夾襖步出監(jiān)牢的瞬間扮叨,已是汗流浹背缤弦。 一陣腳步聲響...
    開封第一講書人閱讀 32,811評論 1 268
  • 我被黑心中介騙來泰國打工, 沒想到剛下飛機就差點兒被人妖公主榨干…… 1. 我叫王不留彻磁,地道東北人碍沐。 一個月前我還...
    沈念sama閱讀 47,685評論 2 368
  • 正文 我出身青樓狸捅,卻偏偏與公主長得像,于是被迫代替她去往敵國和親抢韭。 傳聞我的和親對象是個殘疾皇子薪贫,可洞房花燭夜當晚...
    茶點故事閱讀 44,573評論 2 353

推薦閱讀更多精彩內(nèi)容