from dataclasses import dataclass
from pathlib import Path
from typing import Iterator

import requests
from Bio import SeqIO, Entrez
from Bio.Seq import Seq


Entrez.email = "your@email.example"  # NCBI required


@dataclass
class ExonRegion:
    """遺伝子の1つのエクソンの情報。"""
    gene_name: str
    exon_number: int
    chromosome: str
    strand: str            # "+" または "-"
    start_hg38: int         # 0ベース
    end_hg38: int
    length: int
    sequence: str
    coding_frame: int | None  # UTRの場合はNone、それ以外の場合は0/1/2


def fetch_refseq_mrna(refseq_id: str = "NM_007294") -> str:
    """NCBIからRefSeq mRNA配列を取得（例：BRCA1 = NM_007294）。"""
    with Entrez.efetch(db="nucleotide", id=refseq_id, rettype="fasta", retmode="text") as h:
        record = SeqIO.read(h, "fasta")
    return str(record.seq).upper()


def fetch_gene_exons_ucsc(gene_symbol: str, assembly: str = "hg38") -> list[ExonRegion]:
    """UCSC Genome Browser REST APIを介して遺伝子エクソン座標を取得。
    
    実際には、ローカルでGENCODEアノテーションGTFを解析する方がはるかに安定しています。
    """
    # UCSC Table Browserの代替案：MyGene.info API
    resp = requests.get(
        f"https://mygene.info/v3/query?q=symbol:{gene_symbol}&species=human&fields=exons",
        timeout=30,
    )
    resp.raise_for_status()
    hits = resp.json().get("hits", [])
    if not hits:
        return []
    
    exons_info = hits[0].get("exons", [])
    if not exons_info:
        return []
    
    canonical = exons_info[0]  # カノニカルな転写物
    exon_regions = []
    for i, (start, end) in enumerate(canonical.get("position", [])):
        # UCSCから実際の配列を取得（DAS APIまたはローカルBSgenome）
        seq = _fetch_ucsc_dna(canonical["chr"], start, end, assembly)
        exon_regions.append(ExonRegion(
            gene_name=gene_symbol,
            exon_number=i + 1,
            chromosome=canonical["chr"],
            strand=canonical.get("strand", "+"),
            start_hg38=start, end_hg38=end,
            length=end - start,
            sequence=seq,
            coding_frame=None,  # 詳細な決定にはCDS座標が必要
        ))
    return exon_regions


def _fetch_ucsc_dna(chrom: str, start: int, end: int, assembly: str = "hg38") -> str:
    """UCSC DAS APIを介して特定の領域の配列を返します。"""
    url = f"https://api.genome.ucsc.edu/getData/sequence?genome={assembly};chrom={chrom};start={start};end={end}"
    resp = requests.get(url, timeout=30)
    resp.raise_for_status()
    return resp.json().get("dna", "").upper()

from dataclasses import dataclass
from Bio.Seq import Seq

@dataclass
class GRNACandidate:
    gene_name: str
    exon_number: int
    protospacer: str        # 20 bpのターゲット配列
    pam: str                # 3 bpのPAM
    strand: str             # "+" または "-"
    exon_position: int      # エクソン配列内の開始位置（0ベース）
    genome_position: int    # hg38ゲノム座標（0ベース）
    context_50bp: str       # 15 bp上流 + gRNA + PAM + 15 bp下流 = 53 bp（CRISPRonなどのトレーニング入力）


def scan_pam_ngg(sequence: str, gene_name: str, exon: ExonRegion) -> list[GRNACandidate]:
    """NGG PAM + 20 bpプロトスペーサー候補を両方の鎖から抽出します。"""
    candidates = []
    seq_plus = sequence.upper()
    seq_minus = str(Seq(seq_plus).reverse_complement())
    
    def _scan(strand_seq: str, strand_label: str, seq_len: int) -> None:
        for i in range(20, seq_len - 3):
            pam = strand_seq[i:i+3]
            if pam[1:3] != "GG":
                continue
            protospacer = strand_seq[i-20:i]
            if "N" in protospacer:
                continue
            context_start = max(0, i - 35)
            context_end = min(seq_len, i + 18)
            context = strand_seq[context_start:context_end]
            # ゲノム座標を計算します（鎖に依存した方向の反転）
            if strand_label == "+":
                genome_pos = exon.start_hg38 + (i - 20)
            else:
                genome_pos = exon.end_hg38 - i
            candidates.append(GRNACandidate(
                gene_name=gene_name,
                exon_number=exon.exon_number,
                protospacer=protospacer,
                pam=pam,
                strand=strand_label,
                exon_position=i - 20,
                genome_position=genome_pos,
                context_50bp=context,
            ))
    
    _scan(seq_plus, "+", len(seq_plus))
    _scan(seq_minus, "-", len(seq_minus))
    return candidates


