第4章: ゲノム解析技術

4. ゲノム解析技術

この章を読む理由

ゲノム解析は、FASTQ から BAM / VCF へ至る成果物の流れを最初に固定する章です。 第5章以降の発現解析や臨床応用を理解するためにも、 「どこで品質が落ちるか」「どの成果物を次工程へ渡すか」をここで整理しておく必要があります。

本章で使う共通題材

  • 題材A: SARS-CoV-2 公開データ
    • ラン: SRR11140744
    • 参照配列: MN908947.3
    • 本章では、最小の QC → マッピング → 変異検出の流れを確認するために使います。

学習目標

  • NGS 解析パイプライン(QC→アライメント→後処理→変異検出→アノテーション→要約)と、代表的な入出力(FASTQ/BAM/VCF)を説明できる
  • QC / マッピング / 変異検出で確認すべき品質指標と、その確認手段(FastQC/MultiQC、samtools flagstatbcftools stats 等)を説明できる
  • 参照ゲノム、パラメータ、ツール / DB のバージョン、ログ記録をどう残せば第三者が再実行できるかを整理できる

4.1 配列解析パイプライン

NGS配列解析パイプライン(QC→アライメント→後処理→変異検出→アノテーション→結果要約)
図 4-1: NGS 配列解析パイプラインの概略。入力(FASTQ)から最終成果物(VCF/レポート)までの流れを俯瞰し、各工程の品質指標と代表ツールを把握する。

コマンドライン例(最小ワークフロー): 🧪 概念例(擬似コード)

# 1) 参照インデックス作成(例: BWA)
bwa index ref.fa

# 2) マッピング → SAM
bwa mem -t 8 ref.fa sample_R1.fq.gz sample_R2.fq.gz > sample.sam

# 3) SAM→BAM & ソート
samtools view -bS sample.sam | samtools sort -o sample.sorted.bam
samtools index sample.sorted.bam

# 4) 変異検出(例: bcftools)
bcftools mpileup -Ou -f ref.fa sample.sorted.bam | \
  bcftools call -mv -Ob -o sample.bcf
bcftools view sample.bcf > sample.vcf

# 5) 簡易QC/統計
samtools flagstat sample.sorted.bam > sample.flagstat.txt

QCの基本(FastQC/MultiQC): 🧪 概念例(擬似コード)

