la invisible ciencia básica detrás de las vacunas SARS-CoV-2

Hola,

la pandemia que estamos viviendo, un año después, nos está poniendo a prueba. A pesar de la improvisación de los políticos a escala global, de la desinformación en las redes sociales y los brotes de desconfianza, a pesar de la economía bajo mínimos y del cole en casa, ahora la mayoría tenemos la esperanza de que las vacunas resuelvan el problema. 

En este artículo solamente pretendo recordar el largo camino de la investigación básica que nos ha traído hasta el presente. Estas vacunas son ya grandes hitos de la humanidad, a la altura de la llegada a la luna, pero el camino ha sido largo, de al menos 25 años. Por tanto, nada de milagros, son el fruto de mucho trabajo acumulado que fue explotado con mucho éxito por empresas como BioNTech y Moderna. Como pasó con CRISPR para la edicion de genomas, para que alguien llegara a la cima fueron necesarios muchos pasos previos, muchos de los cuales fuera de contexto serían objeto de "y eso para qué sirve". Aquí enumero los más importantes para las vacunas de ARN, extraídos de este hilo (no están los de los otros tipos de vacunas):

1970 T7 ARN polimerasa: https://nature.com/articles/228227a0 . Esta enzima permite sintetizar moléculas de ARN a medida, como las de las vacunas.

1978 Liposomas para llevar ARN mensajeros (mRNA): https://nature.com/articles/274923a0 . Estos vehículos permiten que el ARN de las vacunas pueda atravesar
la doble capa lipídica de la membrana celular.

1990 Inyecciones de ADN y ARN para expresar genes de manera transitoria en tejidos: https://science.sciencemag.org/content/247/4949/1465.long .  Este trabajo demostró que es posible expresar genes a medida tras ser inyectados en tejidos, de manera que se traducen como proteínas.

2005 Ribonucleótidos modificados no disparan respuestas inmunes: https://cell.com/immunity/fulltext/S1074-7613(05)00211-6 .
Esto permite que el ARN inyectado no desencadene una reacción inmune por si mismo, lo que se pretende es que la reacción la desencadene la proteína codificada por ese ARNm.

2017 Estabilización de proteínas expuestas de los coronavirus MERS-CoV y SARS-CoV: https://pnas.org/content/114/35. Esto permite que la proteína modificada del coronoavirus que expresa el ARNm sea más estable y desencadene una reacción inmune más robusta.

Esta sucesión de descubrimientos y la tecnología actual permitieron que el tiempo de desarrollo de las vacunas haya sido el más corto de la historia:

 

En mi actual institución, el Instituto Europeo de Bioinformática (EBI), también hemos contribuído con el https://www.covid19dataportal.org y el navegador de genomas https://covid-19.ensembl.org

Espero haberos convencido de que la ciencia financiada con fondos públicos y transparente, muchas veces invisible e ingrata, es una parte fundamental e integral del avance de nuestra sociedad, y que de ella beben las empresas que nos venden luego los productos.

Hasta pronto,

Bruno

PD1 Me preguntan si "mucha de esa investigación no la financian las farmacéuticas de forma privada, con escasos beneficios en mucho casos. Esto me recuerda a cuando proponen nacionalizar farmacéuticas". 

 

Mi respuesta es que mi argumento no iba en esa dirección. De hecho creo que esta historia es un buen ejemplo de cómo lo iniciativa privada puede agilizar el desarrollo y las pruebas clínicas de las vacunas llevando a buen puerto en tiempo récord lo que empezó solamente como investigación básica. Mi argumento es que el riesgo es que la opinión pública, y los políticos que deciden, se queden sólo con el final del proceso y no vea la utilidad de la investigación básica. Quién les iba a decir a los autores de los artículos citados que iban a ser instrumentales para dos vacunas en el año 2020?

 

En cuanto a la financiación, como no es mi campo no me atrevo a opinar, pero sí puedo copiar aquí los que dicen al respecto los artículos citados más arriba:

 

1970 "This investigation was supported by a US Public Health Service research grant and training grant from the Institute of General Medical Sciences."

 

1978 "I thank Dr J. R. Tata for the freedom to pursue my own research goals... G.J.D. is in receipt of an EMBO long-term fellowship."

 

1990 "Supported in part by the NIH (grant numbers HD00669-05 and HD03352) and the Lynn F. Taylor Memorial Fund."

 

2005 "This work was supported by National Institutes of Health grants AI060505, AI50484, and DE14825."

 

2017 "This work was supported by Grants P20GM113132 ... and R01AI127521 ..., NIH Contract HHSN261200800001E Agreement 6x142 ..., and intramural funding from National Institute of Allergy and Infectious Diseases to support work at the VRC. Argonne is operated by UChicago Argonne, LLC, for the US Department of Energy (DOE), Office of Biological and Environmental Research under Contract DE-AC02-06CH11357. Use of the Stanford Synchrotron Radiation Lightsource (SSRL), SLAC National Accelerator Laboratory, is supported by the DOE, Office of Science, Office of Basic Energy Sciences under Contract DE-AC02-76SF00515. The SSRL Structural Molecular Biology Program is supported by the DOE Office of Biological and Environmental Research and by the NIH, National Institute of General Medical Sciences (including P41GM103393)."

 

