文章目录

下载
检查blast是否添加至环境路径
使用blast
- 合并fasta文件
- 检查重名
- blast建库
- blast检索目标序列
输出结果
- 提取序列

下载

通过NCBI FTP下载本地blast软件：
https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/
对于Windows用户，直接下载exe可执行文件即可：

然后按照正常软件安装程序安装即可，可自定义安装路径

提示：以下是本篇文章正文内容，下面案例可供参考

检查blast是否添加至环境路径

一般来说新版的blast都会自动添加到环境目录中，这一步一般不用检查，只是为了以防万一，毕竟有些老版本的blast是不能自动添加到环境目录中的，需要手动添加一下

键盘同时按WIN+R，输入sysdm.cpl，

（用的不是自己的电脑，忽略我的英文操作系统）

然后按回车，跳转，选择上方第三个选项卡高级，点下方的环境变量：

检查你安装blast目录下的bin文件夹是否有在用户PATH中：

若环境变量中没有blast目录，则需要手动添加一下：
选中PATH一栏，点击Edit，在弹出窗口中选择NEW，然后在下方栏目中复制粘贴你bin文件的路径，再点OK即可添加

使用blast

首先通过NCBI或Uniprot等数据库下载到自己需要的基因组/蛋白组文件，你的目标序列会在这些基因组/蛋白组文件中进行查询

合并fasta文件

下载完成后，如果你有多个物种的基因组\蛋白组需要检索，则需要将这些物种的基因组fasta文件合并为一个，这里使用了一个简单的python脚本来合并：
使用方法，将需要合并的所有文件放在一个文件夹中，文件夹中不要放置其他文件，同时在文件夹中放入这个python脚本，运行即可。

import osoutname="database.fasta"
with open(outname,'wt') as outfile:for filename in os.listdir():dot=filename.find(".")if (filename[dot+1:]=='txt' or filename[dot+1:]=='fasta') and filename!=outname:species=filename[:dot].replace(" ","_")with open(filename,'rt') as infile:for line in infile:if line[0]=='>':line=line.split('|')outfile.write(">{}_{}\n".format(line[1],species))else:outfile.write(line)print("Finish {}".format(filename))print("Finish all!")

以Uniprot下载得到的蛋白组fasta为例，其序列名都很长，因此我在合并文件时还加入了截短序列名的功能，这里是直接截取了序列名中两个|中间的序列号，同时加入fasta的文件名，因此建议将蛋白组文件名直接命名为其物种名，且尽量只有属名和种加词两个单词
修改前：

修改后：

检查重名

注意：检索的基因组中不能有重复的序列名，否则在建库的过程中会报错
这里同样使用了一个简单的python脚本来检查是否有重复的序列名：

with open("database.fasta",'rt') as infile:namelist=[]repeat=Falsefor line in infile:if line[0]=='>':line=line.strip("")[1:]i=0while i<len(namelist):if namelist[i]==line:repeat=Truei=-1print("序列{}有重复".format(line))breakelif namelist[i]<line:breakelse:i+=1if i>0:namelist.insert(i,line)if not repeat:print("该文件中没有重复的序列")

blast建库

此处建议新建一个空文件夹，将你合并好的蛋白组fasta和要查询的序列fasta放到这个文件夹中，在地址栏输入cmd，回车

在命令行界面中，输入如下命令进行建库

makeblastdb.exe -in database.fasta -parse_seqids -hash_index -dbtype prot

其中database.fasta需要更换为你自己的蛋白组文件名
若你构建的是基因组文件，则需要更换数据类型为核酸：

makeblastdb.exe -in database.fasta -parse_seqids -hash_index -dbtype nucl

建库完成后，会出现一大堆文件：

blast检索目标序列

如果你需要检索的fasta序列已经放置在数据库同一目录下，则可以在上面的cmd中继续输入：

blastp.exe -task blastp -query input.fasta -db database.fasta -out blast_result.txt -outfmt 6 -evalue 1e-5

如果需要比对的是核酸，则使用blastn

blastn.exe -task blastn -query input.fasta -db database.fasta -out blast_result.txt -outfmt 6 -evalue 1e-5

参数说明：
-query 需要检索的序列文件名
-db 需要在哪个数据库中检索
-out 输出文件名
-outfmt 输出格式，一般选6即可，有0-18可选
-evalue E值阈值，高于这个E值的序列不输出到结果中，默认值为10，建议设到10^-5以下

除此之外，还可以加入num_thread等参数，表示同时使用多少个线程来计算。其他更多参数可以使用
-help来查看：

blastp.exe -help

可选输出格式：

输出结果

如果选择的输出格式是6，那么得到的结果会是制表符分隔的表格，使用excel打开即可：
得到的一般会是这样一个有13列的表格

