tooluniverse-protein-sae-variant-interpretation

タンパク質 SAE 変異の解釈

ESMC-6B タンパク質言語モデルから得られる Sparse Autoencoder（SAE）の特徴活性化を、参照配列と変異配列で比較し、単一のミスセンス変異を解釈します。SAE 特徴は、モデルの隠れ状態における解釈可能な潜在次元です。多くの特徴は、生物学的に意味のあるパターン（活性部位、リガンド結合ポケット、PTM シーケオン、構造モチーフ）に対して活性化します。

---

このスキルを使用する場面

次のような場合に使用します。

• 「変異 X（ミスセンス）はなぜ機能喪失を引き起こすのか」と尋ねられ、病原性スコアだけではなく機序に基づく回答が必要な場合
• AlphaMissense / ClinVar で「有害」とされた変異について、どの機能的特徴が損なわれるのか（触媒機能？結合？PTM 部位？構造？）を知りたい場合
• 特定の残基における参照配列と変異配列のタンパク質表現を比較したい場合
• 構造上はわずかな変化（単一のアミノ酸）が大きな機能的影響を及ぼす理由を解釈したい場合

対象外（代わりに他のスキルを使用してください）。

• ACMG による病原性分類 → tooluniverse-variant-interpretation
• 制御領域 / 非コード変異 → tooluniverse-variant-to-mechanism
• 機序を伴わない変異と疾患の関連 → tooluniverse-gene-disease-association
• がんに特化した変異の解釈 → tooluniverse-cancer-variant-interpretation

---

必要な入力

入力	形式	例
タンパク質識別子	UniProt アクセッションまたは HGNC 遺伝子シンボル	P04637 または TP53
変異	1 文字コード：{ref_aa}{position_1idx}{alt_aa}	R175H

任意の入力：

• ウィンドウ半径（デフォルト 8）：解析する変異周辺の残基数
• 参照タンパク質配列（すでに判明している場合は UniProt の検索を省略）

---

前提条件

• 有効な EvolutionaryScale Forge トークンを設定した ESM_API_KEY 環境変数（https://forge.evolutionaryscale.ai）
• SAE をサポートする esm パッケージ：

  pip install 'esm @ git+https://github.com/evolutionaryscale/esm@ee891c52'
  

  PyPI 版の esm には、まだ SAEConfig が含まれていません。上流の機能ブランチからインストールしてください。

ライセンスに関する注意：Forge の SAE 出力には Cambrian Inference Clickthrough License が適用されます。別途商用契約が適用される場合を除き、非商用 / 学術研究での利用に限られます。

---

ワークフロー（5 ステップ）

ステップ 1：遺伝子 → UniProt アクセッションを解決する（必要な場合）

ユーザーから遺伝子シンボルが指定された場合は、レビュー済みのヒトアクセッションを解決します。

UniProt_search(
    query="gene:TP53 AND organism_id:9606 AND reviewed:true",
    fields=["accession", "gene_names", "protein_name"],
)
# → returns accession P04637 as the canonical reviewed human TP53

ステップ 2：標準参照配列を取得する

UniProt_get_sequence_by_accession(accession="P04637")
# → returns the canonical isoform 1 sequence (393 AA for TP53)

ステップ 3：参照残基を検証し、変異配列を作成する

変異文字列（例：R175H）を解析します。

• ref_aa = "R"、position = 175（1 始まり）、alt_aa = "H"
• ref_sequence[174] == "R" であることを確認します（Python は 0 始まり、位置は 1 始まりです）
• 変異配列を作成します：mutant_sequence = ref_sequence[:174] + "H" + ref_sequence[175:]

参照残基が一致しない場合は、明示的なエラーを返してください。誤った位置を黙って変異させてはいけません。

クイックパス（変異解析に推奨）：複合ツール

標準的な変異解釈の用途では、必要な情報量に応じて 2 つの複合ツールのいずれかを使用します。

1 回の呼び出しで最大限の情報を取得 — 変異による破綻、特徴ごとの生物学的カテゴリラベル、機序の要約：

ESM_explain_variant_mechanism(
    sequence=ref_sequence,
    position=175, ref_aa="R", alt_aa="H",
    window=8,
    top_k_features=5,       # describe top 5 lost + top 5 gained
)
# data["mechanism_summary"] e.g. "Disrupted feature categories (lost): catalytic=2, ligand-binding=1"
# data["lost_feature_categories"] / ["gained_feature_categories"] — category counts
# data["top_features_lost"] / ["top_features_gained"] — per-feature delta + category + confidence

