Arhn - архитектура программирования

эффективное разделение файлов fastq по длине

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

@HISEQ:28:H8P69ADXX:1:1101:1462:2036 1:N:0:CTTGTA
NCCATAAAGTAGAAAGCACT
+
#00<FFFFFFFFFIIFIIFF
@HISEQ:28:H8P69ADXX:1:1101:1419:2156 1:N:0:CTTGTA
TGGAGAGAAAGGCAGTTCCTGA
+
BBBFFFFFFFFFFIIIIIIIII
@HISEQ:28:H8P69ADXX:1:1101:1378:2223 1:N:0:CTTGTA
TCCTGTACTGAGCTGCCCCGA
+
BBBFFFFFFFFFFIIIIIIII
@HISEQ:28:H8P69ADXX:1:1101:1585:2081 1:N:0:CTTGTA
AAACCGTTACCATTACTGAGT
+
BBBFFFFFFFFFFIIIIFIII

Прямо сейчас я использую awk для фильтрации последовательностей определенной длины или в определенном диапазоне:

awk 'BEGIN {OFS = "\n"} {header = $0 ; getline seq ; getline qheader ; getline qseq ; if (length(seq) == 22) {print header, seq, qheader, qseq}}'

Если я хочу иметь выходной файл для каждой длины последовательности, я использую цикл for:

for i in {16..33};
awk -v var=$i 'BEGIN {OFS = "\n"} {header = $0 ; getline seq ; getline qheader ; getline qseq ; if (length(seq) == var) {print header, seq, qheader, qseq}}'
done

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

Может ли кто-нибудь помочь мне найти более эффективное решение, чем мой цикл? Если возможно, решение, в котором мне не нужно указывать диапазон, но тот, который проверяет минимальную и максимальную длину и автоматически разбивает их. Я хотел бы сделать это в awk, но я открыт для всего. Спасибо Бенедикт

11.04.2018

  • Напоминает мне об этом: stackoverflow.com/questions/3194349/ 11.04.2018
  • К сожалению, моя проблема немного в другом. Я хочу разбить файл fastq по определенным критериям (длина последовательности), а не только на определенное количество файлов. Если вы разделите по длине, вы получите файлы, которые будут иметь очень разные размеры, поскольку некоторые длины последовательности встречаются редко, а другие очень много. Но спасибо, что нашли время! 11.04.2018

Ответы:


1

что-то вроде этого?

$ awk        '{rec=rec sep $0; sep=ORS} 
       !(NR%4){print rec > fn; rec=sep=""} 
       NR%4==2{fn = length($0)".seq"}' file

сгенерирует эти 3 файла с содержимым

==> 20.seq <==
@HISEQ:28:H8P69ADXX:1:1101:1462:2036 1:N:0:CTTGTA
NCCATAAAGTAGAAAGCACT
+
#00<FFFFFFFFFIIFIIFF

==> 21.seq <==
@HISEQ:28:H8P69ADXX:1:1101:1378:2223 1:N:0:CTTGTA
TCCTGTACTGAGCTGCCCCGA
+
BBBFFFFFFFFFFIIIIIIII
@HISEQ:28:H8P69ADXX:1:1101:1585:2081 1:N:0:CTTGTA
AAACCGTTACCATTACTGAGT
+
BBBFFFFFFFFFFIIIIFIII

==> 22.seq <==
@HISEQ:28:H8P69ADXX:1:1101:1419:2156 1:N:0:CTTGTA
TGGAGAGAAAGGCAGTTCCTGA
+
BBBFFFFFFFFFFIIIIIIIII

поскольку таких выходных файлов будет несколько, нет необходимости их явно закрывать.

Пояснение

{rec=rec sep $0; sep=ORS} построить запись построчно с ORS между строками, с ленивой инициализацией разделителя мы можем устранить оборванный первый разделитель.

!(NR%4) если номер строки кратен 4

{print rec > fn; rec=sep=""} распечатать запись в файл и сбросить запись и разделитель

NR%4==2, если номер строки равен 2 из 4.

{fn = length($0)".seq"} установить имя файла

11.04.2018
  • Большое спасибо и извините, что мне потребовалось некоторое время, прежде чем я смог ответить. Ваш код делает именно то, что я ищу, и менее чем за половину времени по сравнению с моим циклом for. реальный 5m3.676s пользователь 3m25.148s sys 0m8.816s для цикла for и реальный 2m2.489s пользователь 0m32.216s sys 0m1.876s для вашего. 12.04.2018
  • На заметку: не могли бы вы немного объяснить свой код? Я пытаюсь лучше понять awk. Я вижу, что мне нужно изменить, чтобы получить разные имена выходных файлов, но кроме этого я довольно потерян. 12.04.2018
  • Еще раз спасибо за объяснения, я думаю, теперь я понимаю, что вы там делали. 12.04.2018
  • Нет проблем, платите вперед! 12.04.2018
  • Я пытался адаптировать ваш скрипт из файлов fastq в файлы fasta (в основном тот же тип файла, но только с 2 строками в последовательности вместо 4). Я изменил каждое вхождение NR%4 на NR%2, думая, что это будет легко исправить, но получил следующую ошибку: фатальная: выражение для перенаправления `›' имеет нулевое строковое значение. Все остальное то же самое. Любая помощь очень ценится. 18.04.2018
  • строка, в которой установлено имя файла, не работает, вы не можете ожидать, что NR%2==2 будет правдой. измените на NR%2==0, который совпадает с !(NR%2), поэтому вы должны объединить два оператора с одним и тем же условием, очевидно, задав имя файла перед распечаткой. 18.04.2018
  • Новые материалы

    Коллекции публикаций по глубокому обучению
    Последние пару месяцев я создавал коллекции последних академических публикаций по различным подполям глубокого обучения в моем блоге https://amundtveit.com - эта публикация дает обзор 25..

    Представляем: Pepita
    Фреймворк JavaScript с открытым исходным кодом Я знаю, что недостатка в фреймворках JavaScript нет. Но я просто не мог остановиться. Я хотел написать что-то сам, со своими собственными..

    Советы по коду Laravel #2
    1-) Найти // You can specify the columns you need // in when you use the find method on a model User::find(‘id’, [‘email’,’name’]); // You can increment or decrement // a field in..

    Работа с временными рядами спутниковых изображений, часть 3 (аналитика данных)
    Анализ временных рядов спутниковых изображений для данных наблюдений за большой Землей (arXiv) Автор: Рольф Симоэс , Жильберто Камара , Жильберто Кейрос , Фелипе Соуза , Педро Р. Андраде ,..

    3 способа решить квадратное уравнение (3-й мой любимый) -
    1. Методом факторизации — 2. Используя квадратичную формулу — 3. Заполнив квадрат — Давайте поймем это, решив это простое уравнение: Мы пытаемся сделать LHS,..

    Создание VR-миров с A-Frame
    Виртуальная реальность (и дополненная реальность) стали главными модными терминами в образовательных технологиях. С недорогими VR-гарнитурами, такими как Google Cardboard , и использованием..

    Демистификация рекурсии
    КОДЕКС Демистификация рекурсии Упрощенная концепция ошеломляющей О чем весь этот шум? Рекурсия, кажется, единственная тема, от которой у каждого начинающего студента-информатика..