このページでは、次の方法について説明します。
- 入力 VCF ファイルのアノテーション フィールドを解析します。
- オープンソースのアノテーション ツールと Variant Effect Predictor(VEP)などのデータベースを使用して、インポート中の VCF ファイルにアノテーションを付けます。
アノテーション フィールドの解析
Variant Transforms ツールは、VCF 形式のバリアント アノテーションで指定されたアノテーション フィールド形式をサポートしています。VEP はこれと同じ形式を使用します。次の例は、VEP によってアノテーションが付けられたファイルの VCF ヘッダーを示しています。
##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|...>
REF=AAC
と ALT=TAC,A
のサンプル バリアント レコードの対応するアノテーションは、次のように示すことができます。
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
カンマで区切られた複数のアノテーション リストが、バリアント レコードごとに指定されます。各セットの最初のフィールドは、代替対立遺伝子(ALT
)に対応しています。上記の例では、ALT
ごとにサブセットのみを使用しています。実際のファイルには、より多くのフィールドがあり、それぞれの代替対立遺伝子は異なるアノテーションで複数回繰り返すことができます。
Variant Transforms ツールを実行する際に、--annotation_fields
を渡すことができます。このツールは、アノテーション リストを分割し、対応する alternate_bases
レコードの配下に移動します。たとえば、以下の表ではレコード フィールド CSQ
が各 alternate_bases
下に繰り返し見られることを確認できます。
名前 | 型 | モード | 説明 |
---|---|---|---|
alternate_bases.CSQ | RECORD | REPEATED | この代替の CSQ アノテーションのリスト。 |
alternate_bases.CSQ.allele | STRING | NULLABLE | アノテーション フィールドの ALT 部分。 |
alternate_bases.CSQ.Consequence | STRING | Nullable | なし |
alternate_bases.CSQ.IMPACT | STRING | Nullable | なし |
alternate_bases.CSQ.SYMBOL | STRING | Nullable | なし |
alternate_bases.CSQ.Gene | STRING | Nullable | なし |
alternate_bases.CSQ.Feature_type | STRING | Nullable | なし |
alternate_bases.CSQ.Feature | STRING | Nullable | なし |
alternate_bases.CSQ.BIOTYPE | STRING | Nullable | なし |
alternate_bases.CSQ.EXON | STRING | Nullable | なし |
alternate_bases.CSQ.INTRON | STRING | Nullable | なし |
代替対立遺伝子の一致
上記の例では、CSQ
列に表示されている代替対立遺伝子はバリアント レコードの ALT
文字列とは異なります。たとえば、AAC->C
バリアントは「-
」で示され、削除を表します。これは、VEP が --minimal
フラグを持つデータセットで実行されるためです。Variant Transforms は、この minimal
モードでの ALT
フィールドの照合をサポートしています。
ALT
照合の次のモードがサポートされています。
デフォルトでは、
ALT
一致フラグが設定されていない場合は、VCF 形式のバリアント アノテーションが使用されます。上記の例で示すように、
--use_allele_num
フラグを使用すると、Variant Transforms ツールはALLELE_NUM
アノテーション フィールドの値を検索します。これは VEP の機能です。--minimal_vep_alt_matching
フラグを使用して、--minimal
フラグを Variant Transforms ツールでシミュレートできます。ただし、一致があいまいな場合もあります。あいまいな一致は、Dataflow Console でカウントされ、報告されます。作成されたテーブルの新しいフィールドの配下にも追加されます。--minimal
を使用すると(--minimal_vep_alt_matching
ではなく)最適な結果が得られます。
VEP を使用して VCF ファイルに自動的にアノテーションを付ける
通常、アノテーション フィールドを解析して個別のフィールドに分割するには、まず VCF ファイルにアノテーションを追加する必要があります。ただし、Variant Transforms ツールを使用する際、VCF ファイルに事前にアノテーションを追加する必要はありません。
Variant Transforms は、VCF ファイルが BigQuery に読み込まれる際に、自動的にアノテーションを付けることができます。そのためには、Variant Transforms ツールを使用して VEP を実行する必要があります。その場合は、事前に構成された Docker イメージを含む Compute Engine 仮想マシン インスタンスを作成し実行する必要があります。
Variant Transforms ツールを使用した VEP の実行
Variant Transforms ツールを VEP で実行し、BigQuery に読み込まれる際に VCF ファイルに自動的にアノテーションを付けるには、少なくとも次のフラグを使用して Variant Transforms ツールを実行する必要があります。
--run_annotation_pipeline
--max_num_workers MAX_NUM_WORKERS
または--num_workers NUM_WORKERS
を指定します。--annotation_output_dir gs://CLOUD_STORAGE_BUCKET
CLOUD_STORAGE_BUCKET 変数に対して選択したバケットは、プロジェクトで所有している必要があります。
VEP の実行時に、Variant Transforms ツールに渡すことができるフラグの全一覧は、次のとおりです。
フラグ | デフォルト | 説明 |
---|---|---|
--run_annotation_pipeline |
偽 | 自動アノテーションを有効にします。 |
--max_num_workers または --num_workers |
自動的に決定 | それぞれ、Cloud Dataflow が自動スケーリングするジョブのワーカーの最大数、または最初にジョブに割り当てられたワーカーの数を示します。 |
--annotation_output_dir |
なし | VEP 出力ファイル用の Cloud Storage パス。入力ファイルの階層はこの場所に複製されますが、新しいファイルに _vep_output.vcf が追加されます。ディレクトリがすでに存在する場合、Variant Transforms ツールは失敗します。 |
--vep_image_uri |
gcr.io/gcp-variant-annotation/vep_91 |
Cloud Genomics VEP Docker イメージ。Google では、デフォルト画像を保持しています(VEP バージョン 91)。 |
--vep_cache_path |
gs://gcp-variant-annotation-vep-cache/vep_cache_homo_sapiens_GRCh38_91.tar.gz (GRCh38 参照配列に類似したヒトゲノムに推奨) |
VEP キャッシュの圧縮バージョンの Cloud Storage パス。このファイルは、--build_vep_cache.sh スクリプトを使用して作成できます。 |
--vep_info_field |
CSQ_VT |
新規のアノテーションを保持する新しい INFO フィールドの名前。 |
--vep_num_fork |
2 | 1 つのファイルに対して VEP を実行する際に使用される、ローカルのプロセス数。詳細については、VEP ドキュメントの --fork [num_forks] フラグをご覧ください。 |
--shard_variants |
有効 | VEP を使用して VCF ファイルにアノテーションを付ける前に、入力ファイルをより小さい一時 VCF ファイルに分割するかどうかを指定します。入力ファイルが小さい場合(たとえば各 VCF ファイルのバリアントが 50,000 個未満の場合など)、このフラグを True に設定すると、コンピューティング リソースが浪費される可能性があります。 |
--number_of_variants_per_shard |
2,000 | shard_variants が True の場合に、各シャードに書き込まれるバリアントの最大数。サンプル数が多いデータセットの場合は、小さい値を使用することもできます。小さい値を指定すると、パイプラインの終了に時間がかかることがあります。ほとんどの場合に、デフォルト値が推奨されます。 |
トラブルシューティング
VEP を Variant Transforms ツールで実行し、ツールが失敗した場合は、プロセスが実行された仮想マシンのログを示す logs
ディレクトリが作成されます。logs
ディレクトリは、--annotation_output_dir
フラグに使用するパスで指定します。これらのログファイルと、標準の Variant Transforms ツールのログを調べて、問題の原因を特定できます。
Variant Transforms ツールを使用して VEP を実行すると、VEP フラグ --check_ref
が自動的に有効になります。そのため、logs
ディレクトリのファイルには、参照と一致しないバリアントのリストが含まれます。正しい VEP キャッシュを使用していれば、不一致は自動的に破棄され、BigQuery の出力表には表示されないため、不一致は表示されません。
サンプルクエリ
影響の大きいバリアントの抽出
次のサンプルは、ヘモグロビンの β(ベータ)鎖をコードする HBB 遺伝子の影響の大きいバリアントを抽出する方法を示しています。
#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
このクエリを実行すると、次の結果が返されます。
Row | reference_name | start_position | reference_bases | alt | allele | 結果 | 影響 | 記号 | 遺伝子 | Feature_type | 機能 | BIOTYPE | Exon | INTRON |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
1 | 11 | 5246947 | G | GC | C | frameshift_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000335295 | protein_coding | 3/3 | |
2 | 11 | 5246957 | T | C | C | splice_acceptor_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000335295 | proten_coding | 2/2 | |
3 | 11 | 5246957 | T | C | C | splice_acceptor_variant&non_coding_transcript_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000475226 | retained_intron | 1/1 | |
4 | 11 | 5247805 | C | T | T | splice_donor_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000335295 | protein_coding | 2/2 | |
5 | 11 | 5247805 | C | T | T | splice_donor_variant&non_coding_transcript_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000475226 | retained_intron | 1/1 | |
6 | 11 | 5247991 | CAAAG | C | - | frameshift_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000335295 | protein_coding | 2/3 | |
7 | 11 | 5247991 | CAAAG | C | - | frameshift_variant | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000380315 | protein_coding | 4/4 | |
8 | 11 | 5248003 | G | A | A | stop_gained | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000335295 | protein_coding | 2/3 | |
9 | 11 | 5248003 | G | A | A | stop_gained | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000380315 | protein_coding | 4/4 | |
10 | 11 | 5248199 | T | A | A | stop_gained | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000335295 | protein_coding | 1/3 | |
11 | 11 | 5248199 | T | A | A | stop_gained | 高 | HBB | ENSG00000244734 | 文字起こし | ENST00000380315 | protein_coding | 3/4 |
gnomAD データセットには呼び出しまたはサンプルはありませんが、このクエリを拡張して影響の大きいバリアントがあるサンプルを選択できます。
公開 BigQuery データベースとのデータ結合
次のサンプルは、BigQuery で使用可能な pathway データベースなどの他のデータベースとデータを結合する方法を示しています。ISB-CGC では、たとえば次のようなデータベースを使用できます。
次のサンプルは、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
このクエリを実行すると、次の結果が返されます。
Row | reference_name | start_position | reference_bases | alt | 結果 | 影響 | 記号 |
---|---|---|---|---|---|---|---|
1 | 5 | 917449 | G | T | splice_acceptor_variant&non_coding_transcript_variant | 高 | TRIP13 |
2 | 12 | 1022568 | A | C | stop_gained | 高 | RAD52 |
3 | 12 | 1022568 | A | C | stop_gained | 高 | RAD52 |
4 | 12 | 1022568 | A | C | stop_gained | 高 | RAD52 |
5 | 12 | 1023125 | G | A | stop_gained | 高 | RAD52 |
6 | 12 | 1023125 | G | A | stop_gained | 高 | RAD52 |
7 | 12 | 1023125 | G | A | stop_gained | 高 | RAD52 |
8 | 12 | 1023167 | G | A | stop_gained | 高 | RAD52 |
9 | 12 | 1023167 | G | A | stop_gained | 高 | RAD52 |
10 | 12 | 1023167 | G | A | stop_gained | 高 | RAD52 |
11 | 12 | 1023217 | G | T | stop_gained | 高 | RAD52 |
12 | 12 | 1023217 | G | T | stop_gained | 高 | RAD52 |
13 | 12 | 1023217 | G | T | stop_gained | 高 | RAD52 |
14 | 12 | 1025654 | TGA | T | frameshift_variant | 高 | RAD52 |
15 | 12 | 1025699 | T | G | splice_acceptor_variant&non_coding_transcript_variant | 高 | RAD52 |
16 | 12 | 1036004 | GC | G | frameshift_variant | 高 | RAD52 |
17 | 12 | 1038977 | C | CT | frameshift_variant | 高 | RAD52 |