def scan_gene_all_exons(exons: list[ExonRegion]) -> list[GRNACandidate]:
    all_cands = []
    for exon in exons:
        all_cands.extend(scan_pam_ngg(exon.sequence, exons[0].gene_name, exon))
    return all_cands

from Bio.SeqUtils import GC
from Bio.Seq import Seq

def check_self_complementarity(seq: str, min_stem: int = 4) -> bool:
    """自己相補性のチェック：5'末端の4 bpが3'末端の4 bpの逆相補配列である場合、ヘアピンのリスクがあります。"""
    if len(seq) < min_stem * 2:
        return False
    stem_5 = seq[:min_stem]
    stem_3 = seq[-min_stem:]
    return stem_5 == str(Seq(stem_3).reverse_complement())


def filter_grna_physicochemical(
    candidates: list[GRNACandidate],
    min_gc: float = 30.0,
    max_gc: float = 70.0,
    max_poly_t: int = 3,
    prefer_g_start: bool = True,
    check_hairpin: bool = True,
) -> list[GRNACandidate]:
    """物理化学的フィルター。"""
    filtered = []
    for c in candidates:
        # GC含有量
        gc = GC(c.protospacer)
        if not (min_gc <= gc <= max_gc):
            continue
        # ポリ-T
        if "T" * (max_poly_t + 1) in c.protospacer:
            continue
        # G開始（推奨）
        if prefer_g_start and c.protospacer[0] != "G":
            # 完全には除外されませんが、優先順位を下げます（ここで通過します）。
            pass
        # ヘアピン
        if check_hairpin and check_self_complementarity(c.protospacer):
            continue
        filtered.append(c)
    return filtered

import numpy as np
import rpy2.robjects as ro
from rpy2.robjects.packages import importr
from rpy2.robjects.vectors import StrVector


class CrisprScoreEnsemble:
    """RのcrisprScoreから複数のスコアラーをPythonで呼び出します。"""

    def __init__(self):
        # Rパッケージをロードします（BiocManager::installは事前にRコンソールで必要）。
        self.crisprscore = importr("crisprScore")
        self.base = importr("base")

    def crispron(self, contexts: list[str]) -> list[float]:
        """CRISPRonオンターゲットスコア。入力：各項目は30 bpコンテキストです（4 bp上流 + 20 bpプロトスペーサー + 3 bp PAM + 3 bp下流）。"""
        r_input = StrVector(contexts)
        try:
            result = self.crisprscore.getCRISPRonScores(r_input)
            return list(result)
        except Exception as e:
            print(f"CRISPRon failed: {e}")
            return [float("nan")] * len(contexts)

    def deep_hf(self, protospacers_with_pam: list[str], enzyme: str = "WT") -> list[float]:
        """DeepHFスコア。入力：各23 bp（20 bpプロトスペーサー + 3 bp PAM）。enzyme：'WT'、'ESP'、'HF'。"""
        r_input = StrVector(protospacers_with_pam)
        try:
            result = self.crisprscore.getDeepHFScores(r_input, enzyme=enzyme)
            return list(result)
        except Exception as e:
            print(f"DeepHF failed: {e}")
            return [float("nan")] * len(protospacers_with_pam)

    def crispron_be(self, contexts: list[str], editor: str = "ABE8e") -> list[float]:
        """CRISPRon-BEベース編集効率。
        
        エディターオプション：'ABE8e'（アデニン）、'BE4max'（シトシン）など。
        """
        r_input = StrVector(contexts)
        try:
            result = self.crisprscore.getCRISPRonBEScores(r_input, editor=editor)
            return list(result)
        except Exception as e:
            print(f"CRISPRon-BE failed: {e}")
            return [float("nan")] * len(contexts)

    def cfd_off_target(self, protospacer: str, target_dna: str) -> float:
        """切断頻度決定（Doench 2016 [5]）。"""
        try:
            result = self.crisprscore.getCFDScores(
                StrVector([protospacer]),
                StrVector([target_dna]),
            )
            return float(list(result)[0])
        except Exception:
            return 0.0

    def mit_hsu_zhang(self, protospacer: str, target_dna: str) -> float:
        """MIT/Hsu-Zhangオフターゲットスコア（位置依存ミスマッチペナルティ）。"""
        try:
            result = self.crisprscore.getMITScores(
                StrVector([protospacer]),
                StrVector([target_dna]),
            )
            return float(list(result)[0])
        except Exception:
            return 0.0

