Фаста чтение файла Python

Я читаю файл FASTA, который имеет такой формат:

> Г |31563518| исй |NP_852610.1| ассоциированные с микротрубочками белки 1A/1B легкой цепи 3A изоформы b [Homo sapiens]
MKMRFFSSPCGKAAVDPADRCKEVQQIRDQHPSKIPVIIERYKGEKQLPVLDKTKFLVPDHVNMSELVKIIRRRLQLNPTQAFFLLVNQHSMVSVSTPIADIYEQEKDEDGFLYMVYASQETFGF 

Я должен прочитать файл и затем рассчитать расстояние JC (Для пары последовательностей расстояние JC равно -3/4 * ln(1 - 4/3 * p), где p - это доля сайтов, которые отличаются между пара)

Я создал его скелет, но не знаю, как сделать все остальное. После чтения и вычисления расстояния JukesCantor я должен записать его в новый выходной файл, и он должен быть в таблице. Любая помощь, которую я могу получить, очень ценится! спасибо, новый для python и fasta файлов

def readData():
    filename = input("Enter the name of the FASTA file: ")
    infile = open(filename, "r")

def CalculateJC(x,y):
    if x == y:
        return 0
    else:
        return 1 # temporary*

def calcDists(seqs):
    output = []
    for seq1 in seqs:
        newrow = []
        for seq2 in seqs:
            dist = calculateJS(seq1,seq2)
            newrow.append(dist)
        output.append(newrow)
        list(enumerate(seasons))
    return output


def outputDists(distMat):
    pass

def main():
    seqs = readData()
    distMat = calcDists(seqs)
    outputDists(distMat)



if__name__ == "__main__":
    main()

1 ответ

Вы задаете слишком много вопросов одновременно! Сосредоточиться на одном.

Чтение и запись файлов FASTA покрыты в BioPython (как предлагается в комментариях).

Я заметил, что вы еще не вычислили расстояние JC, поэтому, возможно, именно здесь вам нужна помощь. Вот что я придумал:

import itertools, math

def computeJC(seq1, seq2):
    equal = 0
    for base1, base2 in itertools.izip(seq1, seq2):
        equal += (base1 == base2)
    p = equal / float(len(seq1))     
    return -3/4 * math.log(1 - 4/3 * p)  

Трюк itertools.izip объясняется здесь: Как я могу перебирать два списка параллельно Этот фрагмент кода будет работать с любой строкой, и просмотр остановится, когда seq1 или seq2 достигнут конца.

Кто-то еще может придумать "Pythonic one-liner", но сначала попытайтесь понять мой подход. Это позволяет избежать ловушек, в которые попал ваш код: вложенные циклы, ненужное ветвление, рост списка времени выполнения, код спагетти и многие другие. Наслаждайтесь!

Другие вопросы по тегам