问题标签 [sequence-alignment]

对于以下代码，当我执行代码时出现错误，我在下面列出。我想知道是否有人能给我任何关于如何将 CA 原子附加到 tag_atoms/tagged_atoms 列表中的见解，我将使用这些列表进行对齐。并强调我会忽略的代码编写方式中的任何潜在缺陷。我是 python 新手，所以任何见解都会很棒而且很有帮助。

这是下面的错误消息：

我正在编写代码来查找两个序列之间的局部比对。这是我一直在研究的一个最小的工作示例：

以下代码与输出一起运行

但是应该可以得到 '6' 的分数，找到 5 个对齐旁边的 'CC'，如下所示：

关于发生了什么的任何想法？

请帮助我理解这个公式（如果有人想知道，那就是 Needleman-Wunsch 算法），我应该编写一个使用递归的代码，但我不明白怎么做，我已经有了完整的动态编写的版本，所以我知道算法是如何工作的，但我不知道如何用递归/堆栈来实现它。我愿意接受任何建议。

这是公式：

概括：

我从一些关于对齐算法的细节开始，最后，我问了我的问题。如果您了解对齐算法，请通过开头。

考虑我们有两个字符串，例如：

有一些算法，例如：Smith-Waterman或Needleman-Wunsch，它们对齐这两个序列并创建一个矩阵。查看以下部分中的结果：

问题：

我的问题很简单，但也许答案并不像看起来那么容易。我想将一组字符用作单个字符，例如：[A0][C0][A1][B1]。但是在这些算法中，我们必须使用单个字符。我们怎样才能做到这一点？

PS考虑我们有这个序列：#read #write #add #write。然后我将其转换为类似的东西：#read to A .... #write to B.... #add to C. 然后我的序列变为：ABCB。但是我有很多不同的词以 . 开头#。并且 ASCII 表不足以转换所有这些。然后我需要更多的字符。[A0] ... [Z9]唯一的方法是对每个单词使用类似的东西。或使用数字。

PS：此链接中存在 Smith-Waterman 的一些示例代码

PS：还有另一个帖子想要这样的东西，但我想要的是不同的。在这个问题中，我们有一组以 a 开头并以[结尾的字符]。并且不需要使用语义 like eeis equal to i。

我正在寻找一种算法来解决，或者至少是以下问题的正确名称：

我有一组B位串。该算法应该找到一个最小的（定义为“设置最少的位”）位串S，使得：

对于B中的所有b，存在一个移位N（在ℤ</a> 中）使得.(S << N) & b == b

在我看来，这就像某种多序列比对问题。

如果我们可以为B中的每个b _i找到每个N _i ( i = 1 .. | B |)，那么看起来S只是所有 ( b _i >> N _i ) 的按位或。

我的直觉是，第一步是从B中删除每个b ，其中在B中存在另一个位串c和一些移位M使得. 下一步是什么？b & (c << M) == b

我有一个 600 Mb 的 FASTA 文件，其中包含来自 12 个物种的许多对齐块，我想将它们拆分为较小的 FASTA 文件，每个文件包含一个块及其对应的对齐

我有一个看起来像这样的 sed 脚本：

这适用于小规模，但对于 600Mb 的大文件，它需要很长时间，大约 2 天。我不认为这是我正在运行的计算机的问题。

有谁知道如何使它更快？

输入的 Fasta 文件如下所示：

