16S多样性分析docker镜像发布及使用

16S多样性分析docker镜像发布及使用
众所周知微生物多样性分析用到的软件众多,对于初学者而言安装软件是一件非常麻烦的事情,这里我们组学大讲堂发布了专门用于微生物扩增子分析docker镜像,里面内置了几乎所有要用到的关于多样性分析的相关软件,方便大家使用,只需安装docker并下载该镜像即可,避免繁琐的软件安装问题;

扩增子分析流程介绍

目前主流的扩增子分析流程包括,qiime1,qiime2 ,mothur ,usearch等。qiime以其方便好用的流程,被广大研究人员官方使用。逐渐的qiime2会替代qiime1,但是由于数据分析周期较长很多人还在使用qiime1并且也习惯使用qiime1的分析流程。因此我们提供这个ampliseq1镜像仍然安装qiime1,但也会安装:mothur ,usearch 软件方便大家分析数据,qiime2会安装在另一个单独的镜像中,后续发布;
attachments-2020-09-6pz6fYEb5f4f34647c333.png

镜像下载使用方法

1)首先系统中docker安装与使用方法见课程:


attachments-2020-09-q4CKuYCd5f4f347b30229.png


2)镜像搜索、下载、启动使用方法:
> docker search omicsclass   #搜索组学大讲堂提供的所有镜像
NAME                           DESCRIPTION                                     STARS               OFFICIAL            AUTOMATED
omicsclass/gene-family         gene-family analysis docker image               5                  
omicsclass/rnaseq              RNA-seq analysis docker image build by omics…   3                 
omicsclass/gsds-v2             GSDS 2.0 – Gene Structure Display Server        1                  
omicsclass/reseq               whole genome resequence analysis                1                  
omicsclass/biocontainer-base   Biocontainers base Image centos7                1                  
omicsclass/blast-plus          blast+ v2.9.0                                   0                  
omicsclass/isoseq3             isoseq3 v3.3.0 build by omicsclass              0                  
omicsclass/bwa                 BWA v0.7.17 build by omicsclass                 0                  
omicsclass/blastall            legacy blastall v2.2.26                         0                  
omicsclass/sratoolkit          SRAtoolkit v2.10.3 and aspera v3.9.9.177872     0                  
omicsclass/ampliseq-q2         Amplicon sequencing qiime2 v2020.2 image        0                  
omicsclass/ampliseq-q1         Amplicon sequencing qiime1 v1.9.1 image         0                  
omicsclass/samtools            samtools v1.10 build by omicsclass              0                  
omicsclass/bsaseq              NGS Bulk Segregant Analysis image               0                  
omicsclass/gwas                gwas analysis images                            0                  
> docker pull omicsclass/ampliseq-q1  #下载扩增子镜像
>docker run --rm -it -m 4G --cpus 1  -v D:\qiime1-16s:/work omicsclass/ampliseq-q1:latest    #启动并进入镜像

安装软件及版本介绍

OTU分析相关软件
qiime1 v1.9.1
mothur  v.1.25.0
usearch   10.0.240
usearch61 6.1.544
vsearch v2.15.0
序列处理相关分析包
flash v2.15.0 序列合并
bioperl 
biopython
fastqc
multiqc
fastp 

qiime分析相关的包

blast-2.2.22
blat-36
cdhit-3.1
muscle-3.8.31
rdpclassifier-2.2
uclust
功能注释相关分析软件
picrust-1.1.4
bugbase
数据分析与可视化相关的R包:
ggplot2
Ternary
DESeq2
edgeR
ggtree
vegan
pheatmap

高级分析相关包

randomForest  机器学习
scikit-learn  机器学习
circos  圈图绘制
Krona  物种丰度圈图
LefSe 差异比较分析

注意:以上软件包的列表只是部分举例,实际安装的包还有更多。

扩增子镜像使用案例

1.数据准备

fastmap.txt文件, 测序数据文件放在data文件夹中,物种注释数据库文件Greengene silva unite等放在database目录,需要这些测试数据的同学可以关注组学大讲堂公众号,并回复16s,即可得到;准备完成之后,目录结构如下:

