Biopython gives ValueError: Sequences must all be the same length even though sequences are of the same length

Question

我正在尝试通过从我的数据创建一个 .phy 文件来创建一个系统发育树。

我有一个数据框

ndf=

ESV trunc

1 esv1 TACGTAGGTG...

2 esv2 TACGGAGGGT...

3 esv3 TACGGGGGG...

7 esv7 TACGTAGGGT...

我检查了“t运行c”列元素的长度：

length_checker = np.vectorize(len)

arr_len = length_checker(ndf['trunc'])

结果 arr_len 为所有元素提供相同的长度 (=253)。

我将此数据框保存为 .phy 文件，如下所示：

23 253

esv1 TACGTAGGTGGCGAGCGTTATCCGGAATTATTGGGCGTAAAGCGCGCGCAGGCGGTTTCTTAAGTCTGATGTGAAAGCCCACGGCTCAACCGTGGAGGGTCATTGGAAACTGGGAAACTTGAGTGCAGAAGAGGAAAGCGGAATTCCACGTGTAGCGGTGAAATGCGTAGAGATGTGGAGGAACACCAGTGGCGAAGGCGGCTTTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCAAACAGG

esv2 TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG

esv3 TACGGGGGGGGCAAGCGTTGTTCGGAATTACTGGGCGTAAAGGGCGCGTAGGCGGCCAGACCAAGTCGAGTGTGAAATTGCAGGGCTTAACTTTGCAGGGTCGCTCGATACTGGTCGGCTAGAGTGTGGAAGAGGGTACTGGAATTCCCGGTGTAGCGGTGAAATGCGTAGATATCGGGAGGAACACCAGCGGCGAAGGCGGGTACCTGGGCCAACACTGACGCTGAGGCGCGAAAGCTAGGGGAGCAAACAG

这类似于 this tutorial 中使用的文件。

但是，当我运行命令 aln = AlignIO.read('msa.phy', 'phylip')

我收到“ValueError：序列的长度必须相同”

我不知道为什么会出现这个问题或如何解决它。非常感谢任何帮助！

谢谢

Answer 1

一般来说，phylip 是不同程序之间系统发育学中最繁琐的格式。有严格的 phylip 格式和宽松的 phylip 格式等... t 不容易知道正在使用哪个分隔符，一个 space 字符 and/or 一个回车 return.

我认为您似乎在分类单元名称（即序列标签）和序列名称之间留下了 space，即

2. esv2

Phylip 格式正在观察标签和序列数据之间的 space。在本例中，序列长度为 3bp。使用“.”通常也不是一个好主意。整数似乎不表示行号。

另一个问题是您 could/should 尝试将序列与标签保持在同一行并删除回车 return，即。

esv2 TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG

有时一个回车符return确实有效（这可能是宽松的phylip格式），传统格式使用space字符“”。我总是保持统一数量的 spaces 以保持对齐...不确定是否需要这样做。

注意如果您的类群名称超过 10 个字符，您将需要宽松的 phylip 格式，这种格式在任何情况下都是一个好主意。

最终的解决方案是转换为fasta，导入为fasta，然后转换为phylip。如果这一切都失败了...... post 后面还有更多 trouble-shooting

Fasta 格式删除了“23 254”header，然后每个序列看起来像这样，

>esv2
TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG

在“>esv2”和序列之间总是有一个回车return。此外，“>”总是出现在没有任何空间的标签（分类单元名称）的前缀。您可以简单地通过 reg-ex 或 Python 中的 "re" 进行转换。使用 perl one-liner 它将是 s/^([az]+[0-9]+)/>/g 类型代码。我很确定他们会成为一个在线网站来执行此操作。

然后您只需将导入命令中的 "phylip" 替换为 "fasta"。导入后，你要求 BioPython 转换成你想要的任何格式，它应该没有任何问题。

Answer 2

首先，请阅读How to make good reproducible pandas examples. In the future please provide a minimal reproducibl example的答案。

其次，Michael G 绝对正确，phylip 是一种语法非常特殊的格式。

下面的代码将允许您从 Pandas 数据框生成系统发育树。

首先进行一些导入，让我们重新创建您的数据框。

import pandas as pd
from Bio import Phylo
from Bio.Phylo.TreeConstruction import DistanceCalculator
from Bio.Phylo.TreeConstruction import DistanceTreeConstructor
from Bio import AlignIO

