Заменить шаблон из первого столбца в файле2 шаблоном во втором столбце в файле2 в файле1 - PullRequest
3 голосов
/ 17 июня 2020

Я попытался найти решение для следующего: у меня есть файл .gff3, для которого я хочу заменить заголовки генов на упрощенное имя. И исходные заголовки гена, и новое имя гена даны в отдельном файле, с исходным именем в столбце 1 и новым именем в столбце 2. Как я могу использовать sed (я думаю, что sed здесь наиболее подходит) для замены всех вхождений в файле .gff3 с новым сокращенным именем во втором столбце?

Пример строки файла .gff3:

tulip_contig_65_pilon_pilon .   contig  1   93354   .   .   .   ID=tulip_contig_65_pilon_pilon;Name=tulip_contig_65_pilon_pilon
tulip_contig_65_pilon_pilon maker   gene    19497   23038   .   +   .   ID=maker-tulip_contig_65_pilon_pilon-augustus-gene-0.4;Name=maker-tulip_contig_65_pilon_pilon-augustus-gene-0.4
tulip_contig_65_pilon_pilon maker   mRNA    19497   23038   .   +   .   ID=maker-tulip_contig_65_pilon_pilon-augustus-gene-0.4-mRNA-1;Parent=maker-tulip_contig_65_pilon_pilon-augustus-gene-0.4;Name=maker-tulip_contig_65_pilon_pilon-augustus-gene-0.4-mRNA-1;_AED=0.00;_eAED=0.00;_QI=418|1|1|1|0|0|3|2100|206

Пример файла замены строк:

augustus_masked-tulip_contig_306_pilon_pilon-processed-gene-0.1   gene1
maker-tulip_contig_306_pilon_pilon-augustus-gene-0.12 gene2
maker-tulip_contig_65_pilon_pilon-augustus-gene-0.4   gene3

ожидаемый результат:

tulip_contig_65_pilon_pilon   .   contig  1   93354   .   .   .   ID=tulip_contig_65_pilon_pilon;Name=tulip_contig_65_pilon_pilon
tulip_contig_65_pilon_pilon   maker   gene    19497   23038   .   +   .   ID=gene3;Name=gene3
tulip_contig_65_pilon_pilon   maker   mRNA    19497   23038   .   +   .   ID=gene3-mRNA-1;Parent=gene3;Name=gene3-mRNA-1;_AED=0.00;_eAED=0.00;_QI=418|1|1|1|0|0|3|2100|206

Я пробовал использовать:

while read -r pattern replacement; do sed -i "s/$pattern/$replacement/" file.gff3 ; done < rename.txt

Но на основании ответа ниже я сейчас использую AWK. Я использую сценарий ( тот же отступ, что и Эд Мортон , но его копирование здесь немного меняет):

NR==FNR {
 map[$1] = $2
 next 
} 
{
 for (old in map) {
    gsub(old,map[old])
 }
 print 
}

Для запуска я использую:

awk -f tst.awk rename.txt original.gff3 > new.gff3 

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

Длина файла gff составляет 7369803 строки. Файл rename.txt состоит из 18477 строк.

Здесь приветствуются любые советы!

Ответы [ 3 ]

6 голосов
/ 20 июня 2020

Это обеспечивает полное соответствие строки от = до конца -gene=<number> на каждой строке .gff3 и должно быть на несколько порядков быстрее и надежнее, чем то, что мы делали ранее, поскольку оно заменяет только 1 -3 строки, фактически найденные в каждой строке исходного файла .gff3, вместо попытки заменить все 18 000+ строк, которые существуют в файле rename.txt:

$ cat tst.awk
NR==FNR {
    map[$1] = $2
    next
}
{
    head = ""
    tail = $0
    while ( match(tail,/((ID|Parent|Name)=)([^;]+-gene-[0-9]+\.[0-9]+)(.*)/,a) ) {
        old = a[3]
        head = head substr(tail,1,RSTART-1) a[1] (old in map ? map[old] : old)
        tail = a[4]
    }
    print head tail
}

.

$ awk -f tst.awk rename.txt original.gff3
tulip_contig_65_pilon_pilon .   contig  1   93354   .   .   .   ID=tulip_contig_65_pilon_pilon;Name=tulip_contig_65_pilon_pilon
tulip_contig_65_pilon_pilon maker   gene    19497   23038   .   +   .   ID=gene3;Name=gene3
tulip_contig_65_pilon_pilon maker   mRNA    19497   23038   .   +   .   ID=gene3-mRNA-1;Parent=gene3;Name=gene3-mRNA-1;_AED=0.00;_eAED=0.00;_QI=418|1|1|1|0|0|3|2100|206

Он использует GNU awk для третьего аргумента для match () - я предполагаю, что у вас есть GNU awk (или вы можете его установить), так как вы использовали GNU sed.

Итак, match() изолирует строку (который затем сохраняется в old) из текущей строки original.gff3, что может быть в rename.txt (как хранится в map[] в первом блоке), а затем old in map проверка, действительно ли эта строка находится в rename.txt или нет, и, если да, заменяет old соответствующим новым значением из map[]. Это все внутри while, который зацикливается, пока match() продолжает находить новые строки, которые являются кандидатами на замену в текущей строке.

Итак, вместо исходного сценария awk ниже (и сценария sed в вашем вопрос), который повторяется один раз для каждой из 18000+ строк в rename.txt, приведенный выше цикл повторяется только один раз для каждой строки в текущей строке original.gff3, которую, возможно, потребуется заменить, что, согласно вашему опубликованному образцу ввода, составляет только до 3 раз.

