Cette page a été traduite par l'API Cloud Translation.
Switch to English

Analyser et annoter les fichiers VCF

Cette page décrit comment effectuer les opérations suivantes :

  • Analyser les champs d'annotation dans les fichiers VCF d'entrée
  • Annoter des fichiers VCF au cours de leur importation, à l'aide d'outils d'annotation Open Source et de bases de données telles que Variant Effect Predictor (VEP)

Analyser les champs d'annotation

L'outil Variant Transforms accepte un format de champ d'annotation spécifié dans les annotations de variantes au format VCF. VEP utilise le même format. L'exemple suivant illustre un en-tête VCF provenant d'un fichier annoté par 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|...>

Les annotations correspondant à un exemple d'enregistrement de variante où REF=AAC et ALT=TAC,A se présentent comme suit :

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

Plusieurs listes d'annotations séparées par des virgules sont fournies pour chaque enregistrement de variante. Le premier champ de chaque ensemble correspond à un allèle alternatif (ALT). L'exemple ci-dessus n'utilise qu'un sous-ensemble pour chaque allèle ALT. Le fichier réel contient davantage de champs et chaque allèle alternatif peut être répété plusieurs fois avec des annotations différentes.

Lorsque vous exécutez l'outil Variant Transforms, vous pouvez transmettre l'option --annotation_fields qui divise les listes d'annotations et les déplace sous les enregistrements alternate_bases correspondants. Par exemple, vous pouvez voir ci-dessous qu'il existe un champ d'enregistrement répété, CSQ, sous chaque enregistrement alternate_bases :

Nom Type Mode Description
alternate_bases.CSQ RECORD RÉPÉTÉ Liste des annotations CSQ pour cette alternative
alternate_bases.CSQ.allele STRING PEUT ÊTRE VIDE Partie ALT du champ d'annotation
alternate_bases.CSQ.Consequence STRING Peut être vide None
alternate_bases.CSQ.IMPACT STRING Peut être vide None
alternate_bases.CSQ.SYMBOL STRING Peut être vide None
alternate_bases.CSQ.Gene STRING Peut être vide None
alternate_bases.CSQ.Feature_type STRING Peut être vide None
alternate_bases.CSQ.Feature STRING Peut être vide None
alternate_bases.CSQ.BIOTYPE STRING Peut être vide None
alternate_bases.CSQ.EXON STRING Peut être vide None
alternate_bases.CSQ.INTRON STRING Peut être vide Aucune

Correspondance des allèles alternatifs

Dans l'exemple ci-dessus, les allèles alternatifs répertoriés dans la colonne CSQ sont différents des chaînes ALT de l'enregistrement de variante. Par exemple, la variante AAC->C est représentée par un signe -, qui correspond à une suppression. Ce comportement est dû au fait que VEP s'exécute sur l'ensemble de données avec l'option --minimal. L'outil Variant Transforms permet la mise en correspondance des champs ALT dans ce mode minimal.

Les modes acceptés pour la mise en correspondance des champs ALT sont les suivants :

  • Par défaut, si aucune option de correspondance ALT n'est définie, la spécification des annotations de variantes au format VCF est utilisée.

  • Comme le montre l'exemple ci-dessus, lorsque vous utilisez l'option --use_allele_num, l'outil Variant Transforms recherche la valeur du champ d'annotation ALLELE_NUM. Il s'agit d'une fonctionnalité de VEP.

  • Vous pouvez simuler l'option --minimal avec l'outil Variant Transforms à l'aide de l'option --minimal_vep_alt_matching. Toutefois, cette méthode peut entraîner des correspondances ambiguës, faisant l'objet d'un décompte et consignées dans la console Dataflow. Elles sont également ajoutées sous un nouveau champ dans la table créée.

    Pour des résultats optimaux, utilisez --minimal plutôt que --minimal_vep_alt_matching.

Annoter automatiquement des fichiers VCF à l'aide de VEP