data = {'ESV' : ['esv1', 'esv2', 'esv3'],
        'trunc': ['TACGTAGGTGGCGAGCGTTATCCGGAATTATTGGGCGTAAAGCGCGCGCAGGCGGTTTCTTAAGTCTGATGTGAAAGCCCACGGCTCAACCGTGGAGGGTCATTGGAAACTGGGAAACTTGAGTGCAGAAGAGGAAAGCGGAATTCCACGTGTAGCGGTGAAATGCGTAGAGATGTGGAGGAACACCAGTGGCGAAGGCGGCTTTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCAAACAGG',
                 'TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG',
                 'TACGGGGGGGGCAAGCGTTGTTCGGAATTACTGGGCGTAAAGGGCGCGTAGGCGGCCAGACCAAGTCGAGTGTGAAATTGCAGGGCTTAACTTTGCAGGGTCGCTCGATACTGGTCGGCTAGAGTGTGGAAGAGGGTACTGGAATTCCCGGTGTAGCGGTGAAATGCGTAGATATCGGGAGGAACACCAGCGGCGAAGGCGGGTACCTGGGCCAACACTGACGCTGAGGCGCGAAAGCTAGGGGAGCAAACAG']

}

ndf = pd.DataFrame.from_dict(data)
print(ndf)

输出：

    ESV                                              trunc
0  esv1  TACGTAGGTGGCGAGCGTTATCCGGAATTATTGGGCGTAAAGCGCG...
1  esv2  TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTG...
2  esv3  TACGGGGGGGGCAAGCGTTGTTCGGAATTACTGGGCGTAAAGGGCG...

接下来，以正确的格式写入phylip文件。

with open("test.phy", 'w') as f:
    f.write("{:10} {}\n".format(ndf.shape[0], ndf.trunc.str.len()[0]))
    for row in ndf.iterrows():
        f.write("{:10} {}\n".format(*row[1].to_list()))

test.phy 的输出：

         3 253
esv1       TACGTAGGTGGCGAGCGTTATCCGGAATTATTGGGCGTAAAGCGCGCGCAGGCGGTTTCTTAAGTCTGATGTGAAAGCCCACGGCTCAACCGTGGAGGGTCATTGGAAACTGGGAAACTTGAGTGCAGAAGAGGAAAGCGGAATTCCACGTGTAGCGGTGAAATGCGTAGAGATGTGGAGGAACACCAGTGGCGAAGGCGGCTTTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCAAACAGG
esv2       TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG
esv3       TACGGGGGGGGCAAGCGTTGTTCGGAATTACTGGGCGTAAAGGGCGCGTAGGCGGCCAGACCAAGTCGAGTGTGAAATTGCAGGGCTTAACTTTGCAGGGTCGCTCGATACTGGTCGGCTAGAGTGTGGAAGAGGGTACTGGAATTCCCGGTGTAGCGGTGAAATGCGTAGATATCGGGAGGAACACCAGCGGCGAAGGCGGGTACCTGGGCCAACACTGACGCTGAGGCGCGAAAGCTAGGGGAGCAAACAG

现在我们可以开始创建系统发育树了。

# Read the sequences and align
aln = AlignIO.read('test.phy', 'phylip')
print(aln)

输出：

SingleLetterAlphabet() alignment with 3 rows and 253 columns
TACGTAGGTGGCGAGCGTTATCCGGAATTATTGGGCGTAAAGCG...AGG esv1
TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGG...AGG esv2
TACGGGGGGGGCAAGCGTTGTTCGGAATTACTGGGCGTAAAGGG...CAG esv3

计算距离矩阵：

calculator = DistanceCalculator('identity')
dm = calculator.get_distance(aln)
print(dm)

输出：

esv1    0
esv2    0.3003952569169961  0
esv3    0.6086956521739131  0.6245059288537549  0

使用UPGMA算法构建系统发育树并以ascii格式绘制树

constructor = DistanceTreeConstructor()
tree = constructor.upgma(dm)
Phylo.draw_ascii(tree)

输出：

  ________________________________________________________________________ esv3
_|
 |                                     ___________________________________ esv2
 |____________________________________|
                                      |___________________________________ esv1

或者画一个漂亮的树图：

Phylo.draw(tree)

输出：

Biopython gives ValueError: Sequences must all be the same length even though sequences are of the same length

Biopython gives ValueError: Sequences must all be the same length even though sequences are of the same length

python

biopython

dataframe

phylogeny