# Rコンソール内
if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install(c(
    "crisprScore", "crisprDesign",
    "BSgenome.Hsapiens.UCSC.hg38",
    "crisprBase"
))

import requests

def query_crispor_offtargets(
    protospacer: str,
    genome: str = "hg38",
    max_mismatches: int = 4,
    timeout: int = 300,
) -> list[dict]:
    """CRISPORリモートAPIを介してゲノム全体のオフターゲットスキャン。
    
    返されるもの：[{genome_pos、mismatches、cfd_score、gene}、...]
    """
    # 実際のエンドポイント：http://crispor.tefor.net/crispor.py
    # これは概念的なラッパーです。実際のAPIはHTMLを返すため、解析が必要です。
    try:
        resp = requests.get(
            "http://crispor.tefor.net/crispor.py",
            params={
                "seq": protospacer,
                "org": genome,
                "pam": "NGG",
                "showAllOTs": "1",
            },
            timeout=timeout,
        )
        resp.raise_for_status()
        # 実際のCRISPORレスポンスはHTML/TSVです。ここでは例として解析を示します。
        # 実際には、BeautifulSoupまたはローカルのcrispor CLIを使用してTSVを取得します。
        return _parse_crispor_html(resp.text)
    except Exception as e:
        print(f"CRISPORオフターゲットスキャンに失敗しました（{protospacer}）：{e}")
        return []


def _parse_crispor_html(html: str) -> list[dict]:
    """CRISPOR HTMLレスポンスを解析します（概念的なスタブ、実際にはBeautifulSoupを使用します）。"""
    return []  # 本番環境で実装します。


def compute_off_target_summary(off_targets: list[dict]) -> dict:
    """オフターゲットスキャン結果の概要。"""
    if not off_targets:
        return {
            "n_offtargets_mm4": 0,
            "max_cfd_score": 0.0,
            "high_risk_genes": [],
            "genome_safety_score": 1.0,
        }
    max_cfd = max(ot["cfd_score"] for ot in off_targets)
    high_risk = [ot["gene"] for ot in off_targets if ot["cfd_score"] > 0.5 and ot["gene"]]
    # 全体的な安全性：最大CFDを反転します（0 = 危険、1 = 安全）
    safety = 1.0 - max_cfd
    return {
        "n_offtargets_mm4": len(off_targets),
        "max_cfd_score": max_cfd,
        "high_risk_genes": high_risk[:5],
        "genome_safety_score": safety,
    }

import pandas as pd
import numpy as np


DEFAULT_WEIGHTS = {
    "crispron": 0.30,
    "deep_hf": 0.25,
    "crispron_be": 0.10,  # BEシナリオではない場合
    "off_target_safety": 0.35,
}


