微生物比较基因组学¶
**2024.11.9
思维导图¶
flowchart TD
A[Rawdata] --fastp去除adapter,duplication--> B{Cleandata}
B -- Snippy --> C[VCF]
B -- Spades组装 --> D[Scaffolds]
D -- quast,busco -->H[QC data]
C --GWAS--> F[Final Results]
D --prokka注释--> E[基因/功能基因组]
E --kaptive--> G[分型]
E --roary--> I[pan-genome]
I --scoary GWAS-->F
E --pyseer GWAS-->F
G -->F[Final Results]
[
质控QC¶
- fastp(去除低质量序列,去除重复序列)
- Multiqc (整合多个质控结果)
第一次质控¶
# -i 输入序列R1
# -I 输入序列R2
# -o 输出序列R1
# -O 输出序列R2
# --dedup 进行去重复操作
# -w 8 线程数
# -j 输出json格式的报告
# -h 输出网页格式的报告
fastp -i L1EHH1600255--RJ1450.R1.raw.fastq.gz -I L1EHH1600255--RJ1450.R2.raw.fastq.gz -o ./cleandata/L1EHH1600255--RJ1450.R1.clean.fastq.gz -O ./cleandata/L1EHH1600255--RJ1450.R2.clean.fastq.gz --dedup -w 8 -j ./cleandata/L1EHH1600255--RJ1450.fastp.json -h ./cleandata/L1EHH1600255--RJ1450.fastp.html
# 写一个循环批量运行
while read id
do
echo "fastp -i ${id}.R1.raw.fastq.gz -I ${id}.R2.raw.fastq.gz -o ./cleandata/${id}.R1.clean.fastq.gz -O ./cleandata/${id}.R2.clean.fastq.gz --dedup -w 8 -j ./cleandata/${id}.fastp.json -h ./cleandata/${id}.fastp.html "
done<id.txt > 01.fastp.sh
nohup bash 01.fastp.sh &
整合第一次质控结果¶
multiqc ./
第二次质控¶
fastp -i L1EHH1600255--RJ1450.R1.clean.fastq.gz -I L1EHH1600255--RJ1450.R2.clean.fastq.gz -Q -h ./cleandata2/L1EHH1600255--RJ1450.fastp.html -j ./cleandata2/L1EHH1600255--RJ1450.fastp.json -w 8
# 写一个循环批量运行
while read id
do
echo "fastp -i ${id}.R1.clean.fastq.gz -I ${id}.R2.clean.fastq.gz -Q -w 8 -j ./cleandata2/${id}.fastp.json -h ./cleandata2/${id}.fastp.html "
done<id.txt > 02.fastp.stats.sh
nohup bash 02.fastp.stats.sh &
整合第二次质控结果¶
multiqc ./
组装Assembly¶
# --careful 仔细地组装
# -1 质控后的R1
# -2 质控后的R2
# -t 12 线程数为12
# -o 输出目录(指定文件夹名称)
spades.py --careful -1 L1EHH1600255--RJ1450.R1.clean.fastq.gz -2 L1EHH1600255--RJ1450.R2.clean.fastq.gz -t 12 -o 03assembly/RJ1450
# 写一个循环批量运行
while read id1 id2
do
echo "spades.py --careful -1 ${id1}.R1.clean.fastq.gz -2 ${id1}.R2.clean.fastq.gz -t 12 -o 03assembly/${id2}"
done<id12.txt > 03.spades.sh
nohup bash 03.spades.sh&
## !!组装后记得改文件名!!
find -maxdepth 2 -name scaffolds.fasta > id1.txt #找到文件路径
cat id1.txt |cut -d "/" -f 2 > id2.txt #提取文件名
paste id1.txt id2.txt > id12.txt #生成文件路径和文件名的对应表格
while read id1 id2; do cp ${id1} ./${id2}/${id2}.fna; done<id12.txt #生成新的文件
#seqkit replace -p .+ -r "scaffold_{nr}" RJ1450.fna > RJ1450.newname.fasta
while read id; do echo "seqkit replace -p .+ -r '${id}_scaffold_{nr}' ${id}.fasta > ${id}.fna"; done<id2.txt #修改文件内容
组装后的质控¶
-
CheckM 用来评估完整性和污染度
-
Quast 评估基因组组装的连续性和一致性
-
Busco 评估基因组中功能性基因的保守性和完整性
singularity exec ~/download/quast_5.0.2.sif quast ./assmebly/* -r ~/download/ref.fna #会自动生成文件夹
注释Annotation¶
- prokka
singularity exec ~/download/prokka.sif prokka --outdir RJ1450 --prefix RJ1450 --kingdom Bacteria --cpus 8 --locustag RJ1450 --force --evalue 0.001 --proteins ~/download/kp.faa RJ1450.fna
# 自己写一个循环