Mesclar 2 arquivos enormes que correspondem a várias colunas e preservar a ordem (correspondência de impressão e valores não correspondentes) - escala do awk

3

Estou com problemas ao mesclar dados de dois arquivos. São dados genéticos com alelos de cromossomo, posição, referência e alternativos, e eu preciso mesclar os arquivos combinando todas as 4 dessas colunas - com os alelos de referência e alternados de qualquer maneira. Então, eu preciso das colunas $ 1, $ 2, $ 4 e $ 5 ou $ 1, $ 2, $ 5 e $ 4 no arquivo de pesquisa para corresponder exatamente as colunas $ 1, $ 5, $ 6 e $ 7 no arquivo de dados. É realmente importante que eu mantenha a ordem exata no arquivo de dados - não posso classificá-lo (infelizmente não posso usar o join - esta é a resposta sugerida que encontrei em outros exemplos deste tipo de pergunta).

Eu usei awk e fiz o código funcionar para arquivos de amostra com alguns milhares de linhas, mas não aumenta muito para meu grande conjunto de dados (o arquivo de consulta tem > 300 milhões de linhas; o arquivo de dados tem 30 milhões) - presumivelmente como código requer manter na memória 2 grandes matrizes para pesquisa. Alguma sugestão para um código escalável (? Em perl?) Recebido com gratidão!

O formato do arquivo de pesquisa é:

1 10150 rs371194064 C T
1 10165 rs796884232 A AC
1 10177 rs367896724 A AC
1 10177 rs201752861 A C
1 10180 rs201694901 T C
1 10199 rs905327004 A T
1 10231 rs200279319 C A
1 10234 rs145599635 C T
1 10235 rs540431307 T TA
1 10235 rs1035249121 T A
1 10235 rs1035249121 T C
1 10241 rs960927773 T C
1 10247 rs796996180 T C
1 10248 rs148908337 A T
1 10249 rs774211241 AAC A

e o formato do meu arquivo de dados é

1 chr1 chr1:10177 1:10177_A_AC 10177 A AC
1 chr1 chr1:10235 1:10235_T_TA 10235 T TA
1 chr1 chr1:10352 1:10352_T_TA 10352 T TA
1 chr1 chr1:10505 1:10505_A_T 10505 A T
1 chr1 chr1:10506 1:10506_C_G 10506 C G
1 chr1 chr1:10511 1:10511_G_A 10511 G A
1 chr1 chr1:10539 1:10539_C_A 10539 C A
1 chr1 chr1:10542 1:10542_C_T 10542 C T
1 chr1 chr1:10579 1:10579_C_A 10579 C A

A saída deve se parecer com:

1       rs367896724     1:10177_A_AC    10177   A       AC      A       AC
1       rs540431307     1:10235_T_TA    10235   T       TA      T       TA
1       chr1:10352      1:10352_T_TA    10352   T       TA      T       TA
1       chr1:10505      1:10505_A_T     10505   A       T       A       T
1       chr1:10506      1:10506_C_G     10506   C       G       C       G
1       chr1:10511      1:10511_G_A     10511   G       A       G       A
1       chr1:10539      1:10539_C_A     10539   C       A       C       A
1       chr1:10542      1:10542_C_T     10542   C       T       C       T
1       chr1:10579      1:10579_C_A     10579   C       A       C       A

O código awk que consegui trabalhar para o arquivo de amostra é o seguinte:

awk 'BEGIN {OFS = "\t"}
NR==FNR {               #lookup file (323 million rows)
    key = $1 "," $2 "," $4 "," $5
    present[key] = 1
    ID[key] = $3
    Ref[key] = $4
    Alt[key] = $5

    key1 = $1 "," $2 "," $4 "," $5
    present1[key1] = 1
    ID1[key1] = $3
    Ref1[key1] = $4
    Alt1[key1] = $5

    next
}
{                       # my data file (3 million rows)
    key = $1 "," $5 "," $6 "," $7
    key1 = $1 "," $5 "," $7 "," $6
    if (present[key]) print $1, ID[key], $4, $5, $6, $7, Ref[key], Alt[key];
    else if (present1[key1]) print $1, ID1[key1], $4, $5, $6, $7, Ref1[key1], Alt1[key1];
    else              print $1, $3, $4, $5, $6, $7, $6, $7
}' $lookupfile $mydatafile > $outputfile
    
por Katie Fletcher 19.09.2017 / 10:14

1 resposta

0

O problema não é manter esse arquivo na memória, mas varrer a tabela de consulta para cada linha do arquivo de dados. Seu código não mostra, mas nos bastidores você faz 3'000'000 vezes 323'000'000 / 2 = quase meio quadrilhão de comparações, movendo milhares de terabytes pelo seu barramento de memória. Mesmo para memória rápida com 200 GBit / s, isso levará muitas horas.

Portanto, a chave para o problema é como armazenar a tabela de pesquisa. Eu sugiro usar uma árvore binária para reduzir exponencialmente o tempo de execução. Você pode fazer isso em perl ou C ou algum outro idioma, mas neste momento ele ficará fora do tópico aqui.

O conjunto de ferramentas de comando unix não irá ajudá-lo com este problema.

    
por 19.09.2017 / 10:39