[root@aefe86d682b1  13:58:49 /work/amplicon_demo]# tree 
 data
    |-- ERR3975186_1.fastq.gz
    |-- ERR3975186_2.fastq.gz
    |-- ERR3975187_1.fastq.gz
    `-- ERR3975187_2.fastq.gz
 fastmap.txt 

设置一些环境变量方便后续调用:

dbdir=/work/database
workdir=/work/amplicon_demo
datadir=$workdir/data
fastmap=$workdir/fastmap.txt
mkdir /work/tmp
export TMPDIR=/work/tmp  #防止临时目录存储 不够
######################database
#各大数据库地址:
silva_16S_97_seq=$dbdir/SILVA_132_QIIME_release/rep_set/rep_set_16S_only/97/silva_132_97_16S.fna 
silva_16S_97_tax=$dbdir/SILVA_132_QIIME_release/taxonomy/16S_only/97/taxonomy_7_levels.txt 

greengene_16S_97_seq=$dbdir/gg_13_8_otus/rep_set/97_otus.fasta
greengene_16S_97_tax=$dbdir/gg_13_8_otus/taxonomy/97_otu_taxonomy.txt

silva_18S_97_seq=$dbdir/SILVA_132_QIIME_release/rep_set/rep_set_18S_only/97/silva_132_97_18S.fna 
silva_18S_97_tax=$dbdir/SILVA_132_QIIME_release/taxonomy/18S_only/97/taxonomy_7_levels.txt 

unite_ITS_97_seq=$dbdir/unite_ITS8.2/sh_refs_qiime_ver8_97_04.02.2020.fasta
unite_ITS_97_tax=$dbdir/unite_ITS8.2/sh_taxonomy_qiime_ver8_97_04.02.2020.txt

2.合并双端reads,flash

cd $workdir  #回到工作目录
mkdir 1.merge_pe

for i in  `cat $fastmap |grep -v '#'|cut -f 1` ;do 
    echo "RUN CMD: flash $datadir/${i}_1.fastq.gz $datadir/${i}_2.fastq.gz \
        -m 10 -x 0.2 -p 33 -t 1 \
        -o $i -d 1.merge_pe"


    flash $datadir/${i}_1.fastq.gz $datadir/${i}_2.fastq.gz \
        -m 10 -x 0.2 -p 33 -t 1 \
        -o $i -d 1.merge_pe
done
3.对原始数据进行fastqc质控
cd $workdir  #回到工作目录
mkdir 2.fastqc
#fastqc查看数据质量分布等
fastqc -t 2 $workdir/1.merge_pe/*extendedFrags.fastq  -o $workdir/2.fastqc
#质控结果汇总
cd $workdir/2.fastqc
multiqc .

4.数据质控:对原始序列进行去接头,引物,删除低质量的reads等等

cd $workdir  #回到工作目录
mkdir 3.data_qc
cd  3.data_qc
#利用fastp工具去除adapter
#--qualified_quality_phred the quality value that a base is qualified. 
#            Default 15 means phred quality >=Q15 is qualified. (int [=15])
#--unqualified_percent_limit how many percents of bases are allowed to be unqualified
#--n_base_limit if one read's number of N base is >n_base_limit, 
#            then this read/pair is discarded 
#--detect_adapter_for_pe   接头序列未知  可设置软件自动识别常见接头
#
for i in `cat $fastmap |grep -v '#'|cut -f 1`; do 
    echo "RUN CMD: fastp --thread 1 --qualified_quality_phred 10 \
        --unqualified_percent_limit 50 \
        --n_base_limit 10 \
        --length_required 300 \
        --trim_front1 29 \
        --trim_tail1 18 \
        -i $workdir/1.merge_pe/${i}.extendedFrags.fastq \
        -o ${i}.clean_tags.fq.gz \
        --adapter_fasta $workdir/data/illumina_multiplex.fa -h ${i}.html -j ${i}.json"


    fastp --thread 1 --qualified_quality_phred 10 \
        --unqualified_percent_limit 50 \
        --n_base_limit 10 \
        --length_required 300 \
        --trim_front1 29 \
        --trim_tail1 18 \
        -i $workdir/1.merge_pe/${i}.extendedFrags.fastq \
        -o ${i}.clean_tags.fq.gz \
        --detect_adapter_for_pe -h ${i}.html -j ${i}.json
done
去除嵌合体
cd $workdir  #回到工作目录mkdir 4.remove_chimerascd  4.remove_chimeras
#去除嵌合体
for i in `cat $fastmap |grep -v '#'|cut -f 1`; do #相同重复序列合并 vsearch --derep_fulllength $workdir/3.data_qc/${i}.clean_tags.fq.gz \ --sizeout --output ${i}.derep.fa #去嵌合体 vsearch --uchime3_deno ${i}.derep.fa \ --sizein --sizeout \ --nonchimeras ${i}.denovo.nonchimeras.rep.fa #相同序列还原为多个 vsearch --rereplicate ${i}.denovo.nonchimeras.rep.fa --output ${i}.denovo.nonchimeras.fadone
#根据参考序列去除嵌合体for i in `cat $fastmap |grep -v '#'|cut -f 1`; do vsearch --uchime_ref ${i}.denovo.nonchimeras.fa \ --db $dbdir/rdp_gold.fa \ --sizein --sizeout --fasta_width 0 \ --nonchimeras ${i}.ref.nonchimeras.fadone


qiime1分析pick otu 聚类方法
cd $workdir  #回到工作目录
mkdir 5.pick_otu_qiime
cd   5.pick_otu_qiime

#合并fasta文件,并加序列号
for i in `cat $fastmap |grep -v '#'|cut -f 1`do 
    rename_fa_id.pl  -f $workdir/4.remove_chimeras/$i.ref.nonchimeras.fa \
        -n $i -out $i.fa
done

#合并fa文件到qiime.fasta 之后删除所有单个样本的fa文件
cat *fa >qiime.fasta
rm -f *fa

###方法1:pick_de_novo_otus.py
###输出qiime pick otu 参数,更多:http://qiime.org/scripts/pick_otus.html

echo pick_otus:denovo_otu_id_prefix OTU >> otu_params_de_novo.txt
echo pick_otus:similarity 0.97 >> otu_params_de_novo.txt
echo pick_otus:otu_picking_method uclust >> otu_params_de_novo.txt   #sortmerna, mothur, trie, uclust_ref, usearch, usearch_ref, blast, usearch61, usearch61_ref,sumaclust, swarm, prefix_suffix, cdhit, uclust.
echo assign_taxonomy:reference_seqs_fp $silva_16S_97_seq >> otu_params_de_novo.txt
echo assign_taxonomy:id_to_taxonomy_fp $silva_16S_97_tax >> otu_params_de_novo.txt
echo assign_taxonomy:similarity 0.8 >>otu_params_de_novo.txt
echo assign_taxonomy:assignment_method uclust >>otu_params_de_novo.txt    # rdp, blast,rtax, mothur, uclust, sortmerna如果是ITS/18S数据,建议数据库更改为UNITE,方法改为blast。详细使用说明,请读官方文档:http://qiime.org/scripts/assign_taxonomy.html

pick_de_novo_otus.py -i qiime.fasta -f -o pick_de_novo_otus -p otu_params_de_novo.txt  

alpha_diversity 分析

cd $workdir  #回到工作目录
mkdir 8.alpha_diversity
cd   8.alpha_diversity
#alpha多样性指数展示
biom summarize-table  -i $workdir/5.pick_otu_qiime/pick_de_novo_otus/otu_table_clean_rare.biom
echo alpha_diversity:metrics observed_species,PD_whole_tree,shannon,chao1,simpson,goods_coverage > alpha_params.txt
alpha_rarefaction.py -f -i $workdir/5.pick_otu_qiime/pick_de_novo_otus/otu_table_clean.biom -m $fastmap -o ./ -p alpha_params.txt -t $workdir/5.pick_otu_qiime/pick_de_novo_otus/rep_set.tre --retain_intermediate_files --min_rare_depth 40 --max_rare_depth 2032 --num_steps 10

#多样性指数差异比较qiime 自带检验与绘图
compare_alpha_diversity.py \
    -i alpha_div_collated/chao1.txt \
    -o alpha_chao1_stats \
    -m $fastmap \
    -t nonparametric \
    -c city

compare_alpha_diversity.py \
    -i alpha_div_collated/chao1.txt \
    -o alpha_chao1_stats \
    -m $fastmap \
    -t parametric \
    -c city

beta_diversity 分析

cd $workdir  #回到工作目录
mkdir 9.beta_diversity
cd   9.beta_diversity

echo beta_diversity:metrics binary_jaccard,bray_curtis,unweighted_unifrac,weighted_unifrac,binary_euclidean > beta_params.txt

#-e 设置抽平数
beta_diversity_through_plots.py -f -i $workdir/5.pick_otu_qiime/pick_de_novo_otus/otu_table_clean.biom -m $fastmap -o ./ -t $workdir/5.pick_otu_qiime/pick_de_novo_otus/rep_set.tre -e 2844 -p beta_params.txt
#beta多样性adonis检验
compare_categories.py --method adonis -i unweighted_unifrac_dm.txt -m $fastmap -c Treatment -o adonis_out -n 999

部分分析结果展示:


  • 发表于 2020-09-02 13:58
  • 阅读 ( 289 )
  • 分类:宏基因组

4 条评论

请先 登录 后评论
omicsgene
omicsgene

生物信息

367 篇文章

作家榜 »

  1. omicsgene 367 文章
  2. 安生水 217 文章
  3. Daitoue 167 文章
  4. 生物女学霸 120 文章
  5. CORNERSTONE 72 文章
  6. 红橙子 50 文章
  7. 生信老顽童 48 文章
  8. landy 37 文章