0

我有一组 520 个流感序列，我已经对其进行了多序列比对，并计算了成对单位矩阵。如果我想添加另一个序列，我必须重新对齐所有内容，并重新计算整个 PWI 矩阵。是否有任何程序可以用来将这个其他序列“附加”到比对中，并且只计算每个其他序列的 PWI？

一个简单的例子如下。我有一个 2x2 对齐，具有以下分数。

     SeqA SeqB
SeqA 1.00 0.98
SeqB 0.98 1.00

无需重新运行完全对齐，而仅针对所有其他序列运行“SeqC”，我想得到以下矩阵：

     SeqA SeqB SeqC
SeqA 1.00 0.98 0.99
SeqB 0.98 1.00 0.97
SeqC 0.99 0.97 1.00

我正在使用 BioPython 包，Python 是我的首选语言，但如果需要，我也可以使用 Java。

[我在这里声明我是从 BioStars 交叉发布的，以防万一这里有不在 BioStars 上的专家。BioStars 的帖子是： http: //www.biostars.org/p/77607/，但内容完全相同。]

4

1 回答 1

1

如果您的主要问题是重新运行对齐所需的时间（重新计算 PWI 矩阵的计算成本应该很低），那么MUSCLE有能力做您正在寻找的事情，这个过程通常被称为“轮廓-轮廓对齐” ”。

轮廓-轮廓对齐

传递-profile标志时，对齐将“相互重新对齐，保持输入列完整并在需要时插入间隙列。”：

如果您有两个现有的相关序列比对，您可以使用 MUSCLE 的 –profile 选项来比对这两个序列。典型用法是：

   muscle -profile -in1 one.afa -in2 two.afa -out both.afa

在 Biopython 中实现

Biopython有一个MUSCLE 的包装器，但我发现它更容易用于subprocess调用 MUSCLE，然后将结果解析回MultipleSeqAlignment：

# Do profile-profile alignment (one sequence to many aligned)
seq_fn = "influenza_seq.fasta"
aligned_fn = "520_influenza_seqs.afasta"
cmd = ['muscle', '-clwstrict', '-profile', '-in1', seq_fn, '-in2', aligned_fn]
aligner = subprocess.Popen(cmd, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
stdout, stderr = aligner.communicate()

# Get resulting alignment (MultipleSeqAlignment)
alignment =  AlignIO.read(StringIO(stdout), "clustal",
                          alphabet=Alphabet.ProteinAlphabet())
于 2013-08-07T16:50:58.547 回答