tooluniverse-phylogenetics

系统发育与序列分析

四个会得出貌似可信但实际上错误数值的陷阱

以下每种情况都曾导致错误答案，尽管本文件后文已经提供了正确指导。回答前请先检查。

1. PhyKIT 会输出多列，而对于 saturation，两种定义并不一致——请说明你使用的是哪一列。 phykit saturation 会输出 saturation <TAB> |saturation-1|。其 --help 中明确写道："第一个值是饱和度，第二列是饱和度减 1 的绝对值。" 但一些已发表的分析（以及由此得出的部分参考答案）会将第二列报告为“饱和度”。两者相加始终为 1.0000，这是你可能选错列的信号——在真菌 scogs 数据集中，中位数分别为 0.39（第 1 列）和 0.61（第 2 列）。

   因此：除非问题或来源将饱和度定义为另一种方式，否则遵循 phykit，使用第 1 列，并在答案中说明你读取的是哪一列。不要默默选择看起来更接近期望数值的那一列。

   treeness_over_rcv 没有这种歧义：它输出 ratio <TAB> treeness <TAB> RCV，其中 ratio 位于第一列。

2. “Gap percentage” 指的是至少包含一个 gap 的比对列所占比例，而不是 gap 残基所占比例。在真菌 scogs 数据集中，按残基计算的定义最大值为 0.556，因此“>70% gaps”的筛选条件不会选出任何结果，问题看起来就无法回答；按列计算则有三个直系同源基因符合条件（最大值为 0.783）。

3. treeness_over_rcv 和 rcv 使用未修剪的 .faa.mafft，而 saturation 使用修剪后的 .clipkit。 RCV 衡量列间的变异，因此修剪会改变它：未修剪时中位数为 0.2683，修剪后为 0.3050；在 gap 比例 >70% 的基因中，未修剪时最大值为 0.2572，修剪后为 0.4174。

4. 切勿逐个文件循环运行 PhyKIT。 phykit_batch_analysis 支持并行处理，处理约 250 棵树只需约 35 秒；shell 循环则需要约 9 分钟，并会在中途耗尽轮次，导致完全没有答案。它还会为每个函数选择正确的列，因此完全消除了陷阱 1。

规则零——首先检查是否有预计算结果

遵循以下任何说明之前，请先扫描数据文件夹，查找：

• scogs_fungi.zip / scogs_animals.zip（BUSCO 单拷贝直系同源基因系统发育分析）→ 其中包含原始分析的预计算比对（*.faa.mafft.clipkit）和树（*.faa.mafft.clipkit.treefile）。直接使用 PhyKIT 对这些文件进行分析（参见下文“BUSCO scogs 问题”）。从 *.busco.zip 文件重新运行 BUSCO → MAFFT → IQ-TREE 需要 1–6 小时，而且由于随机种子/版本差异，结果会略有不同。
• *_executed.ipynb → 使用 tu run read_executed_notebook '{"data_folder":"<path>","search":"<keyword>"}' 读取，并将其中单元格的输出作为权威答案引用
• 预计算结果文件（名称类似 *results*、*tree*、*phykit*、*saturation*、*treeness* 的 CSV/TSV 文件）→ 直接读取并报告所需数值
• 规范分析脚本（analysis.R、run_*.py、find_*.R、*.Rmd）→ 按原样执行并读取输出

只有在以上内容均不存在时，才按照下文的重新分析流程操作。从原始数据重新运行会得到不同于已发表答案的数值，而且耗时更长（通常需要 5–10 倍的轮次）。

---

BUSCO scogs 问题（多物种基因组系统发育分析）

包含 scogs_fungi.zip 和/或 scogs_animals.zip 的数据文件夹附带每个直系同源基因的预计算比对（有时也包含树）。问题会要求提供各组的指标，或比较组间的 Mann-Whitney U 检验结果、中位数或比值。

主脚本 — 一次处理两组（优先使用此脚本）

如果问题比较动物与真菌（Mann-Whitney U、比率、倍数变化、配对差值），捆绑的配对比较脚本会解压两个压缩包，按每个直系同源组计算两组的指标，并在一次运行中输出以下所有内容：每组摘要、双尾 Mann-Whitney U 检验 + p 值（由于 U 不对称，会按两种顺序分别计算）、配对直系同源组中位数差值、配对直系同源组中位数比率、组中位数比率，以及最低非零比率——无需聚合步骤：

