Update shotter math documentation with Logit Prior algorithm
This commit is contained in:
zly
@@ -682,3 +682,225 @@ flowchart TB
|
||||
- 正文公式保留 score 这一条即可。
|
||||
|
||||
_本文是对 Bttoxin_Shotter v1 实现逻辑的文字化说明，如需与最新代码完全严格对齐，请以代码仓库中的实现为准。_
|
||||
|
||||
---
|
||||
|
||||
## 11. 加入 3 项“基因组背景特征”到菌株×靶标分数（ZWA/Thu/TAA BGC / Mobilome / CRISPR）
|
||||
|
||||
> 你新增的三类指标，本质上都是“**对同一个菌株的全基因组背景做先验修正（strain-level priors）**”。
|
||||
>
|
||||
> - Shotter v1 的 $\text{Score}(\text{strain},\text{order})$ 仅由“毒素命中→靶标分布→noisy-OR”构成（第 3–6 节）。
|
||||
> - 现在我们引入三类**菌株层特征**，用来调整（boost / penalize）每个靶标分数，同时保持：
|
||||
> - 分数仍在 $[0,1]$；
|
||||
> - 解释仍然清晰：毒素命中给“证据”，背景特征给“先验”。
|
||||
>
|
||||
> 三类特征都来自**序列比对/检测的结果**：
|
||||
> - ZWA/Thu/TAA：对已知 BGC/关键酶蛋白做 BLAST/HMM（或用 BGC 预测器）→ 得到“存在/不存在”。
|
||||
> - Mobilome：对转座酶、质粒复制蛋白、噬菌体标志蛋白等做 BLAST/HMM/注释统计 → 得到数量。
|
||||
> - CRISPR/Cas：对 Cas 蛋白做 HMM/BLAST + 阵列检测（repeat/spacer）→ 判定“完整/不完整/不存在”。
|
||||
|
||||
### 11.1 为什么这三项可以影响“杀虫活性潜力”？（三条关键因果链）
|
||||
|
||||
> 你要向专家说明的不是“它们直接提高毒力”，而是：它们改变了菌株在进化与生态上**获得/维持杀虫模块**的概率。
|
||||
|
||||
#### (A) ZWA/Thu/TAA 三类“杀虫素 BGC”：存在即为直接功能模块证据（正向）
|
||||
|
||||
关键因果链（可直接放到讨论里）：
|
||||
|
||||
1. 这些基因簇（BGC）代表一种可表达的生物合成路径（酶系/修饰/转运/自抗等模块齐备）。
|
||||
2. 若用序列比对/结构域检测确认该路径的关键酶与骨架基因存在，则“产物存在的可行性”显著上升。
|
||||
3. 因此，即使 Cry/Vip 等蛋白毒素命中较少，BGC 的存在也代表一种**独立（正交）的杀虫潜力来源**。
|
||||
|
||||
结论：**存在（1）应当加分；不存在（0）不加分。**
|
||||
|
||||
#### (B) 移动元件（转座酶/质粒/噬菌体）数量：mobilome 越丰富，越容易获得/重排毒素模块（总体正向，但要饱和）
|
||||
|
||||
关键因果链：
|
||||
|
||||
1. Bt 及 B. cereus group 中，大量杀虫相关基因（包括 cry 等）常与质粒/移动 DNA 库强相关，且移动元件（IS/转座子等）参与基因重排、模块拼装与在质粒上的迁移。
|
||||
2. mobilome 越丰富，意味着可重排/可迁移的 DNA 元件越多，越容易出现：
|
||||
- 新毒素模块的获得（HGT/共轭/转导等）；
|
||||
- 现有毒素模块的复制、重排、组合与剂量效应；
|
||||
- 数据库未覆盖的新型“other/unknown”杀虫因子的潜力。
|
||||
3. 但 mobilome 指标也会受组装质量、注释阈值影响，因此应做“边际递减/饱和”，避免被噪声拉爆。
|
||||
|
||||
结论：**数量越多总体加分（正向），但用饱和函数。**
|
||||
|
||||
#### (C) CRISPR/Cas 完整度：越完整，越像“限制外源 DNA 的屏障”，从而降低获得毒素/质粒库的先验（负向）
|
||||
|
||||
关键因果链：
|
||||
|
||||
1. CRISPR-Cas 的核心生态功能之一是抵御外源遗传元件（质粒/噬菌体等），本质上会对水平基因转移（HGT）形成选择压力。
|
||||
2. 在 Bt 中，杀虫谱与毒素库的快速扩展常与质粒/移动 DNA 库相关；
|
||||
3. 因此，当 CRISPR/Cas **更完整且更可能功能健全**时，菌株对外源质粒/移动模块的“进入与稳定维持”通常更困难，导致“获取/更新毒素库”的先验下降。
|
||||
4. 反过来，CRISPR 缺失/失活的菌株更可能处于“更开放的 mobilome 交换状态”，更容易累积可迁移毒素模块。
|
||||
|
||||
结论：对“杀虫潜力总评分”应采用 **不存在 > 不完整 > 完整** 的单调顺序（即 CRISPR 越完整越压分）。
|
||||
|
||||
> 注：CRISPR 的实际生态效应可能随环境与 anti-CRISPR 等因素复杂化，但作为“全基因组先验”，上述方向能提供最稳定、最可解释的单调修正。
|
||||
|
||||
### 11.2 三类指标的输入格式（对应你的量化规则）
|
||||
|
||||
你希望的量化方式：
|
||||
|
||||
1. 三种杀虫素（分别为 ZWA、Thu、TAA）生物合成基因簇：存在（1）、不存在（0）。
|
||||
2. 移动元件（转座酶、质粒、噬菌体）：用数量表示。
|
||||
3. CRISPR/Cas 系统：完整（2）、不完整（1）、不存在（0）。
|
||||
|
||||
记作：
|
||||
|
||||
- $b_Z, b_T, b_A \in \{0,1\}$（对应 ZWA/Thu/TAA）
|
||||
- $m \in \mathbb{N}_0$（mobilome 总计数，或分项加和）
|
||||
- $c \in \{0,1,2\}$（CRISPR 状态：0=不存在，1=不完整，2=完整）
|
||||
|
||||
### 11.3 为什么不用把这三类特征直接塞进 noisy-OR？（建模选择）
|
||||
|
||||
noisy-OR 适合合成“多条毒素命中对同一靶标的独立贡献”。
|
||||
而这三类新特征是“菌株整体背景”，并不对应某个具体命中的 $c_i(\text{order})$。
|
||||
|
||||
因此最稳妥的做法是：
|
||||
|
||||
1. 先按 v1 算出 $S_\text{tox}(\text{strain},\text{order})$（即第 6 节的 Score）；
|
||||
2. 再用一个“先验修正项”把它调整为 $S_\text{final}$。
|
||||
|
||||
### 11.4 推荐的并入方式：logit 先验修正（保持 $[0,1]$ 且可解释）
|
||||
|
||||
令：
|
||||
|
||||
$$
|
||||
S_{\text{tox}}(\text{strain},o) = 1 - \prod_{i\in H_{\text{strain}}}\left[1-c_i(o)\right]
|
||||
$$
|
||||
|
||||
我们定义最终分数：
|
||||
|
||||
$$
|
||||
S_{\text{final}}(\text{strain},o)=\sigma\Big(\operatorname{logit}(S_{\text{tox}}(\text{strain},o)+\varepsilon) + \Delta(\text{strain})\Big)
|
||||
$$
|
||||
|
||||
其中：
|
||||
|
||||
- $\sigma(x)=\dfrac{1}{1+e^{-x}}$ 是 sigmoid；
|
||||
- $\operatorname{logit}(p)=\ln\dfrac{p}{1-p}$；
|
||||
- $\varepsilon$ 是极小正数（如 $10^{-6}$），避免 $\logit(0)$ 或 $\logit(1)$。
|
||||
|
||||
把三项特征写进 $\Delta(\text{strain})$：
|
||||
|
||||
$$
|
||||
\Delta(\text{strain}) =
|
||||
\beta_Z\,b_Z + \beta_T\,b_T + \beta_A\,b_A
|
||||
+ \beta_M\,g(m)
|
||||
+ \beta_C\,h(c)
|
||||
$$
|
||||
|
||||
其中 $g(m)$ 与 $h(c)$ 用来保证“饱和”和“单调方向”。
|
||||
|
||||
#### 11.4.1 mobilome 饱和函数 $g(m)$（正向，边际递减）
|
||||
|
||||
推荐二选一：
|
||||
|
||||
- 对数饱和：
|
||||
|
||||
$$
|
||||
g(m)=\ln(1+m)
|
||||
$$
|
||||
|
||||
- 线性截断：
|
||||
|
||||
$$
|
||||
g(m)=\min\left(1,\frac{m}{K}\right)
|
||||
$$
|
||||
|
||||
$K$ 取一个经验阈值（例如 50 或 100），取决于你统计的元件类型与注释粒度。
|
||||
|
||||
#### 11.4.2 CRISPR 单调映射 $h(c)$（“完整压分”）
|
||||
|
||||
为了实现“不存在 > 不完整 > 完整”的顺序，直接把 $c\in\{0,1,2\}$ 映射到 $\{1,0.5,0\}$：
|
||||
|
||||
$$
|
||||
h(c)=1-\frac{c}{2}
|
||||
$$
|
||||
|
||||
解释：
|
||||
|
||||
- 不存在 $c=0\Rightarrow h=1$（最大加成）
|
||||
- 不完整 $c=1\Rightarrow h=0.5$（中等加成）
|
||||
- 完整 $c=2\Rightarrow h=0$（不加成，相当于“最压分”）
|
||||
|
||||
只要取 $\beta_C>0$，就得到你想要的确定性方向：
|
||||
|
||||
> CRISPR 越完整 → $h(c)$ 越小 → $\Delta$ 越小 → $S_{\text{final}}$ 越低
|
||||
|
||||
### 11.5 这三项对评分是“增加还是降低”？（给你一个可写进算法的确定性结论）
|
||||
|
||||
在上述形式下，只要你取参数满足：
|
||||
|
||||
- $\beta_Z,\beta_T,\beta_A > 0$
|
||||
- $\beta_M > 0$
|
||||
- $\beta_C > 0$
|
||||
|
||||
则得到确定性结论：
|
||||
|
||||
- **ZWA/Thu/TAA：存在（1）一定增加分数；不存在（0）不增加。**
|
||||
- **移动元件数量：数量越多一定增加分数（经 $g(m)$ 饱和）。**
|
||||
- **CRISPR：不存在（0）加成最大；不完整（1）次之；完整（2）最小（最压分）。**
|
||||
|
||||
### 11.6 将“三项特征”并入 Shotter 的输出矩阵（菌株×靶标）
|
||||
|
||||
你最终仍输出一个矩阵：
|
||||
|
||||
- `strain_target_scores.tsv`（或新文件名）里，每个单元格从 $S_{\text{tox}}$ 替换为 $S_{\text{final}}$。
|
||||
|
||||
也就是说，Shotter 原本的 per-hit 解释仍然保留（因为 $S_{\text{tox}}$ 没变）；
|
||||
三项新特征只是在菌株层做一个整体“上移/下移”的可解释修正。
|
||||
|
||||
### 11.7 Mermaid：在原流程上加入“背景特征先验”分支（可选）
|
||||
|
||||
> 如果你想把这三项清晰地画进流程图，可以追加一个“Genome context”分支：
|
||||
|
||||
```mermaid
|
||||
flowchart TB
|
||||
G["Genome .fna"] --> ORF[ORF/CDS prediction]
|
||||
ORF --> TR[Translate CDS to proteins]
|
||||
|
||||
TR --> ID["Toxin identification\n(BLAST/HMM)"]
|
||||
ID --> ALL["All_Toxins.txt"]
|
||||
|
||||
BPPRC["BPPRC specificity CSV\n(toxicity-data.csv)"] --> IDX["Build specificity index\nP(order|·), P(species|·)"]
|
||||
|
||||
ALL --> PARSE["Per-hit parsing"]
|
||||
IDX --> WHIT["Compute w_i"]
|
||||
WHIT --> PHIT["Per-hit scoring\nc_i(order)"]
|
||||
PHIT --> COMB["Per-strain noisy-OR\nS_tox(strain,order)"]
|
||||
|
||||
subgraph CONTEXT["Genome context features (BLAST/HMM/detectors)"]
|
||||
BGC["ZWA/Thu/TAA BGC\n(0/1)"]
|
||||
MGE["Mobilome count m"]
|
||||
CRISPR["CRISPR state c\n(0/1/2)"]
|
||||
end
|
||||
|
||||
TR --> BGC
|
||||
TR --> MGE
|
||||
TR --> CRISPR
|
||||
|
||||
COMB --> ADJ["Logit prior adjust\nS_final = sigmoid(logit(S_tox)+Δ)"]
|
||||
CONTEXT --> ADJ
|
||||
|
||||
ADJ --> OUT2["strain_target_scores.tsv\n(strain × order, updated)"]
|
||||
```
|
||||
|
||||
---
|
||||
|
||||
## 12. 参考文献（用于支撑第 11 节的因果链，建议在论文/汇报中引用）
|
||||
|
||||
> 下面列出的是“因果链”的常用支撑点：
|
||||
> - Bt 毒素与质粒/移动 DNA 的关联；
|
||||
> - CRISPR-Cas 对质粒/HGT 的屏障作用；
|
||||
> - Bt 的移动基因库与致病/宿主谱的关系。
|
||||
|
||||
1. Fiedoruk K. et al. *Genetic Environment of cry1 Genes Indicates Their Common Origin.* (2017).（cry 基因位于 Bt 质粒并与特定复制系统/遗传环境相关）
|
||||
2. Lechuga A. et al. *Completed Genomic Sequence of Bacillus thuringiensis…* (2020).（Bt 的昆虫致病能力与质粒携带 cry 基因相关的概述）
|
||||
3. Thomas D.J.I. et al. *Transfer of plasmid pBC16 between Bacillus thuringiensis…* (2002).（Bt 大质粒可经类似共轭机制传播，获得 cry 质粒可扩展生态位/宿主）
|
||||
4. Gillis A. et al. *Role of plasmid plasticity and mobile genetic elements in…* (2018).（综述：Bt 移动基因库与致病/生活史相关）
|
||||
5. Marraffini L.A. & Sontheimer E.J. *CRISPR Interference Limits Horizontal Gene Transfer in Staphylococci.* Science (2008).（经典实验证据：CRISPR 干预能限制共轭/转化等 HGT）
|
||||
6. Wheatley R.M. & MacLean R.C. *CRISPR-Cas systems restrict horizontal gene transfer…* (2020).（综述：CRISPR 对质粒/ICE 等 HGT 的普遍限制作用）
|
||||
|
||||
|
||||
Reference in New Issue
Block a user