Подсчет "диграммы" или пары нуклеотидов в петле for - PullRequest
1 голос
/ 11 марта 2020

Мне нужно написать функцию, которая берет файл fasta и подсчитывает в файле диграммы (AT, CG, TT, CC, et c).

Мой для l oop в настоящее время считывает файл построчно и производит подсчет для этой строки. Затем он перезапускает счет в следующей строке. (Все это организовано в словарь)

Я хочу вести подсчет каждой строки, чтобы получить счетчик для всего файла, а не только для отдельных строк.

Это мой код, который я пытаюсь исправить:

dinucleotides = ['AA','AT','AG','AC',
                 'TA','TT','TG','TC',
                 'GA','GT','GG','GC',
                 'CA','CT','CG','CT']

all_counts = {}

with open('short.fasta', 'r') as dna_file:
    dna_file.readline()

    for line in dna_file:
        my_line = line.strip()

        for pairs in dinucleotides:
            count = my_line.count(pairs)
            all_counts[pairs] = count

Спасибо!

Ответы [ 3 ]

0 голосов
/ 11 марта 2020

Добавьте к последнему счету,

all_counts[pairs] = all_counts.get(pairs, 0) + count
0 голосов
/ 11 марта 2020

Одна идея состоит в том, чтобы инициализировать Python диктовку, отображающую каждую 2грамму в ноль, и увеличивать ее в соответствии с каждой строкой. Здесь я предполагаю, что файл FASQ содержит только базы в «ATG C». Кроме того, повторение каждой возможной пары для каждой строки требует 16 проходов по каждой строке. Этого можно избежать, пропустив каждую строку по очереди и сохранив каждую пару. Возможно, следующим образом:

import random

def random_dnukes(lines=1000, cols=40):
  return [''.join(random.choices('ATGC', k=cols)) for _ in range(lines)]

# e.g.
# ['TGACTCGTCAAAGGTACGTTAATCCTTGGGCAGTTACGGG',
#  'ATTGTTCAATCGAACGTTCGCTACTCGACTCGCGCCCCCT',
#  'TCCCGTGGGACAGGTTCCCAATTGACCGGACGCCGGACTC',
#  'TCGTCGTGCCCCGACATTGCTTCACGGCGGTGCGCGCTGG',
#  'GGTCCGGTCTAGGCGATCCCTAATAGTCAAGCACCGATTA',
#  'CCGAGCCTTGTGTATACTCTGTAAACACTTCTTCCCATAC',
#  'CGGATAGCAGCTAGTGGTTCCCGCAGTACAGGATGACCAA',
#  'CTCGGACGAGAAATCAGGCCAACCTCCACTGGCGACAGAA',
#  'TCTGACCTGCAGTGCAGTCCAGTTATAGTGGAACACCAGC',
#  'GTCAGCCCTTATCCGTTAGCCCAGGTGCCTCAATAGGAGG']

fake_file_iterator = iter(random_dnukes(1000, 40))

from collections import defaultdict
total_counts = defaultdict(int)

for line in fake_file_iterator:
  line = line.strip()
  for i in range(len(line) - 1):
    total_counts[line[i:i+2]] += 1

for k, v in total_counts.items():
  print(k, v)

В результате

GC 2497
CC 2382
CG 2444
GT 2422
TT 2508
TA 2373
AC 2466
GG 2408
TG 2473
CA 2462
AA 2412
CT 2448
AG 2454
GA 2470
TC 2400
AT 2381
0 голосов
/ 11 марта 2020

Вы можете использовать collections.defaultdict с int как default_factory.
И изменить all_counts[pairs] = count на all_counts[pairs] += count.

from collections import defaultdict


dinucleotides = ['AA','AT','AG','AC',
                 'TA','TT','TG','TC',
                 'GA','GT','GG','GC',
                 'CA','CT','CG','CT']

all_counts = defaultdict(int)

with open('short.fasta', 'r') as dna_file:
    dna_file.readline()

    for line in dna_file:
        my_line = line.strip()

        for pairs in dinucleotides:
            count = my_line.count(pairs)
            all_counts[pairs] += count

Или использовать метод dict.setdefault.

...
all_counts = {}
...
            all_counts.setdefault(pairs, 0) += count
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...