rendimiento multihebra de BLAST+ 2.6.0

Hola,
hace unas semanas descubrí gracias a mi colega Pablo Vinuesa que BLAST+ del NCBI iba ya por la versión 2.6.0. Cuando miré el resumen de cambios me llamó la atención que ya desde la versión 2.4 soporta un algoritmo multihebra para la fase reconstrucción hacia atrás del alineamiento, que en la literatura de programación dinámica se llama traceback. Dado que nosotros usamos con mucha frecuencia BLAST quise probar en qué se traducían estas novedades en tiempo de cálculo, dado que ya habíamos observado que BLASTP no paralelizaba bien, razón por la cual desarrollamos split_blast.pl, que recientemente comparamos contra DIAMOND.

El experimento consistió en buscar alineamientos locales de 48.588 secuencias de la variedad de cebada Haruna Nijo entre los 7.927 factores de transcripción de nuestra colección http://floresta.eead.csic.es/footprintdb:

$ ncbi-blast-2.2.30+/bin/makeblastdb -in footprintdb.18042017.tf.fasta -dbtype prot

$ time ~/soft/ncbi-blast-2.2.30+/bin/blastp -query HarunaNijo_proteins.fa \
  -db footprintdb.18042017.tf.fasta -outfmt 6 -max_target_seqs 10 \
  -num_threads 8 > HarunaNijo_proteins.2.2.30.blast

real  53m47.482s
user  122m2.375s
sys 0m9.749s

$ perl split_blast.pl 8 1000 ncbi-blast-2.2.30+/bin/blastp \

  -query HarunaNijo_proteins.fa -db footprintdb.18042017.tf.fasta -outfmt 6 \
  -max_target_seqs 10 -output HarunaNijo_proteins.split.blast

# runtime: 836 wallclock secs ( 0.71 usr  0.20 sys + 6391.38 cusr  5.54 csys = 6397.83 CPU)
# this is ~14m

$ ncbi-blast-2.6.0+/bin/makeblastdb -in footprintdb.18042017.tf.fasta -dbtype prot

$ time ncbi-blast-2.6.0+/bin/blastp -query HarunaNijo_proteins.fa \
  -db footprintdb.18042017.tf.fasta -outfmt 6 -max_target_seqs 10 \
  -num_threads 8 >   HarunaNijo_proteins.2.6.0.blast

real  20m35.969s
user  194m1.715s
sys 2m41.827s

Como podéis ver, al menos para BLASTP esta versión de BLAST+ supone una ganancia clara en procesadores multicore (8 en esta prueba), a costa de un aumento de tamaño del binario, que pasa de 31MB a 38MB, pero sigue siendo más lento que split_blast.pl,
hasta pronto,
Bruno

UniPROBE heterodimers in footprintDB

Hi,
we'd like to let you know that footprintDB, our database of transcription factors (TF), cognate binding sites and interface residues has been updated. Álvaro recently added the latest version of UniPROBE which systematically annotates some TF heterodimers. Their DNA motifs look like this:

Consensus DNA motif recognized by the dimer (MXL1 , MDL1), taken from footprintDB.

The non-redundant collection of motifs has been also updated in RSAT::Plants,
have a good week,
Bruno

ciencia básica = ciencia aplicada (fármacos aprobados)

Buenas,
hoy quisiera comentar un tema recurrente cuando hablamos de ciencia: el de ciencia básica vs ciencia aplicada. En general mi impresión es que la primera se percibe como un esfuerzo romántico que da sentido a la vida de algunos locos, mientras que la segunda es la que vale, puesto que participan ingenieros, y acaba llegando a nuestro smartphone en poco tiempo. Obviamente exagero, pero por ahí van los tiros.

La excusa para sacar este tema hoy es un artículo publicado recientemente en la revista Science (http://science.sciencemag.org/content/356/6333/78.full) donde se analiza cómo se citan los proyectos de Biomedicina financiados por los NIH norteamericanos en patentes. El estudio cubre los años entre 1980 y 2007. Las conclusiones de este trabajo son:

1) que si tenemos en cuenta las citas indirectas, es decir, patentes que citan artículos que a su vez citan proyectos NIH, finannciados con dinero público, hasta un 31% de proyectos en ese periodo son citados en patentes.

2) los proyectos citados en patentes que protegen fármacos aprobados por la FDA son en igual proporción "básicos" y "aplicados", según las definiciones de los autores, que son conscientes de lo resbaladizo de estos términos.

Figura tomada de http://science.sciencemag.org/content/356/6333/78.full

Hasta luego,
Bruno

Apuntes sobre ensamblaje de genomas de plantas

Buenas, ayer asistimos Ernesto Igartua y yo al 6th CNAG Symposium on Genome Research: Agrigenomics, organizado por el Centro Nacional de Análisis Genómico en Barcelona, donde a menudo contratamos servicios de secuenciación.