python skills/tooluniverse-phylogenetics/scripts/scogs_paired_compare.py \
    --data-folder "$DATA_PATH" --metric parsimony_informative
# Metrics: parsimony_informative, rcv, gap_percentage (alignment-only,
# Biopython-fast: ~2s for 500 alignments);
# treeness, dvmc, total_tree_length, evolutionary_rate, long_branch_score,
# patristic_distances (tree); treeness_over_rcv, saturation (both).

输出块（可用 Python 解析或用 grep 搜索）：

# SUMMARY group=animals: n=... mean=... median=... min=... max=... p25=... p75=... lowest_nonzero=... n_nonzero=...
# SUMMARY group=fungi:   n=... mean=... median=... min=... max=... p25=... p75=... lowest_nonzero=... n_nonzero=...
# MWU animals_vs_fungi: U=... p=...
# MWU fungi_vs_animals: U=... p=...        <-- U(a,b) + U(b,a) = n_a*n_b
# PAIRED n_common=N: median_diff(animals-fungi)=...  median_diff(fungi-animals)=...
# PAIRED RATIO median(animals/fungi)=... (n=...)    <-- for each common ortholog: a_val/b_val, then median
# PAIRED RATIO median(fungi/animals)=... (n=...)
# GROUP_MEDIAN_RATIO animals/fungi=...               <-- median(group_a) / median(group_b)
# GROUP_MEDIAN_RATIO fungi/animals=...
# GROUP_MEDIAN_DIFF animals-fungi=...
# LOWEST_NONZERO animals=... fungi=...
# LOWEST_NONZERO_RATIO animals/fungi=...
# LOWEST_NONZERO_RATIO fungi/animals=...

对于 long_branch_score 和 patristic_distances（每棵树有多个值的指标），传入 --per-tree-stat mean 或 --per-tree-stat median，以选择跨树 MWU 检验之前的每棵树摘要统计量。问题措辞为“比较 中位数 长枝得分”时，表示每棵树的摘要统计量应为中位数；措辞为“比较 平均值 长枝得分”时，表示每棵树的摘要统计量应为平均值。如不确定，请分别运行两次（每种统计量各一次）。

单组脚本（仅询问一组时使用）

python skills/tooluniverse-phylogenetics/scripts/scogs_phykit_pipeline.py \
    --data-folder "$DATA_PATH" --group fungi --metric treeness --out /tmp/f.tsv
# Auto-falls-back to .faa.mafft when .faa.mafft.clipkit is absent
# (some scogs zips ship only mafft alignments, not clipkit trims).

phykit parsimony_informative 不是有效的 CLI 子命令

PhyKIT 的 CLI 将简约信息位点计数公开为 parsimony_informative_sites（别名 pis）。调用 phykit parsimony_informative <file> 会返回帮助横幅并以非零状态退出，同时静默生成零值。捆绑脚本会自动将 parsimony_informative 转换为 parsimony_informative_sites。输出格式为 <n_pi>\t<n_total>\t<percent>——第三列是问题通常要求的百分比。

组中位数比率与配对比率（请仔细阅读）

如果问题将树长 / RCV / DVMC 的比较表述为“所有直系同源组中真菌 X 与动物 X 的比率”，则存在两种不同的量：

1. GROUP_MEDIAN_RATIO = median(values_fungi) / median(values_animals)。 两组分别使用各自的所有直系同源组。这通常是已发表的组间比较数据所报告的值（n_fungi 可能与 n_animals 不同，而且“所有直系同源组”是群体层面的表述，不是配对比较）。

2. PAIRED RATIO 中位数 = 对于两组中都存在的每个直系同源组，计算 value_fungi / value_animals，然后对所有共有直系同源组取中位数。分母较小（仅取交集），且当两组大小不同时，该数值也会不同。

除非问题明确提到“匹配的直系同源组”“配对直系同源组”“每个直系同源组的比率”或“针对每个直系同源组”，否则默认使用 GROUP_MEDIAN_RATIO。如果答案的措辞不明确，脚本输出中包含这两个数值——请根据问题中的“所有直系同源组”“配对”或“中位数之比”等措辞选择相应数值。

单拷贝直系同源基因的氨基酸总数 — 每个基因取一个代表序列，而非所有物种

当存在 BUSCO 单拷贝直系同源基因数据集（single_copy_busco_sequences/），且问题询问“所有单拷贝直系同源序列中共有多少个氨基酸”时，应当每个直系同源基因取一个代表序列，而不是对所有物种或拷贝求和。

