0

我有一个 hg19 对齐的 BAM,我希望为其生成一个 DeepVariant VCF。我使用 samtools 提取标题并确保 hg19 参考 FASTA 索引包含相同的重叠群和位置。我最初的目标是在这个 WGS BAM 上只运行一个外显子组模型,使用以下模型和区域:

MODEL=gs://deepvariant/models/DeepVariant/0.7.2/DeepVariant-inception_v3-0.7.2+data-wes_standard

--regions gs://deepvariant/exome-case-study-testdata/refseq.coding_exons.b37.extended50.bed

不幸的是,脚本抗议说 BED 和 BAM / FASTA 参考之间的匹配为 0。我决定运行相同的外显子组模型,但没有指定区域。这是我的脚本:

#!/bin/bash
set -euo pipefail
# Set common settings.
PROJECT_ID=<MY PROJECT>
OUTPUT_BUCKET=gs://<MY BUCKET>
STAGING_FOLDER_NAME=staging
OUTPUT_FILE_NAME=output.vcf
# Model for calling whole genome sequencing data.
MODEL=gs://deepvariant/models/DeepVariant/0.7.2/DeepVariant-inception_v3-0.7.2+data-wes_standard
IMAGE_VERSION=0.7.2
DOCKER_IMAGE=gcr.io/deepvariant-docker/deepvariant:"${IMAGE_VERSION}"
COMMAND="/opt/deepvariant_runner/bin/gcp_deepvariant_runner \
  --project ${PROJECT_ID} \
  --zones us-west1-* \
  --docker_image ${DOCKER_IMAGE} \
  --outfile ${OUTPUT_BUCKET}/${OUTPUT_FILE_NAME} \
  --staging ${OUTPUT_BUCKET}/${STAGING_FOLDER_NAME} \
  --model ${MODEL} \
  --bam gs://my-bucket/wgs_data.bam \
  --ref gs://my-bucket/human_g1k_v37.fa \
  --shards 512 \
  --make_examples_workers 32 \
  --make_examples_cores_per_worker 16 \
  --make_examples_ram_per_worker_gb 60 \
  --make_examples_disk_per_worker_gb 200 \
  --call_variants_workers 32 \
  --call_variants_cores_per_worker 32 \
  --call_variants_ram_per_worker_gb 60 \
  --call_variants_disk_per_worker_gb 50 \
  --gcsfuse"
# Run the pipeline.
gcloud alpha genomics pipelines run \
    --project "${PROJECT_ID}" \
    --service-account-scopes="https://www.googleapis.com/auth/cloud-platform" \
    --logging "${OUTPUT_BUCKET}/${STAGING_FOLDER_NAME}/runner_logs_$(date +%Y%m%d_%H%M%S).log" \
    --regions us-west1 \
    --docker-image gcr.io/cloud-genomics-pipelines/gcp-deepvariant-runner \
    --command-line "${COMMAND}"

BAM 有对应的 BAI,FA 有 FAI 文件。DeepVariant QuickStart 表明这些设置将在 1-2 小时内生成 VCF,但我的管道现在已经运行了 7 多个小时。暂存文件夹现在有一个 call_variants,其中包含 32 个 GZ 文件中的 31 个。Genomics 管道视图显示了 11 个运行 call_variant 的管道,因此我怀疑它正在处理最后一个文件,准备将所有文件合并到一个 VCF 中。

我只是不明白为什么要花这么长时间。我故意排除了抢占式实例,并且文档说外显子组管道应该只需要 20 分钟(WGS 需要 1-2 小时)。为什么会这么慢?

4

1 回答 1

1

您看到的运行时间肯定比 DeepVariant 的预期慢。

开始时的一项观察 - 外显子组捕获床 (refseq.coding_exons.b37.extended50.bed) 和参考 (human_g1k_v37.fa) 的坐标应该匹配。您知道您的 BAM 映射到哪个参考基因组吗?只是为了确认,在您的 FASTA 文件中,第一行应为:>1 且没有“chr”。

使用区域文件时,预期时间应小于 1 小时。

其次,我能否请您按照此页面的说明尝试在单台机器上运行外显子组案例研究:

https://github.com/google/deepvariant/blob/r0.8/docs/deepvariant-exome-case-study.md

运行此程序将有助于确定您看到的问题是否与 DeepVariant 本身有关,或者它是否与与程序分开的 GCP 云实现有关。

谢谢你,安德鲁

于 2019-06-28T15:55:48.823 回答