第3章: データ構造とアルゴリズム

3. データ構造とアルゴリズム

学習目標

章の主要概念（データ構造/アルゴリズムの要点）を説明できる
計算量やメモリ効率を意識した実装上の注意を挙げられる
サンプルコードや実データに対して適切な手法を選択できる

3.0 アルゴリズム設計の理論的基盤

計算複雑性クラス:

P: 多項式時間で解ける決定問題（例: 編集距離がしきい値以下かを判定する問題）
NP: 解の候補を多項式時間で検証できる決定問題
#P: 解の個数を数える問題（RNA構造数え上げ）
PSPACE: 多項式空間（RNA相互作用予測）

近似アルゴリズム理論: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

近似比の定義: ρ = max(OPT/ALG, ALG/OPT)
- PTAS (Polynomial-Time Approximation Scheme): 任意のε > 0に対して(1+ε)近似
- FPTAS: さらに計算時間がεに関して多項式
- APX: 定数近似比を持つ問題クラス

確率的アルゴリズム:

ラスベガス型: 常に正解、実行時間が確率的
モンテカルロ型: 実行時間固定、解の品質が確率的
代表例: MinHashやランダム化探索など、大規模データで近似解を高速に得る手法

文字列アルゴリズムの数学的基礎:

有限オートマトン理論: パターンマッチングの状態機械
Z-algorithm: 線形時間での周期性検出
KMP法: 失敗関数による効率的検索

3.1 配列データ構造

Suffix Tree:

用途: 高速部分文字列検索
構築時間: O(n)
検索時間: O(m) (mは検索文字列長)
メモリ使用量: O(n)

Burrows-Wheeler Transform (BWT):

用途: 配列圧縮と、FM-index などの高速検索構造の基盤
特徴: 後段の圧縮器で高い圧縮効率を得やすい
検索時間: FM-indexと組み合わせると O(m + occ)（occは出現回数）

FM-index: 🔁 疑似コード（実行不可: BWT変換結果の概念図）

入力配列: ATGCATGC（実装側で終端記号 $ を追加）
BWT変換後: CC$GGTTAA
FM-index構築 → 後方検索に利用可能

FM-indexの簡略実装（学習用）: 小規模な文字列で後方検索の流れを確認するための簡略版であり、実運用ではより省メモリな実装が用いられる。 🧪 概念例（Python 3.11+ / 学習用）

class FMIndex:
    """FM-index実装（Burrows-Wheeler変換を使用した高速文字列検索）"""

    def __init__(self, text):
        # 入力 text は終端記号 '$' を含まない前提。実装側で1つだけ追加する。
        self.text = text + '$'
        self.sa = self._build_suffix_array()
        self.bwt = self._build_bwt()
        self.alphabet = sorted(set(self.bwt))
        self.occ = self._build_occ_table()
        self.count = self._build_count_table()

    def _build_suffix_array(self):
        """接尾辞配列の構築"""
        return sorted(range(len(self.text)), key=lambda i: self.text[i:])

    def _build_bwt(self):
        """Burrows-Wheeler変換"""
        return ''.join(self.text[i - 1] if i > 0 else self.text[-1] for i in self.sa)

    def _build_occ_table(self):
        """出現回数テーブルの構築"""
        occ = {char: [0] for char in self.alphabet}
        for char in self.bwt:
            for symbol in self.alphabet:
                occ[symbol].append(occ[symbol][-1] + (symbol == char))
        return occ

    def _build_count_table(self):
        """先頭列で各文字が始まる位置を記録"""
        count = {}
        total = 0
        for char in self.alphabet:
            count[char] = total
            total += self.bwt.count(char)
        return count

    def _get_lf_mapping(self, char, index):
        """LF-mappingで前の文字位置へ写像"""
        return self.count[char] + self.occ[char][index]

    def search(self, pattern):
        """検索パターンの出現位置を返す"""
        if not pattern:
            return []

        top = 0
        bottom = len(self.bwt) - 1

        for char in reversed(pattern):
            if char not in self.occ:
                return []

            top = self._get_lf_mapping(char, top)
            bottom = self._get_lf_mapping(char, bottom + 1) - 1

            if top > bottom:
                return []

        return sorted(self.sa[i] for i in range(top, bottom + 1))