dm3.chr3R(-):17092630-17092781|sequence_index=0|block_index=4|species=dm3|dm3_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droGri2.scaffold_15074(-):2610183-2610334|sequence_index=0|block_index=4|species=droGri2| droGri2_4_0 GGCGGAGATCAAGAATCGTGTTGGGCCGCCGTCGAGCGCCACCGATAACGCTAGCAAAGTGAAAATCGATCAGGGACGCCCAGTGGAAAACAATAGATCTGGTTGCTGCTAAATAA-CTCTGATTGTGAATCATTATTTTATTATACAATTa droMoj3.scaffold_6540(+):33866311-33866462|sequence_index=0|block_index=4|species=droMoj3|droMoj3_4_0 TGCCGAGATTAAGAATCGTGTCGGTCCGCCGTCCAGCGCAACCGACAATGCAAGCAAAGTGAAAATCGATCAGGGACGTCCAGTGGAGAACACCAGATCTGGTTGCTGCTGAATAA-CTCTGATTGTGAATCATTATTTTATTatacaatta droVir3.scaffold_12822(+):1248119-1248270|sequence_index=0|block_index=4|species=droVir3|droVir3_4_0 GGCCGAGATTAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACCGATAATGCTAGCAAAGTGAAAATCGATCAGGGTCGTCCAGTGGAGAACACCAAATCTGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droWil1.scaffold_181130(-):16071336-16071488|sequence_index=0|block_index=4|species=droWil1|droWil1_4_0 GGCCGAGATTAAGAATCGTGTTGGGCCGCCGTCCAGCGCCACTGATAATGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAATACCAAATCCGGTTGCTGCTGAATAAACTCTGATTGTGAATCATTATTTTATTATACAATTA droPer1.super_19(- ):1310088-1310239|sequence_index=0|block_index=4|species=droPer1|droPer1_4_0 GGCTGAGATCAAGAATCGCGTCGGACCGCCGTCCAGCGCCACCGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAAACCCAATTCTGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattattattattatacaatta(-):dp4.chr25593491-5593642|sequence_index=0|block_index=4|species=dp4|dp4_4_0 GGCTGAGATCAAGAATCGCGTCGGACCGCCGTCCAGCGCCACCGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAAGCCCAATTCTGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droAna3.scaffold_13340(-):3754154-3754305|sequence_index=0|block_index=4|species=droAna3|droAna3_4_0 GGCCGAGATCAAGAATCGCGTCGGGCCACCGTCCAGCGCCACCGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAGATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattataaaatta droEre2.scaffold_4770 (+):4567591-4567742|sequence_index=0|block_index=4|species=droEre2|droEre2_4_0 GGCCGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACCGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattatttatttatttattttatacaatta.R5883047-5883198|sequence_index=0|block_index=4|species=droYak2|droYak2_4_0 GGCCGAGATCAAGAATCGCGTCGGGCCGCCATCCAGCGCCACCGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droSec1.super_38(+):36432-36583|sequence_index=0|block_index=4|species=droSec1|droSec1_4_0 GGCGGAGATCAAGAATCGCGTCGGTCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droSim1.chr3R (+):4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattattattattatacaatta4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattattattatatacaatta4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattattattatatacaattasuper_38(+):36432-36583|sequence_index=0|block_index=4|species=droSec1|droSec1_4_0 GGCGGAGATCAAGAATCGCGTCGGTCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droSim1.chr3R(+):4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA -CTCTGATTGTGAATCattattttattatacaattasuper_38(+):36432-36583|sequence_index=0|block_index=4|species=droSec1|droSec1_4_0 GGCGGAGATCAAGAATCGCGTCGGTCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattttattatacaatta droSim1.chr3R(+):4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA -CTCTGATTGTGAATCattattttattatacaatta4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattattattatatacaatta4366350-4366501|sequence_index=0|block_index=4|species=droSim1|droSim1_4_0 GGCGGAGATCAAGAATCGCGTCGGGCCGCCGTCCAGCGCCACTGACAACGCTAGCAAAGTGAAAATCGATCAAGGACGTCCAGTGGAAAACACCAAATCCGGTTGCTGCTGAATAA-CTCTGATTGTGAATCattattattattatatacaatta

