组装细菌基因组
一、从上Genome Announcements网站（https://mra.asm.org/）找一篇细菌基因组文章；找到文章记载的SRA号，如SRR9305630
二、从SRA数据库上用prefetch下载该文件
prefetch SRR9305630
三、Fastq-dump解压
fastq-dump --gzip --split-files SRR9305630.sra
四、Fastqc评价数据质量；
fastqc SRR9305630_1.fastq.gz
fastqc SRR9305630_2.fastq.gz
数据过滤，去接头
java -jar ~/Biosofts/Trimmomatic038/Trimmomatic-0.38/trimmomatic-0.38.jar PE -phred33 SRR9305630_1.fastq.gz SRR9305630_2.fastq.gz ./output_forward_paired.fq.gz ./output_forward_unpaired.fq.gz ./output_reverse_paired.fq.gz ./output_reverse_unpaired.fq.gz ILLUMINACLIP:/home/wangfuhao/Biosofts/Trimmomatic038/Trimmomatic-0.38/adapters/TruSeq2-PE.fa:2:30:10 SLIDINGWINDOW:5:20 LEADING:20 TRAILING:20 MINLEN:75
image.png
五、Spades组装基因组草图
运行 SPAdes 组装细菌基因组
spades.py --careful --pe1-1 SRR9305630_1.fastq.gz --pe1-2 SRR9305630_2.fastq.gz -o ./SPAdesout_7942_new
出现错误，内存不够
用seqtk抽取100000条再试试
先解压
gunzip -c output_forward_paired.fq >output_forward_paired.fq
gunzip -c output_reserve_paired.fq >output_reserve_paired.fq
抽取100000条
seqtk sample -s 60 output_forward_paired.fq 100000 >seqtksample_100000.fq
seqtk sample -s 60 output_reserve_paired.fq 100000 >seqtksample2_100000.fq
再次执行
spades.py --careful --pe1-1 seqtksample_100000.fq --pe1-2 100000 >seqtksample2_100000.fq -o ./SPAdesout_7942_new
有错误，但是算是完成了这一步
六、Quast评价组装的基因组效果
quast.py ~/ncbi/public/sra/SPAdesout_7942_new/contigs.fasta -o ~/SPAdesout_7942_new/quast_out
icarus下载到本地查看