これは変異機序レポートの標準的な選択肢です。ESM_score_variant_sae_disruption と N 回の ESM_describe_sae_feature 呼び出しを連結せずに済みます。差分のみが必要な場合は、include_descriptions=false を設定してラベル付けを省略してください（Forge 呼び出しは 2 回のみです）。

生の差分のみ（カテゴリラベルなし、describe 呼び出しなし — より高速）：

ESM_score_variant_sae_disruption(
    sequence=ref_sequence,
    position=175, ref_aa="R", alt_aa="H",
    window=8, top_k_features=10,
)
# → returns top_features_lost + top_features_gained ranked by |delta|
#   plus ref / mut activation sums per feature

指定された位置の配列と ref_aa が一致しない場合、どちらのツールも明確なエラーを返します（アイソフォームの指定が誤っているか、変異のラベルが間違っています）。以下の詳細な手順は、残基ごとの生の特徴を調べる場合に使用します。

複数の変異を一度に解析する（例：残基 175 の飽和変異 — 19 種類すべての置換）：

alts = "ACDEFGHIKLMNPQRSTVWY".replace("R", "")
variants = [{"position": 175, "ref_aa": "R", "alt_aa": a} for a in alts]
ESM_score_variant_sae_batch(sequence=ref_sequence, variants=variants, top_k_features=5)
# 1 + 19 = 20 Forge calls, not 38

Step 4（長い経路）：参照配列と変異体のSAE特徴を取得する

ref_features = ESM_get_sae_features(
    sequence=ref_sequence,
    position=175,          # 1-indexed mutation position
    window=8,              # +/- 8 residues = 17-residue window
    top_k_per_residue=64,  # full sparsity (k=64 is the SAE's actual k)
)

mut_features = ESM_get_sae_features(
    sequence=mutant_sequence,
    position=175,
    window=8,
    top_k_per_residue=64,
)

各呼び出しは、ウィンドウ内の残基について、{residue_idx_1based, active_features: [{feature_id, activation}]}のリストを返します。通常のレイテンシー：約1～3秒／回（合計約2～6秒）。Forgeのコスト：2クレジット（1回の呼び出しにつき1クレジット）。

Step 5：特徴ごとの活性化の差分を計算する

ウィンドウ全体で、獲得／喪失した特徴を集計します。

# Build per-feature activation arrays across the window
def feature_to_window_sum(features_response):
    sums = {}  # feature_id -> sum of activations across all residues in window
    for residue in features_response["data"]["activations"]:
        for f in residue["active_features"]:
            sums[f["feature_id"]] = sums.get(f["feature_id"], 0.0) + f["activation"]
    return sums

ref_sums = feature_to_window_sum(ref_features)
mut_sums = feature_to_window_sum(mut_features)

# Delta = mut - ref. Positive = gained on mutation. Negative = lost.
all_features = set(ref_sums) | set(mut_sums)
deltas = {
    f: mut_sums.get(f, 0.0) - ref_sums.get(f, 0.0)
    for f in all_features
}

top_lost = sorted(deltas.items(), key=lambda x: x[1])[:10]       # most negative
top_gained = sorted(deltas.items(), key=lambda x: -x[1])[:5]     # most positive

---

解釈表

16,384個のSAE特徴は、UniRef90での活性化パターンに基づき、約6種類の生物学的タイプに分類されています。特徴は、最も優勢なカテゴリに基づいて解釈してください。

カテゴリ	生物学的な意味	検出されることが多い変異の種類
触媒機能	酵素の活性部位またはその近傍の残基で活性化します	酵素活性を損なう変異（キナーゼ、加水分解酵素、転移酵素）
リガンド結合部位	低分子／イオン／ヌクレオチド系リガンドと接触する残基で活性化します	薬物結合、ATP/GTP結合、金属配位、DNA/RNA結合を損なう変異
翻訳後修飾（PTM）	リン酸化部位、糖鎖付加シーケオン、ユビキチン化部位、アセチル化部位で活性化します	リン酸化による制御、N型糖鎖付加、ユビキチンを介した分解を損なう変異
ドメイン／モチーフ	典型的な構造ドメイン（Zn finger、leucine zipper、EF-handなど）で活性化します	ドメイン内の三次構造を損なう変異
構造安定性	局所構造の維持に重要な残基で活性化します	タンパク質を不安定化させる変異 → フォールディングの機能喪失
二次構造／表面	αヘリックス、βストランド、溶媒露出残基で活性化します	特異性が低く、機序を示すシグナルは弱い