fm_index = FMIndex("ATGCATGC")
print(fm_index.search("TGC"))

🧪 概念例（期待される出力）

[1, 5]

実装上の価値: ナイーブな文字列検索はO(nm)の時間計算量を要するが、FM-index などの索引構造により、参照配列に対する反復検索のコストを大幅に削減できる。実際の検索時間やメモリ削減量は、実装、索引構築方法、参照配列、クエリ数、I/O 条件に依存するため、対象データでベンチマークして確認する。

線形参照とパンゲノムグラフの読み替え: BWT / FM-index、k-mer、minimizer は、線形参照上の高速検索や seed-and-extend の考え方を支える。GRCh38 と T2T-CHM13 は座標系やギャップの扱いが異なる線形アセンブリとして選択・記録し、pangenome / パンゲノムは複数の haplotype / assembly を経路として含み得る graph / multipath reference として区別する。HPRC は、単一の線形参照に依存すると reference bias が mapping、variant discovery、association studies に影響し得ると説明している。¹²³

アルゴリズム要素	線形参照での見方	pangenome / graph reference での注意
BWT / FM-index	1本の参照配列や連結配列に対する高速検索	graph へ直接適用するには、経路表現やgraph-aware indexの設計を確認する
k-mer / minimizer	read の候補位置を素早く絞るseed	多様な haplotype を含めると候補が増えるため、重複・反復・population-specific path の扱いを記録する
de Bruijn / overlap graph	assembly や局所構造の表現	variant graph / pangenome graph と目的が異なるため、ノード・エッジの意味を混同しない

この章では graph reference を実装対象にはしない。読者は、線形参照向けの索引構造で学んだ「検索空間を小さくする」考え方が、パンゲノムでは「どの経路・座標系・variant representation を採用したか」という来歴管理と一体になる点を押さえる。実務上は、GRCh38、T2T-CHM13、pangenome のどれを使ったか、または相互変換したかを、第4章の mapping / variant calling、第9章の GWAS / PRS、第10章の reference identifier / checksum と接続して記録する。ただし checksum や GA4GH refget は検索アルゴリズムではなく、参照配列の識別と来歴管理を支える仕組みとして扱う。

3.2 配列アライメントアルゴリズム

動的プログラミングの理論的基盤: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

最適性原理 (Bellman): 最適解は最適部分解から構成される
重複部分問題: F(i,j)の値が複数回使用される
メモ化: 計算済み結果の保存による時間短縮

進化距離の数学的モデル:

Jukes-Cantor モデル: 全置換が等確率（μt = -3/4 ln(1-4p/3)）
Kimura 2-parameter モデル: 転移・転換の区別
HKY85 モデル: 塩基頻度の偏りを考慮
GTR モデル: 6つの置換率パラメータ

統計的有意性評価: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

E-value = Kmn exp(-λS)
- K: 検索空間の補正定数
- m,n: 配列長
- λ: Gumbel分布のパラメータ
- S: アライメントスコア

グローバルアライメント（Needleman-Wunsch法）: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

動的プログラミング:
F(i,j) = max{
    F(i-1,j-1) + s(xi, yj),  // マッチ/ミスマッチ
    F(i-1,j) + gap,          // 挿入
    F(i,j-1) + gap           // 削除
}

初学者向け最小例（Python 3.11+ / 標準ライブラリのみ）: Needleman-Wunsch法と同じ動的計画法の骨格を理解するため、まずは2本の短い配列の編集距離を求める。 🧪 概念例（Python 3.11+）

def edit_distance(seq1: str, seq2: str) -> tuple[int, list[list[int]]]:
    rows = len(seq1) + 1
    cols = len(seq2) + 1
    dp = [[0] * cols for _ in range(rows)]

    for i in range(rows):
        dp[i][0] = i
    for j in range(cols):
        dp[0][j] = j

    for i, base1 in enumerate(seq1, start=1):
        for j, base2 in enumerate(seq2, start=1):
            substitution_cost = 0 if base1 == base2 else 1
            dp[i][j] = min(
                dp[i - 1][j] + 1,                      # 削除
                dp[i][j - 1] + 1,                      # 挿入
                dp[i - 1][j - 1] + substitution_cost, # 置換または一致
            )

    return dp[-1][-1], dp