dm3.chr3R(-):17092781-17092867|sequence_index=0|block_index=5|species=dm3|dm3_5_0 GAGTACGCCGCCCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTTGAGCAGGCCTTCATGACGATGGC droSim1.chr3R(+):4366264-4366350|sequence_index=0|block_index=5|species=droSim1|droSim1_5_0 GAGTACGCCGCCCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTTGAGCAGGCCTTTATGACGATGGC droSec1.super_38(+):36346-36432|sequence_index=0|block_index=5|species=droSec1|droSec1_5_0 GAGTACGCCGCCCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTTGAGCAGGCCTTCATGACGATGGC droYak2.chr3R(-):5883198-5883284|sequence_index=0|block_index=5|species=droYak2|droYak2_5_0 GAGTACGCCGCCCAGTTAGGCATTCCATTCCTTGAAACATCGGCCAAGAGCGCCACCAACGTGGAGCAGGCCTTCATGACGATGGC droEre2.scaffold_4770(+):4567505-4567591|sequence_index=0|block_index=5|species=droEre2|droEre2_5_0 GAGTACGCCGCCCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTGGAGCAGGCCTTCATGACGATGGC droAna3.scaffold_13340(+):20375068-20375148|sequence_index=0|block_index=5|species=droAna3|droAna3_5_0 ------GCCGAAAACTTCGACATGCCCTTCTTCGAGGTCTCTTGCAAGTCAAACATCAATATTGAAGATGCGTTTCTTTCCCTGGC dp4 .chr2(-):5593642-5593728|sequence_index=0|block_index=5|species=dp4|dp4_5_0 GAGTATGCAGCTCAGTTAGGCATTCCATTTCTTGAAACTTCGGCCAAGAGCGCCACGAACGTGGAGCAGGCCTTCATGACGATGGC droPer1.super_19(-):1310239-1310325|sequence_index=0|block_index=5|species=droPer1|droPer1_5_0 GAGTATGCAGCTCAGTTAGGCATTCCATTTCTTGAAACTTCGGCCAAGAGCGCCACGAACGTGGAGCAGGCCTTCATGACGATGGC droWil1 .scaffold_181130(-):16071488-16071574|sequence_index=0|block_index=5|species=droWil1|droWil1_5_0 GAATATGCGGCTCAGTTAGGCATTCCATTCCTTGAAACTTCGGCAAAGAGTGCCACCAATGTGGAGCAGGCCTTTATGACGATGGC droVir3.scaffold_12822(+):1248033-1248119|sequence_index=0|block_index=5|species=droVir3|droVir3_5_0 GAGTACGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTGGAGCAGGCATTTATGACGATGGC droMoj3.scaffold_6540(+): 33866225-33866311|sequence_index=0|block_index=5|species=droMoj3|droMoj3_5_0 GAGTATGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAATGTAGAGCAGGCATTCATGACGATGGC droGri2.scaffold_15074(-):2610334-2610420|sequence_index=0|block_index=5|species=droGri2|droGri2_5_0 GAGTACGCAAATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCGAAGAGTGCCACCAATGTGGAACAGGCATTCATGACGATGGC1248033-1248119|sequence_index=0|block_index=5|species=droVir3|droVir3_5_0 GAGTACGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTGGAGCAGGCATTTATGACGATGGC droMoj3.scaffold_6540(+):33866225-33866311|sequence_index=0|block_index=5|species=droMoj3|droMoj3_5_0 GAGTATGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAATGTAGAGCAGGCATTCATGACGATGGC droGri2.scaffold_15074(-): 2610334-2610420|sequence_index=0|block_index=5|species=droGri2|droGri2_5_0 GAGTACGCAAATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCGAAGAGTGCCACCAATGTGGAACAGGCATTCATGACGATGGC1248033-1248119|sequence_index=0|block_index=5|species=droVir3|droVir3_5_0 GAGTACGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAACGTGGAGCAGGCATTTATGACGATGGC droMoj3.scaffold_6540(+):33866225-33866311|sequence_index=0|block_index=5|species=droMoj3|droMoj3_5_0 GAGTATGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAATGTAGAGCAGGCATTCATGACGATGGC droGri2.scaffold_15074(-): 2610334-2610420|sequence_index=0|block_index=5|species=droGri2|droGri2_5_0 GAGTACGCAAATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCGAAGAGTGCCACCAATGTGGAACAGGCATTCATGACGATGGC33866225-33866311|sequence_index=0|block_index=5|species=droMoj3|droMoj3_5_0 GAGTATGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAATGTAGAGCAGGCATTCATGACGATGGC droGri2.scaffold_15074(-):2610334-2610420|sequence_index=0|block_index=5|species=droGri2|droGri2_5_0 GAGTACGCAAATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCGAAGAGTGCCACCAATGTGGAACAGGCATTCATGACGATGGC33866225-33866311|sequence_index=0|block_index=5|species=droMoj3|droMoj3_5_0 GAGTATGCACATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCCAAGAGCGCCACCAATGTAGAGCAGGCATTCATGACGATGGC droGri2.scaffold_15074(-):2610334-2610420|sequence_index=0|block_index=5|species=droGri2|droGri2_5_0 GAGTACGCAAATCAGTTAGGCATTCCATTCCTTGAAACTTCGGCGAAGAGTGCCACCAATGTGGAACAGGCATTCATGACGATGGC