変異によってSAE特徴が失われた場合、その生物学的機能が損なわれた可能性があります。
SAE特徴が獲得された場合、変異によって非天然のシグナル（多くの場合、無秩序構造／別のフォールド）が生じた可能性があります。

ESM_describe_sae_feature(feature_id=...)を使うと、特定のfeature_idの生物学的カテゴリを調べられます。ある特徴に対する最初の呼び出しは時間がかかります（約30秒、ツールが10種類のタンパク質からなるラベリングパネルでSAEを実行し、UniProtのアノテーションを確認するため、Forgeクレジット約10）。同じ特徴への2回目以降の呼び出しではローカルキャッシュが使われ、即座に結果が返ります。

---

正直な制約事項

1. SAE特徴 ≠ 正解となる機能。SAE特徴は配列パターンから学習されたものです。「Feature 16076はN型糖鎖付加の検出器のように見える」というのは、UniRef90全体でのその特徴の振る舞いに基づく仮説であり、残基ごとの機能アノテーションではありません。証拠として使用し、証明として扱わないでください。
2. 対象は±8残基のウィンドウのみ。このウィンドウは局所的な影響（活性部位の破壊、モチーフの損傷）を捉えます。変異箇所から離れた特徴を変化させる長距離のアロステリック効果、二量体化界面への影響、ドメイン間の再配置は捉えません。
3. 特徴ラベルのカバー範囲は限定的。16,384個の特徴のうち、十分に分類されていないものもあります。活性化が広範囲に及ぶものや、パターンが不明瞭なものがあります。未分類の特徴でも、SAEの差分から変化があったことは分かりますが、生物学的な解釈の確度は低くなります。
4. 単一アイソフォームのみ。UniProtの標準配列を使用します。変異が非標準アイソフォームにある場合、結果が当てはまらない可能性があります。
5. SAEは6bのみ。現在、TUがサポートしているのはesmc-6b-2024-12_k64_codebook16384_layer60 SAEです。より小さい600m SAEも存在しますが、ソースリポジトリのLoFベンチマークとの照合による検証は済んでいません。
6. Forge APIは非商用。EvolutionaryScaleから別途ライセンスを取得しない限り、このスキルの出力を商用目的で使用することはできません。

---

クロス検証パターン（推奨）

SAE に基づく単独の解析だけでは、誤解を招く場合があります。重要な判断を行う際は、構造的・集団的な根拠とクロス検証してください。

レイヤー	ツール	確認できること
集団頻度	gnomad_get_variant / gnomad_search_variants	この変異は希少ですか？（LoF 候補は通常、希少または存在しません）
病原性予測	AlphaMissense_get_variant_score	独立した ML 予測ツールも、これが有害だと判定していますか？
臨床的根拠	ClinVar_search_variants + ClinVar_get_variant_details	この変異はすでに臨床的意義とともにキュレーションされていますか？
構造的コンテキスト	alphafold_get_prediction（残基ごとの pLDDT を返します）	変異した残基は、適切に折りたたまれた領域（pLDDT > 70）にありますか？
機能アノテーション	UniProt_get_disease_variants_by_accession、UniProt_get_function_by_accession	変異した残基は、既知の触媒部位、結合部位、PTM 部位にありますか？

3 つ以上のレイヤーで変異が有害だという結果が一致した場合、SAE 特徴量解析は機序の説明レイヤーになります。「AlphaMissense は有害だと判定し、gnomAD は集団内で見られないことを確認し、SAE は触媒特徴の喪失を示しています → これは触媒機能の LoF 変異です。」

---

レポート形式

ワークフローの完了後、次の形式で要約してください。

Variant: {VARIANT_ID}  (e.g. P04637_R175H = TP53 R175H)
Protein: {gene} ({uniprot_accession})

Top features LOST at position {position} ± {window}:
  - Feature {id}: activation Δ = {delta:.2f}  (category: {category})
    {one-line biological summary}
  ...

Top features GAINED at position {position} ± {window}:
  - Feature {id}: Δ = +{delta:.2f}  (category: {category if known else "uncategorized"})
  ...

Mechanistic interpretation: {2-3 sentences synthesizing the SAE evidence
with cross-validation layers}

Confidence: {high|medium|low}, based on:
  - SAE evidence: {how clear is the dominant category?}
  - Cross-validation: {how many layers agree?}
  - Limitations encountered: {any of the 6 caveats above that apply}