blast 用法汇总
blast 安装：
先在网上找到最新版的blast下载到本地；
this is the download padge: ftp://ftp.ncbi.nih.gov/blast/executables/LATEST/
这里下载的是最近版本：ncbi-blast-2.2.29+-x64-linux.tar.gz
linux下可以通过ftp anonymous 登录，get下载
用filezilla 导入到Linux系统中，找一个目录放下；
先解压：tar zxvf ncbi-blast-2.2.29+-x64-linux.tar.gz
解压后的文件夹ncbi-blast-2.2.29+/bin里面即是可执行的BLAST+程序了
在.barshrc 文件中设置好环境变量即可：
由于2.2.29与之前版本2.2.23有所不同，下面介绍一下2.8.1用法：
blastall
第一步 建库：
/opt/blc/genome/biosoft/blast-2.2.23/bin/formatdb -i ortholog.fa -p T
注：这里用的是安装的绝对路径，
T 表示蛋白库，F表示核酸库；
第二步 运行blastall
eg:
/opt/blc/genome/biosoft/blast-2.2.23/bin/blastall -i query.fasta -d db.seq -o blast.out -p blastn -F F -e 1e-5 -m 8
1．-e参数
-e（value）参数是用来过滤比对较差的结果的，用"-e"参数指定一个实数，blast 会过滤掉期望值大于这个数的比对结果。这样不但简化了结果，还缩短了运行时间和结果占用的空间。
2. -p参数
-p blastp：蛋白序列与蛋白库做比对。
-p blastx：核酸序列对蛋白库的比对。
-p blastn：核酸序列对核酸库的比对。
-p tblastn：蛋白序列对核酸库的比对。
- -F 参数
-F(T/F)参数是用来屏蔽简单重复和低复杂度序列的。如果选“T”，程序在比对过程中会屏蔽query中的简单重复和低复杂度序列；选"F"则不会屏蔽。缺省值为"T"。
4. -m
m8格式如下图：
image.png
举个实例：
image.png
blastn 用法介绍：
第一步：blastn 建index：
makeblastdb -parse_seqids -dbtype nucl -in xxxx.fa
makeblastdb -in /XXX/database/refdata-cellranger-GRCh38-1.2.0/fasta/genome.fa -dbtype nucl -out /database/refdata-cellranger-GRCh38-1.2.0/fasta/genome -parse_seqids
第二步：blastn比对：
常用比对参数：
blastn -db /cygene/database/IMGT_V-QUEST_reference_directory/Homo_sapiens/TR/human_TR_V -query TRAV10.utr.fq2fa.fa -out TRAV10.utr.fq2fa.fa.blastn.out -outfmt 7 -num_threads 20
关于输出格式的参数：
-outfmt <String>
alignment view options:
0 = Pairwise,
1 = Query-anchored showing identities,
2 = Query-anchored no identities,
3 = Flat query-anchored showing identities,
4 = Flat query-anchored no identities,
5 = BLAST XML,
6 = Tabular,
7 = Tabular with comment lines,
8 = Seqalign (Text ASN.1),
9 = Seqalign (Binary ASN.1),
10 = Comma-separated values,
11 = BLAST archive (ASN.1),
12 = Seqalign (JSON),
13 = Multiple-file BLAST JSON,
14 = Multiple-file BLAST XML2,
15 = Single-file BLAST JSON,
16 = Single-file BLAST XML2,
17 = Sequence Alignment/Map (SAM),
18 = Organism Report
关于输出格式，共有19种，我个人常用的是7。
示例如下：会输出每条read的比对信息，如果没有比对上，则 0 hits found
image.png
一般情况下我会设置每条read只输出一种perfect match 的比对结果。在这种条件下统计比对率时就非常简单：需要统计总共有多少条reads，和比对上了多少条reads，然后计算比对率。
le RP1G10E1L2_R1.p2.blastn.out|grep -v "#"|wc 这里统计的是比对上的reads数。
6941468 83297616 628607598
80% pct mapping rate.
le xxx.fa |grep ">" |wc # 这里统计原始fa的总reads数。
比对率就是两数相除的结果。
事实证明这样的方法是错误的！！！
不能用这样的方式计算比对率，因为还有一种情况是一条read比对到了reference 的两个部位。后文给出例子。
线程参数：
-num_threads 20 线程数，根据自身情况设定。
-max_hsps 1 表示只输出一条比对结果，但是有时会输出三条？不懂为什么？
blastn -db /cygene/database/refdata-cellranger-GRCh38-1.2.0/fasta/genome.fa -query $line -out $name.blastn.out -outfmt 7 -num_threads 20 -max_hsps 1
-max_target_seqs 1 只输出一条比对结果，这个参数试了下是没有问题的，只输出一条。（常用）
blastn -db /cygene/database/utr_500bp_index/TRAV.fa -query RP1G10E1L2_R1.p2.fa -out RP1G10E1L2_R1.p2.blastn.out -outfmt 7 -num_threads 20 -max_target_seqs 1
########################
一些经验：
blastn结果中有一种情况是下图这样的，query序列 比对的起始到终止的长度 和 subject序列 比对的起始和终止长度不一致。
image.png
query 比对的长度为76-1=75bp
subject比对的长度为625-551=74bp，那么为什么会造成这样呢？
我将query的序列和subject序列拿出来看了一下。发现query的序列在某个位置上有1bp的插入造成没有对齐。mismatch数量如果不算插入的1bp，那么刚好等于其结果显示的4.
image.png
这个时候如果你需要找到mismatch发生在reference上的位置的时候，那么就不好找了。
blastn对序列reverse的处理：
如果你的序列是reverse的。那么比对一下看看。待补充。
关于设置了-max_target_seqs 1的问题
即使设置了该参数，但是如果数据本身有下面现象也会被输出：
image.png
这种情况就是下图这样的：
image.png
这种情况也会输出两条比对信息。这时统计比对率时就不能单纯用统计行数来表示。