def build_ensemble_ranking(
    candidates: list[GRNACandidate],
    on_scores: dict[str, list[float]],
    off_target_summaries: list[dict],
    weights: dict[str, float] = DEFAULT_WEIGHTS,
) -> pd.DataFrame:
    """統合されたスコアラーランキングDataFrame。"""
    df = pd.DataFrame([{
        "gene": c.gene_name,
        "exon": c.exon_number,
        "protospacer": c.protospacer,
        "pam": c.pam,
        "strand": c.strand,
        "genome_position": c.genome_position,
        "context_50bp": c.context_50bp,
    } for c in candidates])
    
    for name, scores in on_scores.items():
        df[f"on_{name}"] = scores
    
    df["off_n_mm4"] = [s["n_offtargets_mm4"] for s in off_target_summaries]
    df["off_max_cfd"] = [s["max_cfd_score"] for s in off_target_summaries]
    df["off_safety"] = [s["genome_safety_score"] for s in off_target_summaries]
    df["off_high_risk_genes"] = [
        ";".join(s["high_risk_genes"]) for s in off_target_summaries
    ]
    
    # 重み付け平均（正規化後）
    def _normalize(col):
        c = df[col].fillna(df[col].median())
        return (c - c.min()) / (c.max() - c.min() + 1e-8)
    
    final = np.zeros(len(df))
    if "on_crispron" in df.columns:
        final += weights.get("crispron", 0.0) * _normalize("on_crispron")
    if "on_deep_hf" in df.columns:
        final += weights.get("deep_hf", 0.0) * _normalize("on_deep_hf")
    if "on_crispron_be" in df.columns:
        final += weights.get("crispron_be", 0.0) * _normalize("on_crispron_be")
    final += weights["off_target_safety"] * df["off_safety"].fillna(0.5)
    
    df["final_score"] = final
    return df.sort_values("final_score", ascending=False).reset_index(drop=True)

from pathlib import Path
import pandas as pd

def generate_lab_report(
    ranked: pd.DataFrame,
    output_md: Path,
    output_csv: Path,
    top_k: int = 5,
) -> None:
    """ラボで注文可能なレポートを生成します。"""
    top = ranked.head(top_k)
    
    lines = [
        f"# 上位{top_k} gRNA候補の注文レポート",
        "",
        f"ターゲット遺伝子：**{top.iloc[0]['gene']}**",
        f"生成日時：（パイプラインの実行時に記入）",
        "",
        "## 概要テーブル",
        "",
        "| 順位 | エクソン | プロトスペーサー（20 bp） | PAM | 鎖 | オン（CRISPRon） | オン（DeepHF） | オフターゲット安全性 | 最終スコア |",
        "|:----:|:----:|:-------------------:|:---:|:------:|:-------------:|:-----------:|:----------:|:-----:|",
    ]
    for i, row in top.iterrows():
        lines.append(
            f"| {i+1} | {row['exon']} | `{row['protospacer']}` | {row['pam']} | {row['strand']} | "
            f"{row.get('on_crispron', float('nan')):.3f} | {row.get('on_deep_hf', float('nan')):.3f} | "
            f"{row['off_safety']:.3f} | **{row['final_score']:.3f}** |"
        )
    
    lines += [
        "",
        "## 注文配列（5'→3'）",
        "",
        "プロトスペーサー配列の先頭にGがない場合は、先頭にGを追加するか、ExoScribeなどの代替プロモーターを使用します。",
        "",
    ]
    for i, row in top.iterrows():
        spacer = row["protospacer"]
        if spacer[0] != "G":
            forge = f"G{spacer}"
        else:
            forge = spacer
        lines.append(f"### 順位{i+1}")
        lines.append(f"- プロトスペーサー（スパースのみ、20 bp）：`{spacer}`")
        lines.append(f"- 注文用（必要な場合はG接頭辞付き）：`{forge}`（{len(forge)} bp）")
        lines.append(f"- PAM：`{row['pam']}`")
        lines.append(f"- ゲノム位置（hg38）：{row['genome_position']}")
        if row["off_high_risk_genes"]:
            lines.append(f"- ⚠️ 高リスクのオフターゲット遺伝子：{row['off_high_risk_genes']}")
        lines.append("")
    
    lines += [
        "## 事前実験チェックリスト",
        "",
        "- [ ] 注文前にアラインメントを再確認（BLASTnを使用してhg38 / GRCh38.p14に照合）",
        "- [ ] SpCas9-HF1ベクターを準備します（例：pX330またはplentiCRISPR v2）。",
        "- [ ] ターゲット細胞株でマイコプラズマ試験を完了します（K562）。",
        "- [ ] T7E1アッセイまたはアンプリコン深層シーケンス用のプライマーを設計します。",
        "- [ ] 非ターゲットgRNA（例：sgLacZ）を並行して注文します。",
    ]
    
    output_md.write_text("\n".join(lines), encoding="utf-8")
    top.to_csv(output_csv, index=False)

