Tutorial prático do bioBakery — MetaPhlAn 4 (Conda)

Perfil taxonômico com MetaPhlAn 4, a estrela do bioBakery

Fala, pessoal! Bora conhecer o bioBakery Conjunto de ferramentas do Huttenhower Lab (Harvard) e do Segata Lab (Universidade de Trento) para análise de metagenomas. ? Ele não é uma ferramenta só, é uma "padaria" inteira 🥖: o KneadData limpa as reads (inclusive tirando o DNA do hospedeiro), o MetaPhlAn diz quem está na amostra, o HUMAnN diz o que eles estão fazendo (vias metabólicas) e o StrainPhlAn desce até o nível de linhagem.

Neste guia a gente vai focar no MetaPhlAn 4 Metagenomic Phylogenetic Analysis: estima a composição taxonômica usando genes marcadores específicos de cada clado. . A diferença para o Kraken2 é a estratégia: em vez de olhar todos os k-mers de todos os genomas, o MetaPhlAn procura só genes marcadores Genes que existem em um único clado (espécie, gênero…) e em nenhum outro. Se a read bate nele, não tem dúvida de quem é. , milhões deles, organizados em SGBs Species-level Genome Bins: grupos de genomas da mesma espécie, inclusive de espécies que ainda nem têm nome (só conhecidas por metagenomas). . Resultado: menos falsos positivos e uma abundância relativa bem confiável.

Sobre "quem está na festa" e "o que eles estão fazendo", escrevi no Medium: O que existe aqui? e O que eles estão fazendo?.

O que você vai precisar antes:

• Conda — instalação aqui
• WSL/Linux — tutorial aqui
• FASTQs limpos — de preferência depois do trimming (guia do Fastp, de onde este guia continua)
• Máquina parruda — a documentação oficial pede uns 30 GB de RAM e o banco ocupa algumas dezenas de GB em disco

Esse é o guia mais "pesado" da série. Num notebook comum (8–16 GB de RAM) o MetaPhlAn provavelmente não vai rodar. Se você tem acesso a um servidor da universidade ou do laboratório, esse é o lugar certo para ele.

Neste guia, você vai:

• Instalar o MetaPhlAn 4 e baixar o banco de marcadores,
• Gerar o perfil taxonômico de uma amostra e entender a saída,
• Juntar várias amostras numa tabela só,
• Comparar o MetaPhlAn com o Kraken2.

Passo a passo

1. Preparação do ambiente
2. Baixando o banco de marcadores
3. Rodando o MetaPhlAn
4. Entendendo a saída
5. Várias amostras e tabela única
6. MetaPhlAn ou Kraken2?
7. E o HUMAnN?

1. Preparação do ambiente

conda create -n mpa -c conda-forge -c bioconda "metaphlan>=4.2" -y
conda activate mpa

mpa é o nome que a gente deu pro ambiente (vem de MetaPhlAn, e é o mesmo nome que a documentação oficial usa). O Conda já instala junto o bowtie2, que o MetaPhlAn usa para alinhar as reads nos marcadores.

Confira a versão:

metaphlan --version

2. Baixando o banco de marcadores

Crie uma pasta para os bancos (a mesma ideia do guia do Kraken2) e baixe o banco mais recente:

mkdir -p ~/bases/metaphlan
metaphlan --install --db_dir ~/bases/metaphlan --nproc 4

O --install baixa sozinho a versão mais nova do banco (no momento em que escrevi este guia, a vJan25) e já monta o índice do bowtie2. Demora um pouco, então é outra hora boa pra um café ☕.

Guarde a versão do banco! Perfis gerados com bancos diferentes não podem ser comparados nem juntados na mesma tabela. Se for analisar amostras em momentos diferentes, use sempre o mesmo banco. Para conferir qual você tem:

ls ~/bases/metaphlan

3. Rodando o MetaPhlAn

Volte para a pasta do projeto do guia do Fastp e crie uma pasta para os resultados:

cd ~/analise_dados/aula_1
mkdir -p metaphlan

Não fez o guia do Fastp? Use seus próprios FASTQs ou baixe dados públicos com o SRA-Tools. Só troque os nomes dos arquivos.

Agora rode o MetaPhlAn. Repare que o R1 e o R2 vão juntos, separados por vírgula:

metaphlan \
  pos_trim/meta_aula_trim_1_L001_R1_001.fastq.gz,pos_trim/meta_aula_trim_1_L001_R2_001.fastq.gz \
  --input_type fastq \
  --db_dir ~/bases/metaphlan \
  --mapout metaphlan/meta_aula_1.mapout.bz2 \
  --nproc 4 \
  -o metaphlan/meta_aula_1_perfil.txt

Principais parâmetros:

Por que guardar o .mapout.bz2? Porque a parte demorada é o alinhamento. Se depois você quiser rodar de novo com outro parâmetro, dá pra partir dele e o MetaPhlAn termina em segundos:

metaphlan metaphlan/meta_aula_1.mapout.bz2 --input_type mapout --db_dir ~/bases/metaphlan -o metaphlan/meta_aula_1_perfil_v2.txt

Usou Oxford Nanopore ou PacBio? A partir da versão 4.2 o MetaPhlAn aceita reads longas, é só adicionar --long_reads (ele passa a usar o minimap2 no lugar do bowtie2). Ficou em dúvida sobre reads longas x curtas? Escrevi sobre isso no Medium.

4. Entendendo a saída

Dê uma olhada no começo do arquivo:

head -20 metaphlan/meta_aula_1_perfil.txt

As primeiras linhas (com #) dizem a versão do banco e quantas reads foram processadas. Depois vem a tabela:

Quer só as espécies, sem a linhagem gigante? Rode:

grep "s__" metaphlan/meta_aula_1_perfil.txt | grep -v "t__" | \
  awk -F'\t' '{n=split($1,a,"|"); print a[n]"\t"$3}' | sort -t$'\t' -k2,2gr | head -10

Vai aparecer uma linha UNCLASSIFIED. É a fração estimada da amostra que não bateu em nenhum marcador (organismos que não estão no banco, DNA do hospedeiro…). Se ela for muito alta, vale checar se a amostra tem muito DNA humano. O KneadData, do próprio bioBakery, remove essas reads antes do MetaPhlAn.

Atenção: a abundância do MetaPhlAn é de organismos (corrigida pelo tamanho dos genomas), e não de reads. Por isso os números não batem diretamente com os do Kraken2/Bracken, e tá tudo bem!

5. Várias amostras e tabela única

Rode um loop para todas as amostras da pasta pos_trim/:

for r1 in pos_trim/*_R1_001.fastq.gz; do
  r2=${r1/_R1_/_R2_}
  amostra=$(basename "$r1" _R1_001.fastq.gz)

  metaphlan "$r1","$r2" --input_type fastq --db_dir ~/bases/metaphlan \
    --mapout metaphlan/${amostra}.mapout.bz2 --nproc 4 \
    -o metaphlan/${amostra}_perfil.txt
done

Depois junte tudo, com uma coluna por amostra:

merge_metaphlan_tables.py metaphlan/*_perfil.txt > metaphlan/abundancia_todas.txt

E separe só o nível de espécie (mantendo o cabeçalho), pronto para levar pro R ou Python:

grep -E "clade_name|s__" metaphlan/abundancia_todas.txt | grep -v "t__" > metaphlan/especies_todas.txt

6. MetaPhlAn ou Kraken2?

Pergunta clássica! Não existe "o melhor": cada um tem seu ponto forte.

Na dúvida, rode os dois! As espécies que aparecem nos dois com abundância parecida são as que você pode afirmar com mais segurança.

7. E o HUMAnN?

O HUMAnN é o próximo passo natural: ele usa o perfil do MetaPhlAn para descobrir quais genes e vias metabólicas estão presentes na comunidade.

Quando escrevi este guia, o HUMAnN 4 ainda estava em versão alpha e funcionava só com o MetaPhlAn até a versão 4.1.1 (com o banco mpa_vOct22_CHOCOPhlAnSGB_202403), e não com o 4.2 que a gente instalou aqui. Por isso ele vai ganhar um guia próprio, num ambiente separado. Enquanto isso, acompanhe as novidades no fórum do bioBakery.

Conclusão

Você gerou o perfil taxonômico de metagenomas com o MetaPhlAn 4:

• Instalação e download do banco de marcadores
• Perfil de uma amostra e leitura da saída
• Loop para várias amostras e tabela única com merge_metaphlan_tables.py
• E ainda sabe quando usar MetaPhlAn ou Kraken2

Até a próxima!

Desafio!!!

Rode o Kraken2 + Bracken e o MetaPhlAn na mesma amostra e compare o top 10 de espécies. Quantas aparecem nos dois? Quais só aparecem no Kraken2?
Bônus: coloque os dois num pipeline Nextflow e rode tudo com um comando só.
Até a próxima!

Dê um suporte ao meu projeto. Doe um cafézinho ☕.
Pix: biologolee@gmail.com
Bitcoin: bc1qg7qrfhclzt3sm60en53qv8fmwpuacfaxt5v55k

Referências

1. bioBakery — página oficial no GitHub
2. MetaPhlAn — repositório oficial
3. MetaPhlAn 4.2 — wiki oficial
4. Extending and improving metagenomic taxonomic profiling with uncharacterized species using MetaPhlAn 4 (Nature Biotechnology, 2023)
5. Pacote MetaPhlAn no Bioconda
6. Integrating taxonomic, functional, and strain-level profiling of diverse microbial communities with bioBakery 3 (eLife, 2021)
7. Fórum do bioBakery