Analiza y anota archivos VCF

En esta página, se describe cómo realizar las siguientes acciones:

  • Analizar los campos de anotaciones en los archivos VCF de entrada
  • Anotar los archivos VCF mientras se importan, mediante bases de datos y herramientas de anotación de código abierto, como Variant Effect Predictor (VEP)

Analiza campos de anotación

La herramienta Variant Transforms admite un formato de campo de anotación que se especifica en el documento Variant annotations in VCF format (Anotaciones de variantes en formato VCF). VEP utiliza este mismo formato. En el siguiente ejemplo, se muestra un encabezado VCF de un archivo anotado mediante VEP:

##INFO=<ID=CSQ,Number=.,Type=String,Description="Consequence annotations from Ensembl VEP. Format: Allele|Consequence|IMPACT|SYMBOL|Gene|Feature_type|Feature|BIOTYPE|EXON|INTRON|ALLELE_NUM|...>

Las anotaciones correspondientes para un registro de variante de muestra en el que REF=AAC y ALT=TAC,A se ven así:

CSQ=-|downstream_gene_variant|MODIFIER|WASH7P|ENSG00000227232|Transcript|ENST00000423562|unprocessed_pseudogene|||2,...,T|upstream_gene_variant|MODIFIER|DDX11L1|ENSG00000223972|Transcript|ENST00000518655|transcribed_unprocessed_pseudogene|||1

Se proporcionan varias listas de anotaciones separadas por comas para cada registro de variante. El primer campo de cada conjunto corresponde a un alelo alternativo (ALT). En el ejemplo anterior, solo se utiliza un subconjunto para cada ALT. El archivo real tiene más campos, y cada alelo alternativo se puede repetir varias veces con diferentes anotaciones.

Puedes pasar el comando --annotation_fields cuando ejecutes la herramienta Variant Transforms, que divide las listas de anotación y las transfiere a sus registros alternate_bases correspondientes. Por ejemplo, a continuación, puedes ver que hay un campo de registro repetido, CSQ, en cada alternate_bases:

Nombre Tipo Modo Descripción
alternate_bases.CSQ RECORD REPETIDO Lista de anotaciones CSQ para esta alternativa
alternate_bases.CSQ.allele STRING Acepta valor NULL La parte ALT del campo de anotación
alternate_bases.CSQ.Consequence STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.IMPACT STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.SYMBOL STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.Gene STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.Feature_type STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.Feature STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.BIOTYPE STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.EXON STRING Acepta valor NULL Ninguna
alternate_bases.CSQ.INTRON STRING Acepta valor NULL Ninguna

Coincidencia de alelo alternativo

En el ejemplo anterior, los alelos alternativos enumerados en la columna CSQ son distintos de las strings ALT del registro de la variante. Por ejemplo, la variante AAC->C está representada por un -, que implica una eliminación. Esto sucede porque VEP se ejecuta en el conjunto de datos con la marca --minimal. Variant Transforms admite campos ALT coincidentes en este modo minimal.

Se admiten los siguientes modos de coincidencia de ALT:

  • Según la configuración predeterminada, si no se configura una marca de coincidencia ALT, se utiliza la especificación de Anotaciones de variantes en formato VCF.

  • Como se muestra en el ejemplo anterior, cuando se usa la marca --use_allele_num, la herramienta Variant Transforms busca el valor del campo de anotación ALLELE_NUM. Esta es una característica de VEP.

  • Puedes simular la marca --minimal con la herramienta Variant Transforms mediante la marca --minimal_vep_alt_matching. Sin embargo, esto puede generar coincidencias ambiguas. Estas coincidencias ambiguas se consideran y se informan en la Cloud Dataflow Console. También se agregan a un nuevo campo en la tabla creada.

    Para obtener mejores resultados, utiliza --minimal en lugar de --minimal_vep_alt_matching.

Anotación automática de archivos VCF con VEP

Por lo general, si deseas analizar los campos de anotación y dividirlos en campos separados, primero debes anotar previamente tus archivos VCF. Sin embargo, no necesitas hacerlo cuando utilizas la herramienta Variant Transforms.

La herramienta Variant Transforms puede anotar automáticamente los archivos VCF a medida que se cargan en BigQuery. Para hacerlo, debes ejecutar VEP con la herramienta Variant Transforms, lo que requiere crear y ejecutar instancias de máquina virtual de Compute Engine que contengan imágenes de Docker preconfiguradas.

Ejecuta VEP con la herramienta Variant Transforms

