すべてのプロダクト
Search
ドキュメントセンター

Elastic High Performance Computing:BWA、GATK、Samtools を使用した遺伝子シーケンシングの実行

最終更新日:Jan 28, 2026

このトピックでは、E-HPC クラスターを使用して BWA、GATK、Samtools を実行し、遺伝子シーケンシングを行う方法について説明します。

背景情報

生命科学における遺伝子シーケンシング技術の急速な発展により、遺伝子シーケンスデータの量は指数関数的に増加しています。この大量のデータを処理するには、ホモロジー検索、アライメント、変異検出などのタスクのために、大規模なデータ処理と並列計算が必要です。高性能コンピューティング (HPC) は、これらのタスクに必要な強力な計算能力を提供します。HPC は、さまざまなスケジューラを使用して並列性を向上させ、GPU を使用して計算を高速化します。

このトピックでは、典型的な全ゲノムシーケンシング (WGS) ワークフローを例として使用します。GATK ソフトウェアを使用したヒト WGS の一般的なプロセスを示します。実際のバイオインフォマティクス解析では、特定のニーズに基づいてプロセスを調整する必要があります。たとえば、変異検出のための品質管理とフィルタリングを追加する必要がある場合があります。

遺伝子シーケンシングでは、BWA を使用してインデックスを構築し、レコードをアライメントできます。次に、Samtools を使用してアライメントレコードをソートします。最後に、GATK を使用して重複シーケンスを除去し、塩基クオリティスコアを再キャリブレーションし、変異を検出します。

  • BWA (Burrows-Wheeler Alignment Tool) は、DNA シーケンスをヒトゲノムなどのリファレンスゲノムにマッピングするソフトウェアです。

  • GATK (The Genome Analysis Toolkit) は、第二世代の再シーケンシングデータを解析するためのソフトウェアツールキットです。主に、重複シーケンスの除去、塩基クオリティスコアの再キャリブレーション、変異検出に使用されます。

  • Samtools は、SAM および BAM フォーマットを処理するためのツールです。バイナリファイルの表示、ファイル形式の変換、ファイルのソートとマージが可能です。また、SAM ファイルのフラグとタグの情報を使用して、アライメント結果を要約します。

事前準備

  1. E-HPC クラスターを作成します。詳細については、「Standard Edition クラスターの作成」をご参照ください。

    次の表に、クラスター構成の例を示します。

    設定項目

    構成

    シリーズ

    Standard Edition

    デプロイモード

    パブリッククラウドクラスター

    クラスタータイプ

    SLURM

    ノード構成

    次の仕様のマスターノード 1 台、ログインノード 1 台、コンピュートノード 1 台が含まれます。

    • マスターノード:ecs.c7.large インスタンスタイプを使用し、2 vCPU と 4 GiB のメモリを備えています。

    • ログインノード:ecs.c7.large インスタンスタイプを使用し、2 vCPU と 4 GiB のメモリを備えています。

    • コンピュートノード:ecs.g7.16xlarge インスタンスタイプを使用し、64 vCPU と 256 GiB のメモリを備えています。

    クラスターイメージ

    CentOS 7.6 パブリックイメージ

  2. クラスターユーザーを作成します。詳細については、「ユーザー管理」をご参照ください。

    クラスターユーザーは、クラスターにログインしてソフトウェアのコンパイル、ジョブのサブミット、その他の操作を実行するために使用されます。次の例は、このトピックで作成されるユーザーを示しています。

    • ユーザー名:testuser

    • ユーザーグループ:sudo 権限グループ

  3. ソフトウェアをインストールします。

    次の表に、インストールするソフトウェアの詳細を示します。

    ソフトウェア名

    バージョン

    ダウンロード URL

    BWA

    4.2.0.0

    BWA

    GATK

    0.7.17

    GATK

    Samtools

    1.14

    Samtools

手順1:E-HPC クラスターへの接続

  1. E-HPC クラスターに接続します。詳細については、「クラスターへの接続」をご参照ください。

  2. 次のコマンドを実行して、ジョブ実行ディレクトリを作成します。このトピックでは、/home/data/human を例として使用します。

    説明

    データをダウンロードして処理するための十分なディスク領域を確保するために、クラスターにマウントされている /home 共有ストレージディレクトリにジョブ実行ディレクトリを作成してください。

    sudo mkdir -p /home/data/human

