Existen variadas herramientas para realizar una "predicción funcional" a partir de amplicones r16S o ITS. En esta página dejaré instrucciones y comandos para las herramientas que uso habitualmente.
Para continuar requiero que tengas OTUs o ASVs con su respectiva identificación taxonómica y una tabla de conteo del número de reads por para cada OTU/ASV en cada muestra.
SITIO EN CONSTRUCCIÓN
1. Picrust2 (sitio oficial)
Esta herramienta sólo está diseñada para amplicones r16S (procariontes). Existe un pipeline alternativo para realizar predicción desde ITS (fungi) pero es incipiente y por ahora las asignaciones no son muy confiables.
Picrust utiliza las secuencias de los OTUs/ASVs para, mediante homología, encontrar los mejores match en una base de datos de r16S que, a su vez, posee información de funcionalidad de muchas especies de bacterias, tanto en metacyc_pathways (default), COG, EC y KO. Utilza la tabla de conteo para tratar de relacionar la abundancia de cada categoría funcional en cada muestra. Picrust realiza una normalización interna de la abundancia de cada función respecto del número de copias conocidas del gen r16S para cada especie intentificada.
1.1 Instalación
conda create -n picrust241 -c bioconda -c conda-forge picrust2
1.2 Archivos necesarios (input)
a) Archivo fasta con las secuencias "representantes", ya sean OTUs o ASVs. [archivo.fasta]
b) Tabla del conteo de OTUs o ASVs en cada muestra (de preferencia debe ser una tabla rarefaccionada). Esta tabla debe estar en un formato de columnas delimitadas por tabulaciones y no requiere la taxonomía. El encabezado no debe tener el símbolo #. [tabla.tsv]
1.3 Ejecución
picrust2_pipeline.py -s archivo.fasta -i table.txt -o directorio_de_salida -p 16 --in_traits EC,KO,GO,GENE_NAMES --remove_intermediate --stratified
-o : directorio donde quedarán los archivos de salida
-p : número de procesadores a utilizar
--in_traits : tipos de bases de datos desde las cuales obtendremos información. Se generará un resultado por cada base de datos. Metacyc es el resultado por default.
--remove_intermediate : eliminar archivos intermediarios
--stratified : entrega el nombre de los OTUs o ASVs que aportan a cada función asignada.
Versión con datos qiime2 usando datos rarefaccionados según ESTE PROTOCOLO (coloco la lista completa de comandos).
Es extremadamente importante usar datos rarefaccionados pues el análisis arrojará resultados que dependen de la profundidad de secuenciación. Las librerías que tengan mayor profundidad tendrán mayor representación. Para evitarlo, debemos "normalizar" la profundidad de secuenciación, por ejemplo, usando rarefacción.
cd core-metrics-phylogenetic-metrics_r7k
#obtenemos la tabla rarefaccionada en formato biom
qiime tools export --input-path rarefied_table.qza --output-path ../exported_r7k
#transformamos la tabla de formato biom a formato texto y limpiamos los encabezados
cd ../exported_r7k
biom convert -i feature-table.biom -o feature-table.tsv --header-key taxonomy --output-metadata-id "taxonomy" --to-tsv
sed 's/^#OTU ID/id/' feature-table.tsv
sed -i '1d' feature-table.tsv
cd ..
# Ejecutamos picrust2
picrust2_pipeline.py -s representative_sequences.fasta -i exported_r7k/feature-table.tsv -o picrust2 -p 64 --in_traits EC,KO,GO,GENE_NAMES --remove_intermediate --stratified
1.4 Post-proceso de resultados
1.4.1 Traducción de códigos de genes y/o rutas a nombres o descripciones.
Picrust por default entrega como resultado una tabla por cada base datos consultada, Dichas tablas sólo contienen el código interno de la base de datos, no la descripción. Para asignar la descripción (y que sea legible para humanos) entonces tendremos que usar estos comandos.
cd directorio_de_salida
add_descriptions.py -i pathways_out/path_abun_unstrat.tsv.gz -m METACYC -o pathways_out/path_abun_unstrat_descrip.tsv.gz
add_descriptions.py -i KO_metagenome_out/pred_metagenome_unstrat.tsv.gz -m KO -o KO_metagenome_out/pred_metagenome_unstrat_descrip.tsv.gz
add_descriptions.py -i EC_metagenome_out/pred_metagenome_unstrat.tsv.gz -m EC -o EC_metagenome_out/pred_metagenome_unstrat_descrip.tsv.gz
add_descriptions.py -i COG_metagenome_out/pred_metagenome_unstrat.tsv.gz -m COG -o COG_metagenome_out/pred_metagenome_unstrat_descrip.tsv.gz
1.4.2 Resultados Metacyc pathways (carpeta pathways_out)
La base de datos Metacyc_pathways tiene distintos niveles de jerarquía. El resultado "default" de Picrust entrega todas las Pathways sin un orden de jerarquía. Para simplificar el análisis vamos a procesar el output de Picrust para obtener categorías de 1ra y 2da jerarquía. Para ello necesitamos el "mapa" de las categorías, es decir, una archivo que agrupe categorías y subcategorías. Los mapas que yo utilizo los encontré en este sitio web.
Archivos:
metacyc_pathways_info_prokaryotes_top_level.tsv
metacyc_pathways_info_prokaryotes_sec_level.tsv
add_descriptions.py -i pathways_out/path_abun_unstrat.tsv.gz -o pathways_out/path_abun_unstrat_descrip_top.tsv.gz --custom_map_tabl metacyc_pathways_info_prokaryotes_top_level.tsv
add_descriptions.py -i pathways_out/path_abun_unstrat.tsv.gz -o pathways_out/path_abun_unstrat_descrip_sec.tsv.gz --custom_map_tabl metacyc_pathways_info_prokaryotes_sec_level.tsv
El resultado que más he usado es el de Metacyc pues contiene información fácil de interpretar y agrupada por jerarquía. Serían los archivos path_abun_unstrat_descrip_top.tsv.gz y path_abun_unstrat_descrip_sec.tsv.gz. Estos archivos los descomprimo y manipulo en Excel. Mi recomendación es transformar la abundancia absoluta en abundancia relativa por muestra (columna) y luego utilizar la tabla de datos para producir un colorido heatmap en R o en Morpheus.
1.4.3 Resutlados KEGG Onthology (carpeta KO_metagenome_out)
La tabla de resultados KEGG contiene una lista de los genes que podrían estar presentes en las comunidades estudiadas. Recordemos que, a partir de una homología del fragmento del gen r16S hemos entrado "especies similares", las cuales sabemos por bibliografía que su genoma contiene ciertos genes. Dada la homología del r16S, por transitividad "asumimos" que nuestras bacterias (OTUs/ASVs) contienen los mismos genes (predicción indirecta de funcionalidad). La lista de genes es interesante pero no nos permite ver "procesos" a menos que hagamos un trabajo detallado de unir rutas metabólicas con los genes encontrados. Este trabajo "manual" no es necesario pues KEGG tiene categorizados los genes en categorias BRITE, es decir, podemos asociar cada gen encontrado con alguna o varias categorías funcionales. Lamentablemente Picrust2 no posee una forma de transformar nuestra tabla de genes en categorías (Picrust1 si lo hacía) y hay una razón que se explica en este link. De todas formas, en el mismo link ofrecen una forma de realizar la conversión usando un script de R. A continuación describo los pasos en base al script de R que he modificado.
1.4.3.1 Localizar el archivo pred_metagenome_unstrat.tsv.gz (carpeta KO_metagenome_out), descomprimirlo y renombrarlo a ko.pred_metagenome_unstrat.tsv
cd KO_metagenome_out
gunzip -c pred_metagenome_unstrat.tsv.gz > ko.pred_metagenome_unstrat.tsv
sed -i 's/^ko://' ko.pred_metagenome_unstrat.tsv
1.4.3.2 Descargar mi script de R en la misma carpeta donde está el archivo ko.pred_metagenome_unstrat.tsv
wget https://github.com/jomaldon/scripts_bioinfo/raw/master/picrust1_categorize_by_func.R
1.4.3.3 Descargar la base de datos que relaciona códigos KEGG con categorías funcionales
wget https://github.com/jomaldon/scripts_bioinfo/raw/master/supporting_files/picrust1_KO_BRITE_map.tsv
1.4.3.4 Ejecutar el script de R
En un sistema linux con R instalado (puede ser mediante conda) se puede ejecutar la siguiente línea.
Rscript picrust1_categorize_by_func.R
Otra opción es usar RStudio o cualquier entorno de ejecución de R compatible.
1.4.3.5 Análisis de resultados.
EL resultado se encuentra en los archivos ko_L1_sorted.txt, ko_L2_sorted.txt y ko_L3_sorted.txt los cuales representan las categorías funcionales hasta nivel1, hasta nivel2 y hasta nivel3 respectivamente. Estos archivos los manipulo en Excel. Mi recomendación es transformar la abundancia absoluta en abundancia relativa por muestra (columna) y luego utilizar la tabla de datos para producir un colorido heatmap en R o en Morpheus.
2. Faprotax (sitio oficial)
Esta herramienta está diseñada sólo para amplicones r16S (procariontes).
Faprotax utiliza sólo la información de taxonomía como "clave" para contrastar contra una base de datos de bacterias que, a su vez, posee un registro de funciones conocidas (evidencia). Las categorías funcionales usadas en Faprotax son mas fáciles de interpretar pues apuntan a "procesos biológicos" conocidos. Faprotax utilza la tabla de conteo para tratar de relacionar la abundancia de cada categoría funcional en cada muestra.
2.1 Instalación
Descargar la última versión desde este link.
Descomprimir el archivo .zip descargado.
2.2 Archivos necesarios (input)
a) Tabla del conteo de OTUs o ASVs en cada muestra (de preferencia debe ser una tabla rarefaccionada). Esta tabla debe estar en un formato de columnas delimitadas por tabulaciones y requiere la información de taxonomía la cual debe estar en la última columna. El encabezado no debe tener el símbolo #. [tabla.con_taxonomia.tsv]
2.3 Ejecución
~/bin/instaladores/FAPROTAX_1.2.4/collapse_table.py -i tabla.con_taxonomia.tsv -o tabla.con_taxonomia.faprotax.tsv -g ~/bin/instaladores/FAPROTAX_1.2.4/FAPROTAX.txt -d "taxonomy" -v 1>1a2.txt 2>2a2.txt
-o : nombre del archivo donde se guardará el resultado del análisis
-g : base de datos de relación "taxonomía -> función" la cual viene incluída en Faprotax
-d : nombre de la columna en la cual se encuentra la taxonomía
-v : entrega en pantalla todos los detalles de la ejecución del análisis
2.4 Post-proceso de resultados
El archivo tabla.con_taxonomia.faprotax.tsv no requiere un postproceso
2.5 Análisis de resultados
El archivo tabla.con_taxonomia.faprotax.tsv se puede revisar directamente ya sea en R o en Excel.
3. FunGuild (sitio oficial)
Esta herramienta está diseñada sólo para Hongos (ITS o r18S).
Se utiliza sólo la información de taxonomía como "clave" para contrastar contra una base de datos de hongos que, a su vez, posee un registro de estilos de vida y alimentación de especies de hongos conocidos (evidencia).
3.1 Instalación
git clone https://github.com/UMNFuN/FUNGuild
cd FUNGuild/
python FUNGuild.py taxa -otu example/otu_table.txt -format tsv -column taxonomy -classifier unite
python FUNGuild.py guild -taxa example/otu_table.taxa.txt
3.2 Archivos necesarios (input)
a) Tabla del conteo de OTUs o ASVs en cada muestra (de preferencia debe ser una tabla rarefaccionada). Esta tabla debe estar en un formato de columnas delimitadas por tabulaciones y requiere la información de taxonomía la cual debe estar en la última columna. El encabezado no debe tener el símbolo #. [tabla.con_taxonomia.tsv]
3.3 Ejecución
python FUNGuild/FUNGuild.py taxa -otu tabla.con_taxonomia.tsv -format tsv -column taxonomy -classifier unite
python FUNGuild/FUNGuild.py guild -taxa tabla.con_taxonomia.tsv
Esto generará el archivo tabla.con_taxonomia.taxa.guilds.txt que contiene toda la info de resultados.
3.4 Post-proceso de resultados
El archivo tabla.con_taxonomia.taxa.guilds.txt no requiere un postproceso
3.5 Análisis de resultados
El archivo tabla.con_taxonomia.taxa.guilds.txt se puede revisar directamente ya sea en R o en Excel.