tooluniverse-regulatory-variant-analysis

説明ではなく計算を実行する

解析に計算（統計、データ処理、スコアリング、情報の補完）が必要な場合は、Bash 経由で Python コードを記述して実行します。実行する内容を説明するのではなく、実際に実行して結果を報告します。ToolUniverse ツールでデータを取得し、Python（pandas、scipy、statsmodels、matplotlib）で分析します。

Regulatory Variant Analysis スキル

体系的な制御バリアント解釈：GWAS から形質との関連を特定し、eQTL の効果をマッピングし、クロマチンの状況を注釈し、制御エレメントとの重なりを評価し、AlphaGenome／AlphaGenome Atlas で配列レベルの直接的な効果を予測して、非コードバリアントの機能的影響を証拠の強さに応じて評価します。

使用する場面

• 「rs12913832 にはどのような GWAS 関連がありますか？」
• 「脳組織における APOE 遺伝子座の eQTL を見つけてください」
• 「このバリアント領域と重なる制御エレメントは何ですか？」
• 「GWAS で 2 型糖尿病と関連する SNP を教えてください」
• 「このイントロンバリアントは活性エンハンサー内にありますか？」
• 「rs429358 の RegulomeDB スコアはいくつですか？」
• 「BRCA1 プロモーター領域の ENCODE ヒストンマークを見つけてください」
• 「『血圧』の形質オントロジー用語を EFO ID にマッピングしてください」
• 「このバリアントには GWAS／eQTL／ENCODE のヒットがまったくありません。配列モデルは何らかの影響を予測しますか？」
• 「この制御領域の候補をスキャンして、どの 1 塩基置換が最も重要か特定してください」

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

• コーディングバリアントの病原性評価 → tooluniverse-variant-interpretation を使用
• 臨床バリアントの包括的な分類（ACMG） → tooluniverse-variant-interpretation を使用
• 遺伝子と疾患の関連（バリアント固有ではないもの） → tooluniverse-gene-disease-association を使用
• 薬理ゲノミクスのバリアント注釈 → tooluniverse-pharmacogenomics を使用
• エピゲノミクスデータの処理（BED／narrowPeak ファイル） → tooluniverse-epigenomics を使用

---

非コード変異の影響を推論する

非コード変異を評価する際は、5つの問いに沿って根拠を集めます。最初の4つはアノテーションに基づく問い（この集団において、この位置についてすでに何がわかっているか）で、5つ目はモデルに基づく問い（その配列自体が何を予測するか。そこにまだアノテーションがなくても予測します）です。それぞれ異なる問いに答えるため、互いの代わりにするのではなく、相互に裏付けとなるように使います。

1. 変異は調節エレメント内にありますか？ RegulomeDBを使って、変異が転写因子結合部位、クロマチンアクセシビリティピーク、または既知の調節アノテーションと重なるかを評価します。RegulomeDBスコアが低い（カテゴリー1a～2a）場合、その位置が機能的に活性である強い根拠となります。ENCODEのヒストンマークで確認します。H3K27acは活性エンハンサーと活性プロモーターを示し、H3K4me1のみの場合は準備状態のエンハンサーを示します。H3K4me3は活性プロモーターを示し、H3K27me3はサイレンシングされた領域を示します。

2. 転写因子結合部位を変化させますか？ RegulomeDBの転写因子結合に関する根拠と、ENCODEの転写因子ChIP-seq実験を確認します。変異が転写因子のフットプリント内にあり、コンセンサスモチーフを破壊する場合、特にその転写因子が疾患組織で重要だとわかっているときは、機序に基づく介入可能な知見となります。

3. 遺伝子との関連を示すeQTLの根拠はありますか？ GTExを検索し、その変異（またはLDが強い変異）が近傍遺伝子の発現を、組織特異的または広範に調節しているかを判断します。組織特異的なeQTLは細胞型特異的な調節を示唆し、広範に見られるeQTLは中核的な調節エレメントを示唆します。解釈には、NESの方向（正＝代替アレルが発現を増加させる、負＝減少させる）と効果量が重要です。

