2.3 Blast Practice
徐英泽 2024011267
返回引导页
0. 对于序列MSTRSVSSSSYRRMFGGPGTASRPSSSRSYVTTSTRTYSLGSALRPSTSRSLYASSPGGVYATRSSAVRL
1. 请使用网页版的 blastp, 将上面的蛋白序列只与 mouse protein database 进行比对， 设置输出结果最多保留10个， E 值最大为 0.5。将操作过程和结果截图，并解释一下 E value和 P value 的实际意义
-
E value：在数据库搜索的背景下，E-value 是与得分 S 相关联的预期匹配次数，表示在数据库搜索中，预期会随机出现多少个得分 ≥ S 的不同比对。与比对分数成指数关系，范围是0到无穷大。E-value 回答的是”预期会随机出现多少个这样的匹配？”
-
P value：在序列比对的背景下，在泊松分布下，发现分数大于等于本条比对分数的概率。P-value 是与得分 S 相关联的概率，表示随机获得一个至少等于 S 的得分的概率。P-value 回答的是”有多大的概率会随机出现这样好的匹配？”
2. 请使用 Bash 脚本编程：将上面的蛋白序列随机打乱生成10个， 然后对这10个序列两两之间进行 blast 比对，输出并解释结果。（请上传bash脚本，注意做好重要code的注释；同时上传一个结果文件用来示例程序输出的结果以及对这些结果的解释。）
#!/bin/bash
#初始蛋白序列
origin_protein="MSTRSVSSSSYRRMFGGPGTASRPSSSRSYVTTSTRTYSLGSALRPSTSRSLYASSPGGVYATRSSAVRL"
#定义函数shuffle
shuffle_seq(){
local seq="$1" #原始序列
echo "$seq" | grep -o . | shuf | tr -d '\n' #分行，打乱，重排
}
#数组储存
sequences=()
sequence_names=()
for i in {1..10}; do
#调用函数，打乱序列
shuffled=$(shuffle_seq "$origin_protein")
sequences+=("$shuffled")
sequence_names+=("shuffled_$i")
done
#结果输出文件
results_file="/home/test/blast/output/result_blast.txt"
> "results_file"
#循环比对
for i in {1..10}; do
for j in $(seq $((i+1)) 10);do
seq_i="${sequences[$((i-1))]}"
seq_j="${sequences[$((j-1))]}"
name_i="${sequence_names[$((i-1))]}"
name_j="${sequence_names[$((j-1))]}"
echo "正在比对: $name_i vs $name_j"
blastp -query <(echo -e ">query\n$seq_i") \
-subject <(echo -e">subject\n$seq_j")\
-outfmt 6\
-out "$results_file.tmp"
#将临时结果导入到目标文件中
if [ -s "$results_file.tmp" ]; then
echo "# $seq_i vs $seq_j" >> "$results_file"
cat "$results_file.tmp" >> "$results_file"
echo "" >> "$results_file"
fi
rm -f "$results_file.tmp"
done
done
echo 0
结果：results_file.txt 文件：
# STVRSSTGSVTSRSAPTRSLPMAAGRGSLRSGATTPMGSPRRTRYFSGRYLSSAYYSSSTSLYVSSSRVS vs RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV
query unnamed 54.545 11 5 0 3 13 14 24 3.9 10.4
query unnamed 40.000 10 6 0 50 59 41 50 4.7 10.4
# STVRSSTGSVTSRSAPTRSLPMAAGRGSLRSGATTPMGSPRRTRYFSGRYLSSAYYSSSTSLYVSSSRVS vs GSFAMSGTYGSSTTRYGSPARTSSSYTRSRSRPVLLSRSSSRTYPSTSLGVSMSAVRGSRYATRVSASLP
query unnamed 71.429 7 2 0 34 40 22 28 1.6 11.5
# STVRSSTGSVTSRSAPTRSLPMAAGRGSLRSGATTPMGSPRRTRYFSGRYLSSAYYSSSTSLYVSSSRVS vs GSSASASSYGASRSPLGYRSGSRRTPTSTTYFTSSLRGSAYRRVSMPSAVRVSSSTRLSRTVMTLYSPSG
query unnamed 41.176 17 10 0 39 55 34 50 0.27 13.9
# STVRSSTGSVTSRSAPTRSLPMAAGRGSLRSGATTPMGSPRRTRYFSGRYLSSAYYSSSTSLYVSSSRVS vs TGMRVMRSTTTTGSSGLSPSSFPSYLARTGRRSSPVSSSRSATGSSSYRTRSSYPAVVRRGALSYLSYAS
query unnamed 66.667 6 2 0 34 39 18 23 9.5 9.2
# STVRSSTGSVTSRSAPTRSLPMAAGRGSLRSGATTPMGSPRRTRYFSGRYLSSAYYSSSTSLYVSSSRVS vs SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV
query unnamed 80.000 5 1 0 36 40 37 41 0.43 13.1
query unnamed 58.333 12 5 0 8 19 64 75 4.0 10.4
# STVRSSTGSVTSRSAPTRSLPMAAGRGSLRSGATTPMGSPRRTRYFSGRYLSSAYYSSSTSLYVSSSRVS vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 64.286 14 3 1 25 38 43 54 0.008 17.7
# RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV vs GSFAMSGTYGSSTTRYGSPARTSSSYTRSRSRPVLLSRSSSRTYPSTSLGVSMSAVRGSRYATRVSASLP
query unnamed 45.455 11 6 0 22 32 8 18 7.6 9.6
# RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV vs GSSASASSYGASRSPLGYRSGSRRTPTSTTYFTSSLRGSAYRRVSMPSAVRVSSSTRLSRTVMTLYSPSG
query unnamed 40.000 25 12 1 17 38 26 50 2.2 11.2
# RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV vs SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV
query unnamed 52.941 17 5 1 26 42 54 67 0.44 13.1
# RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV vs RSSRTSSGVFRMYSAGRSSVYTSGRYTLRLPSGVSASLASPYMRARPRSSTSGRSTATVTSPTYSSSSGL
query unnamed 40.625 32 16 1 34 62 42 73 0.17 14.2
query unnamed 51.852 27 12 1 5 30 52 78 0.31 13.5
# RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 85.714 7 1 0 63 69 21 27 0.12 14.6
query unnamed 50.000 12 6 0 21 32 26 37 4.0 10.4
# RFARVRSRGSASTSRRGYRTSTPSSSSSSGSYMGLSAYSSAPGRPMSGTVTRYLSLVRSSTYLSPTSATV vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 44.444 18 8 1 31 48 18 33 1.6 11.5
# TSRGRRSFGPVTSVMYRTALRRYSPTSSSGSSARPYASMYVSVAYSSSSRSSPGTLSATTLLTRGSGSSR vs GSSASASSYGASRSPLGYRSGSRRTPTSTTYFTSSLRGSAYRRVSMPSAVRVSSSTRLSRTVMTLYSPSG
query unnamed 33.333 48 31 1 16 63 27 73 0.61 12.7
# TSRGRRSFGPVTSVMYRTALRRYSPTSSSGSSARPYASMYVSVAYSSSSRSSPGTLSATTLLTRGSGSSR vs SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV
query unnamed 40.000 25 15 0 18 42 45 69 0.27 13.9
# TSRGRRSFGPVTSVMYRTALRRYSPTSSSGSSARPYASMYVSVAYSSSSRSSPGTLSATTLLTRGSGSSR vs RSSRTSSGVFRMYSAGRSSVYTSGRYTLRLPSGVSASLASPYMRARPRSSTSGRSTATVTSPTYSSSSGL
query unnamed 100.000 5 0 0 45 49 73 77 2.9 10.8
# TSRGRRSFGPVTSVMYRTALRRYSPTSSSGSSARPYASMYVSVAYSSSSRSSPGTLSATTLLTRGSGSSR vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 38.889 18 11 0 11 28 60 77 1.1 11.9
# TSRGRRSFGPVTSVMYRTALRRYSPTSSSGSSARPYASMYVSVAYSSSSRSSPGTLSATTLLTRGSGSSR vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 42.105 19 8 1 11 29 16 31 0.41 13.1
query unnamed 40.741 27 15 1 28 54 8 33 2.0 11.2
# GSFAMSGTYGSSTTRYGSPARTSSSYTRSRSRPVLLSRSSSRTYPSTSLGVSMSAVRGSRYATRVSASLP vs SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV
query unnamed 83.333 6 1 0 59 64 48 53 0.21 13.9
# GSFAMSGTYGSSTTRYGSPARTSSSYTRSRSRPVLLSRSSSRTYPSTSLGVSMSAVRGSRYATRVSASLP vs RSSRTSSGVFRMYSAGRSSVYTSGRYTLRLPSGVSASLASPYMRARPRSSTSGRSTATVTSPTYSSSSGL
query unnamed 100.000 5 0 0 65 69 43 47 1.1 11.9
query unnamed 50.000 20 10 0 6 25 58 77 1.5 11.5
# GSSASASSYGASRSPLGYRSGSRRTPTSTTYFTSSLRGSAYRRVSMPSAVRVSSSTRLSRTVMTLYSPSG vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 55.000 20 5 1 2 17 35 54 0.38 13.1
query unnamed 50.000 10 5 0 59 68 14 23 6.0 10.0
# GSSASASSYGASRSPLGYRSGSRRTPTSTTYFTSSLRGSAYRRVSMPSAVRVSSSTRLSRTVMTLYSPSG vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 46.667 15 8 0 42 56 59 73 0.20 13.9
# TGMRVMRSTTTTGSSGLSPSSFPSYLARTGRRSSPVSSSRSATGSSSYRTRSSYPAVVRRGALSYLSYAS vs SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV
query unnamed 83.333 6 1 0 48 53 50 55 1.7 11.5
# TGMRVMRSTTTTGSSGLSPSSFPSYLARTGRRSSPVSSSRSATGSSSYRTRSSYPAVVRRGALSYLSYAS vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 56.250 16 6 1 41 56 25 39 2.6 11.2
query unnamed 80.000 5 1 0 45 49 62 66 5.9 10.0
# TGMRVMRSTTTTGSSGLSPSSFPSYLARTGRRSSPVSSSRSATGSSSYRTRSSYPAVVRRGALSYLSYAS vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 77.778 9 2 0 39 47 28 36 0.20 13.9
query unnamed 80.000 5 1 0 11 15 31 35 2.5 11.2
# SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV vs RSSRTSSGVFRMYSAGRSSVYTSGRYTLRLPSGVSASLASPYMRARPRSSTSGRSTATVTSPTYSSSSGL
query unnamed 57.143 14 6 0 12 25 20 33 0.16 14.2
# SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 66.667 6 2 0 5 10 29 34 1.4 11.9
query unnamed 33.333 21 14 0 15 35 11 31 5.1 10.0
# SGRSSVTFSRRRTYSSGSRSSLTSGSTPAGSPALSTVYTRYATRSSSGAMYMAYGSVPLSSLRTRSPSRV vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 66.667 6 2 0 44 49 29 34 1.6 11.5
# RSSRTSSGVFRMYSAGRSSVYTSGRYTLRLPSGVSASLASPYMRARPRSSTSGRSTATVTSPTYSSSSGL vs SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR
query unnamed 50.000 10 5 0 53 62 43 52 2.0 11.2
# RSSRTSSGVFRMYSAGRSSVYTSGRYTLRLPSGVSASLASPYMRARPRSSTSGRSTATVTSPTYSSSSGL vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 53.333 15 7 0 38 52 61 75 0.053 15.4
# SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 45.455 11 6 0 8 18 21 31 3.8 10.4
query unnamed 40.000 5 3 0 31 35 55 59 9.7 9.2
结果分析：
# SSPGSRTYARTVSPTSATPSLAFSRSSYSASLYGRGSSGGTSPLGRVLRSASSMSYRVSRTMYTVTSSRR vs GSRRRSMVSYMYRVTYSPSRSTTGASSTLAVTGYSLGSSSSSGPSPSFGRRLASPYVTRRASSTASTSLR
query unnamed 45.455 11 6 0 8 18 21 31 3.8 10.4
query unnamed 40.000 5 3 0 31 35 55 59 9.7 9.2
以该结果为例，在这两个序列中，比对到了两对相似片段，产生两行输出。 每行输出的前两项为query和subject的名称，第三项为一致的百分比，第四、五、六项分别为比对长度、错配（mismatches）个数和间隙个数（gap openings）。 第7-10项分别为比对区域在query中的起点、终点和在subject中的起点、终点。第11项为E-value，用于评估比对结果可靠性。最后一项为bit score，即归一化处理后的比对得分。
3. 解释blast 中除了动态规划还利用了什么方法来提高速度，为什么可以提高速度。
-
在运行搜索前，BLAST会对整个数据库（或查询序列，取决于具体实现版本，如blastn是对数据库建索引）建立索引。BLAST不直接对整条序列进行矩阵计算，而是先在查询序列和数据库序列之间寻找精确匹配的短片段。 BLAST通过构建哈希表（Hash Table）或后缀数组索引，可以在预期时间内找到所有匹配的“种子”，将问题从“全矩阵计算”转化为“极少数的局部区域计算”。
-
BLAST在比对过程中设置了许多经验性的阈值。 传统的动态规划必须计算出最优解才能停止。BLAST允许在发现当前局部比对无法达到显著性（如期望值大于0.001）时，提前终止该区域的延伸。以牺牲灵敏度为代价，换取速度提升。
-
双命中（Two-Hit）策略。它要求在同一条对角线上，两个非重叠的种子距离在 A 个碱基以内。只有当两个种子同时出现时，才开始进行延伸。 随机匹配（噪声）往往是孤立的、单点的。要求两个独立的种子在附近出现，极大地降低了假阳性延伸的概率。
4. 我们常见的PAM250有如下图所示的两种（一种对称、一种不对称），请阅读一下 “Symmetry of the PAM matrices” ,然后总结和解释一下这两种（对称和不对称）PAM250不一样的原因及其在应用上的不同。
-
在 PAM 的构建过程中，替换概率矩阵一开始不一定是对称的，因为从氨基酸 A 突变为 B 与从 B 突变为 A 在进化上未必对称（比如考虑突变率差异、密码子偏好、选择压力等）。但在实际用于序列比对打分时，通常会对矩阵进行对称化处理，因为比对算法（如 Smith-Waterman、BLAST）使用的打分矩阵要求对称，才能保证比对结果与序列顺序无关。因此常见的 PAM250 打分矩阵是对称的。
- 而不对称的PAM250则会反映突变的真实情况。
- 对称的PAM250：广泛用于蛋白质序列比对（BLAST、FASTA、全局/局部比对），用于评估同源性。
- 不对称的PAM250：可能用于系统发育分析中构建有根树，或用于模拟定向进化过程。