# 入力FASTQの品質評価(並列実行例)
mkdir -p qc/fastqc
fastqc -t 8 -o qc/fastqc data/*.fq.gz

# レポート統合
multiqc -o qc qc/fastqc

前処理フェーズ: 🧪 概念例(擬似コード)

def preprocess_sequence(raw_sequence):
    # 品質管理
    quality_filtered = quality_control(raw_sequence)
    
    # アダプター除去
    trimmed = adapter_trimming(quality_filtered)
    
    # 重複除去
    deduplicated = remove_duplicates(trimmed)
    
    return deduplicated

実装例(要環境調整): このコードは構造の例です。依存ライブラリ(BioPython / pysam 等)や入出力、例外設計は環境に合わせて調整してください(CIの実行対象ではありません)。 🧪 概念例(擬似コード)

import logging
import os
from typing import List, Dict, Optional, Tuple
from dataclasses import dataclass
from Bio import SeqIO
import pysam
import numpy as np

# ロギング設定
logging.basicConfig(
    level=logging.INFO,
    format='%(asctime)s - %(name)s - %(levelname)s - %(message)s'
)
logger = logging.getLogger(__name__)

@dataclass
class QualityMetrics:
    """品質管理メトリクス"""
    total_reads: int
    passed_reads: int
    mean_quality: float
    q20_percentage: float
    q30_percentage: float
    gc_content: float

class SequenceAnalysisPipeline:
    """エラーハンドリングを含む配列解析パイプライン"""
    
    def __init__(self, min_quality: int = 20, min_length: int = 50):
        self.min_quality = min_quality
        self.min_length = min_length
        
    def quality_control(self, fastq_file: str) -> Tuple[str, QualityMetrics]:
        """
        品質管理とフィルタリング
        
        Args:
            fastq_file: 入力FASTQファイル
            
        Returns:
            tuple: (フィルタ済みファイル名, 品質メトリクス)
        """
        try:
            output_file = fastq_file.replace('.fastq', '_filtered.fastq')
            
            total_reads = 0
            passed_reads = 0
            quality_scores = []
            gc_counts = []
            
            with open(output_file, 'w') as out_handle:
                for record in SeqIO.parse(fastq_file, 'fastq'):
                    total_reads += 1
                    
                    # 品質スコアの計算
                    quals = record.letter_annotations['phred_quality']
                    mean_qual = np.mean(quals)
                    
                    # 長さチェック
                    if len(record.seq) < self.min_length:
                        logger.debug(f"Read {record.id} too short: {len(record.seq)}")
                        continue
                    
                    # 品質チェック
                    if mean_qual < self.min_quality:
                        logger.debug(f"Read {record.id} low quality: {mean_qual:.2f}")
                        continue
                    
                    # GC含量計算
                    gc_count = record.seq.count('G') + record.seq.count('C')
                    gc_content = gc_count / len(record.seq)
                    
                    # フィルタを通過
                    passed_reads += 1
                    quality_scores.extend(quals)
                    gc_counts.append(gc_content)
                    
                    # 出力
                    SeqIO.write(record, out_handle, 'fastq')
            
            # メトリクス計算
            if passed_reads == 0:
                raise ValueError("No reads passed quality filters")
            
            quality_array = np.array(quality_scores)
            metrics = QualityMetrics(
                total_reads=total_reads,
                passed_reads=passed_reads,
                mean_quality=np.mean(quality_array),
                q20_percentage=np.sum(quality_array >= 20) / len(quality_array) * 100,
                q30_percentage=np.sum(quality_array >= 30) / len(quality_array) * 100,
                gc_content=np.mean(gc_counts)
            )
            
            logger.info(f"Quality control completed: {passed_reads}/{total_reads} reads passed")
            return output_file, metrics
            
        except FileNotFoundError:
            logger.error(f"Input file not found: {fastq_file}")
            raise
        except Exception as e:
            logger.error(f"Quality control failed: {e}")
            raise

マッピングフェーズ:

  • BWA: Burrows-Wheeler Alignerによる高速マッピング
  • SAMtools: マッピング結果の操作・統計
  • GATK: 変異検出パイプライン

後処理フェーズ:

  • 変異アノテーション: 機能予測・病原性評価
  • 品質管理: 偽陽性率の制御
  • 統計解析: 有意性検定・多重検定補正

参照ゲノムの前提と最新動向(Telomere-to-Telomere (T2T)/パンゲノム):

  • 線形参照(例: GRCh38)が依然標準だが、ギャップのないTelomere-to-Telomere (T2T) 参照(例: T2T-CHM13)やパンゲノム参照が登場している。
  • パンゲノムは複数個体の多様性を表現するため、グラフ参照や多系統の表現が用いられる。
  • 実務では、目的(臨床/研究/集団解析)に合わせて参照を選び、解析結果の互換性と再現性を確保する。

解析パイプラインの3系統(概要):

  • short-read: QC → 線形参照へのマッピング → 変異コーリング → 注釈付与
  • long-read: QC → ロングリードマッピング/アセンブリ → 構造変異検出 → ハプロタイプ解析
  • pangenome/graph: グラフ参照の構築 → グラフアラインメント → 変異検出 → 注釈付与

臨床診断への直接貢献: 標準化されたパイプラインにより、ゲノムデータから臨床的に有用な情報をより安定して抽出できる。品質管理プロセスにより、偽陽性・偽陰性を最小化し、診断精度を向上させる。自動化により、解析時間を短縮し、緊急症例での迅速診断を可能にする。標準化により、施設間での結果の一貫性を確保し、医療の質を向上させる。

4.2 変異検出と機能解析

体細胞変異検出: 🧪 概念例(擬似コード)

class SomaticVariantCaller:
    """腫瘍-正常ペアからの体細胞変異検出"""
    
    def __init__(self, tumor_bam, normal_bam, reference):
        self.tumor_bam = tumor_bam
        self.normal_bam = normal_bam
        self.reference = reference
    
    def call_variants(self, min_vaf=0.05, min_depth=20):
        """
        体細胞変異の検出
        
        Args:
            min_vaf: 最小変異アレル頻度
            min_depth: 最小リード深度
        """
        # MuTect2スタイルの変異検出
        variants = []
        
        # 腫瘍特異的変異の統計的検定
        # Fisher's exact testによる有意性評価
        
        return variants

構造変異(SV)検出:

  • 欠失・挿入・逆位・転座の検出
  • Split-read/Paired-end情報の統合
  • CNV(コピー数変異)解析

深層学習ベースの変異コーリング:

  • 体細胞変異の文脈では、DeepSomaticやNeuSomaticなどの深層学習ベースcallerが提案されている。
  • 生殖細胞系列では、DeepVariant等が画像化したリード情報を用いて変異を推定する。
  • これらのツールは参照ゲノムのバージョン、学習データの前提、カバレッジ条件に依存するため、適用範囲の確認が必要。

変異の機能的影響予測:

  • SIFT/PolyPhen-2による病原性予測
  • VEP(Variant Effect Predictor)
  • スプライシング影響予測

4.3 ロングリードシークエンシング

PacBio/Oxford Nanopore技術:

  • リード長: 10kb-100kb+
  • エラー率: 5〜15%(生データ)
  • リアルタイムシークエンシング

ハイブリッド解析(short+long):

  • short-readの精度とlong-readの連続性を組み合わせ、構造変異や反復領域の解析精度を高める。
  • 解析設計としては、short-readを基準としつつ、long-readで難所を補完する形が一般的。

アセンブリアルゴリズム: 🧪 概念例(擬似コード)

class LongReadAssembler:
    """ロングリードアセンブリ"""
    
    def overlap_layout_consensus(self, reads):
        # 1. オーバーラップグラフ構築
        overlap_graph = self.build_overlap_graph(reads)
        
        # 2. レイアウト(グラフ簡約化)
        layout = self.simplify_graph(overlap_graph)
        
        # 3. コンセンサス配列生成
        consensus = self.generate_consensus(layout)
        
        return consensus

ハプロタイプ分離アセンブリ:

  • 父母由来の染色体を分離
  • HiFiリードによる高精度化
  • 構造変異の正確な検出

🎯 認定試験ポイント

重要概念チェックリスト

ゲノム解析基礎 ⭐⭐⭐

  • ゲノムサイズ・遺伝子数・非コード領域の割合をヒトで把握している
  • SNP・Indel・CNV・構造変異の定義と検出手法を理解している
  • 遺伝率とは何か、その推定方法を説明できる
  • GWAS(ゲノムワイド関連解析)の基本原理と限界を理解している

配列解析技術 ⭐⭐⭐

  • Sanger法と次世代シークエンシング技術の特徴比較ができる
  • リードマッピング・アセンブリの基本戦略を理解している
  • ペアエンドリード・ロングリードの利点と使い分けを把握している
  • カバレッジ・精度・N50等の品質指標の意味を説明できる

変異解析 ⭐⭐

  • SNP・Indelコーリングのアルゴリズムと品質フィルタリング
  • アレル頻度・ハーディー・ワインベルグ平衡の計算
  • 連鎖不平衡(LD)とハプロタイプの概念
  • 集団遺伝学パラメータ(FST、Tajima’s D等)の解釈

アセンブリ手法 ⭐⭐

  • グリーディアルゴリズムとグラフベースアセンブリの違い
  • De Bruijnグラフの構築とパス探索
  • アセンブリ品質評価(N50、コンティグサイズ分布等)
  • ギャップクローニングとスキャフォールディング

典型的な出題パターン

【遺伝学計算】 🧪 概念例(擬似コード)

問題例: 日本人集団でのSNPアレル頻度がA=0.7、a=0.3の場合、
ハーディー・ワインベルグ平衡における各遺伝型頻度は?

解答: AA = 0.7² = 0.49 (49%)
      Aa = 2×0.7×0.3 = 0.42 (42%)
      aa = 0.3² = 0.09 (9%)

【技術選択】 🧪 概念例(擬似コード)

問題例: 新規ゲノムのde novoアセンブリを行う際、
ショートリードとロングリードの使い分けを説明せよ。

解答: ショートリード(Illumina): 高精度、反復配列に課題
      ロングリード(PacBio/ONT): 反復配列解決、精度向上が課題
      → ハイブリッドアプローチで両者の利点を活用

【統計・品質管理】 🧪 概念例(擬似コード)

問題例: GWAS解析で多重検定補正が必要な理由と、
一般的な補正手法を2つ挙げよ。

解答: 理由: 数百万SNPを同時検定するため偽陽性率が急増
      手法: 1) Bonferroni補正(保守的)
           2) FDR補正(Benjamini-Hochberg法)

【データ解釈】 🧪 概念例(擬似コード)

問題例: RNA-seqデータで発現変動遺伝子(DEG)を検出した後、
生物学的意味を調べるための解析手法を3つ挙げよ。

解答: 1) Gene Ontology (GO) enrichment解析
      2) KEGGパスウェイ解析
      3) Gene Set Enrichment Analysis (GSEA)

関連する付録・章

まとめ

  • 第4章で固定すべきなのは、FASTQ / BAM / VCF の役割と、各工程で確認する QC 指標です。
  • 最初から大規模データを扱うのではなく、SRR11140744 のような小規模公開データで最小パイプラインを再現すると理解しやすくなります。
  • 後続章でも、ここで作った成果物とログの残し方が再現性の基盤になります。

Source notes / 次の一歩

次章への橋渡し

ゲノム解析で「配列から変異を取り出す」流れを固定したら、次は 「その変化が発現量や機能にどう現れるか」を追う必要があります。 第5章では、同じ再現性の考え方を保ちながら RNA-seq と発現解析へ進みます。

最小入出力(期待成果物/期待ログ)

  • 入力: 参照ゲノム(FASTA)とリード(FASTQ)(小規模データで可)
  • 出力(期待成果物): ソート済みBAM+インデックス、変異(VCF/BCF)+インデックス、QCレポート(FastQC/MultiQC)
  • 期待ログ(例): マッピング率・read数・変異数などのサマリ(例: samtools flagstat, bcftools stats

前へ: データ構造とアルゴリズム 目次 次へ: トランスクリプトーム解析

演習

  1. SRR11140744MN908947.3 を用い、QC→アライメント→BAM 整形→変異検出までの最小パイプラインを実行し、成果物(BAM / VCF / QC レポート)を作成せよ。
  2. トリミング閾値、最小品質、最小深度を 1 つずつ変え、マッピング率・変異数・主要 QC 指標の変化を比較せよ。
  3. 使用したコマンド、パラメータ、ツール / 参照 DB のバージョンを README 形式で残し、第三者が再実行できるよう整理せよ。

具体課題例

  • samtools flagstatbcftools stats を使い、最小限の品質サマリを 1 ページで説明する。
  • マッピングツールまたは変異検出ツールを 1 つ差し替え、結果差分と採用判断を文章でまとめる。
  • 第10章で再利用できるよう、出力ファイル名・サンプル ID・参照配列 ID の命名規則を決める。