4. 形質との関連を示すGWASの根拠はありますか？ GWAS CatalogでrsIDまたは周辺の遺伝子座を検索します。関連する形質でゲノムワイド有意な関連（p < 5×10⁻⁸）があれば、その変異の生物学的重要性を裏付ける根拠となります。複数のGWAS研究による遺伝子座から遺伝子へのマッピングについては、OpenTargetsと照合します。

5. 配列ベースのモデルは何を直接予測しますか？ AlphaGenome_atlas_lookup_variantを実行すると、ほぼ即座にAVI_SCOREを取得できます（AlphaGenomeとAlphaMissenseを統合した影響推定値です。事前計算済みのため、以下のアノテーション工程を始める前でも、候補となるすべてのSNVについて最初に低コストで確認できます）。他の根拠がないのにAVI_SCOREが高い場合は、まさに「明らかな役割がわからない変異」のケースです。これは、アノテーションデータベースの更新が追いついていないことを意味し、変異が不活性であることを意味するわけではありません。AVI_SCOREが高い場合、または位置がAtlasでは扱えないindel／合成配列の場合は、ライブモデルに切り替えます。トラックごとの影響の内訳にはAlphaGenome_score_variantを、周辺ウィンドウをスキャンして、どの塩基変化がどの生物学的指標（発現、スプライシング、クロマチン）に影響するかを特定するにはAlphaGenome_score_ism_variantsを使います。5つの問いのうち、相関的な回答ではなく、「この置換によって、それまで存在しなかった転写因子結合モチーフが作られる」のような機序に基づく回答が得られるのは、これだけです。

根拠の統合：複数の層からなる根拠を組み立てます。GWASの有意性、eQTLの根拠、RegulomeDBスコア1a～2a、関連組織における活性クロマチン（H3K27ac）がそろった変異は、調節への影響が高い信頼度で示されます。さらに、AVI_SCOREまたはISMの結果が一致すれば、相関的なアノテーションを追加したものではなく、独立した因果メカニズムの根拠であるため、さらに確かな裏付けとなります。根拠が2～3種類収束している場合（例：eQTLと活性エンハンサー）は、中程度の信頼度です。根拠が1種類のみの場合、または活性ではなく準備状態の調節環境にある変異の場合は、信頼度は低くなります。ただし、「根拠なし」と判断する前に、以下のAlphaGenomeのみのケースを確認してください。

---

ワークフローの概要

Input (rsID, genomic coordinates, trait/disease, gene)
  |
  v
Phase 0: Variant/Trait Resolution
  Resolve rsIDs, map trait names to EFO/MONDO IDs via OLS
  |
  v
Phase 0.5: Sequence-Based Triage (AlphaGenome Atlas)
  Instant AVI_SCORE per candidate SNV -- cheap, run before the annotation phases
  |
  v
Phase 1: GWAS Association Lookup
  GWAS Catalog associations, p-values, effect sizes, study metadata
  |
  v
Phase 2: eQTL Analysis
  GTEx tissue-specific eQTLs, target gene identification
  |
  v
Phase 3: Regulatory Element Annotation
  ENCODE histone marks, RegulomeDB scores, chromatin state
  |
  v
Phase 4: OpenTargets GWAS Integration
  OpenTargets GWAS study aggregation, locus-to-gene mapping
  |
  v
Phase 4.5: Sequence-Based Deep-Dive (if Phases 1-4 are silent, or the Phase 0.5
  AVI_SCORE was high and needs a mechanism) -- live AlphaGenome scan
  |
  v
Phase 5: Functional Impact Synthesis
  Integrate all evidence, assign regulatory impact level
  |
  v
Phase 6: Report
  Evidence-graded regulatory variant report

---

フェーズ 0：バリアント／形質の解決

GWAS クエリの前に、ols_search_terms を使って形質名をオントロジー ID に解決します。GWAS 形質には ontology="efo" を指定してください。OpenTargets では MONDO ID が優先されます（例：2 型糖尿病には EFO_0001360 ではなく MONDO_0005148）。初期の影響注釈と最寄りの遺伝子の特定には EnsemblVEP_annotate_rsid を使います（パラメーターは rsid ではなく variant_id です）。

---

フェーズ 0.5：配列に基づくトリアージ（AlphaGenome Atlas）