En règle générale, si vous souhaitez analyser les champs d'annotation et les diviser en champs distincts, vous devez commencer par préannoter vos fichiers VCF. Cependant, cette opération n'est pas nécessaire lorsque vous utilisez l'outil Variant Transforms.

Variant Transforms peut annoter automatiquement les fichiers VCF lors de leur chargement dans BigQuery. Pour ce faire, vous devez exécuter VEP avec l'outil Variant Transforms, ce qui nécessite de créer et d'exécuter des instances de VM Compute Engine contenant des images Docker préconfigurées.

Exécuter VEP avec l'outil Variant Transforms

Pour exécuter l'outil Variant Transforms avec VEP et annoter automatiquement les fichiers VCF lors de leur chargement dans BigQuery, vous devez exécuter l'outil Variant Transforms avec au moins les indicateurs suivants :

  • --run_annotation_pipeline
  • --max_num_workers MAX_NUM_WORKERS ou --num_workers NUM_WORKERS

  • --annotation_output_dir gs://CLOUD_STORAGE_BUCKET

    Le bucket que vous choisissez pour la variable CLOUD_STORAGE_BUCKET doit appartenir à votre projet.

Vous trouverez ci-dessous la liste complète des indicateurs que vous pouvez transmettre à l'outil Variant Transforms lors de l'exécution de VEP :

Indicateur Valeur par défaut Description
--run_annotation_pipeline False Active l'annotation automatique.
--max_num_workers ou --num_workers Défini automatiquement Respectivement, le nombre maximal de nœuds de calcul associés à la tâche et autorisés pour l'autoscaling de Cloud Dataflow, ou le nombre de nœuds de calcul initialement affectés à la tâche.
--annotation_output_dir Aucune Une chemin d'accès Cloud Storage pour les fichiers de sortie de VEP. La hiérarchie des fichiers d'entrée est dupliquée à cet emplacement, mais les nouveaux fichiers sont suivis de _vep_output.vcf. Si le répertoire existe déjà, l'outil Variant Transforms échoue.
--vep_image_uri gcr.io/gcp-variant-annotation/vep_91 L'image Docker VEP Cloud Genomics. Google assure la maintenance de l'image par défaut (VEP version 91).
--vep_cache_path gs://gcp-variant-annotation-vep-cache/vep_cache_homo_sapiens_GRCh38_91.tar.gz (recommandé pour un génome humain aligné sur une séquence de référence GRCh38) Un chemin d'accès Cloud Storage pour la version compressée du cache VEP. Vous pouvez créer ce fichier à l'aide du script --build_vep_cache.sh.
--vep_info_field CSQ_VT Le nom d'un nouveau champ INFO contenant les nouvelles annotations.
--vep_num_fork 2 Le nombre de processus locaux utilisés lors de l'exécution de VEP sur un fichier individuel. Pour en savoir plus, référez-vous à l'option --fork [num_forks] dans la documentation VEP.
--shard_variants Vrai Détermine si les fichiers d'entrée sont segmentés en fichiers VCF temporaires plus petits avant annotation des fichiers VCF par VEP. Si les fichiers d'entrée sont petits (par exemple, si chaque fichier VCF contient moins de 50 000 variantes), définir cette option sur True risque de gaspiller des ressources de calcul.
--number_of_variants_per_shard 2 000 Si la valeur de shard_variants est True, le nombre maximal de variantes est écrit pour chaque segment. Si votre ensemble de données contient de nombreux échantillons, il peut être judicieux d'utiliser une valeur plus faible. Une valeur plus faible peut toutefois conduire à allonger le délai de traitement du pipeline. La valeur par défaut est recommandée dans la plupart des cas.

Dépannage

Si vous exécutez VEP avec l'outil Variant Transforms et que celui-ci échoue, un répertoire de journaux logs est créé. Il fournit les journaux des VM sur lesquelles le processus a été exécuté. Le répertoire logs est spécifié dans le chemin renseigné pour l'option --annotation_output_dir. Vous pouvez examiner ces fichiers journaux en plus des journaux standards de l'outil Variant Transforms afin de déterminer la cause du problème.