手順2:データのダウンロードと前処理

  1. 次のコマンドを実行して、data_pre.slurm という名前のジョブスクリプトファイルを作成して開きます。

    sudo vim data_pre.slurm

    スクリプトは次のとおりです。

    #!/bin/bash -l
    #SBATCH -J data_pre
    #SBATCH --partition comp
    #SBATCH -N 1
    #SBATCH -n 1
    #SBATCH --output=data_pre.log
    
    cd /home/data/human/
    
    # ヒトのリファレンスゲノムと 2 つのヒト遺伝子シーケンシングサンプルをダウンロードして解凍します。
    wget https://public-ehpc-package.oss-cn-hangzhou.aliyuncs.com/lifescience/b37_human_g1k_v37.fasta
    wget https://public-ehpc-package.oss-cn-hangzhou.aliyuncs.com/lifescience/b37_1000G_phase1.indels.b37.vcf
    wget https://public-ehpc-package.oss-cn-hangzhou.aliyuncs.com/lifescience/b37_Mills_and_1000G_gold_standard.indels.b37.vcf
    wget https://public-ehpc-package.oss-cn-hangzhou.aliyuncs.com/lifescience/b37_dbsnp_138.b37.vcf.zip
    unzip b37_dbsnp_138.b37.vcf.zip
    
    wget https://public-ehpc-package.oss-cn-hangzhou.aliyuncs.com/lifescience/gatk-examples_example1_NA20318_ERR250971_1.filt.fastq.gz
    gunzip gatk-examples_example1_NA20318_ERR250971_1.filt.fastq.gz
    
    wget https://public-ehpc-package.oss-cn-hangzhou.aliyuncs.com/lifescience/gatk-examples_example1_NA20318_ERR250971_2.filt.fastq.gz
    gunzip gatk-examples_example1_NA20318_ERR250971_2.filt.fastq.gz
    
    # bgzip を使用してファイルを圧縮します。
    bgzip -f b37_1000G_phase1.indels.b37.vcf
    bgzip -f b37_Mills_and_1000G_gold_standard.indels.b37.vcf
    bgzip -f b37_dbsnp_138.b37.vcf
    
    # tabix を使用してインデックスを構築します。
    tabix -p vcf  b37_1000G_phase1.indels.b37.vcf.gz
    tabix -p vcf  b37_Mills_and_1000G_gold_standard.indels.b37.vcf.gz
    tabix -p vcf  b37_dbsnp_138.b37.vcf.gz
  2. 次のコマンドを実行してジョブをサブミットします。

    sbatch data_pre.slurm

    次の出力は、ジョブが正常にサブミットされたことを示しています。

    image

  3. 次のコマンドを実行してジョブステータスを確認します。

    squeue
  4. ジョブが完了したら、次のコマンドを実行して data_pre.log ファイルで結果を表示します。

    cat data_pre.log

手順3:リファレンスゲノムインデックスの構築

  1. 次のコマンドを実行して、bwa_index.slurm という名前のジョブスクリプトファイルを作成して開きます。

    sudo vim bwa_index.slurm

    スクリプトは次のとおりです。

    #!/bin/bash -l
    #SBATCH -J bwa_index
    #SBATCH --partition comp
    #SBATCH -N 1
    #SBATCH -n 1
    #SBATCH --output=bwa_index.log
    
    cd /home/data/human/
    
    # BWA アライメントのためにリファレンスゲノムのインデックスデータベースを構築します。
    time bwa index b37_human_g1k_v37.fasta
    
    # FASTA シーケンスフォーマットのインデックスを作成します。
    samtools faidx b37_human_g1k_v37.fasta
  2. 次のコマンドを実行してジョブをサブミットします。

    sbatch --dependency=afterok:jobid1 bwa_index.slurm
  3. 次のコマンドを実行してジョブステータスを確認します。

    squeue
  4. ジョブが完了したら、次のコマンドを実行して bwa_index.log ファイルで結果を表示します。

    cat bwa_index.log