ALPHA_GENOME_API_KEY が必要です。 設定されていない場合、AlphaGenome ツールはツールセットに一切表示されません。その場合はこのフェーズとフェーズ 4.5 を完全にスキップし、フェーズ 1〜4 に進んでください。これらだけで完全に自立したパイプラインを構成します（AlphaGenome のサポートが追加される前は、このスキルもそのように動作しており、以下の内容は実行されていることを前提としていません）。

特定のバリアントが 1 つ（染色体、位置、参照塩基／代替塩基）わかっている場合は、以下の注釈フェーズに時間をかける前に AlphaGenome_atlas_lookup_variant を実行します。これは、実質的に人間のほぼすべての SNV を網羅する事前計算データベースを参照するため、最後の手段ではなく、標準の最初のステップとして実行できるほど低コストです。返される AVI_SCORE は、AlphaGenome の調節予測と AlphaMissense のコーディング影響モデルを統合した単一の数値であり、コーディングバリアントと非コーディングバリアントのどちらにも同じように適用できます。

既知のバリアントが 1 つではなく候補のリストがある場合（たとえば、フェーズ 1 の gwas_get_variants_for_trait がある形質について数十個の SNP を返した場合や、ファインマッピングを行う場合）は、候補ごとに atlas_lookup_variant を繰り返すのではなく、対象領域に対して AlphaGenome_atlas_scan_interval を使います。1 回の呼び出しで、最大 10 kb のウィンドウ内にあるすべての SNV を AVI_SCORE で順位付けします。これは、高コストの追加調査を行う前に予測された分子影響で候補をグループ化する方法そのもので、大規模コホート解析では、候補をすべて同じ調査価値として扱うよりも、形質との関連を大幅に多く発見できました。

どちらの結果も、最終判断ではなく事前情報として扱ってください。注釈の根拠が弱いのに AVI_SCORE が高い場合は、「ここをさらに詳しく調べる必要があります。データベースがまだ追いついていません」という意味です（フェーズ 4.5 を参照）。GWAS／eQTL の根拠が強い一方で AVI_SCORE が低い場合は、因果バリアントが実際には LD 内の別のバリアントではないか、再検討する価値があります。Atlas が対象とするのは一塩基置換のみです。インデルの場合は、AlphaGenome_score_variant（ライブモデル）に直接進んでください。

---

フェーズ 1：GWAS 関連の検索

主要なツールは gwas_search_associations です。disease_trait（自由テキスト）、efo_id（精度を優先する場合に推奨）、rs_id、および p_value 閾値を受け付けます。ゲノム全体で有意な関連には p_value=5e-8 を使います。遺伝子座レベルで調べる場合、gwas_get_variants_for_trait を使うと、ある形質に関連するすべての SNP を取得できます。gwas_get_snps_for_gene を使うと、特定の遺伝子にマッピングされた GWAS カタログ掲載の SNP を検索できます。

推論のヒント：自由テキストの形質で検索して GWAS Catalog が結果を返さない場合は、efo_id パラメーターに切り替えてください。カタログは統制語彙を使用しているため、自由テキストによる照合の精度は高くありません。

---

フェーズ 2：eQTL 解析

GTEx_query_eqtl は遺伝子シンボル（GENCODE ID に自動解決）または Ensembl 遺伝子 ID を受け付けます。組織ごとの SNP と遺伝子の関連、および各組織の NES（正規化効果量）と p 値を返します。

結果を解釈する際は、eQTL の影響が疾患に最も関連する組織で生じているかを確認してください。神経変性疾患のバリアントの場合、脳に特異的な eQTL のほうが、広範な組織で見られるものより有力です。eQTL の根拠を重視する前に、GTEx_get_median_gene_expression を使い、標的遺伝子が関連する組織で実際に発現していることを確認してください。

注記：GTEx API は v8 のデータを使用します。クエリによっては、gtex_v10 エンドポイントから結果が返らないことがあります。

---

フェーズ 3：調節エレメントのアノテーション

RegulomeDB_query_variant（パラメーター：rsid）は、調節スコアと特徴アノテーションを返します。カテゴリー1a～2aのスコアは、強い調節エビデンス（eQTLとの重複＋TF結合＋クロマチンアクセシビリティ）を示します。スコア3a～6は、段階的に弱いエビデンスを示します。

