Удаление дубликатов PCR из файла Fastq, содержащего уникальные молекулярные идентификаторы

Я пытаюсь отредактировать файл Fastq, содержащий геномные данные и уникальные молекулярные идентификаторы, фланкирующие каждую последовательность.

Пример первых двух прочитанных приведен ниже:

1 @HISEQ:230:C6G45ANXX:3:1101:1395:2141 1:N:0:ACAGTGGTTGAACCTT 2 TGACGGCACTTTCTCTTCCCAACCACGTGGCTGCAGACTTCTTGCTCTCAAGTTGTCCTGACATGCTCTGAGAGCACACACAACATACATACAACACCTGGATCTGTGAATTAATTACTGCCTAGG 3 + 4 BB//<<BFBFFF<FFFFBBB<<<F/FBBB<FF/B<FFFFFFFFFFFFFFBFFFBFB/FBFFB//F//B<FFF</</BF<BBBFFFFF//B<FBFF/77F/B/BF7/FF/<BF/7FFFFBBF//B7B 5 @HISEQ:230:C6G45ANXX:3:1101:1498:2162 1:N:0:ACAGTGGTTGAACCTT 6 TGACGGCACTTTCTCTTCCCAACCACGTGGCTGCAGACTTCTTGCTCTCAAGTTGTCCTGACATGCTCTGAGAGCACACACAACATACATACAACACCTGGATCTGTGAATTAATTACTGCCTAGG 7 + 8 BBB<B<F<FFFFFFFBFFFFFFBFFFFBFF/F<FFFFBBFFFFFFFFFFBFB/BFFFFFFFFFFFBFFB/<<<FFFFFFFFFFFFFFBFFFF################################## 

Эти строки объясняются ниже:

 1 Information 2 Sequence 3 + 4 Quality Scoring 5 Information 6 Sequence 7 + 8 Quality Scoring 

Мне нужен выходной файл, в котором все точные повторы заданной последовательности (и ее соответствующей информации) были удалены. То есть, мне нужно удалить те блоки из 4 строк, в которых второй уже появился в файле.

Таким образом, в приведенном выше примере, поскольку последовательность соответствует строкам 2 и 6, выходной файл должен содержать строки 1,2,3 и 4, но не 5,6,7 и 8.

Результат выходного файла:

 1 @HISEQ:230:C6G45ANXX:3:1101:1395:2141 1:N:0:ACAGTGGTTGAACCTT 2 TGACGGCACTTTCTCTTCCCAACCACGTGGCTGCAGACTTCTTGCTCTCAAGTTGTCCTGACATGCTCTGAGAGCACACACAACATACATACAACACCTGGATCTGTGAATTAATTACTGCCTAGG 3 + 4 BB//<<BFBFFF<FFFFBBB<<<F/FBBB<FF/B<FFFFFFFFFFFFFFBFFFBFB/FBFFB//F//B<FFF</</BF<BBBFFFFF//B<FBFF/77F/B/BF7/FF/<BF/7FFFFBBF//B7B 

2 Solutions collect form web for “Удаление дубликатов PCR из файла Fastq, содержащего уникальные молекулярные идентификаторы”

Кажется, это идеальный случай, когда мы дважды прокручиваем файл: сначала для вычисления дубликатов, а затем для печати соответствующих строк:

 awk 'FNR==NR { if (FNR%4==2) { a[$2]++ if (a[$2]>1) b[int(FNR/4)]=1 } next} b[int(FNR/4)]==0' file file 

Ключевым моментом здесь является воспроизведение с файлом 4K + 2 в файле и отслеживание того, какие из них появились до сих пор. Если они это сделают, мы сохраняем K (от 4K+2 ), так что в следующем цикле файла мы избегаем тех строк, которые находятся в форме 4K+0/1/2/3 .

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

Контрольная работа

 $ awk 'FNR==NR {if (FNR%4==2) {a[$2]++; if (a[$2]>1) b[int(FNR/4)]=1} next} b[int(FNR/4)]==0' aa @HISEQ:230:C6G45ANXX:3:1101:1395:2141 1:N:0:ACAGTGGTTGAACCTT TGACGGCACTTTCTCTTCCCAACCACGTGGCTGCAGACTTCTTGCTCTCAAGTTGTCCTGACATGCTCTGAGAGCACACACAACATACATACAACACCTGGATCTGTGAATTAATTACTGCCTAGG + BBB<B<F<FFFFFFFBFFFFFFBFFFFBFF/F<FFFFBBFFFFFFFFFFBFB/BFFFFFFFFFFFBFFB/<<<FFFFFFFFFFFFFFBFFFF################################## 

Я думаю, что вы отключены по одной строке в заказе FASTQ. В вашем примере:

 1 Junk 2 Junk 3 Junk 4 Information +-> 5 Sequence | these four lines constitute a single record 6 + | 7 Quality Scoring +-> 9 Information 10 Sequence 11 + 

Таким образом, строки 1-3 (Junk) на самом деле являются верхними 3 строками предыдущей записи, а строки 9-11 – это верхние 3 следующей записи.

В любом случае, я предлагаю вам использовать SeqIO BioPython для анализа вашего файла FASTQ и удаления дубликатов.

http://biopython.org/wiki/SeqIO

Один из основных способов сделать это:

 from Bio import SeqIO from Bio.SeqUtils.CheckSum import seguid detected = [] unique = [] for rec in SeqIO.parse(open('inputfile.fastq', 'rU'), 'fastq'): cksum = seguid(rec.seq) if cksum not in detected: unique.append(rec) detected.append(cksum) SeqIO.write(unique, open('deduplicated.fastq','w'), 'fastq') 

Это считывает каждую запись и вычисляет контрольную сумму последовательности для хранения в списке. Каждая последующая запись добавляется только к списку для вывода («unique»), если в его последовательности нет контрольной суммы, уже обнаруженной.

  • Каковы различия между Perl, Python, AWK и sed?
  • Как использовать регулярное выражение python для замены с помощью захваченной группы?
  • Как случайным образом удалить несколько строк из большого файла?
  • Как конкатенировать идентификатор, указанный в двух строках?
  • Рекурсивно перемещать каталог и заменять вызовы функций
  • Создание нескольких CSV-файлов из данных в файле csv
  • Как получить плоский XML, чтобы внешние объекты были объединены на верхний уровень
  • извлекать каждый файл последовательности в виде отдельного файла
  • Python - лучший язык программирования в мире.