这几列对应的表头是：
query id：查询序列ID标识
subject id：比对上的目标序列ID标识
identity：序列比对的一致性百分比
alignment length：符合比对的比对区域的长度
mismatches：比对区域的错配数
gap openings：比对区域的gap数目
query start：比对区域在查询序列(Query id)上的起始位点
query end：比对区域在查询序列(Query id)上的终止位点
subject start：比对区域在目标序列(Subject id)上的起始位点
subject end：比对区域在目标序列(Subject id)上的终止位点
e-value：表明在随机的情况下，其它序列与目标序列相似度要大于S值的可能性，e-value的分值越低越好
bit score：比对结果的bit score值

后续就可以根据这些blast结果从中找到自己感兴趣的序列来进行分析

提取序列

由于需要根据序列号从基因组/蛋白组文件中逐个找回序列的工作十分繁琐，这里同样使用一个简单的python脚本来进行批量查找：

import osprint("请在本文件夹下input.txt文件中输入需要查找的序列,按回车分隔")
print("待查找的序列需与该程序置于同一文件夹中,且后缀名只能为.txt或.fasta")
try:with open('input.txt', 'rt', encoding='utf-8') as input_sequence_file:search_name=[n.strip() for n in input_sequence_file.read().strip().split("\n")]
except IOError:print("本文件夹下没有input.txt文件")outfilename="result.fasta"
with open(outfilename,'wt') as outfile:for search_file_name in os.listdir():filetype=search_file_name[search_file_name.find(".")+1:]if (filetype!='txt' and filetype!='fasta') or search_file_name=='input.txt' or search_file_name==outfilename:continuewith open(search_file_name,'rt') as search_file:key=''for line in search_file:if line[0]=='>':if len(key)>0:outfile.write(">{}\n{}\n".format(key,seqdict[key]))key=''if line[1:].strip() in search_name:key=line[1:].strip()seqdict={key:''}elif len(line.strip())>0 and len(key)>0:seqdict[key]+=line.strip()

使用方法：将要查找的序列名放在input.txt中，每行一个；然后将数据库fasta文件和input.txt放在同一个文件夹目录下，不要放置其他多余文件，将该脚本放置在该文件夹中，运行，得到的结果会输出在output.txt中

【实用教程】本地blast使用及简单python脚本辅助相关推荐

简单python脚本实例-python脚本例子
广告关闭 2017年12月,云+社区对外发布,从最开始的技术博客到现在拥有多个社区产品.未来,我们一起乘风破浪,创造无限可能. 简单python脚本:注:1.python 变量不需要$2.for里面的 ...
【实用教程】一文学会安装Python环境，Python语言给我整起来
文章目录 Python 环境搭建 Python下载 Python安装 Unix & Linux 平台安装 Python: Window 平台安装 Python: MAC 平台安装 Python ...
简单python脚本实例-Python简单实现阴阳师挂机脚本
参考大佬们的博客,凑出了阴阳师的简单挂机脚本,对各个博主进行感谢,鞠躬.jpg 实现原理: 对各个按钮等需要点击的地方进行截图,然后对整个屏幕截屏,通过模板匹配得到截图在截屏中的坐标位置,进行点击实 ...
【python教程入门学习】通过运行python脚本来更改Windows背景
通过运行python脚本来更改Windows背景在我们开始之前,一定要注意这篇文章只针对Windows用户!对于那些使用Windows的人来说,这是一个有趣的想法. 如果您想使用python更改桌面 ...
简单python脚本实例-简单了解python调用其他脚本方法实例
1.用python调用python脚本 #!/usr/local/bin/python3.7 import time import os count = 0 str = ('python b.py') ...
csv转为矩阵 python_超简单python脚本实现从resfinder+blast下机数据到csv抗性基因矩阵...
今天还是在做张阿姨的项目,还是那100多个菌,她说要blast一下抗性基因,那就比咯,最出名与权威的的细菌抗性基因库自然肯定是resfinder了. for i in $(cat list); do ...
简单python脚本实例-对Python实现简单的API接口实例讲解
get方法代码实现 # coding:utf-8 import json from urlparse import parse_qs from wsgiref.simple_server impor ...
简单python脚本实例-终于晓得python入门脚本实例
Python面向对象跟Java的面向对象大同小异,这篇我们简单介绍一下Python面向对象.以下是小编为你整理的python入门脚本实例类(class)是通俗的说就是事物的属性和行为的抽象.下面我们 ...
简单python脚本实例-30个Python 小例子，帮你快速上手Python
可能很多人学编程有个误区,总想着从最基础的原理开始,看了一大堆书一堆视频还是不会,兴趣也早没了.而自己动手写代码是最好的开始方式,不管TM三七二十一,直接就是干,无论是看书还是看视频,确保身边有台电脑 ...
简单python脚本实例-超实用的 30 段 Python 案例
Python是目前最流行的语言之一,它在数据科学.机器学习.web开发.脚本编写.自动化方面被许多人广泛使用. 它的简单和易用性造就了它如此流行的原因. 如果你正在阅读本文,那么你或多或少已经使用过P ...

【实用教程】本地blast使用及简单python脚本辅助