def full_pipeline(
    gene_symbol: str,
    refseq_id: str,
    output_prefix: str,
    editor: str | None = None,  # None：ヌクレアーゼKO、"ABE8e"/"BE4max"：ベース編集
    top_k: int = 5,
) -> pd.DataFrame:
    print(f"[1/7] {gene_symbol}エクソン情報取得（RefSeq {refseq_id}）")
    # 実際には、GENCODE GTFからのエクソン座標をローカルで解析することをお勧めします。
    exons = fetch_gene_exons_ucsc(gene_symbol)
    if not exons:
        # フォールバック：mRNA全体を1つの「エクソン」として扱います（単純なデモンストレーション）。
        seq = fetch_refseq_mrna(refseq_id)
        exons = [ExonRegion(
            gene_name=gene_symbol, exon_number=1, chromosome="",
            strand="+", start_hg38=0, end_hg38=len(seq),
            length=len(seq), sequence=seq, coding_frame=0,
        )]
    print(f"      エクソン数：{len(exons)}")
    
    print("[2/7] PAM（NGG）スキャン")
    all_cands = scan_gene_all_exons(exons)
    print(f"      生の候補数：{len(all_cands)}")
    
    print("[3/7] 物理化学的フィルター")
    filtered = filter_grna_physicochemical(all_cands)
    print(f"      フィルター通過数：{len(filtered)}")
    
    print("[4/7] オンターゲットスコアラーアンサンブル")
    scorer = CrisprScoreEnsemble()
    contexts = [c.context_50bp for c in filtered]
    proto_pam = [c.protospacer + c.pam for c in filtered]
    on_scores = {
        "crispron": scorer.crispron(contexts),
        "deep_hf": scorer.deep_hf(proto_pam, enzyme="HF"),  # SpCas9-HF1
    }
    if editor:
        on_scores["crispron_be"] = scorer.crispron_be(contexts, editor=editor)
    
    print("[5/7] オフターゲットスキャン（CRISPORリモート）")
    off_summaries = []
    for i, c in enumerate(filtered):
        if i % 50 == 0:
            print(f"      進捗：{i}/{len(filtered)}")
        offs = query_crispor_offtargets(c.protospacer)
        off_summaries.append(compute_off_target_summary(offs))
    
    print("[6/7] アンサンブルランキング")
    ranked = build_ensemble_ranking(filtered, on_scores, off_summaries)
    ranked.to_csv(f"{output_prefix}_all_candidates.csv", index=False)
    
    print(f"[7/7] 上位{top_k}注文レポート")
    generate_lab_report(
        ranked,
        Path(f"{output_prefix}_top{top_k}.md"),
        Path(f"{output_prefix}_top{top_k}.csv"),
        top_k=top_k,
    )
    print(f"      完了：{output_prefix}_top{top_k}.md")
    return ranked


# 実行例（BRCA1 KOシナリオ）
# ranked = full_pipeline(
#     gene_symbol="BRCA1",
#     refseq_id="NM_007294",
#     output_prefix="brca1_k562",
#     editor=None,  # SpCas9-HF1ヌクレアーゼKO
#     top_k=5,
# )

[1/7] BRCA1エクソン情報取得（RefSeq NM_007294）
      エクソン数：24
[2/7] PAM（NGG）スキャン
      生の候補数：423
[3/7] 物理化学的フィルター
      フィルター通過数：197
[4/7] オンターゲットスコアラーアンサンブル
[5/7] オフターゲットスキャン（CRISPORリモート）
      進捗：0/197
      進捗：50/197
      ...
[6/7] アンサンブルランキング
[7/7] 上位5注文レポート
      完了：brca1_k562_top5.md

## パフォーマンス・コスト・既知の失敗事例

