mVISTA:在线程序展示叶绿体基因组相似性小实例

import os
import sys
import argparse
parser = argparse.ArgumentParser(description='This script was used to download gb or fasta file of cp genome from NCBI nucleotides database')
parser.add_argument('-f','--format',help='Please input file format you want to download',required=True)
parser.add_argument('-a','--accession',help='file name contain accession number of cp genome you want to download',required=True)
args = parser.parse_args()


fr = open(args.accession,'r')
acc = []
for line in fr:
    acc.append(line.strip())
print("The number of sequence will be downloaded is: ",len(acc))

from Bio import Entrez
from Bio import SeqIO
Entrez.email = "" #在这里添加邮箱

fold_name = "downladed_seq_"+args.format
os.mkdir(fold_name)
os.chdir(fold_name)

for line in acc:
    print(line + " will be downloaded!")
    fw = open(line+"."+args.format,'w')
    cp = Entrez.efetch(db='nucleotide',id=[line],rettype=args.format)
    seq = SeqIO.read(cp,args.format)
    SeqIO.write(seq,fw,args.format)
    fw.close()
    print(line + " is downloaded!")

print(len(os.listdir()),"sequence have been downloaded")

#下载fasta格式
python download_gb_or_fa_from_NCBI_cp_genome_database_1.py -f fasta -a accession_numbers.txt
# 下载genbank格式
python download_gb_or_fa_from_NCBI_cp_genome_database_1.py -f gb -a accession_numbers.txt
# 序列号放到文本文件中每行一个

import argparse

from Bio import SeqIO

parser = argparse.ArgumentParser(description = 'This script was used to get mVISTA annotation file used to cp genome comparative analysis')
parser.add_argument('-i','--genbank',help="Please input genBank format file", required = True)
args = parser.parse_args()

fw = open(args.genbank + "_mVISTA_annotation","w")

for rec in SeqIO.parse(args.genbank,"gb"):
    for feature in rec.features:
        if feature.type == "gene":
            for part in feature.location.parts:
                if part.strand == -1:
                    start_location = str(part.start)
                    end_location = str(part.end)
                    gene_name = feature.qualifiers["gene"][0]
                    fw.write("<\t%s\t%s\t%s\n"%(start_location,end_location,gene_name))
                else:
                    start_location = str(part.start)
                    end_location = str(part.end)
                    gene_name = feature.qualifiers["gene"][0]
                    fw.write("<\t%s\t%s\t%s\n"%(start_location,end_location,gene_name))
        elif feature.type == "CDS":
            for part in feature.location.parts:
                if part.strand == -1:
                    start_location = str(part.start)
                    end_location = str(part.end)
                    fw.write("%s\t%s\t%s\n"%(start_location,end_location,"exon"))
                else:
                    start_location = str(part.start)
                    end_location = str(part.end)
                    fw.write("%s\t%s\t%s\n"%(start_location,end_location,"exon"))
        elif feature.type == "tRNA" or feature.type == "rRNA":
            for part in feature.location.parts:
                start_location = str(feature.location.start)
                end_location = str(feature.location.end)
                fw.write("%s\t%s\t%s\n"%(start_location,end_location,"utr"))
        else:
            print("%s %s"%(feature.type,"is not needed"))

print("The process is done!")

fw.close()

# 使用方法
python get_mVISTA_annotation_file_from_genbank_1.py -i NC_031163.gb