1.上Genome Announcements网站（https://mra.asm.org/）找一篇细菌基因组文章；找到文章记载的SRA号；
1
2
- 由于步骤都是一样的，我们只对菌株4041进行组装。
2.从SRA数据库上用prefetch下载该文件；
- 代码：
prefetch SRR5513009
prefetch -
因为中间有失去连接，我们不敢保证下载的序列是否完整，上ftp官网看一下大小，都是600M，应该下载完整了。
下载的序列
3. Fastq-dump解压，解压为gz文件，可以节省空间。因为需要点时间，我们让它在后台运行。
fastq-dump --gzip --split-files ~/ncbi/public/sra/SRR5513009.sra &
5
4. Fastqc质控
wwwww77@wwwww77-VirtualBox:~/assembly$ fastqc SRR5513009_1.fastq.gz
Started analysis of SRR5513009_1.fastq.gz
Approx 5% complete for SRR5513009_1.fastq.gz
Approx 10% complete for SRR5513009_1.fastq.gz
Approx 15% complete for SRR5513009_1.fastq.gz
Approx 20% complete for SRR5513009_1.fastq.gz
Approx 25% complete for SRR5513009_1.fastq.gz
Approx 30% complete for SRR5513009_1.fastq.gz
Approx 35% complete for SRR5513009_1.fastq.gz
Approx 40% complete for SRR5513009_1.fastq.gz
Approx 45% complete for SRR5513009_1.fastq.gz
Approx 50% complete for SRR5513009_1.fastq.gz
Approx 55% complete for SRR5513009_1.fastq.gz
Approx 60% complete for SRR5513009_1.fastq.gz
Approx 65% complete for SRR5513009_1.fastq.gz
Approx 70% complete for SRR5513009_1.fastq.gz
Approx 75% complete for SRR5513009_1.fastq.gz
Approx 80% complete for SRR5513009_1.fastq.gz
Approx 85% complete for SRR5513009_1.fastq.gz
Approx 90% complete for SRR5513009_1.fastq.gz
Approx 95% complete for SRR5513009_1.fastq.gz
Analysis complete for SRR5513009_1.fastq.gz
wwwww77@wwwww77-VirtualBox:~/assembly$ fastqc SRR5513009_2.fastq.gz
Started analysis of SRR5513009_2.fastq.gz
Approx 5% complete for SRR5513009_2.fastq.gz
Approx 10% complete for SRR5513009_2.fastq.gz
Approx 15% complete for SRR5513009_2.fastq.gz
Approx 20% complete for SRR5513009_2.fastq.gz
Approx 25% complete for SRR5513009_2.fastq.gz
Approx 30% complete for SRR5513009_2.fastq.gz
Approx 35% complete for SRR5513009_2.fastq.gz
Approx 40% complete for SRR5513009_2.fastq.gz
Approx 45% complete for SRR5513009_2.fastq.gz
Approx 50% complete for SRR5513009_2.fastq.gz
Approx 55% complete for SRR5513009_2.fastq.gz
Approx 60% complete for SRR5513009_2.fastq.gz
Approx 65% complete for SRR5513009_2.fastq.gz
Approx 70% complete for SRR5513009_2.fastq.gz
Approx 75% complete for SRR5513009_2.fastq.gz
Approx 80% complete for SRR5513009_2.fastq.gz
Approx 85% complete for SRR5513009_2.fastq.gz
Approx 90% complete for SRR5513009_2.fastq.gz
Approx 95% complete for SRR5513009_2.fastq.gz
Analysis complete for SRR5513009_2.fastq.gz
- 我们可以下载html文件到Windows端看一下结果。
-
从中我们可以知道输入文本的reads的数量是5843752，测序长度是35-151，GC含量是67%，有点高，但由于由于二代测序GC偏好性高，且深度越高，GC含量会越高。
-
从Per base sequence quality来看我们的reads大部分都在绿色区域，说明质量比较高。
-
碱基总体质量值也都在高质量区域。
fastqc.html
SRR5513009_1.fastq.gz
SRR5513009_1.fastq.gz
SRR5513009_1.fastq.gz
SRR5513009_2.fastq.gz
SRR5513009_2.fastq.gz
SRR5513009_2.fastq.gz
5.Trimmomatic去接头：
-
由文章可知这些数据是由illumina平台测序得到的，我们用Trimmomatic去除接头，因为这个软件其实就是专为illumina平台数据而设计的。
mkdir trim_out
java -jar ~/Biosofts/Trimmomatic038/Trimmomatic-0.38/trimmomatic-0.38.jar PE -phred33 SRR5513009_1.fastq.gz SRR5513009_2.fastq.gz ./trim_out/output_forward_paired.fq.gz ./trim_out/output_forward_unpaired.fq.gz ./trim_out/output_reverse_paired.fq.gz ./trim_out/output_reverse_unpaired.fq.gz ILLUMINACLIP:/home/wwwww77/Biosofts/Trimmomatic038/Trimmomatic-0.38/adapters/TruSeq3-PE.fa:2:30:10 SLIDINGWINDOW:4:15 LEADING:5 TRAILING:5 MINLEN:50
Trimmomatic
trim_out
6.再次FastQC对过滤后的数据进行质量测评
wwwww77@wwwww77-VirtualBox:~/assembly$ fastqc trim_out/output_forward_paired.fq.gz
Started analysis of output_forward_paired.fq.gz
Approx 5% complete for output_forward_paired.fq.gz
Approx 10% complete for output_forward_paired.fq.gz
Approx 15% complete for output_forward_paired.fq.gz
Approx 20% complete for output_forward_paired.fq.gz
Approx 25% complete for output_forward_paired.fq.gz
Approx 30% complete for output_forward_paired.fq.gz
Approx 35% complete for output_forward_paired.fq.gz
Approx 40% complete for output_forward_paired.fq.gz
Approx 45% complete for output_forward_paired.fq.gz
Approx 50% complete for output_forward_paired.fq.gz
Approx 55% complete for output_forward_paired.fq.gz
Approx 60% complete for output_forward_paired.fq.gz
Approx 65% complete for output_forward_paired.fq.gz
Approx 70% complete for output_forward_paired.fq.gz
Approx 75% complete for output_forward_paired.fq.gz
Approx 80% complete for output_forward_paired.fq.gz
Approx 85% complete for output_forward_paired.fq.gz
Approx 90% complete for output_forward_paired.fq.gz
Approx 95% complete for output_forward_paired.fq.gz
Analysis complete for output_forward_paired.fq.gz
wwwww77@wwwww77-VirtualBox:~/assembly$ fastqc trim_out/output_reverse_paired.fq.gz
Started analysis of output_reverse_paired.fq.gz
Approx 5% complete for output_reverse_paired.fq.gz
Approx 10% complete for output_reverse_paired.fq.gz
Approx 15% complete for output_reverse_paired.fq.gz
Approx 20% complete for output_reverse_paired.fq.gz
Approx 25% complete for output_reverse_paired.fq.gz
Approx 30% complete for output_reverse_paired.fq.gz
Approx 35% complete for output_reverse_paired.fq.gz
Approx 40% complete for output_reverse_paired.fq.gz
Approx 45% complete for output_reverse_paired.fq.gz
Approx 50% complete for output_reverse_paired.fq.gz
Approx 55% complete for output_reverse_paired.fq.gz
Approx 60% complete for output_reverse_paired.fq.gz
Approx 65% complete for output_reverse_paired.fq.gz
Approx 70% complete for output_reverse_paired.fq.gz
Approx 75% complete for output_reverse_paired.fq.gz
Approx 80% complete for output_reverse_paired.fq.gz
Approx 85% complete for output_reverse_paired.fq.gz
Approx 90% complete for output_reverse_paired.fq.gz
Approx 95% complete for output_reverse_paired.fq.gz
Analysis complete for output_reverse_paired.fq.gz
-
为了方便看过滤后数据的质量对比，我们用MultiQC把结果整合成一个HTLM网页交互式报告。
过滤后的正反序列质量报告居然一样，multiqc直接把它们识别为一个报告文件了。
multiqc *.zip
multiqc -
把multiqc_report.html用WinSCP下载到本地查看
发现其实过滤效果并不十分明显，其中reads重复率降低了一点，还有就是SRR5513009_2过滤前的每条reads各位置N碱基含量比例高了一点点，但其实也是处于高质量区域。
另外图三也说明了原序列基本没什么接头污染。
1
2
3
7.Spades组装基因组草图：
- 原文有提到参数要求,文件是paired-end reads,要选用--careful来减少错误和插入缺失
Genome assemblies were produced with SPAdes genome assembler version 3.10 (14), set in “paired-end assembly, careful mode,”
wwwww77@wwwww77-VirtualBox:~/assembly/trim_out$ spades.py --careful --pe1-1 output_forward_paired.fq.gz --pe1-2 output_reverse_paired.fq.gz -o ./SPAdes_out
- 出现报错了，我上网查了查SPAdes的err code ：255是由于RAM不够造成的。我们关闭虚拟机，把虚拟机的内存大小调大一点，我调到了5058MB。重启后再执行这个语句
报错
内存调整
-
SPAdes组装完成
image.png
8.Quast评价组装的基因组效果
- 可以自定义参数，skip contigs shorter than 200 bp
wwwww77@wwwww77-VirtualBox:~/assembly/trim_out$ quast.py SPAdes_out/contigs.fasta --min-contig 200 -o SPAdes_out/quast_out
quast
- quast执行完成后结果有很多，我们可以直接查看report.txt；
我们也可查看其中的网页版报告，一般看icarus.html，其为导航页面，更便于查看更多结果。
quast结果·
-
一般contigs/scaffolds序列总数越少、序列总长度合理、N50等值越高长，组装结果越好
从report.txt可以看到Arthrobacter sp. 4041总基因组长度为3912868 bp，GC含量为67.65%，N50值为536987bp 。
report.txt
image.png