### パフォーマンスの参考 (公開されたベンチマークの引用)

| スコアラー | ベンチマークデータセット | スピアマンのρ | 備考 | 出典 |
|---|---|:----------:|---|---|
| Rule 4 (Doench 2014) | 独自のベンチマーク | 0.42 | 手動で設計されたベースライン | Doench et al., Nat Biotechnol 2014 |
| DeepCRISPR | Kim 2017データセット | 0.65 | 1D CNN | Chuai et al., Genome Biol 2018 [1] |
| CRISPRon (Xiang 2021) | Kim 2019 | 0.72 | コンテキストを考慮したCNN | Xiang et al., Nat Commun 2021 [2] |
| DeepHF (WT SpCas9) | Wang 2019 | 0.73 | RNN + アテンション | Wang et al., Nat Commun 2019 [3] |
| DeepHF (SpCas9-HF1) | Wang 2019 | 0.75 | HF1バリアントに特化 | Wang et al. 2019 [3] |
| CRISPRon-BE (ABE8e) | ABEmaxデータセット | 0.68 | 塩基編集 | Kim et al. 2022 [4] |
| BE-Hive | BE4/ABEデータセット | 0.75 | 編集ウィンドウ予測に優れている | Arbab et al. 2020 |
| アンサンブル (本稿で近似) | ベンチマークの再現 | 0.77〜0.82 (推定) | 複数のスコアラーの平均 | コミュニティベンチマーク [7] |

### 推定される学習モデルの再現にかかるコスト

- APIコストはかからない（完全にローカルまたは無料のCRISPORリモート）。
- 小規模な消費向けGPUで500個のgRNAをスコアリングする場合、約1〜3分。
- CPUで代替する場合、5〜15分。
- CRISPORリモートAPI：推奨は1秒あたり1リクエスト未満、500個のgRNAのスクリーニングには約10〜20分。
- ローカルのcrispritz代替案：CPUで30分 + 20GBのゲノムインデックスのダウンロード。

### 5つの既知の失敗事例（コミュニティ/論文から収集）

1. **rpy2環境の競合（R 4.3 + Bioconductor 3.18 + rpy2 3.5）**
   症状：Rライブラリのロードに失敗する（例：`libR.so symbol lookup error`、`PermissionError: could not load package`）。
   原因：condaとシステムRの間でRライブラリのパスが二重に登録されている（特にmacOS/Linuxの混合使用の場合に多い）。
   解決策：（a）Rをconda-forgeを通じて、conda環境内で一貫してインストールする（`conda install -c conda-forge r-base bioconductor-crisprscore`）、（b）`R_HOME`と`LD_LIBRARY_PATH`の環境変数を明示的に設定する、（c）Dockerコンテナで分離する（rocker/tidyverse + crisprScoreのインストール）、（d）rpy2をバイパスして、Rscriptを直接subprocess経由で呼び出す。
   出典：rpy2 GitHub Issues [8]、BioconductorサポートフォーラムのcrisprScoreのスレッド。

2. **計算されたCFDスコアと実際の実験におけるオフターゲットの結果の乖離**
   症状：CFDで安全（低いスコア）と判断されたgRNAが、実際の実験で予期しないオフターゲット編集を引き起こす（GUIDE-seq、CIRCLE-seq）。
   原因：CFDは、学習されたミスマッチ位置の重みを基にしている（Doench 2016 [5]）が、クロマチン状態、DNAメチル化、ヌクレオソームの位置を反映していない。
   解決策：（a）複数のオフターゲットアルゴリズム（CFD + MIT/Hsu-Zhang + Elevation）を組み合わせる、（b）TSSやエンハンサー領域に対して、別途ペナルティを適用する（ENCODEアノテーションを参照）、（c）GUIDE-seqやDigenome-seqなどの測定されたオフターゲットデータで検証する、（d）最終的な選択の前に、上位3つの候補を並行して測定する。
   出典：Tsai et al. "GUIDE-seq" Nat Biotechnol 2015 [9]、Doench 2016 CFD論文 [5]、CRISPORドキュメント [6]。

