[toc]

寫(xiě)在前面

kraken基于mini數(shù)據(jù)庫(kù)。并且這個(gè)序列也比較少，所以，很快就能完成

繼續(xù)處理

膠水操作：提取序列名稱zcat ./trimmomatic/SUBERR793599_forward_paired.fq.gz | awk '(NR%4 == 1){{print $0}}' | cut -d' ' -f 1 | cut -c 2- > ./trimmomatic/SUBERR793599_1.names.txt

# ${base}

for i in ./rawdata_0.01/*_1.fq.gz

do

base=$(basename $i _1.fq.gz)

echo $base

zcat ./trimmomatic/${base}_forward_paired.fq.gz | awk '(NR%4 == 1){{print $0}}' | cut -d' ' -f 1 | cut -c 2- > ./trimmomatic/${base}_1.names.txt

done

膠水操作：合并雙端數(shù)據(jù)

合并測(cè)序文件：這部分?jǐn)?shù)據(jù)存在于./unpack/文件夾，是未經(jīng)過(guò)質(zhì)控和barcode去除的最原始的數(shù)據(jù)。# 用于豐度統(tǒng)計(jì) 加州大學(xué)戴維斯分校的成果

pip install khmer

interleave-reads.py -o ./unpack/SUBERR793599.fastq ./unpack/SUBERR793599_1.fastq.gz ./unpack/SUBERR793599_2.fastq.gz

# ${base}

for i in ./rawdata_0.01/*_1.fq.gz

do

base=$(basename $i _1.fq.gz)

echo $base

interleave-reads.py -o ./unpack/${base}.fastq ./unpack/${base}_1.fastq.gz ./unpack/${base}_2.fastq.gz

done

膠水操作：合并全部樣本的雙端序列cat ./unpack/SUBERR793599.fastq > ./interleaved_merged/all.fastq

# `tail -n+2 ./map.txt |cut -f 1|sed 's/^//unpack\//;s/$/_.fastq/'| tr '\n' ','|sed 's/,$//'`

這個(gè)到這里就先停下了，又開(kāi)始了另外一條注釋線。

kraken注釋物種信息

kraken物種注釋

kraken也是基于read的注釋，這里呢，我們使用triming文件夾中的數(shù)據(jù)進(jìn)行注釋

安裝krakenconda install kraken -c bioconda

kraken 注釋流程mkdir kraken

# sudo apt install pigz

pigz ./trimming/*_1.fq

pigz ./trimming/*_2.fq

kraken --db ~/db/kraken/minikraken_20171019_8GB --threads 6 --fastq-input ./trimming/SUBERR793599_1.fq.gz ./trimming/SUBERR793599_2.fq.gz --paired --check-names --gzip-compressed --output ./kraken/SUBERR793599.kraken 2> ./kraken/SUBERR793599.kraken.log

kraken-translate --db ~/db/kraken/minikraken_20171019_8GB --mpa-format ./kraken/SUBERR793599.kraken > ./kraken/SUBERR793599.kraken.taxonomy

kraken-mpa-report --db ~/db/kraken/minikraken_20171019_8GB ./kraken/SUBERR793599.kraken > ./kraken/SUBERR793599.kraken.report

批量運(yùn)行# ${base}

for i in ./rawdata_0.01/*_1.fq.gz

do

base=$(basename $i _1.fq.gz)

echo $base

kraken --db ~/db/kraken/minikraken_20171019_8GB --threads 6 --fastq-input ./trimming/${base}_1.fq.gz ./trimming/${base}_2.fq.gz --paired --check-names --gzip-compressed --output ./kraken/${base}.kraken 2> ./kraken/${base}.kraken.log

kraken-translate --db ~/db/kraken/minikraken_20171019_8GB --mpa-format ./kraken/${base}.kraken > ./kraken/${base}.kraken.taxonomy

kraken-mpa-report --db ~/db/kraken/minikraken_20171019_8GB ./kraken/${base}.kraken > ./kraken/${base}.kraken.report

done

后續(xù)轉(zhuǎn)化物種注釋文件：kraken2phyloseq 這部分?jǐn)?shù)據(jù)清晰沒(méi)什么卵用不過(guò)咱們也不需要，直接在R進(jìn)行然后可視化就可以了。# echo -e \"#OTU ID\tSUBERR793599\tConsensus Lineage\" > ./kraken/SUBERR793599.kraken.qiime.taxonomy

# awk 'BEGIN {{OFS=\"\t\"}} {{ print $1,1,$2 }}' ./kraken/SUBERR793599.kraken.taxonomy | sed 's/|/;/g' >> ./kraken/SUBERR793599.kraken.qiime.taxonomy

作者寫(xiě)了這幾行，證明大工作走完了，開(kāi)始作圖了== #TODO: ==

== # Make plots for raw reads tax assignment with kraken/diamond ==

== # Convert to a dummy otutable to load with phyloseq with a fake count column ==

== # Add qiime header (todo) ==

# awk 'BEGIN {OFS="\t"} { print $1,1,$2 }' 1511KMI-0007/kraken/MPR1-1n-EST.kraken.taxonomy | sed 's/|/;/g' > test2.qiime

# make plots with gpplot

到這兒，咱們解決了整個(gè)代碼的500行內(nèi)容了。

配置數(shù)據(jù)庫(kù)

這將下載NCBI分類信息，以及RefSeq中細(xì)菌，古細(xì)菌和病毒域的完整基因組。下載所有這些數(shù)據(jù)之后，構(gòu)建過(guò)程開(kāi)始；這是最耗時(shí)的步驟。如果您具有多個(gè)處理核心，則可以使用多個(gè)線程來(lái)運(yùn)行此過(guò)程，例如：kraken-build --standard --threads 24 --db ~/db/kraken

請(qǐng)注意，如果任何步驟（包括初始下載）失敗，則構(gòu)建過(guò)程將中止。但是，kraken-build將在整個(gè)安裝過(guò)程中產(chǎn)生檢查點(diǎn)，并且如果您嘗試在部分構(gòu)建的數(shù)據(jù)庫(kù)上再次運(yùn)行同一命令，則會(huì)在最后一個(gè)未完成的步驟中重新啟動(dòng)構(gòu)建。

構(gòu)建數(shù)據(jù)庫(kù)之后，要除去任何不必要的文件（包括不再需要的庫(kù)文件），請(qǐng)運(yùn)行以下命令：kraken-build --db ~/db/kraken --clean

數(shù)據(jù)庫(kù)完成后應(yīng)該有這幾個(gè)文件database.kdb：包含k -mer到分類單元的映射

database.idx：在database.kdb中包含最小化程序偏移量位置

taxonomy/nodes.dmp：分類樹(shù)結(jié)構(gòu)+等級(jí)

taxonomy/names.dmp：分類名稱

如果要構(gòu)建自定義數(shù)據(jù)庫(kù)，請(qǐng)參考：http://ccb./software/kraken/MANUAL.html#standard-kraken-database

這里考慮到大家的內(nèi)存

您可以獲取現(xiàn)有的Kraken數(shù)據(jù)庫(kù)并從中創(chuàng)建較小的MiniKraken數(shù)據(jù)庫(kù)。使用此選項(xiàng)將刪除除指定數(shù)量的k -mer / taxon對(duì)，以創(chuàng)建一個(gè)新的較小的數(shù)據(jù)庫(kù)。例如：kraken-build --shrink 10000 --db ~/db/minikraken --new-db minikraken

根際互作生物學(xué)研究室 簡(jiǎn)介

根際互作生物學(xué)研究室是沈其榮教授土壤微生物與有機(jī)肥團(tuán)隊(duì)下的一個(gè)關(guān)注于根際互作的研究小組。本小組由袁軍副教授帶領(lǐng)，主要關(guān)注：1.植物和微生物互作在抗病過(guò)程中的作用；2 環(huán)境微生物大數(shù)據(jù)整合研究；3 環(huán)境代謝組及其與微生物過(guò)程研究體系開(kāi)發(fā)和應(yīng)用。團(tuán)隊(duì)在過(guò)去三年中在 isme J， Microbiome， PCE，SBB，Horticulture Research等期刊上發(fā)表了多篇文章。歡迎關(guān)注 微生信生物 公眾號(hào)對(duì)本研究小組進(jìn)行了解。

團(tuán)隊(duì)工作及其成果 （點(diǎn)擊查看）

團(tuán)隊(duì)關(guān)注

團(tuán)隊(duì)文章成果

團(tuán)隊(duì)成果-EasyStat專題

ggClusterNet專題

袁老師小小組

了解 交流 合作

團(tuán)隊(duì)成員郵箱 袁軍：

junyuan@njau.edu.cn；

文濤：

2018203048@njau.edu.cn

團(tuán)隊(duì)公眾號(hào)：