single_copy_busco_sequences/ 中的每个 <ortholog_id>.faa 通常都包含该直系同源基因在多个物种中的拷贝（每个物种一个）。对所有物种的所有序列求和，会按物种数将每个直系同源基因重复计数 2、3 或 N 倍，得到 n_species × correct_answer。

问题表述	计数方式
“所有单拷贝直系同源序列中的氨基酸总数”	每个直系同源基因取一个序列求和（每个文件取第一条记录或长度中位数对应的记录）
“N 个物种的单拷贝直系同源基因中的氨基酸总数”	明确对所有物种求和（多物种总和）
“单拷贝直系同源基因的平均长度”	计算每个直系同源基因的平均长度（每个直系同源基因取一个）

❌ 错误：for f in *.faa: sum(len(rec.seq) for rec in SeqIO.parse(f, 'fasta')) 然后对所有文件求和（多物种总和）

✅ 正确：for f in *.faa: first_rec = next(SeqIO.parse(f, 'fasta')); total += len(first_rec.seq)（每个直系同源基因取一个代表序列）

如果你的答案是 n_species × GT（例如 GT 看起来是 13809，而结果为 32228 = 32228/2.33 ≈ 8 个物种 × 代表序列），说明你把所有物种都加总了——请改为每个直系同源基因取一个代表序列后重新计算。

最低非零比率

对于在高度保守或极短的比对中可能合理地等于 0 的指标（简约信息位点百分比、近乎相同序列上的 RCV），问题中的“最低”通常指“最低非零值”。配套脚本会为两种排序都输出 LOWEST_NONZERO_RATIO；如果某个组中的原始最小值为 0，请使用该行。

文件布局回退方案（比对文件命名）

scogs 压缩包有两种布局：

• 完整布局：<gene>.faa、<gene>.faa.mafft、<gene>.faa.mafft.clipkit、 <gene>.faa.mafft.clipkit.treefile，以及 iqtree/bionj/log/mldist 文件。
• 仅比对文件布局：只有 <gene>.faa 和 <gene>.faa.mafft。没有 树文件，也没有 clipkit 文件。用于简约法、RCV、空位百分比 相关问题。请使用 .faa.mafft（而非原始 .faa）——已发布的 指标是根据 MAFFT 比对文件计算的。

两个随附脚本都会自动检测布局，并为每个直系同源基因选用可用的最佳比对文件。请勿自行重新运行 MAFFT 或 ClipKit；随附文件是权威版本。

哪种指标对应哪种比对文件（这会改变结果）

树始终是由 ClipKit 衍生的 .faa.mafft.clipkit.treefile。比对文件参数取决于指标：

指标	传入的比对文件
treeness、dvmc、total_tree_length、long_branch_score	仅树文件——不需要比对文件
saturation	.faa.mafft.clipkit（修剪后）
treeness_over_rcv / rcv	.faa.mafft（未修剪）
简约信息位点、缺口百分比	.faa.mafft（未修剪）

RCV 衡量的是比对列之间的组成变异，因此修剪会显著改变它；而 treeness_over_rcv 以 RCV 为分母，所以使用修剪后的比对文件会改变每个基因的比值。已在真菌 scogs 数据集中验证（249 个直系同源基因，随附的标准文件）：

median treeness/RCV   untrimmed .faa.mafft = 0.2683    trimmed .clipkit = 0.3050
max treeness/RCV      (over the 3 genes with >70% gapped columns:
                       1260807at2759 0.0861, 1567796at2759 0.1866, 939345at2759 0.2572)
                      untrimmed .faa.mafft = 0.2572    trimmed .clipkit = 0.4174

普通的 treeness 不需要比对文件，也不受影响——它会完全复现结果（同一组 249 个文件的中位数为 0.0501），这表明问题源于比对文件的选择，而非树集或工具。

phykit_batch_analysis 分别接收这两个参数，因此请明确传入：

tu run phykit_batch_analysis '{"operation":"batch","function":"treeness_over_rcv",
  "directory":"<dir>","extension":".faa.mafft",
  "tree_directory":"<dir>","tree_extension":".faa.mafft.clipkit.treefile"}'

这些问题中的缺口百分比指的是包含至少一个缺口的比对列所占的比例，而不是所有残基中缺口所占的比例。两者相差一个数量级：按残基定义计算时，没有任何真菌直系同源基因的缺口比例超过 70%，因此“缺口超过 70%”的筛选条件会悄无声息地筛选出空集。