seq1 = "GATTACA"
seq2 = "GACTATA"

distance, matrix = edit_distance(seq1, seq2)
print(f"{seq1} と {seq2} の編集距離: {distance}")
print("DP行列の最終行:", matrix[-1])

🧪 概念例（期待される出力）

GATTACA と GACTATA の編集距離: 2
DP行列の最終行: [7, 6, 5, 4, 4, 3, 3, 2]

編集距離が小さいほど、2配列を一致させるための挿入・削除・置換の回数が少ないことを意味する。まずは短い配列でDP行列の更新規則を確認し、その後にスコア行列やギャップペナルティを導入すると、Needleman-Wunsch法やSmith-Waterman法の理解につながる。

ローカルアライメント（Smith-Waterman法）:

グローバル法を局所最適化に拡張
0以下の値をリセット → 局所的類似領域を検出

高速近似法（BLAST）:

シード検索: 短い完全一致領域を特定
拡張: シード周辺の類似性評価
統計的有意性検定

スコア行列の理論的基礎: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

PAM行列 (Point Accepted Mutation):
- 進化距離1 PAMは1%のアミノ酸変化
- PAM250: 平均2.5回の置換（遠縁関係）

BLOSUM行列 (BLOcks SUbstitution Matrix):
- 保存ブロック内の実際の置換頻度から導出
- BLOSUM62: 62%以上の同一性を持つ配列から構築

ギャップペナルティの最適化:

線形ギャップ: gap penalty = g × length
アフィンギャップ: penalty = open + extend × (length-1)
凸ギャップ: 生物学的により現実的なモデル

応用における重要性: 配列アライメントは生物学的類似性を定量化するための基礎技術である。進化関係、機能候補、variant / 変異の注釈候補を整理する前処理として広く使われる。厳密解法は計算量がO(nm)で大規模データに適用しにくいが、BLAST などの近似解法により、既知データベースとの比較候補を短時間で得られる場合がある。結果は候補提示であり、機能や疾患関連性の最終判断には追加検証と専門家レビューが必要である。

3.3 グラフアルゴリズム

グラフ理論の基礎概念: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

グラフG = (V, E): 頂点集合Vと辺集合E
- 次数分布: P(k) = ネットワーク中の次数kのノード割合
- クラスタ係数: C = 3×三角形数 / 連結三組数
- 最短経路長: 全ノード対間の最短距離の平均
- 中心性: 重要ノードの定量的指標

生物学的ネットワークの位相的特性:

スケールフリー性: P(k) ∝ k^(-γ), γ ≈ 2-3
スモールワールド性: 高クラスタ係数 + 短い特性経路長
階層構造: モジュール内高密度、モジュール間低密度
堅牢性: ランダム攻撃に強く、標的攻撃に脆弱

ネットワーク中心性指標: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

次数中心性: CD(v) = deg(v) / (n-1)
近接中心性: CC(v) = (n-1) / Σ d(v,u)
媒介中心性: CB(v) = Σ σst(v) / σst
固有ベクトル中心性: Av = λv （主固有ベクトル）
PageRank: PR(v) = (1-d)/n + d Σ PR(u)/deg(u)

コミュニティ検出アルゴリズム:

Modularity最適化: Q = 1/2m Σ [Aij - kikj/2m]δ(ci,cj)
Louvain法: 階層的クラスタリング
Label Propagation: ラベル伝播による高速クラスタリング
スペクトラル法: ラプラシアン行列の固有ベクトル

タンパク質相互作用ネットワーク:

ノード: タンパク質
エッジ: 相互作用関係
問題: 最短経路、クラスタリング、中心性解析

代謝パスウェイ解析:

有向グラフ: 代謝反応の流れ
重み: 反応速度定数
最適化: フラックスバランス解析（FBA）

フラックスバランス解析の数学的定式化: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

線形計画問題:
max cᵀv
subject to: Sv = 0 (定常状態制約)
           lb ≤ v ≤ ub (反応速度制約)

S: 化学量論行列
v: フラックスベクトル
c: 目的関数係数

ネットワーク解析の計算複雑性:

最短経路: Dijkstra O((V+E)logV), Floyd-Warshall O(V³)
最大フロー: Ford-Fulkerson O(E f* ), Dinic O(V²E)
最小カット: Karger O(V²logV), Stoer-Wagner O(VE+V²logV)
TSP (Traveling Salesman): NP困難、近似比3/2