L'option VEP --check_ref est automatiquement activée lorsque vous exécutez VEP avec l'outil Variant Transforms. Par conséquent, les fichiers du répertoire logs contiennent une liste des variantes qui ne correspondent pas à la référence. Si vous utilisez le cache VEP adéquat, aucune incohérence ne doit apparaître, car celles-ci sont automatiquement supprimées et n'apparaissent pas dans la table de résultats BigQuery.

Exemples de requêtes

Extraire des variantes à fort impact

L'exemple suivant montre comment extraire des variantes à fort impact au niveau du gène HBB, qui code la chaîne bêta de l'hémoglobine :

#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

L'exécution de la requête renvoie les éléments suivants :

Ligne reference_name start_position reference_bases alt allele Effet IMPACT SYMBOLE Gêne Feature_type Fonctionnalité BIOTYPE EXON INTRON
1 11 5246947 G GC C frameshift_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000335295 protein_coding 3/3
2 11 5246957 T C C splice_acceptor_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000335295 protein_coding 2/2
3 11 5246957 T C C splice_acceptor_variant&non_coding_transcript_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000475226 retained_intron 1/1
4 11 5247805 C T T splice_donor_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000335295 protein_coding 2/2
5 11 5247805 C T T splice_donor_variant&non_coding_transcript_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000475226 retained_intron 1/1
6 11 5247991 CAAAG C - frameshift_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000335295 protein_coding 2/3
7 11 5247991 CAAAG C - frameshift_variant ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000380315 protein_coding 4/4
8 11 5248003 G A A stop_gained ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000335295 protein_coding 2/3
9 11 5248003 G A A stop_gained ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000380315 protein_coding 4/4
10 11 5248199 T A A stop_gained ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000335295 protein_coding 1/3
11 11 5248199 T A A stop_gained ÉLEVÉ HBB ENSG00000244734 Transcription ENST00000380315 protein_coding 3/4

L'ensemble de données gnomAD ne comporte aucun appel ni échantillon, mais vous pouvez étendre cette requête afin de sélectionner des échantillons présentant des variantes à fort impact.

Effectuer une jointure des données avec des bases de données BigQuery publiques

L'exemple suivant montre comment effectuer une jointure des données avec d'autres bases de données disponibles dans BigQuery, telles que des bases de données de voies de signalisation. Voici quelques-unes des bases de données disponibles via l'ISB-CGC :

L'exemple suivant montre comment extraire des variantes à fort impact dans le cadre de la voie de signalisation associée à la réparation des cassures double brin, tel que définie par le processus biologique GO (Gene Ontology) 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

L'exécution de la requête renvoie les éléments suivants :

Ligne reference_name start_position reference_bases alt Effet IMPACT SYMBOLE
1 5 917449 G T splice_acceptor_variant&non_coding_transcript_variant ÉLEVÉ TRIP13
2 12 1022568 A C stop_gained ÉLEVÉ RAD52
3 12 1022568 A C stop_gained ÉLEVÉ RAD52
4 12 1022568 A C stop_gained ÉLEVÉ RAD52
5 12 1023125 G A stop_gained ÉLEVÉ RAD52
6 12 1023125 G A stop_gained ÉLEVÉ RAD52
7 12 1023125 G A stop_gained ÉLEVÉ RAD52
8 12 1023167 G A stop_gained ÉLEVÉ RAD52
9 12 1023167 G A stop_gained ÉLEVÉ RAD52
10 12 1023167 G A stop_gained ÉLEVÉ RAD52
11 12 1023217 G T stop_gained ÉLEVÉ RAD52
12 12 1023217 G T stop_gained ÉLEVÉ RAD52
13 12 1023217 G T stop_gained ÉLEVÉ RAD52
14 12 1025654 TGA T frameshift_variant ÉLEVÉ RAD52
15 12 1025699 T G splice_acceptor_variant&non_coding_transcript_variant ÉLEVÉ RAD52
16 12 1036004 GC G frameshift_variant ÉLEVÉ RAD52
17 12 1038977 C CT frameshift_variant ÉLEVÉ RAD52