ENCODE_search_histone_experiments は histone_mark（例：「H3K27ac」）と biosample_term_name（組織または細胞株の名前。疾患名ではありません。ENCODEでは「liver」や「breast epithelium」のような生物学的サンプル名を使用します）を受け取ります。TF結合データを検索する際は、assay_title="TF ChIP-seq"（単なる「ChIP-seq」ではありません）を使用してください。

推論のヒント：RegulomeDBは、ENCODE、Roadmap、その他のデータを集約しています。ENCODEに特定のバイオサンプルがなくても、RegulomeDBには関連する細胞種から集約されたエビデンスがある場合があります。

バリアントがGRCh38座標でわかっている場合、FAVOR_annotate_variant(variant="chr-pos-ref-alt") は調節アノテーションブロック（保存性、頻度、CADDも含む）を1回の呼び出しで返します。RegulomeDB／ENCODEを詳しく調べる前に、その位置がアノテーション済みの調節エレメントに含まれるかどうかをすばやく確認するために使用してください。

---

フェーズ 4：OpenTargets GWASの統合

OpenTargets_search_gwas_studies_by_disease は、MONDO IDの配列を diseaseIds として受け取ります。複数のGWAS研究から遺伝子座と遺伝子の関連（L2G）スコアを提供します。このスコアは、単純な近接性にとどまらず、共局在、eQTL、クロマチンデータを取り入れています。まず OpenTargets_multi_entity_search_by_query_string または OpenTargets_get_disease_id_description_by_name を使用して、疾患名をMONDO／EFO IDに解決してください。

---

フェーズ 4.5：配列ベースの詳細分析（条件付き）

次のいずれかの条件に当てはまる場合に、このフェーズを実行します：（a）フェーズ1～4で結果が得られない、または弱い場合（典型的な「明確な役割がない」バリアント）、または（b）フェーズ0.5のAVI_SCOREが高く、レポートに数値だけでなく実際のメカニズムが必要な場合。

• AlphaGenome_score_variant：重要なモダリティ（発現、スプライシング、アクセシビリティ、特定のヒストンマーク）がすでにわかっている場合に、バリアントのトラックごとの参照配列と代替配列の影響を評価する推奨ツールです。
• AlphaGenome_score_ism_variants：バリアント周辺の短いウィンドウ（≤500 bp）で考えられるすべての置換をスキャンし、影響の大きさで順位付けします。どの塩基が重要なのかはまだわからず、領域内の何かが影響していることだけがわかっている場合に使用します。たとえば、隣接する位置で同じ破壊パターンが見られるかどうかを調べ、バリアントがモチーフ内にあることを確認できます。
• AlphaGenome_predict_interval：単一バリアントによる差分だけでなく、より広い調節環境が必要な場合に、より広範なウィンドウ（最大1 Mb）全体のトラックプロファイルを予測します。

2つまたは3つのトラックで一致する結果（例：エンハンサー活性を示すヒストンマークと発現の両方が同じ方向に変化する場合）は、単一トラックでのヒットより強いエビデンスです。これは、フェーズ5ですでにアノテーションのエビデンスに適用している考え方と同じです。

AlphaGenomeはToolUniverseにある唯一の配列予測モデルではありませんが、最も幅広く、このフェーズでデフォルトとして使用されるモデルです。利用できない場合（キーがない場合）、または特定のケースで別のモデルが適している場合（ホスト型APIではなく公開済みでセルフホスト可能なモデルを使う場合や、幅広いマルチモダリティの結果ではなくクロマチンアクセシビリティを特に調べる場合）は、tooluniverse-regulatory-genomics を参照してください。Enformer、Borzoi、ChromBPNet、Evo 2が紹介されており、異なるトレードオフで同種の分析を行えます。

---

フェーズ 5：機能的影響の統合

エビデンスを収集したら、各レイヤーを踏まえて判断します。