Allí presentamos nuestro trabajo con cebada, junto a otros colegas que trabajan en ganadería, piscicultura y agricultura y utilizan herramientas de la genómica contemporánea.

Como curiosidades me apunté que André Eggen, de Illumina, mencionó que comparando razas bovinas habían imputado SNPs mezclando genotipos de baja densidad (chips de ~10K SNPs), con genomas completos, alcanzando millones de SNPs. Por cierto, habían usado el software propietario DeNovoMAGIC para ensamblar genomas bovinos.

Otra cosa fue que los peces que estudian Franscesc Piferrer y su grupo tienen un mecanismo de metilación en función de la temperatura para controlar la producción de hormonas sexuales, algo que me recordó mucho a la memoria de vernalización en las plantas.

Pero además de estas charlas, y de visitar las salas de secuenciación y de servidores del CNAG, tuvimos dos sesiones casi seguidas donde repasamos los últimos métodos de ensamblaje y validación de genomas de plantas de la mano de Tyler Alioto y Gareth Linsmith. Éstas son mis notas:

Detección de contaminantes en las lecturas/reads: 
kraken : https://ccb.jhu.edu/software/kraken

Ensamblajes híbridos y diploides, combinando lecturas cortas y largas y estrategias más complejas para genomas de individuos heterocigotos.

• reads cortos, generalmente Illumina, de entre 100 y 300b, para alcanzar profunidades de al menos 30X en cada tipo de librería: 
• paired-end (PE) con insertos de por ejemplo 400 y 730pb 
• mate-pair (MP) con insertos de 4 y 8Kb para superar la longitud de la mayoría de secuencias repetidas

• reads largos, generalmente PacBio o de Oxford Nanopore. EN CNAG usan secuenciadores minIon para producir lecturas de 11.5Kb de media, alcanzando longitudes máximas > 100kb. Gareth comentó que en manzano necesitaron 60x, y eso que era material doble haploide. Este tipo de reads requieren consensos calculados con software como Sparc, Racon o Nanopolish.

En cuanto a ensambladores, Tyler destacó DISCOVAR de novo y Platanus, más adecuado para individuos con moderadas tasas de sitios heterocigotos. Pero advirtió del efecto negativo que tiene la heterocigosis sobre N50. En cambio, Gareth mencionó que primero ensambla las lecturas cortas con SOAPdenovo sin resolver las burbujas de Bruijn para luego luego combinar los reads largos con DBG2OLC y CANU.

Estrategias complementarias de ensamblaje
Datos de RNAseq para scaffolding con AGOUTI y Rascaf.

Pools de fósmidos como los empleados en el genoma de la ostra.
Mapas ópticos con enzimas nickasas que cortan cada 10Kb, con Bionano.
Dovetail genomics, aproximación basada en Hi-C.

Herramientas para corregir y finalizar genomas
PILON : https://github.com/broadinstitute/pilon/wiki
BESST : https://github.com/ksahlin/BESST

Estrategias para evaluar y validar genomas
Aparte del criterio clásico de sintenia respecto a especies cercanas, ambos mencionaron los problemas de evaluar un ensamblaje solamente por su N50 sin mirar por ejemplo los genes core anotados, por ejemplo con BUSCO, el sucesor de CEGMA. Gareth mencionó ALE para calcular la verosimilitud de un ensamblaje dadas las librerías de secuencias y KAT para comparar los k-meros originales de los reads con los del ensamblaje, que deberían coincidir, o para determinar la fracción de sitios heterocigotos:

Frecuencias de k-meros de los genotipos B73 y Mo17 de maíz, tomada de http://www.nature.com/articles/srep42444.

Casi se me olvida mencionar la comparación entre el mapa físico y el genético como criterio de calidad, muy útil en el genoma de manzano o en el de la cebada:

Comparación entre las posiciones de marcadores en una población de mapeo en cebada y sus posiciones en los mapas físico IBSC y POPSEQ de cebada, tomada de http://link.springer.com/article/10.1007%2Fs11032-015-0253-1.

Hasta  pronto,
Bruno

Tutorial: pan-genome analysis with GET_HOMOLOGUES

Hi,
a new tutorial on the analysis pan-genomes using GET_HOMOLOGUES and GET_HOMOLOGUES-EST is now available. After a short introduction, where the main concepts are illustrated, the remaining sections cover the installation and typical operations required to analyze and annotate genomes and transcriptomes from a pan-genome perspective, in which individuals or species contribute genetic material to a pool.

The examples include both bacterial sequences in GenBank format and plant transcripts. This tutorial has been created for a two-day workshop to be held at BIOS (Manizales, Colombia) next week, with title "From genomes to pangenomes: understanding variation among individuals and species":

The tutorial can be found at: http://digital.csic.es/handle/10261/146411 

Code, sample datasets and documentation are available at:
https://github.com/eead-csic-compbio/get_homologues

Suggestions and error reports are welcome,
Bruno