最近在摸索著學(xué)習(xí)宏基因組的數(shù)據(jù)分析摩钙,記錄一下,方便以后看说墨,哪里有錯誤希望看到的人可以指出來,如果有幸?guī)偷搅四悴园兀覍荛_心尼斧。
第二篇:宏基因組學(xué)習(xí)記錄-基因預(yù)測
一、質(zhì)控
- 原始數(shù)據(jù)質(zhì)量信息
這次用的fastqc试吁,當(dāng)然還有其他軟件突颊,比如fastp進(jìn)行質(zhì)控,還可以生成報告潘悼,也很方便
# conda 無腦安裝
conda install -c bioconda fastqc
# fastqc
fastqc -t 10 -o fastqc sample1_1.fq sample1_2.fq ...
- 數(shù)據(jù)質(zhì)控
# 下載直接解壓使用,網(wǎng)頁:http://www.usadellab.org/cms/?page=trimmomatic
wetget http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/Trimmomatic-0.39.zip
# 解壓
unzip Trimmomatic-0.39.zip
# 質(zhì)控 shell寫的小循環(huán)
ls raw_data/*fq.gz | while read id
do
echo $id
java -jar /softs_path/Trimmomatic-0.39/trimmomatic-0.39.jar PE -threads 10 ${id%%.*}.R1.fq.gz ${id%%.*}.R2.fq.gz ${id%%.*}.R1.trim.fq ${id%%.*}.R1.unpaired.fq ${id%%.*}.R2.trim.fq ${id%%.*}.R2.unpaired.fq ILLUMINACLIP:~/softs/Trimmomatic-0.39/adapters/TruSeq2-PE.fa:2:30:10 LEADING:10 TRAILING:10 SLIDINGWINDOW:4:20 MINLEN:70;done
3.去宿主
# 下載人基因組文件eg:hg38
wget ftp://ftp.ensembl.org/pub/release-101/fasta/homo_sapiens/dna/Homo_sapiens.GRCh38.dna_rm.primary_assembly.fa.gz
# conda 無腦安裝 bowtie2
conda install -c bioconda bowtie2
# 構(gòu)建index
bowtie2-build genome.fasta genome
#比對去宿主
ls *trim.fq | while read id
do
echo $id
bowtie2 -p 10 -x /index_path/genome -1 ${id%%.*}.R1.trim.fq -2 ${id%%.*}.R2.trim.fq -S ${id%%.*}.sam 2> ${id%%.*}.bowtie2.log --un-conc-gz ${id%%.*} ;done
# 生成文件sample.1 ,sample.2是去除宿主的序列,重命名
mv sample.1 sample.R1.fq.gz
mv sample.2 sample.R2.fq.gz
最后可以用fastqc看一下clean_reads質(zhì)量情況
二爬橡、組裝
組裝軟件很多,網(wǎng)上也有很多大神總結(jié)的各種軟件的優(yōu)缺點(diǎn)治唤,自行查閱,各取所需
我暫時選擇了megahit這個軟件糙申,資源消耗較少宾添,速度較快,準(zhǔn)確性也滿足基本需求了柜裸。
- 軟件安裝
#當(dāng)然是選擇conda了
conda install -c bioconda megahit
# 查看軟件使用方法
megahit -h
2.合并數(shù)據(jù)
可以分開單組組裝缕陕,也可以合并組裝,此次選擇合并組裝
cat *_R1.fq.gz > all_reads_R1.fq.gz
cat *_R2.fq.gz > all_reads_R2.fq.gz
3.組裝
megahit -1 all_reads_R1.fq.gz -2 all_reads_R2.fq.gz -o assembly/ --out-prefix assembly -t 10 --min-contig-len 300
#簡單介紹一下參數(shù)
-1: R1 reads
-2: R2 reads
-o: 結(jié)果目錄
--out-prefix: 輸出結(jié)果前綴
-t: 線程數(shù)
--min-contig-len: 最小組裝長度
4.組裝質(zhì)量評估
#下載疙挺,解壓扛邑,安裝
wget https://nchc.dl.sourceforge.net/project/quast/quast-5.0.2.tar.gz
tar -xzf quast-5.0.2.tar.gz
cd quast-5.0.2/
python setup.py install_full
#conda 安裝
conda install -c bioconda quast
quast.py assembly.contigs.fasta -o quast
碼字不易,轉(zhuǎn)載請注明出處铐然,謝謝~