BLAST中，与HSP和query_coverage(qcov)相关的关键参数
2018-11-17 13:40阅读：
首先，在BLAST的输出结果中，输入的序列被称为query，database的序列被称为subject（通过makeblastdb建库的序列）
通过BLAST得到的结果是局部的比对（local alignment），而不是全局的比对（global
alignment）。所谓局部，就是一条序列的一部分。如果序列A的一部分（e.g.
几十个碱基，而整条序列有几百个碱基），与序列B的一部分（同样几十个碱基）比对起来，满足一些预先设定的条件（e.g.
identity，e-value等），BLAST就会将其作为一个HSP输出。BLAST结果中的每一行就是一个HSP。
HSP是high scoring
pair的首字母缩写。局部比对的特点就是一条query和一条subject会产生多个HSP。如下面的BLAST结果：
从左到右的14列数据依次是: query id, subject id, identity, alignment length,
query coverage per subject, query coverage per hsp, query coverage
per uniq subject, mismatches, gap opens, query start, query end,
subject start, subject end, evalue.
Query 'scaffold249' 与 subject 'GGLTR3G2'
产生了两个HSP：
第一个HSP是query的第239~519个碱基，与subject的第289~603个碱基；
第二个HSP是query的第28~119个碱基，与subject的第21~112个碱基。
那么在理解了局部比对与BLAST输出结果之间的关系以后，就不难理解BLAST中与HSP和query
coverage相关的参数了。
-max_target_seqs
规定了BLAST结果中每一条query对应的subject的最大数量。注意不是HSP的最大数量。例如，若将其设置为10，那么对每一条query，最多有10条subject与之对应，但可能产生不止10个HSP，因为每一条subject都至少产生一个HSP。上例中的一条subject就产生了两个HSP。BLAST结果中每一行就是一个HSP（或者叫做hit），所以每一query在BLAST结果中都不止10行。
-max_hsps
规定了BLAST结果中每一条query对应的HSP的最大数量。注意与-max_target_seqs区分。由于BLAST结果中每一行就是一个HSP，所以如果将其设置为10，那么每一query在BLAST结果中都只有10行。
此外，在BLAST的标准输出（outfmt=6/7）中，是没有query coverage的相关信息的。但由于是local
alignment，query
coverage的信息对于判断比对的可信度又非常重要，因此我们常常要通过设定outfmt参数的方式来增加输出query
coverage的信息。
设定的方式也非常简单，只需要在outfmt 6 或 7 后指定相应的域就可以了。如：
-outfmt '7 qseqid sseqid pident length qcovs qcovhsp qcovus
mismatch gapopen qstart qend sstart send evalue
bitscore'
指定后，输出结果中就会包含对应的域。
与query coverage有关的域有三个，它们又分别具有不同的含义：
-qcovs, 即 query coverage per subject
-qcovhsp, 即 query coverage per hsp
-qcovus, 即 query coverage per uniq subject
前面我们已经知道了，一条query和一条subject会产生多个不同的HSP，每一个HSP都对应着query/subject不同的部分。qcovs的计算公式是用这对query/subject所有的HSP的query
length之和除以query序列的总长度。query length的计算方法是query end - query start +
1. 也因此同一对query/subject的所有HSP的qcovs都是一样的。
在上述scaffold249/GGLTR3G2的例子中，第10列是query
start，11列是query end。所以两个HSP的query length分别是
519-239+1=281，119-28+1=92，它们的和是373，而query序列的总长度为533，所以qcov=373/533≈70%,
即是第5列的内容。可以看到一对query/subject的所有HSP的qcovs都一样。
而qcovhsp是将每一个HSP都分别计算，用该HSP的query
length除以query序列的总长度。上例中，两个HSP的qcovhsp分别是281/533≈53%, 92/533≈17%,
即是第6列的内容。
最后，qcovus计算所采用的也是一对query/subject所有HSP query
length之和除以query序列的总长度，但其query length之和的计算方法，不再是简单地将所有HSP的query
length相加，而是要去除掉不同HSP的query
length之间的重叠部分。如果不同的HSP之间没有重叠，那么qcovus=qcovs, 如上例第7列