Para ejecutar la herramienta Variant Transforms con VEP y anotar automáticamente archivos VCF a medida que se cargan en BigQuery, debes ejecutar la herramienta Variant Transforms, al menos, con las siguientes marcas:

  • --run_annotation_pipeline
  • --max_num_workers MAX_NUM_WORKERS o --num_workers NUM_WORKERS
  • --annotation_output_dir gs://CLOUD_STORAGE_BUCKET

    El depósito que selecciones para la variable CLOUD_STORAGE_BUCKET debe ser propiedad de tu proyecto.

A continuación, se proporciona una lista completa de las marcas que puedes pasar a la herramienta Variant Transforms cuando ejecutas VEP:

Marca Valor predeterminado Descripción
--run_annotation_pipeline Falso Habilita la anotación automática.
--max_num_workers o --num_workers Determinado automáticamente La cantidad máxima de trabajadores en el trabajo a la que Cloud Dataflow ajustará su escala en forma automática o la cantidad de trabajadores a quienes se les asignó el trabajo desde un principio, respectivamente.
--annotation_output_dir Ninguno Una ruta de Cloud Storage para los archivos de salida de VEP. La jerarquía de los archivos de entrada se replica en esta ubicación, pero los archivos nuevos se agregan con _vep_output.vcf. Si el directorio ya existe, la herramienta Variant Transforms fallará.
--vep_image_uri gcr.io/gcp-variant-annotation/vep_91 La imagen de Docker de VEP para Cloud Genomics. Google mantiene la imagen predeterminada (versión 91 de VEP).
--vep_cache_path gs://gcp-variant-annotation-vep-cache/vep_cache_homo_sapiens_GRCh38_91.tar.gz (recomendado para un genoma humano alineado con una secuencia de referencia GRCh38) Una ruta de Cloud Storage para la versión comprimida de la caché de VEP. Puedes crear este archivo con la secuencia de comandos --build_vep_cache.sh.
--vep_info_field CSQ_VT El nombre de un nuevo campo INFO que contiene las anotaciones nuevas.
--vep_num_fork 2 La cantidad de procesos locales que se utilizan cuando se ejecuta VEP en un solo archivo. Consulta la marca --fork [num_forks] en la documentación de VEP para obtener más información.
--shard_variants True Determina si los archivos de entrada se fragmentan en archivos VCF temporales más pequeños antes de anotar los archivos VCF con VEP. Si los archivos de entrada son pequeños (por ejemplo, si cada archivo VCF contiene menos de 50,000 variantes), establecer esta marca en True puede desperdiciar recursos de procesamiento.
--number_of_variants_per_shard 2,000 Si shard_variants es True, la cantidad máxima de variantes se escribirá en cada fragmento. Si tienes un conjunto de datos con muchas muestras, te recomendamos usar un valor menor. Es posible que la canalización demore más tiempo en completarse cuando se proporciona un valor menor. En la mayoría de los casos, se recomienda el valor predeterminado.

Soluciona problemas

Si ejecutas VEP con la herramienta Variant Transforms y la herramienta falla, se crea un directorio de logs que proporciona registros de las máquinas virtuales donde se ejecutó el proceso. El directorio logs se especifica en la ruta utilizada para la marca --annotation_output_dir. Puedes inspeccionar estos archivos de registro, además de los registros estándares de la herramienta Variant Transforms, para determinar la causa del problema.

La marca --check_ref de VEP se habilita automáticamente cuando ejecutas VEP con la herramienta Variant Transforms. Como resultado, los archivos en el directorio de logs contendrán una lista de variantes que no coinciden con la referencia. Si estás utilizando la caché de VEP correcta, no deberían aparecer valores no coincidentes, porque estos se eliminan de forma automática y no aparecen en la tabla de salida de BigQuery.

Consultas de muestra

Extrae variantes de alto impacto

En el siguiente ejemplo, se muestra cómo extraer variantes de alto impacto en el gen HBB, que codifica la cadena beta de la hemoglobina:

#standardSQL
SELECT reference_name, start_position, reference_bases, ALT.alt, CSQ.*
FROM vcf_imports_external.gnomad_genomes_chr_hg19 AS T, T.alternate_bases AS ALT, ALT.CSQ AS CSQ
WHERE CSQ.SYMBOL = "HBB" AND CSQ.IMPACT = "HIGH"
ORDER BY start_position

Una vez que hayas ejecutado la consulta, se mostrará el siguiente resultado:

Fila reference_name start_position reference_bases alt alelo Consecuencia IMPACTO SÍMBOLO Gen Feature_type Función BIOTIPO EXON INTRON
1 11 5246947 G GC C frameshift_variant ALTO HBB ENSG00000244734 Transcripción ENST00000335295 protein_coding 3/3
2 11 5246957 T C C splice_acceptor_variant ALTO HBB ENSG00000244734 Transcripción ENST00000335295 proten_coding 2/2
3 11 5246957 T C C splice_acceptor_variant&non_coding_transcript_variant ALTO HBB ENSG00000244734 Transcripción ENST00000475226 retained_intron 1/1
4 11 5247805 C T T splice_donor_variant ALTO HBB ENSG00000244734 Transcripción ENST00000335295 protein_coding 2/2
5 11 5247805 C T T splice_donor_variant&non_coding_transcript_variant ALTO HBB ENSG00000244734 Transcripción ENST00000475226 retained_intron 1/1
6 11 5247991 CAAAG C - frameshift_variant ALTO HBB ENSG00000244734 Transcripción ENST00000335295 protein_coding 2/3
7 11 5247991 CAAAG C - frameshift_variant ALTO HBB ENSG00000244734 Transcripción ENST00000380315 protein_coding 4/4
8 11 5248003 G A A stop_gained ALTO HBB ENSG00000244734 Transcripción ENST00000335295 protein_coding 2/3
9 11 5248003 G A A stop_gained ALTO HBB ENSG00000244734 Transcripción ENST00000380315 protein_coding 4/4
10 11 5248199 T A A stop_gained ALTO HBB ENSG00000244734 Transcripción ENST00000335295 protein_coding 1/3
11 11 5248199 T A A stop_gained ALTO HBB ENSG00000244734 Transcripción ENST00000380315 protein_coding 3/4

El conjunto de datos gnomAD no tiene llamadas ni muestras, pero puedes ampliar esta consulta para seleccionar muestras que tengan variantes de alto impacto.

Une datos con bases de datos públicas de BigQuery

En el siguiente ejemplo, se muestra cómo unir datos con otras bases de datos, como bases de datos de rutas, que están disponibles en BigQuery. Algunas de las bases de datos disponibles a través del ISB-CGC son las siguientes:

En el siguiente ejemplo, se muestra cómo extraer variantes de alto impacto en la ruta de reparación por recombinación, tal como se define en los procesos biológicos de GO (ontología génica) (GO:0006302):

#standardSQL
SELECT
  reference_name, start_position, reference_bases, ALT.alt,
  CSQ.Consequence, CSQ.Impact, CSQ.SYMBOL
FROM
  `vcf_imports_external.gnomad_genomes_chr_hg19` AS T,
  T.alternate_bases AS ALT, ALT.CSQ AS CSQ
WHERE
  # Note: Matching based on symbol is "best effort" as the names may not be
  # standardized across sources.
  CSQ.SYMBOL IN (SELECT DB_Object_Symbol
                 FROM `isb-cgc.genome_reference.GO_Annotations`
                 WHERE GO_ID = 'GO:0006302')
  AND CSQ.IMPACT = "HIGH"
ORDER BY
  start_position

Una vez que hayas ejecutado la consulta, se mostrará el siguiente resultado:

Fila reference_name start_position reference_bases alt Consecuencia IMPACTO SÍMBOLO
1 5 917449 G T splice_acceptor_variant&non_coding_transcript_variant ALTO TRIP13
2 12 1022568 A C stop_gained ALTO RAD52
3 12 1022568 A C stop_gained ALTO RAD52
4 12 1022568 A C stop_gained ALTO RAD52
5 12 1023125 G A stop_gained ALTO RAD52
6 12 1023125 G A stop_gained ALTO RAD52
7 12 1023125 G A stop_gained ALTO RAD52
8 12 1023167 G A stop_gained ALTO RAD52
9 12 1023167 G A stop_gained ALTO RAD52
10 12 1023167 G A stop_gained ALTO RAD52
11 12 1023217 G T stop_gained ALTO RAD52
12 12 1023217 G T stop_gained ALTO RAD52
13 12 1023217 G T stop_gained ALTO RAD52
14 12 1025654 TGA T frameshift_variant ALTO RAD52
15 12 1025699 T G splice_acceptor_variant&non_coding_transcript_variant ALTO RAD52
16 12 1036004 GC G frameshift_variant ALTO RAD52
17 12 1038977 C CT frameshift_variant ALTO RAD52