extrai entradas fasta da lista usando enquanto lido

Question

Eu tenho 28 arquivos, cada um com ~ 14.000 "entradas". Uma única entrada consiste em um cabeçalho, denotado por > string, uma nova linha e, em seguida, uma sequência que é uma string. Cada entrada tem sequência / sequência de comprimento variável. Em todos os 28 arquivos, há cabeçalhos de entrada idênticos, mas a sequência de cada entrada é variável.

Por exemplo, um arquivo CR1_ref.fasta seria parecido com

>FBgn0080937
ATGGATAAAAGGCTCAGCGATAGTCCCGGAGATTGTCGCGTAACCAGATCCAGCATGACGCCCACCCTCCGCTTGGAGCACAGTCCCCGGCGGCAACAACAGCAACAACA
>FBgn0076379
ATGCTGCGCACCCTTTTCGCCGTGCGTGGTCAGTGCCAGCAGCTGCTGAGGAGAACATTCACCCCCCATTGCAGTGGCCAACGA
>FBgn0070974
ATGCAGACGCGTCCGAGCAGTGAACCGCAGCGCGCCAAGGAGCAACTCCTGCGGGAGCTGCCGCCGCAGAAATGCTCCAGCGCCACGCTGGCCAAGAAGGTGCTGTCGCAGAGCCCGCCGGCAGCCCCGCCGCCCACACCGGCCACAATTGTGCCGCTCACTGCGGTGCCCGTCATCCAGCTGACGCCTCCGTCGCACTCCGGCGACACGCCGCAAAAGCCAGCACCTCCGGCGCCGCCGCCGCC

O objetivo geral é criar ~ 14.000 novos arquivos. Onde cada arquivo é a entrada associada a um ID / cabeçalho específico em todos os 28 arquivos.

Para extrair uma única entrada de um único arquivo, posso usar o seguinte comando

sed -n '/^>FBgn0080937$/{p;n;p;}' CR1_ref.fasta

Para extrair essa entrada em todos os 28 arquivos, cada um termina em ref.fasta, eu posso fazer

for i in *ref.fasta; do sed -n '/^>FBgn0080937$/{p;n;p}' $i; done > FBgn0080937.fasta

Eu tenho um arquivo de texto separado que tem 14.000 linhas cada linha correspondente a um cabeçalho para uma entrada chamada gene.txt. As primeiras linhas deste arquivo parecem

FBgn0080937
FBgn0076379
FBgn0070974
FBgn0081668
FBgn0076576
FBgn0076572
FBgn0079684
FBgn0070907
FBgn0080226
FBgn0072746

Gostaria de ler este arquivo criando um novo arquivo de texto por ID de cabeçalho. Abaixo de $ F está extraindo entradas para um cabeçalho particular (FBgn *) e armazenando isto em um novo arquivo. Eu estou usando o comando de substituição para renomear as seqüências com base no arquivo ref.fasta de onde elas vêm.

while read -r line;
do F=$line
for i in *ref.fasta
do sed -n "/^>$F$/{s/FB.*/$i/;p;n;p;}" $i > $line.fasta
done
done < "gene.txt"

Atualmente, esse script cria 14.000 arquivos, mas cada arquivo tem apenas uma única sequência.

>Z9_ref.fasta
ATGCAGACGCGTCCGAGCAGTGAACCGCAGCGCGCCAAGGAGCAAC

Estou esperando 28 sequências de uma sequência por arquivo * ref.fasta. O comando sed está exibindo a última entrada. A saída esperada seria

    >CR1_ref.fasta
    ATGCAGACGCGTCCGAGCAGTGAACC
    >FH2_ref.fasta
    AGCAGTGAACCGCAGCGCGCCAAGGAGCAAC
    >MSH10_ref.fasta
    CGCGTCCGAGCAGTGAACCGCAGCGCGCCAAGGAGCAAC
    >Z9_ref.fasta
    ATGCAGACGCGTCCGAGCAGTGAACCGCAGCGCGCCAAGGAGCAAC

    

Answer 1

O shell não é realmente adequado para esse tipo de análise. Você pode ver em seu próprio código que você leu a totalidade de cada arquivo uma vez para cada dos nomes dos genes lidos no arquivo gene.txt .

O seguinte comando awk único faria a mesma coisa mais rápido.

awk -F '>' '
    FNR == NR           { genes[$1]; next }
    /^>/ && $2 in genes { if (out != "") close(out);
                          out = $2 ".fa"
                          split(FILENAME, a, "_")
                          $0 = ">" a[1] "_" $2 }
    out != ""           { print >>out }' genes.txt *_ref.fasta

Primeiro, lê o arquivo genes.txt e cria uma matriz associativa chamada genes com os nomes dos genes como chaves.

Quando chega aos arquivos Fasta (o código assume que todos eles são chamados de XXX_ref.fasta ), quando lemos um cabeçalho Fasta, e o gene no cabeçalho é uma chave na nossa lista genes , então criamos um nome de arquivo de saída do nome do gene como genename.fa e reescrevemos o cabeçalho para incluir a parte do nome do arquivo atual antes do sublinhado.

Se o cabeçalho original em XXX_ref.fasta for

>genename

então isso seria transformado em

>XXX_genename

A última parte do script awk envia todas as linhas para o arquivo de saída apropriado.

Testar isso com os dados que você forneceu gera três arquivos:

$ ls *.fa
FBgn0070974.fa FBgn0076379.fa FBgn0080937.fa

$ cat FBgn0076379.fa
>CR1_FBgn0076379
ATGCTGCGCACCCTTTTCGCCGTGCGTGGTCAGTGCCAGCAGCTGCTGAGGAGAACATTCACCCCCCATTGCAGTGGCCAACGA