Análises taxonômicas
INTRO
Primeiramente temos de obter os arquivos de taxonomia. Isto pode ser feito através do RDP (preferencial para sequências 16S) ou do BLAST (para sequências 18S ou funcionais).
Tabela de conteúdo |
Através do RDP
Como mencionado acima, se o arquivo de entrada for sequencias de 16S então o melhor caminho é alinhar pelo RDP.
De forma local:
- RDP Classifier
- Navegue até a pasta do RDP classifier e:
java -Xmx1g -jar rdp_classifier-2.3.jar -q A07.fasta -o A07out.txt -f fixrank
- RDP Multiclassifer
- Navegue até a pasta do RDP multiclassifier e:
java -Xmx1g -jar multiclassifier.jar --assign_outfile=A07.class --hier_outfile=A07.hier --format=fixRank A07.fasta
No entanto ao utilizar qualquer versão local do RDP é necessário formataro arquivo de saída:
Desta forma, o arquivo de saída:
gi|156447538|gb|EU074232.1| Bacteria domain 1.0 "Actinobacteria" phylum 1.0 Actinobacteria class 1.0 Actinomycetales orde0.41r 0.99 Propionibacteriaceae family 0.66 Friedmanniella genus 0.42 gi|156447537|gb|EU074231.1| Bacteria domain 1.0 "Actinobacteria" phylum 1.0 Actinobacteria class 1.0 Actinomycetales order 1.0 Propionibacteriaceae family 0.79 Friedmanniella genus
Deve ser formatdo para percentuais e ";", assim:
gi|156447538|gb|EU074232.1|; ; Bacteria; 100%; "Actinobacteria"; 100%; Actinobacteria; 100%; Actinomycetales; 99%; Propionibacteriaceae; 66%; Friedmanniella; gi|156447537|gb|EU074231.1|; ; Bacteria; 100%; "Actinobacteria"; 100%; Actinobacteria; 100%; Actinomycetales; 100%; Propionibacteriaceae; 79%; Friedmanniella;
Este arquivo em percentual é que será importado pelo MEGAN.
- Ao realizar a classificação pelo site do RDP (ver abaixo) o arquivo final já sairá em %.
O Maurício Cantão, LNCC, escreu um script em PERL bastante útil, para a reformatação dos arquivos de classificação produzidos pelo RDP classifier local (rdp2megan.pl).
Com o PERL instalado basta abrir um console e executar
rdp2megan.pl arquivo.rdpout
No servidor online
No site do RDP faça login e submeta as sequências desejadas.
Após o alinhamento clique no sinal + e vá em Classifier, indique as sequências selecionadas e espere.
Quando aparecerem os resultados, clique em [show assignment detail] e depois em "download as text file".
ATENÇÃO!
No momento há grandes diferenças entre os alinhamentos LOCAL' (Standalone) e ONLINE do RDP Classifier. Recomendo (por enquanto) a utlização do método ONLINE por resultar em classificações mais detalhadas.
Importando para o MEGAN
Agora que já temos o arquivo de classificação gerado pelo RDP e formatado em % podemos importá-lo para o MEGAN.
- Para isto vá em "File>Import" ou pressione CRTL+SHIFT+O
- A seguir, indique o diretório dos arquivos requisitados
Seja qual a forma de produção do arquivo de classificação do RDP (seja local ou online) o campo "Max number of matches per read" (dentro das opção de importação) deve ser 1.
Pois, segundo o Maurício Cantão:
- Quando você importa o resultado do RDP, então use 1 para o Max number... pois cada read só terá um hit, e os demais parâmetros ficam em 1 para o Min support, 80 para o percentual de confiança e 100 para Top percent.
- Quando os resultados são importados de uma placa de shotgun metagenômico, então é feito um BLASTx ou n contra os reads e os valores dos parâmetros para importação mudam. Para uma corrida metagenômica de 454, utilizo os valores: Max number of matches... 5 (ele vai buscar o ancestral comum dos 5 primeiros hits/organismos, e atribuir na árvore gerada), Min support: 5, um taxon somente será válido se 5 reads tiverem resultados iguais (isso evita falsos positivos), caso uma corrida tenha poucos reads (menos de 200 mil), eu diminuo para 3. Para o top percent utilizo 10 (o score do segundo hit tem que ser no mínimo 10% do score do primeiro, idem para os demais hits). E em score eu coloco 60, para não utilizar qualquer resultado no blast; o valor 60 já elimina resultados ao acaso.
Através do BLAST
Esta opção é indicada para sequências taxonômicas eucarióticas ou eucarióticas e procarióticas funcionais.
De forma local
Segue um exemplo genérico:
blastall -i seq.fasta -o seq.fasta.nr.out -p blastx -d nr -e 1e-03 -b 10 -v 10 -F F
É possível especificar a formatação do arquivo de saída utilizando o parâmetro '-m' para o Blast legacy:
pairwise (flag:-m 0) master-slave showing identities (flag:-m 1) master-slave no identities (flag:-m 2) flat master-slave, show identities (flag:-m 3) flat master-slave, no identities (flag:-m 4) master-slave no identities and blunt ends (flag:-m 5) flat master-slave, no identities and blunt ends (flag:-m 6) XML Blast output (flag:-m 7) tabular (flag:-m 8) tabular with comment lines (flag:-m 9)
Para o Blast+, utiliza-se o parâmetro '-outfmt':
0 = pairwise, 1 = query-anchored showing identities, 2 = query-anchored no identities, 3 = flat query-anchored, show identities, 4 = flat query-anchored, no identities, 5 = XML Blast output, 6 = tabular, 7 = tabular with comment lines, 8 = Text ASN.1, 9 = Binary ASN.1, 10 = Comma-separated values, 11 = BLAST archive format (ASN.1)
Caso não seja especificado o formato de saída, este assumirá a formatação 'pairwise' (que é a mesma fornecida pelo Blast online no servidor do NCBI)
Caso a saída selecionada seja tabular (ótimizada para abrir no excel, por exemplo), ao importar este arquivo para o MEGAN é necessário que sejam indicadas nas configurações de importação o arquivo que converte as GIs e um arquivo de sinônimos do SILVA database, ambos disponíveis no site de download do MEGAN.
Abaixo seguem alguns exemplos do Blast em execução
blastall -p blastx (Blast legacy, tal como no linux)
blastall -i mydata\A07.fasta -o mydata\A07.fasta.nr.out -p blastx -d c:\blast\db\nr -e 1e-03 -b 10 -v 10 -F F
blastx (Blast+)
blastx -query mydata\A07.fasta -out mydata\A07.blastout.txt -db nr -evalue 1000 -num_alignments 10 -num_descriptions 10 -seg no
blastn (Blast+)
blastn -query mydata\A07.fasta -out mydata\A07.blastout.nt.txt -db nt -evalue 1000 -num_alignments 10 -num_descriptions 10 -dust no