• 影響大：GWAS でゲノムワイド有意であり、意味のある NES を持つ eQTL があり、RegulomeDB スコアが ≤ 2 で、関連組織に活性クロマチン（H3K27ac）が存在します。複数の独立したエビデンスが同じ座位と遺伝子に収束しています。一致する AlphaGenome/AVI_SCORE の結果があれば、これに機序に関するエビデンスが加わりますが、アノテーションのエビデンスがすでに収束している場合、信頼度を高く判断するうえで必須ではありません。
• 中程度の影響：2～3 種類のエビデンス（例：eQTL と活性エンハンサーの重複、または GWAS で有意かつ RegulomeDB ≤ 3）があるものの、完全には収束していません。
• 影響小：エビデンスが 1 種類のみ、または機能データを伴わない計算によるアノテーション（VEP の結果カテゴリ）のみです。
• 配列モデルで支持されるが、アノテーションでは検出されない：GWAS/eQTL/ENCODE/RegulomeDB のエビデンスはありませんが、フェーズ 4.5 で高い AVI_SCORE や、明確な複数トラックの機序（例：転写因子モチーフの新規形成、スプライスジャンクションの破壊）が示されています。これは「エビデンスなし」ではなく、独立した段階として扱います。既存のデータベースが追いつけないほど新規性が高い、または希少なため、変異に機能がある可能性は十分にあります。何も起きていないという意味ではありません。信頼度は正直に示してください。これは、独立した十分に検証されたモデルによる有力な仮説であり、上記の段階に見られるような複数の情報源の収束ではありません。
• エビデンスなし：どの情報源にも調節アノテーションがなく、フェーズ 4.5（実行した場合）でも AVI_SCORE が低い、または強い配列効果が示されていません。変異は機能のない領域にある可能性があります。

---

フォールバック戦略

• GWAS Catalog が空の結果を返す場合：フリーテキストの disease_trait から efo_id に切り替えます。形質名をより広い語にします。
• 遺伝子の GTEx eQTL が空の場合：遺伝子記号のスペルを確認します。Ensembl ID を試し、size パラメーターを増やします。
• RegulomeDB がデータを返さない場合：ENCODE を直接検索するか、調節および保存性のアノテーションを得るために（GRCh38 座標で）FAVOR_annotate_variant を実行します。利用可能なデータにその変異の調節アノテーションがない可能性があります。
• OpenTargets GWAS が None を返す場合：MONDO/EFO ID の形式を確認します。正しい ID を確認するため、まず OpenTargets_multi_entity_search_by_query_string を試します。
• ENCODE で組織が見つからない場合：ENCODE では特定のバイオサンプル名が使われています。RegulomeDB は多くの細胞種のデータを集約しているため、不足を補える可能性があります。
• フェーズ 1～4 がすべて空または弱い結果の場合：「エビデンスなし」と結論付けるのはまだ早いです。フェーズ 4.5 を実行してください。AlphaGenome は配列から直接予測するため、その変異（または近傍の領域）がこれまで研究されている必要はありません。データベースに何もない場合でも機能する唯一のエビデンス源です。
• 変異が一塩基置換ではなくインデルの場合：AlphaGenome_atlas_lookup_variant/AlphaGenome_atlas_scan_interval は SNV のみを対象とします。代わりにライブの AlphaGenome_score_variant を使用してください。
• 使用できるツールセットに AlphaGenome ツールがない場合：ALPHA_GENOME_API_KEY が設定されていません。フェーズ 0.5 と 4.5 はスキップし、フェーズ 1～4 を独立したパイプラインとして実行してください。提供されない可能性のあるキーを待って解析を止めてはいけません。このケースで配列モデルによる予測が引き続き重要であれば、tooluniverse-regulatory-genomics には、この特定のキーを必要としないセルフホスト可能な代替手段（Enformer、Borzoi、ChromBPNet）があります。

---

ワークフロー例

GWAS 変異の機能アノテーション（rs429358 / APOE）

