Removendo novos caracteres de linha ao redor de padrões

Question

Eu tenho um arquivo fasta de milhões de seqüências pareadas, que se parece com isso:

>7001289F:56:HKH3FBCXX:2:1101:1692:2074 1:N:0:CGATGT
GAGCAGAGGCACCGCTGAGCAGACAGCGAGCGAGTGAAGGGGTCAGGGGCCAGTCAGCAATCTCGTGTAGAAAGAATCACGGTCGAGCGGTGCACGCATG
>NNNNN
GACACCTTCATTTCCACTTTATTGAGCAGCGGCGCATGCGTGCACCGCTCGACCGTGATTCTTTCTACACGAGATTGCTGACTGGCCCCTGACCCCTTCA
>7001289F:56:HKH3FBCXX:2:1101:1522:2186 1:N:0:CGATGT
GTAGATGATGAATACAGCTGTTGCTGCAGCAACTGGTGCTGAGTAAGCAACTGCGATCCATGGACGCATACCTAAACGGAAAGATAATTCCCAC
>NNNNN
GTGGGAATTATCTTTCCGTTTAGGTATGCGTCCATGGATCGCAGTTGCTTACTCAGCACCAGTTGCTGCAGCAACAGCTGTATTCATCATCTAC

Eu preciso formatá-lo como abaixo:

>7001289F:56:HKH3FBCXX:2:1101:1692:2074 1:N:0:CGATGT
GAGCAGAGGCACCGCTGAGCAGACAGCGAGCGAGTGAAGGGGTCAGGGGCCAGTCAGCAATCTCGTGTAGAAAGAATCACGGTCGAGCGGTGCACGCATGNNNNNGACACCTTCATTTCCACTTTATTGAGCAGCGGCGCATGCGTGCACCGCTCGACCGTGATTCTTTCTACACGAGATTGCTGACTGGCCCCTGACCCCTTCA
>7001289F:56:HKH3FBCXX:2:1101:1522:2186 1:N:0:CGATGT
GTAGATGATGAATACAGCTGTTGCTGCAGCAACTGGTGCTGAGTAAGCAACTGCGATCCATGGACGCATACCTAAACGGAAAGATAATTCCCACNNNNNGTGGGAATTATCTTTCCGTTTAGGTATGCGTCCATGGATCGCAGTTGCTTACTCAGCACCAGTTGCTGCAGCAACAGCTGTATTCATCATCTAC

Basicamente, os cabeçalhos complexos representam as leituras de avanço da sequência de DNA e o cabeçalho imediatamente a seguir representa a leitura inversa correspondente com os cabeçalhos NNNNN. Eu preciso acrescentar essas leituras reversas para as leituras de encaminhamento separadas com apenas NNNNN, mas estou lutando para remover os novos caracteres de linha com sed. Alguém pode esclarecer isso, por favor?

    

Answer 1

Se o seu arquivo é pequeno o suficiente para caber na memória, você pode fazer:

perl -00pe 's/\n>(NNNNN)\n/$1/g' file

Como seu arquivo quase certamente será grande demais para ser carregado em sua RAM, você pode usar isso:

$ perl -pe '$c++; if($c==2){chomp}
                  elsif($c==3){s/[>\n]//g;}
                  elsif($c==4){$c=0}' file.fa
>7001289F:56:HKH3FBCXX:2:1101:1692:2074 1:N:0:CGATGT
GAGCAGAGGCACCGCTGAGCAGACAGCGAGCGAGTGAAGGGGTCAGGGGCCAGTCAGCAATCTCGTGTAGAAAGAATCACGGTCGAGCGGTGCACGCATGNNNNNGACACCTTCATTTCCACTTTATTGAGCAGCGGCGCATGCGTGCACCGCTCGACCGTGATTCTTTCTACACGAGATTGCTGACTGGCCCCTGACCCCTTCA
>7001289F:56:HKH3FBCXX:2:1101:1522:2186 1:N:0:CGATGT
GTAGATGATGAATACAGCTGTTGCTGCAGCAACTGGTGCTGAGTAAGCAACTGCGATCCATGGACGCATACCTAAACGGAAAGATAATTCCCACNNNNNGTGGGAATTATCTTTCCGTTTAGGTATGCGTCCATGGATCGCAGTTGCTTACTCAGCACCAGTTGCTGCAGCAACAGCTGTATTCATCATCTAC

Explicação

• perl -pe '...' file.fa : para cada linha do arquivo de entrada file.fa , execute o script fornecido por -e e -p rint.
• $c++ : incrementa a variável $c em um para cada linha.
• if($c==2){chomp} : se o valor atual de $c for 2 , remova a nova linha do final da linha. Isso corresponde às linhas de seqüências avançadas.
• elsif($c==3){s/[>\n]//g;} : se $c for 3 , as linhas >NNNNN , remova os caracteres > e nova linha.
• elsif($c==4){$c=0}' : se $c for 4 , defina-o novamente para 0 novamente.

Observe que isso pressupõe leituras pareadas. Ele falhará se você não tiver exatamente uma leitura reversa para cada leitura direta no arquivo. Também assume que você tenha suas seqüências em uma única linha. Arquivos Fasta geralmente têm várias linhas por sequência, o padrão é cortar 60 caracteres. Isso mudou nos últimos anos, mas o formato ainda permite sequências de várias linhas.