Generación de árboles filogenéticos UPGMA a partir de un código de proteína, en Python

En este post se presenta una aplicación para la generación de árboles filogenéticos UPGMA a partir de un código de proteína en formato FASTA.

La aplicación se puede descargar de:

https://sites.google.com/site/elviajedelnavegante/codigo/FILO_TREE.zip?attredirects=0&d=1

Los árboles filogenéticos son útiles para encontrar los familiares más cercanos a nuestro organismo de estudio. Se ha implementado mediante la creación de un script Python, que utiliza la biblioteca externa BioPython para la creación de árboles filogenéticos. La aplicación permite obtener secuencias similares a partir de una secuencia de proteína dada. Estas secuencias obtenidas por Blast son procesadas con un algoritmo de alineamiento múltiple (Clustal), cuyo resultado sirve para crear la matriz de distancias y el árbol filogenético UPGMA. La aplicación es multiplataforma, para sistemas Microsoft Windows y Linux.

La configuración permite seleccionar cuántas secuencias como máximo se quiere utilizar y el umbral de similitud para seleccionar las secuencias. El script se puede lanzar en modo atendido (menú de aplicación) y desatendido (para usar en procesos batch de procesamiento secuencial). Así tendremos las siguientes opciones de inicio de aplicación:

python filo_tree.py

Muestra ayuda de funcionamiento de la aplicación.

python filo_tree.py -g

Lanza el siguiente menú de aplicación:

Las opciones de aplicación son las siguientes:

• CARGAR PROTEÍNA: Carga datos de proteína en formato FASTA a partir de un fichero.
• VISUALIZAR PROTEÍNA: Visualiza la proteína cargada en formato FASTA.
• BUSCAR SECUENCIAS SIMILARES: A partir de la proteína cargada se realiza una búsqueda de secuencias similares mediante el algoritmo blastp en una base de datos online. En esta opción se escoge tanto el número de secuencias a devolver como la base de datos a utilizar. Se genera un fichero XML con dichas secuencias.
• VISUALIZAR SECUENCIAS SIMILARES: Se visualizan por pantalla las secuencias generadas por blastp que están guardadas en el fichero XML.
• GENERAR ÁRBOL FILOGENÉTICO UPGMA: A partir de las secuencias guardadas en el fichero XML se generan las secuencias de alineamiento múltiple, con el fin de crear la matriz de distancias y poder crear el árbol filogenético UPGMA. En esta opción se escoge el umbral de similitud.
• ACERCA DE: Créditos de la aplicación.
• SALIR: Salir de la aplicación.

python filo_tree.py -p FICHERO.FASTA -nseq N_SEQ -e UMBRAL_SIMILITUD -bd N_BD

Lanza el script filo_tree.py, buscando secuencias similares a la secuencia proteica definida en formato FASTA en el fichero FICHERO.FASTA. El número de secuencias a buscar está definido por N_SEQ y el umbral de similitud por UMBRAL_SIMILITUD. Además se tiene que especificar la base de datos a la que conectarse para realizar la búsqueda, donde N_BD es un número entero correspondiente con:

1. nr: Non-redundant GenBank CDS translations + RefSeq + PDB + SwissProt + PIR + PRF, excluding those in PAT, TSA, and evr_nr.
2. refseq_protein: Protein sequences from NCBI Reference Sequence project.
3. swissprot: Las major release of the UniProtKB/SWISS-PROT protein sequence database (no incremental updates).
4. pat: Proteins from the Patent division of GenBank.
5. pdb: Protein sequences from the 3-dimensional structure records from the Protein Data Bank.
6. env_nr: Protein sequences translated from the CDS annotation of metagenomic nicleotide sequences.
7. tsa_nr: Protein sequences translated from CDSs annotated on trascriptome shotgun assemblies

El script realiza alineamiento múltiple de las secuencias obtenidas mediante Clustal, generando la matriz de distancias y árbol filogenético UPGMA.

Ejemplo: Búsqueda de proteína contenida en H0XFT0.fasta en la base de datos de swissprot, recuperando 50 secuencias con e-value 1.

python filo_tree.py -p H0XFT0.fasta -nseq 50 -e 1 -bd 3

Devuelve lo siguiente:

Y el siguiente fichero que contiene la matriz de distancias y el árbol filogenético:

https://sites.google.com/site/elviajedelnavegante/documentacion/TREE_UPGMA_20150523121542810823.TXT?attredirects=0&d=1

Esquema general de la aplicación

El esquema general de la aplicación consiste en:

1. Pasar como argumento (parámetro) un fichero con una secuencia de proteína en formato FASTA, al algoritmo BLAST (implementado en una función dentro del paquete BioPython). Dicha función busca un número de secuencias determinado y las alinea por pares a partir de la secuencia FASTA dada. Es decir, alinea N secuencias a partir de una (no es alineamiento múltiple, es alineamiento por pares). Esta información (la secuencia FASTA original, la secuencia similar y el alineamiento, por cada secuencia similar encontrada) se guarda en un fichero XML.
2. A partir del fichero BLAST_FILE.XML se obtienen las secuencias similares (no los alineamientos),y se guardan en un fichero en formato FASTA (FASTA_MULTIPLE.FASTA). Por tanto ahora tenemos las N secuencias similares en un formato adecuado, listas para un alineamiento múltiple.
3. Se realiza un alineamiento múltiple a partir del algoritmo de Clustal por medio de la aplicación externa clustalw2, que se puede obtener de http://www.clustal.org/. El script Python busca la aplicación en el PATH y lanza el programa (algoritmo), guardando la salida (el alineamiento múltiple) en el fichero FASTA_MULTIPLE.aln.
4. El fichero de alineamiento múltiple es el parámetro de entrada para el cálculo de la matriz de distancias (implementado en la clase DistanceCalculator de BioPython). Los datos de la matriz de distancia se vuelcan o a un fichero de salida TREE_UPGMA_fecha_hora.TXT o por pantalla. 
   
5. La variable DISTANCE MATRIX (matriz_distancia) es a su vez el parámetro de entrada del método upgma del objeto producto de la instanciación de la clase DistanceTreeConstructor para la construcción del árbol UPGMA:
   
6. Una vez calculado el árbol filogenético UPGMA ya se pueden obtener todos los datos, volcando a fichero (TREE_UPGMA_fecha_hora.TXT) o pantalla los datos propios del árbol, la representación del árbol en ASCII (mediante la función Phylo.draw_ascii de BioPython) y la generación de un fichero gráfico PNG (TREE_UPGMA_fecha_hora.PNG) a partir de Phylo.draw_graphviz de BioPython. Opcionalmente se puede visualizar por pantalla mediante Phylo.Show() de BioPython.

Parámetro de entrada - Fichero FASTA 

La aplicación tiene que tener como parámetro de entrada un fichero con un código de proteína en formato FASTA. Por ejemplo el fichero H0XFT0.fasta (adjunto con el script a modo de ejemplo):

Este fichero se puede obtener desde el portal web de UNIPROT http://www.uniprot.org/.

Saludos.