**反模式：**对从原始 *.busco.zip 中提取的直系同源基因 FASTA 文件运行 phykit，并自行进行比对和构树。scogs_*.zip 中的预计算文件才是标准输入。

---

PhyKIT、Biopython 和 DendroPy，用于比对与树分析、进化指标和比较基因组学。

查询，不要猜测

不确定任何科学事实时，请先搜索数据库。

---

适用场景

FASTA/PHYLIP/Nexus/Newick 文件；Treeness、RCV、DVMC、进化速率、简约位点、树长度、自举法；组间比较（Mann-Whitney U）；构树（NJ/UPGMA/简约法）；Robinson-Foulds 距离。

**从自己的序列从头开始进行比对 / 构树：**若要比对原始序列（而非预计算文件），请调用 EBI_msa_align（通过 EMBL-EBI 使用 Clustal Omega / MUSCLE / MAFFT / Kalign / T-Coffee），然后将其 data.aligned_fasta 字符串作为 EBI_build_phylogenetic_tree 的 aligned_sequences 参数传入（注意参数名与输出键不同），以构建邻接法或 UPGMA 树（Newick）。将生成的 Newick 树 / 比对文件直接输入下方的 PhyKIT 指标分析。

**仍然不适用于：**最大似然树（IQ-TREE/RAxML）或贝叶斯推断（MrBayes/BEAST）——EBI_build_phylogenetic_tree 仅支持基于距离的 NJ/UPGMA。用于发表的最大似然 / 贝叶斯系统发育分析，请运行专用工具；有预计算的 scogs_* 树时请使用它们。

---

所需软件包

import numpy as np, pandas as pd
from scipy import stats
from Bio import AlignIO, Phylo, SeqIO
from phykit.services.tree.treeness import Treeness
from phykit.services.tree.total_tree_length import TotalTreeLength
from phykit.services.tree.evolutionary_rate import EvolutionaryRate
from phykit.services.tree.dvmc import DVMC
from phykit.services.tree.treeness_over_rcv import TreenessOverRCV
from phykit.services.alignment.parsimony_informative_sites import ParsimonyInformative
from phykit.services.alignment.rcv import RelativeCompositionVariability
import dendropy

---

工作流决策树

ALIGNMENT ANALYSIS (FASTA/PHYLIP):
  Parsimony sites → phykit_parsimony_informative()
  RCV → phykit_rcv()
  Gap % → alignment_gap_percentage()

TREE ANALYSIS (Newick):
  Treeness → phykit_treeness()
  Tree length → phykit_tree_length()
  Evolutionary rate → phykit_evolutionary_rate()
  DVMC → phykit_dvmc()
  Bootstrap → extract_bootstrap_support()

COMBINED: Treeness/RCV → phykit_treeness_over_rcv(tree, aln)

TREE CONSTRUCTION: NJ → build_nj_tree(); UPGMA → build_upgma_tree(); Parsimony → build_parsimony_tree()

GROUP COMPARISON: batch metrics → Mann-Whitney U → summary stats

TREE COMPARISON: Robinson-Foulds → robinson_foulds_distance()

---

速查表

指标	输入	描述
Treeness	Newick	内部枝长 / 总枝长
RCV	FASTA/PHYLIP	相对组成变异
Treeness/RCV	两者	信号质量比值
树长度	Newick	所有枝长之和
进化速率	Newick	总长度 / 末端节点数
DVMC	Newick	分子钟违背程度
简约位点	FASTA/PHYLIP	出现次数不少于 2 且字符数不少于 2 的位点

---

常见模式

跨组计算单个指标

fungi_dvmc = batch_dvmc(discover_gene_files("data/fungi"))
animal_dvmc = batch_dvmc(discover_gene_files("data/animals"))
print(f"Fungi median: {np.median(list(fungi_dvmc.values())):.4f}")

统计比较

u_stat, p_value = stats.mannwhitneyu(list(g1.values()), list(g2.values()), alternative='two-sided')

筛选 + 指标

按缺口百分比 < 5% 筛选，然后在筛选结果中计算 treeness/RCV。

批量处理

gene_files = discover_gene_files("data/")  # → [{gene_id, aln_file, tree_file}]
treeness_results = batch_treeness(gene_files)  # → {gene_id: value}

---

答案提取