PD2 Parece ser que Uğur Şahin, uno de los dos fundadores de BioNTech, dirige actualmente un proyecto ERC financiado con dinero público de la UE: https://twitter.com/ERC_Research/status/1372962936856190982

contenido GC de un fichero FASTA

Hola,

una pregunta habitual cuando analizas un fichero de nucleótidos, por ejemplo un ensamblaje de un genoma, es qué porcentaje GC tiene.  Asumiendo que el fichero está en formato FASTA, podemos obtener fácilmente ese valor con un mini-programa (one-liner) escrito en lenguaje perl. Por ejemplo, para el genoma comprimido de Brachypodium distachyon obtenido de Ensembl Plants, podríamos obtenerlo así:

zcat Brachypodium_distachyon.Brachypodium_distachyon_v3.0.dna.toplevel.fa.gz | \
   perl -lne 'if(!/^>/){ $SQ=uc($_); while($SQ =~ /([ACTG])/g){ $stat{$1}++; $tot++ } } 
   END{ printf("%%GC=%1.1f\n",100*($stat{"G"}+$stat{"C"})/$tot);  
      foreach $nt (keys(%stat)){ print "$nt\t$stat{$nt}" } }'

%GC=46.4
A	72549289
T	72561114
C	62839311
G	62789747

Si quieres calcular el %GC solamente para ciertas regiones del genoma entonces puedes codificarlas en un fichero BED y usar bedtools nuc, como se explica en https://www.biostars.org/p/47047

Hasta pronto,

Bruno

mapeando genes con liftoff

Hola,

esta semana en Ensembl tuvimos que probar a mapear genes de un ensamblaje genómico de arroz sobre otro ensamblaje más reciente. Esto se llama lift-over en la literatura. Para ello probamos un software que se llama Liftoff, publicado en Bioinformatics, con código fuente en https://github.com/agshumate/Liftoff

Como resume la figura para el tránscrito humano ENST00000598723.5, Liftoff calcula por medio de minimap2 alineamientos parciales al genoma de referencia y luego calcula el grafo más corto que conecte los alineamientos:

Resumo ahora cómo instalé este software escrito en python y cómo lo probé.

## pyenv 

# set up a pyenv folder
export PYENV_ROOT="$HOME/.pyenv"
 
curl https://pyenv.run | bash
 
pyenv install 3.7.9
 
# create dedicated env
pyenv virtualenv 3.7.9 liftoff

# install Liftoff inside that environment
pyenv shell 3.7.9 liftoff
pip install --upgrade pip
pip install Liftoff

# add the following to your .bashrc:
export PATH="$PYENV_ROOT/bin:$PATH"
eval "$(pyenv init -)"
eval "$(pyenv virtualenv-init -)"

## Test run:

pyenv shell 3.7.9 liftoff
liftoff -g old.gff3 -o old.new.liftoff.gff -p 4 new.genome.fna old.genome.fasta

 

El fichero de salida es un GFF con contenido como éste:

1       Liftoff gene    2903    10817   .       +       .       ID=LOC_Os01g01010;Name=LOC_Os01g01010;Note=TBC domain containing protein, expressed;coverage=1.0;sequence_ID=1.0;extra_copy_number=0;copy_num_ID=LOC_Os01g01010_0
1       Liftoff mRNA    2903    10817   .       +       .       ID=LOC_Os01g01010.1;Name=LOC_Os01g01010.1;Parent=LOC_Os01g01010;extra_copy_number=0
1       Liftoff exon    2903    3268    .       +       .       ID=LOC_Os01g01010.1:exon_1;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    3354    3616    .       +       .       ID=LOC_Os01g01010.1:exon_2;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    4357    4455    .       +       .       ID=LOC_Os01g01010.1:exon_3;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    5457    5560    .       +       .       ID=LOC_Os01g01010.1:exon_4;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    7136    7944    .       +       .       ID=LOC_Os01g01010.1:exon_5;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    8028    8150    .       +       .       ID=LOC_Os01g01010.1:exon_6;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    8232    8320    .       +       .       ID=LOC_Os01g01010.1:exon_7;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    8408    8608    .       +       .       ID=LOC_Os01g01010.1:exon_8;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    9210    9617    .       +       .       ID=LOC_Os01g01010.1:exon_9;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    10104   10187   .       +       .       ID=LOC_Os01g01010.1:exon_10;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    10274   10430   .       +       .       ID=LOC_Os01g01010.1:exon_11;Parent=LOC_Os01g01010.1;extra_copy_number=0
1       Liftoff exon    10504   10817   .       +       .       ID=LOC_Os01g01010.1:exon_12;Parent=LOC_Os01g01010.1;extra_copy_number=0