Исходный ответ, основанный на простом ускорении вашей оболочки l oop вызов sed:

Что-то вроде этого вам нужно:

$ cat tst.awk
NR==FNR {
    map[$1] = $2
    next
}
{
    for (old in map) {
        gsub(old,map[old])
    }
    print
}

.

$ awk -f tst.awk repl.txt foo.gff3
tulip_contig_65_pilon_pilon .   contig  1   93354   .   .   .   ID=tulip_contig_65_pilon_pilon;Name=tulip_contig_65_pilon_pilon
tulip_contig_65_pilon_pilon maker   gene    19497   23038   .   +   .   ID=gene3;Name=gene3
tulip_contig_65_pilon_pilon maker   mRNA    19497   23038   .   +   .   ID=gene3-mRNA-1;Parent=gene3;Name=gene3-mRNA-1;_AED=0.00;_eAED=0.00;_QI=418|1|1|1|0|0|3|2100|206

Существуют некоторые решения относительно сопоставления строк и регулярных выражений и полного и частичного сопоставления, которые также применяются к вашей оболочке + sed l oop, поэтому подумайте о своих полных требованиях и предоставьте образцы ввода / вывода для тестирования, а затем мы можем обработать это для подходит, если это не совсем то, что вам нужно. Прямо сейчас он выполняет частичное сопоставление регулярных выражений, как это делает команда sed в вашем вопросе.

1 голос
/ 26 июня 2020

Количество попыток сценария оболочки обычно составляет O (n × m), что для 7369803 × 18477 огромно.

Я не собираюсь писать код, но для чего-то такого большого я бы напишите специальную c программу на таком языке, как C или python:

  • Добавьте номер строки в начало каждой строки файла данных.
  • Сортировка файла данных по полю ключевого слова.
  • Сортировка файла переименования по полю ключевого слова.
  • Параллельный просмотр двух отсортированных файлов (программа «ad_ho c»):
    • если ключ карты сортируется до ключа данных, прочитать следующую запись карты.
    • иначе:
      • если ключи совпадают, сделайте замену.
      • записать запись данных
      • прочитать следующую запись данных
  • Сортировать результат по номеру строки.
  • Удалить номера строк.

Вычислительная сложность будет полностью зависеть от операций сортировки O (n × log (n)), которые будут невероятно быстрее, чем любые шарнир с вложенными петлями. Сама программа ad_ho c является линейной, O (n + m), просто считывая два файла и обрабатывая каждую строку один раз, без каких-либо циклов.

#!/bin/bash
#Assumes "%" does not occur in either file

ad_hoc \
<( sed <map  -e 's/\([^ ]*\)  *\([^ ]*\)$/\1%\2/' \
 | sort -t '%' -k 1,1 ) \
<( sed  <input  '/./=' \
 | sed '/./N; s/\n/ /' \
 | sed -e 's/ID=\([^;]*\);/ID=%\1%;/' \
 | sort -t '%' -k 2,2 ) \
| sort -t ' ' -k 1,1n \
| sed -e 's/^[0-9]* *//'

Все, что вам нужно, это «ad_ho c "программа.

1 голос
/ 23 июня 2020

С perl, очень похоже в logi c на awk ответ

perl -ne 'if(!$#ARGV){ @n = split; $h{$n[0]}=$n[1] } else{ 
          s/(ID|Name|Parent)=\K[^;]+-gene-\d+\.\d+/exists $h{$&} ? $h{$&} : $&/ge;
          print }' rename.txt original.gff3
  • if(!$#ARGV) для обработки только первого файла
  • @n = split чтобы получить ключ и значение для переименования карты
    • $h{$n[0]}=$n[1] карта
  • (ID|Name|Parent)=\K для отметки начального местоположения, не будет частью соответствующей части
  • [^;]+-gene-\d+\.\d+ пожалуйста, проверьте, все ли ваши гены соответствуют этому формату, он не может содержать символ ; и всегда заканчивается -, за которым следует gene, за которым следуют цифры, за которыми следует . с последующими цифрами, иначе это не сработает
  • exists $h{$&} проверьте, существует ли совпадающая часть в качестве ключа на карте, и если да, используйте значение с карты, иначе не меняйте имя гена

Я провел некоторый тест скорости, но он может существенно отличаться для вашего реального случая использования. Edit - только что понял, что созданный ниже t2 бесполезен, так как он будет иметь только 3 разных сопоставления.

$ perl -0777 -ne 'print $_ x 1000000' original.gff3 > t1
$ perl -0777 -ne 'print $_ x 10000' rename.txt > t2

$ time perl -ne 'if(!$#ARGV){ @n = split; $h{$n[0]}=$n[1] } else{ 
            s/(ID|Name|Parent)=\K[^;]+-gene-\d+\.\d+/exists $h{$&} ? $h{$&} : $&/ge;
            print }' t2 t1 > m1 
real    0m17.473s

$ time awk -f tst.awk t2 t1 > m2    
real    2m53.598s

$ # assuming your input is ASCII
$ time LC_ALL=C awk -f tst.awk t2 t1 > m3
real    1m53.976s

$ du -h t1 t2 m1
591M    t1
1.9M    t2
371M    m1

$ diff -sq m1 m2
Files m1 and m2 are identical
$ diff -sq m1 m3
Files m1 and m3 are identical
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...