手順4:遺伝子シーケンシング解析の実行

  1. 次のコマンドを実行して、wgs_main_process.slurm という名前のジョブスクリプトファイルを作成して開きます。

    sudo vim wgs_main_process.slurm

    スクリプトは次のとおりです。

    重要

    6 行目の #SBATCH -w compute000 で、compute000 をご利用のコンピュートノードの名前に置き換えてください。

    #!/bin/bash -l
    #SBATCH -J wgs_main_process
    #SBATCH --partition comp
    #SBATCH -N 1
    #SBATCH -n 64
    #SBATCH -w compute000
    #SBATCH --output=wgs_main_process.log
    
    cd /home/data/human/
    
    # サンプルシーケンシングデータのリードをヒトのリファレンスゲノムにアライメントし、出力ファイルを BAM フォーマットに変換してディスク領域を節約します。
    time bwa mem -t 64 \
    -R '@RG\tID:ehpc\tPL:illumina\tLB:library\tSM:b37' \
    b37_human_g1k_v37.fasta \
    gatk-examples_example1_NA20318_ERR250971_1.filt.fastq \
    gatk-examples_example1_NA20318_ERR250971_2.filt.fastq | samtools view -S -b - > ERR250971.bam
    
    # アライメントレコードを昇順でソートします。
    time samtools sort -@ 64 -O bam -o ERR250971.sorted.bam ERR250971.bam
    
    # DNA シーケンスの PCR 増幅によって生成された重複リードを除去します。
    time gatk MarkDuplicates \
    -I ERR250971.sorted.bam -O ERR250971.markdup.bam \
    -M ERR250971.sorted.markdup_metrics.txt
    
    # markdup シーケンシング BAM ファイルのインデックスを構築します。
    samtools index ERR250971.markdup.bam
    
    # リファレンスゲノムの dict ファイルを生成します。
    time gatk CreateSequenceDictionary  \
    -R b37_human_g1k_v37.fasta \
    -O b37_human_g1k_v37.dict
    
    # 既存の SNP および indel データベースを使用してモデルを構築し、再キャリブレーションテーブルを作成します。
    time gatk BaseRecalibrator \
    -R b37_human_g1k_v37.fasta \
    -I ERR250971.markdup.bam  \
    --known-sites b37_1000G_phase1.indels.b37.vcf.gz \
    --known-sites b37_Mills_and_1000G_gold_standard.indels.b37.vcf.gz \
    --known-sites b37_dbsnp_138.b37.vcf.gz -O ERR250971.BQSR.table
    
    # 元の塩基クオリティスコアを調整します。
    time gatk ApplyBQSR \
    --bqsr-recal-file ERR250971.BQSR.table \
    -R b37_human_g1k_v37.fasta \
    -I ERR250971.markdup.bam \
    -O ERR250971.BQSR.bam
    
    # BQSR シーケンシング BAM ファイルのインデックスを構築します。
    time samtools index ERR250971.BQSR.bam
    
    # 変異検出結果を VCF ファイルに出力します。
    time gatk HaplotypeCaller \
    -R b37_human_g1k_v37.fasta  \
    -I ERR250971.BQSR.bam \
    -O ERR250971.HC.vcf.gz
  2. 次のコマンドを実行してジョブをサブミットします。

    sbatch --dependency=afterok:jobid2 wgs_main_process.slurm
  3. 次のコマンドを実行してジョブステータスを確認します。

    squeue
  4. ジョブが完了したら、次のコマンドを実行して wgs_main_process.log ファイルで結果を表示します。

    cat wgs_main_process.log

シーケンシング時間

このトピックで使用されるコンピュートノードは ecs.g7.16xlarge です。これはデモンストレーションのみを目的としています。全ゲノムシーケンシングプロセス全体には約 6~7 時間かかります。次の表に、参考として各手順のおおよその時間を示します。これらの時間は、最適なパフォーマンスのベンチマークではありません。

プロセス

機能

実行時間 (分)

bwa index

インデックスの構築

54

bwa mem

アライメント

25

samtools sort

ソート

2

gatk MarkDuplicates

重複の除去

13

gatk CreateSequenceDictionary

辞書の作成

0.15

gatk BaseRecalibrator

キャリブレーションテーブルの構築

40

gatk ApplyBQSR

再キャリブレーション

21

gatk HaplotypeCaller

変異の検出

211