Hasta pronto, 

Bruno

eliminar redundancia en grandes ficheros de nucleótidos (linclust)

Hola,

recientemente me vi en la necesidad de eliminar redundancia de un fichero con millones de secuencias de nucleótidos. En mi caso se trataba de secuencias repetidas del genoma de cebada (n=4.638.834), y lo intenté con CD-HIT-EST, una herramienta que he probado muchas veces, y que tenía por muy eficiente.

Sin embargo, para esta tarea, a pesar de reservar más de 10GB de RAM, no terminá en varias horas usando 20 cores.

Por tanto, me puse a buscar opciones:

• dnaclust, usa demasiada RAM
• SigClust, muy rápido, pero usa el algoritmos K-medias y por tanto le tienes que decir cuántos clusters quieres, que en mi caso es lo que quiero averiguar
• swarm, es demasiado fino separando, solamente es eficiente para eliminar copias idénticas

Mi colega Pablo Vinuesa me habló de otras opciones que se ocupan de calcular distancias entre genomas o metagenomas, un problema relacionado:

• mash
• dRep
  

Después de mucho buscar encontré linclust, cuyo código fuente está en https://github.com/soedinglab/MMseqs2. Este programa se describió originalmente para secuencias de péptidos (artículo) pero los autores, entre ellos J Soding, añadieron después la posibilidad de agregar nucleótidos. En mis manos este es la mejor opción ahora mismo.

A continuación os muestra un banco de pruebas con secuencias repetidas de la planta Arabidopsis thaliana (n=66.752), analizadas con el binario más rápido distribuido por los autores (AVX2).

El programa lo ejecuté de esta manera:

$ command time -v mmseqs/bin/mmseqs easy-linclust \

  arabidopsis_thaliana.repeats.nondeg.fasta --threads 4 \

  arabidopsis_thaliana.48.repeats.nr0.95 ./ --min-seq-id 0.95

Verás que guarda datos temporales en el directorio actual (./) que luego debes eliminar.

Lo que observé en mis pruebas con cebada y A. thaliana es que a partir de un umbral de identidad del 70% el programa converge:

umbral	RAM (kbytes)	segundos	secuencias
0.99	311164	4.31	58285
0.95	303764	4.22	54353
0.90	307804	4.10	51359
0.80	306648	4.05	48722
0.70	304736	3.85	47770
0.50	307140	3.34	46433

Hasta pronto,

Bruno

Consultas GraphQL al Protein Data Bank

Hola, 

durante el mantenimiento anual de footprintDB he descubierto que el Protein Data Bank ha cerrado su interfaz de web servicios (WS) https://www.rcsb.org/pages/webservices/rest-fetch .  Leyendo en https://data.rcsb.org/index.html#data-api veo que ahora se pueden hacer consultas REST o GraphQL, y opté por la segunda por aprender un poco.

Hay documentación en https://data.rcsb.org/migration-guide.html#legacy-fetch-api para migrar las consultas WS antiguas a GraphQL. Las nuevas consultas tienen este aspecto:

query={
 entry(entry_id:"9ANT"{
  polymer_entities{
   rcsb_polymer_entity{
    pdbx_description
   }
   rcsb_entity_source_organism{
    scientific_name
   }
   rcsb_polymer_entity_container_identifiers{
    entry_id
    auth_asym_ids
    reference_sequence_identifiers{
     database_accession,
     database_name
    }
   }
  }
  struct{
   title
  }
  rcsb_primary_citation{
   pdbx_database_id_PubMed
  }
 }
}

Y se pueden convertir a una URL como https://data.rcsb.org/graphql?query={entry(entry_id:%229ANT%22){polymer_entities{rcsb_polymer_entity{pdbx_description}rcsb_entity_source_organism{scientific_name},rcsb_polymer_entity_container_identifiers{entry_id,auth_asym_ids,reference_sequence_identifiers{database_accession,database_name}}}struct{title}rcsb_primary_citation{pdbx_database_id_PubMed}}}

 

Si haces click verás que en Chrome o Firefox no obtienes el resultado esperado, con un error:

Invalid character found in the request target. The valid characters are defined in RFC 7230 and RFC 3986

Creo que se debe a un problema del servidor Tomcat del PDB. Sin embargo, sí funciona con wget y así obtienes la salida JSON que puedes procesar con tu lenguaje favorito:

wget 'https://data.rcsb.org/graphql?query={entry(entry_id:%229ANT%22){polymer_entities{rcsb_polymer_entity{pdbx_description}rcsb_entity_source_organism{scientific_name},rcsb_polymer_entity_container_identifiers{entry_id,auth_asym_ids,reference_sequence_identifiers{database_accession,database_name}}}struct{title}rcsb_primary_citation{pdbx_database_id_PubMed}}}' -O-

Hasta pronto,

Bruno