Фильтровать последовательности с более чем 8 одинаковыми последовательными нуклеотидами в файле fastq? - PullRequest
Новая
       61

Фильтровать последовательности с более чем 8 одинаковыми последовательными нуклеотидами в файле fastq?

3 голосов
/ 03 ноября 2019

Я хочу отфильтровать мои последовательности, содержащие более 8 одинаковых последовательных нуклеотидов, таких как "GGGGGGGG", "CCCCCCCC" и т. Д. В моих fastq файлах.

Как мне это сделать?

Ответы [ 3 ]

1 голос
/ 05 ноября 2019

Я решил использовать приведенные ниже коды в R и решил свою проблему. 

library(ShortRead)

fq <- FastqFile("/Users/path/to/file")
reads_fq <- readFastq(fq)


trimmed_fq <- reads_fq[grep("GGGGGGGG|TTTTTTTTT|AAAAAAAAA|CCCCCCCCC", 
                       sread(reads_fq), invert = TRUE)]

writeFastq(trimmed_fq, "new_name_for_fq.fastq", compress = FALSE)
1 голос
/ 03 ноября 2019
  1. Быстрый и неправильный способ, который может быть достаточно близок: grep -E -B1 -A2 'A{8}|C{8}|G{8}|T{8}' yourfile.fastq. Это пропустит блоки, где 8-мер разделен на две строки (например, первая строка заканчивается AAAA, а вторая начинается с AAAA). Также предполагается, что выходные данные имеют блоки по 4 строки в каждом.
  2. Правильный способ: напишите небольшую программу (на Python или на ваш язык), которая буферизует один блок FASTQ (например, 4 строки) и проверяет, чтоконкатенация последовательности предыдущего (буферизованного) блока и последовательности текущего блока не имеет 8-членного числа, как указано выше. Если это так, выведите буферизованный блок.
0 голосов
/ 04 ноября 2019

Вы можете использовать пакет Python биотит для него (https://www.biotite -python.org ).

Допустим, у вас есть следующий файл FASTQ: 

@Read:01
CCCAAGGGCCCCCCCCCACTGCGATCACCTGGTTGCTGCCGGGAAAGGAGACCCAGGAGGTGAAACGGACTGGTGAATTG
CGGGGGTAGATATGGCGGGTGACACAAAAACATATAATCGGGCC
+
.+.+:'-FEAC-4'4CA-3-5#/4+?*G@?,<)<E&5(*82C9FH4G315F*DF8-4%F"9?H5535F7%?7@+6!FDC&
+4=4+,@2A)8!1B#,HA18)1*D1A-.HGAED%?-G10'6>:2
@Read:02
AACACTACTTCGCTGTCGCCAAAGGTTGGTGTAGGTCGGACTTCGAATTATCGATACTAGTTAGTAGTACGTCGCGTGGC
GTCAGCTCGTATGCTCTCAGAACAGGGAGAACTAGCACCGTAAGTAACCTAGCTCCCAAC
+
6%9,@'4A0&%.19,1E)E?!9/$.#?(!H2?+E"")?6:=F&FE91-*&',,;;$&?@2A"F.$1)%'"CB?5$<.F/$
7055E>#+/650B6H<8+A%$!A=0>?'@",8:#5%18&+3>'8:28+:5F0);E9<=,+

Это скрипт, который должен выполнять работу: 

import biotite.sequence.io.fastq as fastq
import biotite.sequence as seq


# 'GGGGGGGG', 'CCCCCCCC', etc.
consecutive_nucs = [seq.NucleotideSequence(nuc * 8) for nuc in "ACGT"]

fastq_file = fastq.FastqFile("Sanger")
fastq_file.read("example.fastq")
# Iterate over sequence entries in file
for header in fastq_file:
    sequence = fastq_file.get_sequence(header)
    # Iterative over each of the consecutive sequences
    for consecutive_nuc in consecutive_nucs:
        # Find all indices, where a match was found
        matches = seq.find_subsequence(sequence, consecutive_nuc)
        if len(matches) > 0:
            # If any match was found report it
            print(
                f"Found '{consecutive_nuc}' "
                f"in sequence '{header}' at position {matches[0]}"
            )

Это вывод: 

Found 'CCCCCCCC' in sequence 'Read:01' at pos 8
...