1. NT数据库下载（跳过，这里用的是下好的）
2. NT库比对 (XML格式输出)
   1. 将组装后的fa文件和nt库进行比对，生成 XML 格式输出的 BLAST 结果
   2. query，目标组装的序列
   3. out，输出文件表格（？）
   4. outfmt5，列表的一种格式来着
   5. db，nt数据库，这里是在nt数据库的上两个层级的目录下挂在了
   6. num_threads，现成数量
   7. perc_identity，设定最低的百分比相似度为 50%，即只有相似度大于等于 50% 的比对结果才会被保留
   8. subject_besthit，保留最佳比对结果
   9. 拷贝一下建库产生的文件，因为nt库下面文件太多，不方便操作

cp nt.n* /home/asm_first/ctg/

1. 下载ncbi的脚本包，运行比对命令

/home/ncbi-blast-2.15.0+/bin/blastn \
  -query pulchrum_asm_third.hic.p_ctg.fa \
  -out pulchrum_asm_third.hic.p_ctg.fa.blast_nt.xml \
  -outfmt 5 \
  -db /home/server/nt_makeblastdb/decompress_nt/nt \
  -num_threads 30 \
  -evalue 1e-5 \
  -perc_identity 50 \
  -subject_besthit

1. 使用 blast_xml_parse.py 脚本解析生成的 XML 文件，输出 CSV 格式的比对结果。

nohup python blast_xml_parse.py \
  -i pulchrum_asm_third.hic.p_ctg.fa.blast_nt.xml \
  -q pulchrum_asm_third.hic.p_ctg.fa &

1. 筛选保留每个查询序列的最优5个结果
   1. BLAST 默认按比对质量由高到低输出结果，因此第一条记录即为最优比对。
   2. 基因组组装去污染

awk 'BEGIN {OFS="\t"} {if (count[$1] < 5) {print $0; count[$1]++}}' blast_nt.csv > top5_blast_nt.csv

1. 提取ribosomal（核糖体）、mitochondrion（线粒体）、chloroplast（叶绿体）和细菌（bacteria）所在行的内容

grep 'ribosomal' top5_blast_nt.csv | cut -f 1 | sort | uniq > ribosomal.id
grep 'mitochondrion' top5_blast_nt.csv | cut -f 1 | sort | uniq > mitochondrion.id
grep 'chloroplast' top5_blast_nt.csv | cut -f 1 | sort | uniq > chloroplast.id
grep 'bacteria' top5_blast_nt.csv | cut -f 1 | sort | uniq > bacteria.id

1. 将比对到的核糖体和线粒体的序列id合并到delete.id，作为要去掉的序列

cat mitochondrion.id ribosomal.id > delete.id

1. 删除delete.id文档中的序列

seqkit grep -f delete.id -v pulchrum_asm_third.hic.p_ctg.fa > pulchrum_asm_third.hic.p_ctg_filtered.fa

1. 统计fa的bp数量

seqkit stats *.fa