Bio python дает ValueError: все последовательности должны иметь одинаковую длину, даже если последовательности имеют одинаковую длину - PullRequest
Новая
       68

Bio python дает ValueError: все последовательности должны иметь одинаковую длину, даже если последовательности имеют одинаковую длину

1 голос
/ 18 февраля 2020

Я пытаюсь создать дерево филогенети c, создавая файл .phy из моих данных.

У меня есть фрейм данных

ndf =

ESV trun c

1 esv1 TACGTAGGTG ...

2 esv2 TACGGAGGGT ...

3 esv3 TACGGGGGG ...

7 esv7 TACGTAGGGT ...

Я проверил длину элементов столбца " trun c ":

length_checker = np.vectorize(len)

arr_len = length_checker(ndf['trunc'])

Полученный arr_len дает одинаковую длину (= 253) для всех элементов.

1028 * Я сохранил этот dataframe как .phy файл, который выглядит следующим образом: * * * 1029 1030

23 253

1033

esv1 TACGTAGGTGGCGAGCGTTATCCGGAATTATTGGGCGTAAAGCGCGCGCAGGCGGTTTCTTAAGTCTGATGTGAAAGCCCACGGCTCAACCGTGGAGGGTCATTGGAAACTGGGAAACTTGAGTGCAGAAGAGGAAAGCGGAATTCCACGTGTAGCGGTGAAATGCGTAGAGATGTGGAGGAACACCAGTGGCGAAGGCGGCTTTCTGGTCTGTAACTGACGCTGAGGCGCGAAAGCGTGGGGAGCAAACAGG

esv2 TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG

1042 *

esv3 TACGGGGGGGGCAAGCGTTGTTCGGAATTACTGGGCGTAAAGGGCGCGTAGGCGGCCAGACCAAGTCGAGTGTGAAATTGCAGGGCTTAACTTTGCAGGGTCGCTCGATACTGGTCGGCTAGAGTGTGGAAGAGGGTACTGGAATTCCCGGTGTAGCGGTGAAATGCGTAGATATCGGGAGGAACACCAGCGGCGAAGGCGGGTACCTGGGCCAACACTGACGCTGAGGCGCGAAAGCTAGGGGAGCAAACAG

Это похоже на файл, используемый в этого урока .

Однако, когда я запускаю команду aln = AlignIO.read('msa.phy', 'phylip')

, я получаю «ValueError: все последовательности должны быть одинаковой длины»

Я не знаю, почему я получаю это или как это исправить. Любая помощь с благодарностью!

Спасибо

Ответы [ 2 ]

2 голосов
/ 19 февраля 2020

Как правило, phylip - самый простой формат в филогенетике между различными программами. Существует строгий формат phylip и упрощенный формат phylip и т. Д. c ... t нелегко узнать, какой разделитель используется, пробел и / или возврат каретки.

Я думаю, что вы появляетесь чтобы оставить пробел между именем таксона (то есть меткой последовательности) и именем последовательности, а именно: 

2. esv2

Формат Phylip отслеживает пространство между меткой и данными последовательности. В этом примере последовательность будет иметь длину 3 б.п. Использование "." как правило, не очень хорошая идея, а также. Целое число, по-видимому, не обозначает номер строки.

Другая проблема заключается в том, что вы можете / должны попытаться сохранить последовательность в той же строке, что и метка, и удалить возврат каретки, а именно. 

esv2 TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG

Иногда возврат каретки работает (это может быть упрощенный формат phylip), в традиционном формате используется пробел "". Я всегда сохранял одинаковое количество пробелов, чтобы сохранить выравнивание ... не уверен, если это необходимо.

Примечание если имя таксона превышает 10 символов, вам понадобится смягченный формат phylip, и это Формат в любом случае, как правило, хорошая идея.

Окончательное решение, в котором все остальное терпит неудачу - это преобразовать в fasta, импортировать как fasta, а затем преобразовать в phylip. Если все это не сработает ... отправьте обратно, есть еще проблемы устранения

Формат Fasta удаляет заголовок "23 254", а затем каждая последовательность выглядит следующим образом: 

>esv2
TACGGAGGGTGCAAGCGTTATCCGGATTCACTGGGTTTAAAGGGTGCGTAGGTGGGTTGGTAAGTCAGTGGTGAAATCTCCGAGCTTAACTTGGAAACTGCCATTGATACTATTAATCTTGAATATTGTGGAGGTTAGCGGAATATGTCATGTAGCGGTGAAATGCTTAGAGATGACATAGAACACCAATTGCGAAGGCAGCTGGCTACACATATATTGACACTGAGGCACGAAAGCGTGGGGATCAAACAGG

Там всегда есть возврат каретки между "> esv2" и последовательностью. Кроме того, «>» всегда присутствует в качестве префикса метки (названия таксона) без какой-либо точки. Вы можете просто конвертировать через reg-ex или "re" в Python. При использовании однострочного perl это будет код типа s/^([az]+[0-9]+)/>$1/g. Я почти уверен, что это будет онлайн-сайт, который будет делать это.

Затем вы просто замените "phylip" на "fasta" в вашей команде импорта. После импорта вы просите Bio Python конвертировать в любой формат, который вам нужен, и с ним не должно возникнуть никаких проблем.

1 голос
/ 19 февраля 2020

Сначала прочитайте ответ на Как сделать хорошие воспроизводимые pandas примеры . В будущем, пожалуйста, предоставьте минимальный воспроизводимый пример .

Во-вторых, Michael G абсолютно прав, что phylip - это формат, который очень своеобразен по своему синтаксису.

Приведенный ниже код позволит вам сгенерировать филологетическое дерево c из вашего 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...

Далее запишите файл в правильном формате. 

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

Теперь мы можем Начнем с создания нашего дерева филогенетов c. 

# 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

Постройте дерево филогенет c, используя алгоритм UPGMA, и нарисуйте дерево в ascii 

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

Вывод: 

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

Или составьте хороший график из дерево: 

Phylo.draw(tree)

Выход:

enter image description here

...