システム医学への貢献: 生物学的システムはネットワーク構造を持つため、グラフ理論は疾患関連経路、候補遺伝子、薬剤標的候補を整理する研究支援に使える。中心性やモジュール解析は優先順位付けの補助になるが、治療効果や副作用を単独で保証するものではない。臨床応用には、実験的検証、専門家レビュー、施設ごとの品質管理と倫理・規制手続きが必要である。

3.4 シングルセルデータ構造

スパース行列表現: SciPy が必要な、細胞×遺伝子発現行列の保存形式を示す概念例です（CIの実行対象ではありません）。検証状況: CIおよび現在のローカル環境（Linux / Python 3.12.11）では SciPy 未導入のため、構文確認のみを対象とする。代替確認方法: SciPy 導入済み環境で python -c "import scipy.sparse" が成功することを確認してから実行する。 ⚠️ 環境依存（SciPy依存の概念例）

# 細胞×遺伝子発現行列の効率的表現
import scipy.sparse as sp

class SingleCellDataStructure:
    def __init__(self, expression_matrix):
        # CSR形式でのスパース行列保存
        self.data = sp.csr_matrix(expression_matrix)
        self.n_cells, self.n_genes = self.data.shape

k-NN グラフ構造:

細胞間類似性の効率的計算
Approximate Nearest Neighbor (ANN)アルゴリズム
LSH (Locality Sensitive Hashing)による高速化

軌跡推定のための構造:

最小全域木（MST）
主曲線（Principal Curves）
拡散マップ（Diffusion Maps）

3.5 集団遺伝学のアルゴリズム

ハプロタイプフェージング:

隠れマルコフモデル（HMM）によるアプローチ
Li-Stephensモデル
SHAPEIT, Eagle等の実装

連鎖不平衡（LD）計算: 🔁 疑似コード（実行不可: 数式/定義・処理手順の説明）

r² = (pAB - pA×pB)² / (pA×pa×pB×pb)
D' = D / Dmax

集団構造解析:

主成分分析（PCA）による祖先成分推定
ADMIXTURE/STRUCTUREアルゴリズム
f統計量（FST, f3, f4）

Source notes / 次の一歩

2026年版出典監査メモ: GRCh38 / T2T / pangenome / reference bias（HPRC、GRC、NCBI、HPRC Data Use の確認日）
第4章: ゲノム解析技術（参照配列選択、mapping、variant calling、annotation への接続）
第9章: 集団ゲノミクス（GWAS / PRS の ancestry、population structure、reference bias への接続）
第10章: データベース技術（accession、version、checksum、reference identifier の管理）

最小入出力（期待成果物/期待ログ）

入力: 小規模な配列データ（FASTA等）と、簡単な行列/グラフ（トイデータ）
出力（期待成果物）: 代表アルゴリズムの実装（例: 編集距離/アライメント/グラフ探索）と、計算量の見積もり・実測の比較メモ
期待ログ（例）: 入力サイズ・実行時間・メモリの計測結果が残る（再現できる形で記録）

前へ: 計算インフラストラクチャ

次へ: ゲノム解析技術

演習

章内のサンプルコードをベースに、入力規模を10倍にした際の実行時間/メモリを測定し、ボトルネックを特定せよ。
同じ課題を別手法（データ構造/アルゴリズム）で解き、計算量/実行時間の差を比較して考察せよ。

具体課題例

公開データを用いた再現（SRA/GEO/ArrayExpressから実在アクセッションを選定し）、前処理→主解析→結果要約まで実施。
代替ツールの比較（例: ツールA vs ツールB）。処理時間/メモリ/精度など評価指標を定義し、比較表を作成。
成果物一式（レポート、使用コマンド/パラメータ、ツール/ライブラリのバージョン、入力/出力、実行ログ、MultiQC等のレポート、図表）を添付。

Human Pangenome Reference Consortium, Human Pangenome Reference Consortium（参照日: 2026-05-13 JST） ↩
Human Pangenome Reference Consortium, HPRC Data Use: Best Practices（参照日: 2026-05-13 JST） ↩
Genome Reference Consortium, GRC FAQ（参照日: 2026-05-13 JST） ↩