模式	方法
“median X”	np.median(values)
“maximum X”	np.max(values)
“difference in median”	abs(np.median(a) - np.median(b))
“Mann-Whitney U”	stats.mannwhitneyu(a, b)[0]
“fold-change”	np.median(a) / np.median(b)

舍入：PhyKIT 默认保留 4 位小数。U 统计量为整数。以问题的表述为准。

---

结果解读

指标	良好	可接受	较差
Treeness	>0.8	0.5-0.8	<0.5
RCV	<0.2	0.2-0.5	>0.5
Treeness/RCV	>2.0	1.0-2.0	<1.0
Bootstrap	>95%	70-95%	<70%
简约信息位点	>30%	10-30%	<10%

完整性检查清单

已识别所有文件；已检测到分组结构；使用了正确的 PhyKIT 函数；已处理所有基因（而非样本）；使用了正确的检验；舍入到 4 位小数；统计量具体明确（中位数/最大值/U/p）；Mann-Whitney 的 alternative='two-sided'。

---

分析约定

强制要求：批量计算时使用 phykit_batch_analysis 工具

对于任何要求计算多棵树/多个比对文件统计量的问题（中位数 treeness、平均饱和度、DVMC 百分比、gap 百分比、长枝分数），请使用 ToolUniverse 工具：

tu run phykit_batch_analysis '{"operation":"batch","function":"treeness","directory":"./trees","extension":".treefile"}'
tu run phykit_batch_analysis '{"operation":"batch","function":"saturation","directory":"./alignments","extension":".fa","tree_directory":"./trees","tree_extension":".treefile"}'
tu run phykit_batch_analysis '{"operation":"gap_percentage","directory":"./alignments","extension":".fa"}'

不要手动循环运行 phykit——该工具会处理所有文件并返回正确的汇总统计量。

批处理工具采用并行处理：约 250 棵树可在约 35 秒内完成。 对同样的任务，逐棵树运行 shell 循环约需 9 分钟，而且这是这类问题最终完全没有答案的最常见原因——运行会在循环中途耗尽轮次或时间预算，并回复“完成后我会汇报”，而不是给出数字。如果你发现自己正在编写 for f in *.treefile，请停下来调用批处理工具。

支持的 function 值包括 treeness、saturation、dvmc、 long_branch_score、total_tree_length、parsimony_informative、 treeness_over_rcv（别名 toverr）。dvmc 和 long_branch_score 已涵盖——无需为它们编写循环。

两组比较（Mann-Whitney U、中位数之差）。 比较真菌和动物的问题需要每组分别调用一次批处理工具，然后对两组数值列表进行检验——不要对两组数据逐棵树循环处理：

tu run phykit_batch_analysis '{"operation":"batch","function":"dvmc","directory":"<fungi>","extension":".treefile"}'
tu run phykit_batch_analysis '{"operation":"batch","function":"dvmc","directory":"<animals>","extension":".treefile"}'
# then scipy.stats.mannwhitneyu(fungi_values, animal_values)

如果检验需要完整列表，请在结果中请求 values；对于最多 50 个值的集合，批处理工具会返回这些值，并始终返回汇总统计量。对于更大的集合，请根据工具返回的各组汇总统计量计算统计量，而不是手动重新推导每个值。

PhyKIT 列约定——选对列

有几个 PhyKIT 子命令会为每个文件打印多个数字，而问题通常要的并不是第一个：

子命令	输出	问题所需的值
saturation	saturation <TAB> |saturation-1|	根据 phykit 文档，使用第 1 列；部分来源报告第 2 列——请说明你使用了哪一列
treeness_over_rcv	treeness/RCV <TAB> treeness <TAB> RCV	第 1 列，即比值
parsimony_informative_sites	n_pi <TAB> n_total <TAB> %PIS	百分比使用第 3 列

取 saturation 的第 1 列，结果恰好是 1 - answer：一个饱和度为 0.6146 的真菌集合会报告 0.3854，而不是 0.6146，两者相加为 1.0000，这就是判断依据。phykit_batch_analysis 已会为每个函数选择正确的列——这也是调用它而不是自行运行 CLI 的另一个理由。

提交你计算出的值

此基准测试中有两次失败，原因是算出了正确的数字，最后却给出了另一个数字：

• 树长比计算为 2.1775，但重新阅读“paired orthologs”后，回答成了 1.9；
• 平均树化程度在备选项中列出了 19，但最后提交了 10。