3. **BE（塩基編集）の編集ウィンドウの誤解釈 + バイスタンダー効果**
   症状：CRISPRon-BEでスコアリングされた塩基エディターが、意図された編集ターゲットに加えて、他のヌクレオチドも編集する（バイスタンダー）。
   原因：BEは、プロトスペーサー内の特定のウィンドウ（例：4〜8塩基）を編集すると想定されているが、実際には5〜10塩基まで拡張し、アデニンデアミナーゼ（ABE）は特定の配列コンテキストにバイアスがかかっている。
   解決策：（a）複数のモデル（BE-Hive、BE-DICT、CRISPRon-BE）の編集ウィンドウのアノテーションを組み合わせる、（b）実験設計において、バイスタンダー候補を別途マークする、（c）Prime editing（PE）スコアラー（PRIDICT）を使用して、正確な編集を行うことを検討する、（d）コドン単位でシミュレーションを行い、編集ウィンドウ内のターゲットのC/Aの位置が不要なコドンの変化を引き起こさないか確認する。
   出典：Anzalone et al. "Prime editing" Nature 2019 [10]、Arbab et al. "BE-Hive" Cell 2020。

4. **CRISPORリモートAPIのレート制限 + 応答遅延**
   症状：500個のgRNAを連続してスキャンすると、半数以上が失敗またはタイムアウトする。
   原因：CRISPORリモートサーバーは、CPUの制約を受けている。これはコミュニティリソースである。
   解決策：（a）ローカルのcrispritz + BWA + ゲノムインデックスを使用してバッチスキャンする（20GBのゲノムインデックスのダウンロードが必要）、（b）CRISPORをローカルにインストールする（Pythonスクリプトをローカルにインストールする）、（c）5秒以上の間隔でリクエストする、（d）上位20個の候補のみをCRISPORリモートに送信する（100個以上の候補をローカルで処理する）。
   出典：CRISPORの公式ドキュメント「バッチ使用」セクション [6]。

5. **BRCA1のオルタナティブスプライシングによるエクソン座標の不一致**
   症状：BRCA1には複数のスプライスバリアントがある。NM_007294がカノニカルであるが、他の転写物（NM_007297など）は異なるエクソン数/座標を持つ。
   原因：MyGene.infoは1つのカノニカル転写物のみを返す。ターゲットの組織/細胞株で実際に発現しているアイソフォームが異なる場合がある。
   解決策：（a）GENCODE GTFから複数の転写物の座標を解析する、（b）ターゲットの組織/細胞株の主要なアイソフォームを測定された発現データ（GTEx、Human Protein Atlas）を使用して確認する、（c）複数のアイソフォームに共通のエクソン（構成的エクソン）を優先する。
   出典：GENCODEアノテーションドキュメント [11]、UCSC Genome BrowserのBRCA1トラック。

## 拡張アイデア

- **Prime Editing (PE) スコアリング拡張機能**: PE は複雑な 5' pegRNA デザインを必要としますが、高精度な編集に最適です。最新のスコアリングツール（例: PRIDICT、DeepPE）を組み込みます。
- **CRISPRi/CRISPRa (活性化/干渉)**: KO の代わりに、遺伝子発現の調節を行います。`getCrispraiScores` を crisprScore で使用し、dCas9-VP64（活性化）または dCas9-KRAB（干渉）を適用します。
- **マルチ gRNA (多重) 組み合わせの最適化**: 複数の遺伝子を同時にノックアウトする場合、相互のオフターゲット干渉を最小限に抑えます（組み合わせ最適化、整数計画法）。
- **患者特異的な遺伝的背景**: 特定の SNPs を持つ患者では、オフターゲットの位置が異なります（個別化されたオフターゲット、gnomAD や UK Biobank のアノテーションを参照）。
- **ヌクレオソーム配置の統合**: MNase-seq データを使用して、ヌクレオソームで覆われた領域に対して個別のペナルティを適用します。
- **T7 in vitro 合成 vs オリゴプール**: gRNA の注文形式（シングル vs プール）に応じて、上位 5 件と上位 200 件に対するスコアリング戦略が異なります。

## 次のトピック