我确信这是一件容易做的事情，但我的生物信息学经验非常有限。

我有许多 -100,000- FASTA 文件，其中包含相同 12 个物种的不同基因的比对。每个文件看起来像这样：

并非所有文件都以相同的方式排序，我希望它们都是。如果这更容易，它们可以按字母顺序排序，只要所有文件均等排序，它们的排序方式并不重要。按字母顺序应该是这样的：

任何自动执行此操作的脚本将不胜感激。

编辑：我一直在使用使用 sed 的 shell 脚本，该脚本有效但有问题。当文件数量不是很大时，它可以工作，但在这种特殊情况下，它会创建具有不同名称的重复文件。脚本内容如下：

未排序的文件称为 bloque.fasta.trim，此脚本会创建一堆名为 b7.fasta.trim$ 的文件，这些文件应该会创建一个 b7。每个块的文件。文件。问题是有时它会复制一个文件但名称不同。我相信大多数情况下都会有一种更简单的方法，不会出现重复错误。

我有大约 13K 序列和 120 个碱基，我想比较它们以找到诸如保守区域、它们之间的平均差异或非常不同的异常值之类的东西。

问题是，有了这么多的序列，我尝试过的事情是不可行的。

那么有没有人做过类似的事情并且可以给我一些提示如何实现它？或者也许只是我应该寻找的一些提示？

我已经看过这里和其他论坛，但找不到我的问题的答案。

我想为目标富集测序方法设计诱饵，并使用 MarkerMiner 搜索来自四个不同基因组的直系同源基因座的输出，其中
拟南芥作为参考。这些输出比对是每个拟南芥注释基因的单独 Fasta 文件，我的数据集中的序列与之对齐。
我已经运行了一个脚本来过滤掉那些被我的四个输入数据集中的至少两个支持为直系同源的基因座。

然而，现在，我被难住了。
我的比对很差，因为输入数据主要是 RNAseq，而参考也包含内含子。所以它看起来像这样：

AT01G1234567 ATCGATCGATGCGCGCTAGCTGAATCGATCGGATCGCGGTAGCTGGAGCTAGSTCGGATCGC MyData1
CGATGCGCGC-----------CGGATCGCGG---------------CGGATCGC
MyData2
CGCTGCGCGC------------GGATAGCGG--- ------------CGGATCCC

为了有效地设计诱饵，我现在需要从文件中提取所有对齐的部分，以便最终得到单独的文件；或在文件中单独对齐；对于在 MyData 和参考序列之间对齐的部分，排除了所有有间隙的部分。这些 fasta 文件大约有 1300 个，因此无法手动进行。我在 python 和 Linux 命令行工具方面有一些编程经验，但是我完全不知道如何去做。我会很感激一个提示，关于我可以使用什么样的工具或者我需要想出什么样的算法。谢谢你。干杯

我想知道如何量化 Needleman-Wunsch 算法的结果（通常用于比对核苷酸/蛋白质序列）。

考虑一些固定的评分方案和两个不同长度S1和的序列S2。假设我们通过蛮力计算所有可能的对齐方式，S1并且得分最高的对齐方式有一个 score 。当然，这比 Needleman-Wunsch 方法的复杂性要高得多。S2x

当使用 Needleman-Wunsch 算法查找序列比对时，假设它有一个 score y。

考虑r是通过 Needleman-Wunsch 为两个随机序列R1和生成的分数R2。

x相比如何y？是否y总是大于r已知同源性的两个序列？

总的来说，我确实知道我们使用 Needleman-Wunsch 算法来显着加快序列比对（相对于蛮力方法），但不了解随之而来的准确性成本（如果有的话）。我尝试阅读了原始论文（Needleman & Wunsch，1970），但仍然留下了这个问题。