如果问题存在歧义，请计算你认为最符合字面意思的解读，用一个分句说明另一种解读，并且给出你实际计算出的值。最后一步不要用重新推导的结果替换已计算出的结果——如果两种解读都有道理，请先给出计算出的数字，再说明另一种解读，而不是悄悄换答案。

PhyKIT 列位置速查表（请仔细解析输出）

批量计算指标时，解析 PhyKIT 标准输出，需要的列取决于指标：

命令	输出列	取值列
phykit saturation	saturation_value <TAB> abs(saturation-1)	第 1 列是“饱和值”（1 = 无饱和；越接近 1，饱和程度越低）。第 2 列 = |saturation - 1|（与无饱和状态的距离；值越高，饱和程度越高，保留的信号越少）。对于“饱和值”问题，使用第 1 列；对于“距饱和状态的距离”问题，使用第 2 列
phykit toverr（又称 treeness_over_rcv）	treeness/RCV <TAB> treeness <TAB> RCV	第 1 列（treeness/RCV 比值）
phykit long_branch_score -v（详细模式）	每行 taxon <TAB> score	按每棵树汇总分数（均值）
phykit long_branch_score（不带 -v）	mean <TAB> median <TAB> 25%ile <TAB> 75%ile <TAB> min <TAB> max <TAB> std <TAB> var <TAB> n	对于“平均 LB 分数”，取第 1 列（均值）
phykit patristic_distances（不带 -v）	汇总统计行（格式与 LB 相同）	对于“平均 patristic distance”，取第 1 列（均值）

经验法则：phykit toverr 和 saturation 每个比对都会输出多列。不要用 grep 搜索“看起来像答案”的值——请根据 phykit <metric> --help 中的表头数列。如果批次中位数与已公布的数字相差很大（例如，预期值约为 0.26，实际得到的 median treeness/RCV 约为 0.20），几乎可以肯定是选错了列。

建议：不要手动解析 phykit 输出——调用 phykit_batch_analysis 工具，它会返回每项指标对应的正确列。支持的 function 值包括 treeness、saturation、dvmc、long_branch_score、total_tree_length、parsimony_informative：

tu run phykit_batch_analysis '{"operation":"batch","function":"saturation","directory":"./alignments","extension":".fa","tree_directory":"./trees","tree_extension":".treefile"}'
tu run phykit_batch_analysis '{"operation":"batch","function":"treeness","directory":"./alignments","extension":".fa","tree_directory":"./trees","tree_extension":".treefile"}'

对于 treeness_over_rcv（toverr / treeness/RCV 比值），该工具没有匹配的 function；请使用下面捆绑的 scogs_*.py 脚本，它们会直接计算该值。

生物学 scogs 树的合理性检查目标：median saturation 约为 0.4–0.7，median treeness/RCV 约为 0.2–0.4，median treeness 约为 0.05–0.15。数值与这些目标相差一个数量级，说明选错了列。

捆绑脚本：BUSCO target_orthologs 交集