- トピック 07 `drug-target-gnn`: GNN を使用した薬剤-標的予測により、gRNA でノックアウトする候補標的を事前に絞り込みます。
- トピック 10 `med-llm-reproduction`: 複数の LLM に gRNA スコアリングのベンチマークを再現させ、評価を行います。
- トピック 12 `single-cell-perturbation`: 実際にノックアウトする前に、gRNA 候補に対する細胞応答を in silico で予測します（Perturb-seq シミュレーション）。
- トピック 14 `bio-mcp-agent`: gRNA デザインを MCP ツールとして公開し、「BRCA1 エキソン 3〜5 から KO gRNA を抽出する」などの自律的な実行を可能にします。

## 参考文献

1. Chuai G, Ma H, Yan J, et al. "DeepCRISPR: optimized CRISPR guide RNA design by deep learning." Genome Biology 2018. `https://genomebiology.biomedcentral.com/articles/10.1186/s13059-018-1459-4`
2. Xiang X, Corsi GI, Anthon C, et al. "Enhancing CRISPR-Cas9 gRNA efficiency prediction by deep learning (CRISPRon)." Nature Communications 2021. `https://www.nature.com/articles/s41467-021-23576-0`
3. Wang D, Zhang C, Wang B, et al. "Optimized CRISPR guide RNA design for two high-fidelity Cas9 variants by deep learning (DeepHF)." Nature Communications 2019. `https://www.nature.com/articles/s41467-019-12281-8`
4. Kim HK, Yu G, Park J, et al. "Predicting the efficiency of prime editing guide RNAs in human cells (PRIDICT); base editing companion (CRISPRon-BE)." Nature Biotechnology 2023 (関連部分). `https://www.nature.com/articles/s41587-022-01613-7`
5. Doench JG, Fusi N, Sullender M, et al. "Optimized sgRNA design to maximize activity and minimize off-target effects of CRISPR-Cas9 (CFD original paper)." Nature Biotechnology 2016. `https://www.nature.com/articles/nbt.3437`
6. CRISPOR ウェブツールとドキュメント: `http://crispor.tefor.net/` / `https://github.com/maximilianh/crisporWebsite`
7. crisprScore Bioconductor パッケージ: `https://bioconductor.org/packages/release/bioc/html/crisprScore.html`
8. rpy2 GitHub Issues: `https://github.com/rpy2/rpy2/issues`
9. Tsai SQ, Zheng Z, Nguyen NT, et al. "GUIDE-seq enables genome-wide profiling of off-target cleavage by CRISPR-Cas nucleases." Nature Biotechnology 2015. `https://www.nature.com/articles/nbt.3117`
10. Anzalone AV, Randolph PB, Davis JR, et al. "Search-and-replace genome editing without double-strand breaks or donor DNA (Prime Editing)." Nature 2019. `https://www.nature.com/articles/s41586-019-1711-4`
11. GENCODE ヒトゲノムアノテーション: `https://www.gencodegenes.org/human/`
12. Kim HK et al. "Predicting the efficiency of prime editing guide RNAs in human cells (PRIDICT)." Nat Biotechnol 2023.
13. Broad Institute CRISPResso2 (編集結果の分析): `https://github.com/pinellolab/CRISPResso2`
14. Bioconductor BSgenome.Hsapiens.UCSC.hg38: `https://bioconductor.org/packages/release/data/annotation/html/BSgenome.Hsapiens.UCSC.hg38.html`
15. RTH-tools CRISPRon-BE GitHub: `https://github.com/RTH-tools/crispron-BE`
16. BE-Hive: Arbab M et al. "Determinants of Base Editing Outcomes from Target Library Analysis and Machine Learning." Cell 2020. `https://github.com/maxwshen/be_predict_efficiency`
17. Addgene sgRNA 注文標準プロトコル: `https://www.addgene.org/guides/crispr/`
18. crispritz (ローカルオフターゲットスキャン): `https://github.com/pinellolab/CRISPRitz`
19. NCBI RefSeq NM_007294 (BRCA1 カノニカル mRNA): `https://www.ncbi.nlm.nih.gov/nuccore/NM_007294`
20. MyGene.info API: `https://mygene.info/v3/api`