CRISPR spacers数据库；CRT和PILER-CR用于MAGs的spacers搜索

CAAS_IFR_zp

已于 2025-03-15 15:23:47 修改

阅读量538

点赞数 3

文章标签：数据库

于 2025-02-16 20:11:33 首次发布

本文链接：https://blog.youkuaiyun.com/m0_53945548/article/details/145626188

版权

iPHoP：病毒宿主预测-优快云博客

之前介绍了这个方法来预测病毒宿主，今天来介绍另一种比较用的多的方法CRISPR比对

CRISPR spacers数据库

Dash 在这可以下载作者搜集的spacers用于后期比对

CRT和PILER-CR

使用 CRT 和 PILERCR 识别 CRISPR 间隔区，合并冗余 CRISPR 阵列，并格式化输出

用的是别人写好的代码（好用就行，这两软件太老了，别人帮忙下好，配置好了）

Metagenomic compendium of 189,680 DNA viruses from the human gut microbiome | Nature Microbiology

MGV/crispr_spacers at master · snayfach/MGV · GitHub

#运行
software="/home/zhongpei/hard_disk_sda2/zhongpei/Software/MGV/crispr_spacers"
out_dir="/home/zhongpei/hard_disk_sda2/zhongpei/database/SPACER_rumen_MAGs"
gunzip -k *
dir_name=$(basename $(pwd))
mkdir ${out_dir}/${dir_name}
for i in *.fa
do
	num=${i%%.fa}
	dir=$(pwd)
	cd ${software}
	identify_crispr.py -i ${dir}/${num}.fa -o ${out_dir}/${dir_name}/${num}_spacer
	merge_crispr.py ${out_dir}/${dir_name}/${num}_spacer/crt ${out_dir}/${dir_name}/${num}_spacer/pilercr ${out_dir}/${dir_name}/${num}_spacer/merged
	cd ${dir}
done
rm *.fa

提取merge.spacer

#!/usr/bin/env python
#########################################################
# filter blastn crisper
# written by PeiZhong in IFR of CAAS

import argparse
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
from Bio import SeqIO

parser = argparse.ArgumentParser(description='Filter BLASTN CRISPR')
parser.add_argument('--in_txt', required=True, help='Input file containing multiple spacers')
parser.add_argument('--out_fa', required=True, help='Output FASTA of these spacers')

args = parser.parse_args()
in_txt = args.in_txt
out_fa = args.out_fa

with open(in_txt, "r") as fa:
    count = 0
    fasta_records = []
    for i in fa.readlines():
        if "spacer_seq" not in i:
            spacer = i.strip("\n").split("\t")
            contig_id = spacer[0] + f"_SPACER{count}"
            spacer_seq = spacer[5].strip()
            if spacer_seq:  # Skip if spacer_seq is empty
                record = SeqRecord(Seq(spacer_seq), id=contig_id, description="")
                fasta_records.append(record)
                count += 1

with open(out_fa, "w") as fasta:
    SeqIO.write(fasta_records, fasta, 'fasta')

print("FASTA file has been generated.")

spacers比对

我们现在有了两组spacers，一组是NAR文章整理的，一组是宿主MAG提取的。

使用Blastn进行比对（viral contigs和spacers）

参考文献的Viral host prediction部分：A metagenomic catalog of the early-life human gut virome | Nature Communications

blastn -query ${fa}/fetal_5.0_95.fa -db ${database} -evalue 0.0000001 -gapopen 10 -gapextend 2 -reward 1 -penalty -1 -word_size 5 -perc_identity 100 -max_target_seqs 10000 -out ${out}/fetal_vOTU_spacerDB.txt -outfmt 6 -num_threads 150