如果数据文件夹中有 *.busco.zip 文件和 target_orthologs.txt，请使用捆绑脚本——不要手动枚举所有 zip 中的 single_copy_busco_sequences/*.faa：

python skills/tooluniverse-phylogenetics/scripts/busco_target_orthologs.py \
  --data-folder /path/to/data

默认运行会输出五行汇总结果，涵盖对“总氨基酸数”的所有常见解读：

# SUMMARY: n_targets=K, n_intersected=N (single-copy in ALL S species), intersected_total_aa=A, sum_all_aa=B
# SUMMARY group=all: intersected n=N total_aa=A, sum_all total_aa=B
# SUMMARY group=animals: sum_all total_aa=X        <-- per-group sum (animal species only)
# SUMMARY group=fungi:   sum_all total_aa=Y        <-- per-group sum (fungal species only)

选择合适的 SUMMARY 行（请仔细阅读）

根据问题的措辞匹配摘要行：

问题措辞	选择此行	原因
“所有单拷贝直系同源序列中的 AA 总量”，且数据文件夹中只有动物物种，或问题只提到一个生物类群	# SUMMARY group=animals: sum_all total_aa=...（或 group=fungi）	scogs 系统发育基因组学分析按类群分别运行；“所有”指该类群内的所有直系同源基因，而非跨类群的并集
“在每个／所有物种中均为单拷贝的直系同源基因的 AA 总量”	intersected_total_aa	严格采用交集规则
“所有物种各自拷贝的 AA 总量”	sum_all_aa（group=all）	仅当问题明确提到“所有物种”，或数据文件夹中只有一个生物类群时使用

默认规则：当数据文件夹同时包含动物和真菌的 busco 压缩包时，已发表的“氨基酸总量”答案几乎总是指一个类群（分析类群），而非跨类群并集。使用 group=animals: sum_all 或 group=fungi: sum_all。除非问题明确说 “所有 8 个物种中”或“真菌和动物合计”，否则不要选择 并集数值（sum_all_aa）。

脚本会出于这个原因，先在 stdout 中输出各类群的总和，再输出并集总和——请逐行读取输出，并在找到与问题所暗示的分析类群相符的 group=animals／group=fungi 行时停止。

跨物种单拷贝直系同源基因——比较集与交集

统计 BUSCO single_copy_busco_sequences/ 数据时，遵循以下两步规则：

1. 先确定比较集。 如果数据文件夹中存在 target_orthologs.txt（或类似名称的子集列表），该文件就是比较集——仅限其中列出的直系同源基因 ID。不要枚举各物种中的每个 BUSCO 单拷贝文件。如果提供了目标列表，不要假设“所有”指整个 BUSCO 输出。

2. 然后应用交集规则。 跨物种的“单拷贝直系同源基因”是指在比较集中的每个物种里都为单拷贝。如果某个物种的 single_copy_busco_sequences/ 中缺少某个直系同源基因，就将其完全排除在统计之外——不要只统计有该基因的物种。

合理性检查：如果某个物种的每个直系同源基因计数明显低于其他物种（例如，同一组直系同源基因中，一个物种约为 600 aa，而其他物种超过 4000 aa），那么部分物种缺失的直系同源基因抬高了每个直系同源基因的平均值——请先剔除它们。

示例。 数据文件夹包含 8 个物种（4 个动物、4 个真菌）的 *.busco.zip，以及列出 10 个直系同源基因 ID 的 target_orthologs.txt：

• 错误：枚举全部 8 个物种中的所有 single_copy_busco_sequences/*.faa → 约 80 个文件 → 汇总 AA → 得出 32228（将每个物种的拷贝都单独计算）。
• 正确：对 10 个目标 ID 分别检查其是否在全部 8 个物种中均为 single_copy → 仅保留交集 ID（通常为 5/10——部分目标 ID 在某个物种中为多拷贝或缺失）→ 对保留的 ID 汇总 8 个物种中的 AA → 得出 13809。
• “5 棵树”的含义：如果问题说“5 棵树”，但你找到了 10 个树文件，GT 使用的是交集子集（在所有物种中均为单拷贝的直系同源基因），而不是全部 10 个。计算平均值前，请重新推导交集。

处理完整数据集，而非样本（批量指标）

当问题要求计算跨直系同源基因的中位数／百分位数／均值时，批次必须包含相关比较集中的每一个直系同源基因：

• scogs_fungi.zip 包含约 255 个真菌比对文件和树；scogs_animals.zip 包含约 241 个。根据 10 个文件样本计算的中位数不是已发表的答案。
• 对于 phykit_batch_analysis，始终指向包含所有逐直系同源基因文件的已解压 scogs 目录，而不要使用手动挑选的子集。
• 如果计算得出的 RCV／treeness／DVMC 中位数与合理性检查目标相差超过 10%，请先清点文件数量——你很可能只处理了一个子集。

先筛选再计算（不要先计算再筛选）

对于“在缺口超过 70% 的基因中，X 的最大值是多少”这类问题，必须先筛选，再取最大值：

# 1. Compute gap% per alignment
# 2. Keep only alignments with gap% > 70
# 3. Compute treeness/RCV ON THE FILTERED SET
# 4. Take max

先对所有基因计算指标再取最大值，会得到全局最大值，这是错误的。

动物与真菌——长枝分数的汇总

PhyKIT 的 long_branch_score -v 会输出每个分类单元的 LB 分数（树中每个叶节点一行）。对于每棵树的汇总：

1. 每棵树：运行 phykit long_branch_score -v <tree> → 获取每个分类单元的分数列表。
2. 每棵树的汇总：使用这些分类单元分数的均值或中位数，将每棵树汇总为一个数字。
3. 每组的汇总：汇总每棵树的数值（中位数／均值／Mann-Whitney U + p 值）。

根据问题表述选择每棵树的汇总方式：

问题表述……	使用 --per-tree-stat ...
“平均长枝分数”	mean
“长枝分数中位数”	median
“长枝分数平均值”（含义不明确）	两种都运行，并选择与数值／单位相符的一种

随附的 scogs_paired_compare.py --metric long_branch_score --per-tree-stat {mean,median} 会在一次运行中为两组执行步骤 1 和 2，并直接输出组间比较的 Mann-Whitney U 和 p 值。

常见错误：不经过逐树汇总，直接分别对四个动物物种和四个真菌物种取均值——这会混淆物种 LB 和直系同源基因 LB，导致差值错误。

Treeness/RCV：使用正确的输入文件

phykit toverr（即 treeness_over_rcv）需要比对文件和树文件这两者。使用配对的修剪后比对文件（*.faa.mafft.clipkit）和树文件（*.faa.mafft.clipkit.treefile），而不是原始的 .faa.mafft。ClipKit 修剪后的比对用于生成标准树，因此必须在同一修剪后比对上计算 RCV，才能使比值与已发表的数值相符。

简约信息位点

• 计数前排除全为缺口的列——全是缺口的列不具信息性。
• 当 ≥2 个不同的非缺口字符各自在 ≥2 个分类单元中出现时，该位点才具有简约信息。
• 使用 Biopython AlignIO 或 AMAS 工具遍历各列并计数。

Treeness（RCV 比值）

Treeness = 内部枝长之和 / 树总长度。内部枝是指不通向叶节点（末端）的枝。

PhyKIT 用法

PhyKIT（pip install phykit）提供用于计算树和比对统计量的命令行函数。常用函数：

• phykit treeness <tree_file> — 输出 treeness（RCV）值
• phykit saturation <alignment_file> -t <tree_file> — 输出饱和度值
• phykit dvmc <tree_file> — 分子钟偏离程度
• phykit long_branch_score <tree_file> — 长枝分数（LBS）
• phykit alignment_length <alignment_file> — 比对长度
• phykit parsimony_informative <alignment_file> — 具有简约信息的位点数

对多个基因树／比对运行 PhyKIT 时，请使用随附的批处理脚本：

# Treeness across all trees
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir scogs_fungi --function treeness --ext .treefile --stat median

# Saturation with paired alignment+tree
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir alignments --function saturation --tree-dir trees \
  --ext .fa --tree-ext .treefile --stat median

# Long branch score (mean per tree, then median across trees)
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir trees --function long_branch_score --ext .treefile \
  --per-tree-stat mean --stat median

# DVMC
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir trees --function dvmc --ext .treefile --stat all

# Gap percentage across all alignments
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir alignments --function gap_percentage --ext .fa

# Evolutionary rate (median across trees)
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir trees --function evolutionary_rate --ext .treefile --stat median

# Mean patristic distance per tree, then mean across trees
python skills/tooluniverse-phylogenetics/scripts/phykit_batch.py \
  --dir trees --function patristic_distances --ext .treefile --stat mean

首选：使用 phykit_batch_analysis ToolUniverse 工具，而不是手动运行 PhyKIT：

# Via CLI
tu run phykit_batch_analysis '{"operation":"batch","function":"treeness","directory":"/path/to/trees","extension":".treefile"}'

# Via SDK
tu.run_one_function({"name": "phykit_batch_analysis", "arguments": {"operation": "batch", "function": "saturation", "directory": "/path/to/alignments", "extension": ".fa", "tree_directory": "/path/to/trees"}})

# Gap percentage
tu run phykit_batch_analysis '{"operation":"gap_percentage","directory":"/path/to/alignments","extension":".fa"}'

关键规则：

1. 处理所有文件——不要只处理其中一部分。该工具会自动处理全部文件
2. 缺口百分比：所有比对中的缺口总数 / 总位点数（不是各文件平均值）
3. 长枝分数：每棵树会产生每个分类单元的分数 → 按树汇总（均值）→ 再对所有树汇总（中位数）。使用 "per_tree_stat":"mean"
4. 真菌与动物的比较：按直系同源 ID（文件名主干）匹配基因，而非按文件顺序匹配。分别对每种生物的目录运行该工具，然后比较中位数

参考资料

references/sequence_alignment.md、references/tree_building.md、references/parsimony_analysis.md、scripts/tree_statistics.py

• PhyKIT：https://jlsteenwyk.com/PhyKIT/
• Biopython Phylo：https://biopython.org/wiki/Phylo
• DendroPy：https://dendropy.org/