Step 1: gwas_search_associations(rs_id="rs429358")
  -> All trait associations (Alzheimer's disease, LDL cholesterol, etc.)

Step 2: GTEx_query_eqtl(gene_symbol="APOE")
  -> Tissue-specific eQTL evidence; note effect in brain vs liver

Step 3: RegulomeDB_query_variant(rsid="rs429358")
  -> Regulatory score and TF binding annotations

Step 4: ENCODE_search_histone_experiments(histone_mark="H3K27ac", biosample_term_name="brain")
  -> Active enhancer context near the variant

Step 5: Synthesize: does GWAS significance + eQTL + active chromatin converge on one gene?

非コード変異の評価（イントロン／UTR 変異）

Step 1: EnsemblVEP_annotate_rsid(variant_id="rs12345678")
  -> Confirm non-coding consequence, identify nearest gene

Step 2: RegulomeDB_query_variant(rsid="rs12345678")
  -> Is this position in a regulatory context?

Step 3: gwas_search_associations(rs_id="rs12345678")
  -> Any GWAS associations in relevant traits?

Step 4: GTEx_query_eqtl(gene_symbol=nearest_gene)
  -> Does this variant or nearby variants modulate expression?

Step 5: ENCODE_search_histone_experiments(histone_mark="H3K27ac", biosample_term_name=relevant_tissue)
  -> Active chromatin confirmation

Step 6: Classify impact based on convergence of evidence lines

配列予測によるアノテーション未検出変異の解明

GREGoR Consortium が AlphaGenome Atlas を用いて、深部イントロンの DNM1 変異による未解決のてんかん症例を解決した実例などに見られるパターンです。エクソンにヒットせず、GWAS の関連もなく、アノテーションの各フェーズでも何も見つかりません。配列モデルに直接尋ねるまでは、何もわからないのです。

Step 1: AlphaGenome_atlas_lookup_variant(chromosome=..., position=..., reference_bases=..., alternate_bases=...)
  -> High AVI_SCORE despite the variant being deep intronic

Step 2: gwas_search_associations, GTEx_query_eqtl, RegulomeDB_query_variant, ENCODE_search_histone_experiments
  -> All empty or weak -- nothing here for standard annotation to catch

Step 3: AlphaGenome_score_ism_variants(chromosome=..., start=..., end=..., output_type="SPLICE_SITES")
  -> Scan the surrounding window; identify that this substitution creates a cryptic splice site

Step 4: Synthesize: high AVI_SCORE + a concrete, mechanistic splice-site prediction, despite zero
  annotation-database hits -> classify as "sequence-model-supported, annotation-silent", not
  "no evidence". Recommend experimental validation (e.g., a minigene splicing assay) rather than
  treating the model output as a final answer.

GWAS 遺伝子座における多数の候補の優先順位付け

GWAS のヒットが単一のバリアントだけに関与することはほとんどありません。LD により、通常、1 つの遺伝子座には数十個の候補が存在します。すべての候補に完全なアノテーションパイプライン（フェーズ 1〜4）を実行するとコストがかかるため、まず候補を絞り込みます。

Step 1: gwas_get_variants_for_trait(trait/efo_id=...)
  -> Dozens of candidate SNPs at the associated locus

Step 2: AlphaGenome_atlas_scan_interval(chromosome=..., start=..., end=...)
  -> One call, AVI_SCORE for every SNV across the locus (<=10 kb per call; tile larger loci)

Step 3: Rank candidates by AVI_SCORE; take the top handful (not just the GWAS lead SNP --
  a variant in tight LD with a stronger predicted effect is often the better causal candidate)

Step 4: Run Phases 1-4 (GWAS, eQTL, RegulomeDB/ENCODE, OpenTargets) only on those top candidates,
  not the full list -- the expensive annotation calls are now spent where they're likely to pay off

---

制限事項

• GWAS Catalog が対象とするのは公表済みの GWAS のみであり、未公表の研究は含まれません。
• GTEx の eQTL データは v8 のものです。v10 のエンドポイントでは空の結果が返される場合があります。
• RegulomeDB のアノテーションは、特定の細胞型で利用可能な ENCODE/Roadmap データに依存します。
• eQTL 解析で特定できるのは相関であり、因果関係ではありません。因果バリアントの特定にはファインマッピングが必要です。
• RegulomeDB のスコアはヒューリスティックです。スコアが 1a でも、機能的影響が保証されるわけではありません。
• GWAS の関連は集団レベルのものです。個々のバリアントの影響は遺伝的背景によって異なります。
• AlphaGenome/AVI_SCORE の予測はモデルによる推定であり、実験による測定値ではありません。アノテーションの根拠を裏付けたり、検証に値する仮説を生み出したりできますが、高いスコアだけでは機能的影響が確認されたことにはなりません。AlphaGenome_score_ism_variants の対象範囲は最大 500 bp、AlphaGenome_atlas_scan_interval は呼び出しごとに最大 10,000 bp です。どちらを呼び出す場合も、事前に対象領域を絞り込んでください。