BioGeoBEARS 生物地理分析

概览

BioGeoBEARS (BioGeography with Bayesian and Likelihood Evolutionary Analysis in R Scripts) 执行系统发育树上祖先地理范围的概率推断。此技能通过以下方式帮助设置完整的生物地理分析：

验证和重新格式化输入文件（系统发育树和地理分布数据）
生成有组织的分析文件夹结构
创建定制的RMarkdown分析脚本
指导用户选择参数和模型选择
生产出适合出版的可视化结果

何时使用此技能

当用户请求：

“在我的系统发育树上分析生物地理”
“重建我物种的祖先范围”
“运行BioGeoBEARS分析”
“我的祖先占据了哪些区域？”
“测试生物地理模型（DEC, DIVALIKE, BAYAREALIKE）”

此技能在用户提及系统发育生物地理学、祖先区域重建或提供树+分布数据时触发。

所需输入

用户必须提供：

系统发育树（Newick格式，.nwk, .tre, 或 .tree文件）
- 必须是有根的
- 树尖标签将与地理文件匹配
- 需要分支长度
地理分布数据（任何表格格式）
- 物种名称（与树尖匹配）
- 不同地理区域的存在/缺失数据
- 可以是CSV, TSV, Excel, 或已经处于PHYLIP格式

工作流程

第1步：收集信息

当用户请求BioGeoBEARS分析时，询问：

输入文件路径：
- “你的系统发育树文件路径是什么？”
- “你的地理分布文件路径是什么？”
分析参数（如果未指定）：
- 最大范围大小（一个物种可以同时占据多少区域？）
- 比较哪些模型（默认：全部六个 - DEC, DEC+J, DIVALIKE, DIVALIKE+J, BAYAREALIKE, BAYAREALIKE+J）
- 输出目录名称（默认：“biogeobears_analysis”）

使用AskUserQuestion工具高效收集此信息：

示例问题：
- “最大范围大小” - 基于区域数量的选项（例如，对于4个区域：“全部4个区域”, “3个区域”, “2个区域”）
- “比较模型” - 选项：“全部6个模型（推荐）”, “仅基础模型（DEC, DIVALIKE, BAYAREALIKE）”, “仅+J模型”, “自定义选择”
- “可视化类型” - 选项：“饼图（显示概率）”, “文本标签（显示最可能的状态）”, “两者”

第2步：验证和准备输入文件

验证树文件

使用Read工具检查树文件：

# 在R中，基本验证：
library(ape)
tr <- read.tree("path/to/tree.nwk")
print(paste("Tips:", length(tr$tip.label)))
print(paste("Rooted:", is.rooted(tr)))
print(tr$tip.label)  # 检查物种名称

验证：

文件可以解析为Newick
树是有根的（如果不是，询问用户使用哪个外群）
记录树尖标签以验证地理文件

验证和重新格式化地理文件

使用scripts/validate_geography_file.py验证或重新格式化地理文件。

如果文件已经是PHYLIP格式（以数字开头）：

python scripts/validate_geography_file.py path/to/geography.txt --validate --tree path/to/tree.nwk

这检查：

正确的制表符分隔符
物种名称与树尖匹配
二进制代码长度正确
物种名称或二进制代码中没有空格

如果文件是CSV/TSV格式（需要重新格式化）：

python scripts/validate_geography_file.py path/to/distribution.csv --reformat -o geography.data --delimiter ","

或对于制表符分隔：

python scripts/validate_geography_file.py path/to/distribution.txt --reformat -o geography.data --delimiter tab

该脚本将：

从标题行检测区域名称
将存在/缺失数据转换为二进制（处理“1”, “present”, “TRUE”等）
从物种名称中移除空格（用下划线替换）
创建正确格式化的PHYLIP文件

在继续之前始终验证重新格式化的文件：

python scripts/validate_geography_file.py geography.data --validate --tree path/to/tree.nwk

第3步：设置分析文件夹结构

为分析创建一个有组织的目录：

biogeobears_analysis/
├── input/
│   ├── tree.nwk                 # 原始或复制的树
│   ├── geography.data            # 验证/重新格式化的地理文件
│   └── original_data/            # 原始输入文件
│       ├── original_tree.nwk
│       └── original_distribution.csv
├── scripts/
│   └── run_biogeobears.Rmd       # 生成的RMarkdown脚本
├── results/                      # 由分析创建（输出目录）
│   ├── [MODEL]_result.Rdata      # 保存的模型结果
│   └── plots/                    # 可视化输出
│       ├── [MODEL]_pie.pdf
│       └── [MODEL]_text.pdf
└── README.md                     # 分析文档

以编程方式创建此结构：

mkdir -p biogeobears_analysis/input/original_data
mkdir -p biogeobears_analysis/scripts
mkdir -p biogeobears_analysis/results/plots

# 复制文件
cp path/to/tree.nwk biogeobears_analysis/input/
cp geography.data biogeobears_analysis/input/
cp original_files biogeobears_analysis/input/original_data/

第4步：生成RMarkdown分析脚本

使用scripts/biogeobears_analysis_template.Rmd模板，并根据用户参数进行定制。

复制和定制模板：

cp scripts/biogeobears_analysis_template.Rmd biogeobears_analysis/scripts/run_biogeobears.Rmd

创建参数文件或修改Rmd中的YAML标题以使用用户的特定设置：

示例定制通过R代码：

# 通过程序编辑YAML参数或在渲染时提供参数
rmarkdown::render(
  "biogeobears_analysis/scripts/run_biogeobears.Rmd",
  params = list(
    tree_file = "../input/tree.nwk",
    geog_file = "../input/geography.data",
    max_range_size = 4,
    models = "DEC,DEC+J,DIVALIKE,DIVALIKE+J,BAYAREALIKE,BAYAREALIKE+J",
    output_dir = "../results"
  ),
  output_file = "../results/biogeobears_report.html"
)

或创建运行脚本：

# biogeobears_analysis/run_analysis.sh
#!/bin/bash
cd "$(dirname "$0")/scripts"

R -e "rmarkdown::render('run_biogeobears.Rmd', params = list(
  tree_file = '../input/tree.nwk',
  geog_file = '../input/geography.data',
  max_range_size = 4,
  models = 'DEC,DEC+J,DIVALIKE,DIVALIKE+J,BAYAREALIKE,BAYAREALIKE+J',
  output_dir = '../results'
), output_file = '../results/biogeobears_report.html')"

第5步：创建README文档

在分析目录中生成README.md，解释：

有哪些文件
如何运行分析
使用了哪些参数
如何解释结果

示例：

# BioGeoBEARS分析

## 概览

[NUMBER]物种在[NUMBER]地理区域的生物地理分析。

## 输入数据

- **树**：`input/tree.nwk` ([NUMBER]树尖)
- **地理**：`input/geography.data` ([NUMBER]物种 × [NUMBER]区域)
- **区域**：[A, B, C, ...]

## 参数

- 最大范围大小：[NUMBER]
- 测试模型：[LIST]

## 运行分析

### 选项1：直接使用RMarkdown

```r
library(rmarkdown)
render("scripts/run_biogeobears.Rmd",
       output_file = "../results/biogeobears_report.html")

选项2：使用运行脚本

bash run_analysis.sh

输出

结果将保存在results/：

biogeobears_report.html - 完整的分析报告和可视化
[MODEL]_result.Rdata - 每个模型保存的R对象
plots/[MODEL]_pie.pdf - 祖先范围重建（饼图）
plots/[MODEL]_text.pdf - 祖先范围重建（文本标签）

解释结果

HTML报告包括：

模型比较 - AIC分数，AIC权重，最佳拟合模型
参数估计 - 扩散（d），灭绝（e），创始事件（j）率
似然比测试 - 嵌套模型的统计比较
祖先范围图 - 在系统发育树上的可视化
会话信息 - R包版本以重现

模型描述

DEC：扩散-灭绝-谱系起源（通用）
DIVALIKE：强调 vicariance
BAYAREALIKE：强调同域物种形成
+J：添加创始事件物种形成参数

见references/biogeobears_details.md了解详细的模型描述。

安装要求

# 安装BioGeoBEARS
install.packages("rexpokit")
install.packages("cladoRcpp")
library(devtools)
devtools::install_github(repo="nmatzke/BioGeoBEARS")

# 其他包
install.packages(c("ape", "rmarkdown", "knitr", "kableExtra"))

第6步：向用户提供使用说明

设置分析后，向用户提供清晰的说明：

分析设置完成！

目录结构创建于：biogeobears_analysis/

📁 文件创建：
   ✓ input/tree.nwk - 系统发育树（[N]树尖）
   ✓ input/geography.data - 地理分布数据（验证）
   ✓ scripts/run_biogeobears.Rmd - RMarkdown分析脚本
   ✓ README.md - 文档和说明
   ✓ run_analysis.sh - 方便运行分析的脚本

📋 下一步：

1. 查看README.md了解分析细节

2. 如果尚未安装BioGeoBEARS，请安装：
   ```r
   install.packages("rexpokit")
   install.packages("cladoRcpp")
   library(devtools)
   devtools::install_github(repo="nmatzke/BioGeoBEARS")

运行分析：

cd biogeobears_analysis
bash run_analysis.sh

或在R中：

setwd("biogeobears_analysis")
rmarkdown::render("scripts/run_biogeobears.Rmd",
                  output_file = "../results/biogeobears_report.html")

查看结果：
- 在网页浏览器中打开results/biogeobears_report.html
- 检查results/plots/中的PDF可视化

⏱️ 预计运行时间：[基于树大小的估计]

小树（<50树尖）：5-15分钟
中等树（50-100树尖）：15-60分钟
大树（>100树尖）：1-4小时

💡 HTML报告包括模型比较、参数估计和您的系统发育树上祖先范围的可视化。


## 分析参数指导

当用户询问参数指导时，参考`references/biogeobears_details.md`并提供建议：

### 最大范围大小

**询问**：“您组中的物种最多可以实际占据多少个区域？”

常见方法：
- **保守**：区域数量 - 1（防止不切实际的泛世界祖先范围）
- **宽容**：所有区域（如果生物学上合理）
- **数据驱动**：现存物种中观察到的最大值

**影响**：更大的值会指数级增加计算时间

### 模型选择

**默认建议**：运行全部6个模型进行全面比较

- DEC, DIVALIKE, BAYAREALIKE（基础模型）
- DEC+J, DIVALIKE+J, BAYAREALIKE+J（+J变体）

**理由**：
- 模型比较是推断的关键
- +J参数通常很重要
- 计算成本增加很小

如果计算是问题，建议从DEC和DEC+J开始。

### 可视化选项

**饼图**（`plotwhat = "pie"`）：
- 显示所有可能状态的概率分布
- 更适合传达不确定性
- 许多区域时可能会混乱

**文本标签**（`plotwhat = "text"`）：
- 只显示最大似然状态
- 清洁，易于阅读
- 不显示不确定性

**建议**：在分析中生成两者（模板自动执行此操作）

## 常见问题和故障排除

### 物种名称不匹配

**症状**：关于树中物种不在地理文件中（反之亦然）的错误

**解决方案**：使用带有`--tree`选项的验证脚本识别不匹配，然后要么：
1. 编辑地理文件以匹配树尖标签
2. 编辑树尖标签以匹配地理文件
3. 删除不在两者中的物种

### 树没有根

**症状**：关于未根树的错误

**解决方案**：
```r
library(ape)
tr <- read.tree("tree.nwk")
tr <- root(tr, outgroup = "outgroup_species_name")
write.tree(tr, "tree_rooted.nwk")

询问用户使用哪个物种作为外群。

地理文件格式错误

症状：关于制表符、空格或二进制代码的验证错误

解决方案：使用重新格式化选项：

python scripts/validate_geography_file.py input.csv --reformat -o geography.data

优化未能收敛

症状：参数估计中的NA值或非常负的对数似然值

可能的原因：

树和地理数据不匹配
所有物种在同一区域（没有变化）
不切实际的max_range_size

解决方案：检查输入数据质量，尝试更简单的模型（仅DEC）

运行时间非常慢

原因：

区域数量大（>6-7区域变慢）
最大范围大小大
许多树尖（>200）

解决方案：

减少max_range_size
如适当，合并地理区域
使用force_sparse = TRUE在运行对象中
在HPC集群上运行

资源

此技能包括：

scripts/

validate_geography_file.py - 验证和重新格式化地理文件
- 检查PHYLIP格式合规性
- 与树尖标签验证
- 从CSV/TSV重新格式化为PHYLIP
- 使用方法：python validate_geography_file.py --help
biogeobears_analysis_template.Rmd - 完整的RMarkdown分析模板
- DEC, DIVALIKE, BAYAREALIKE（带/不带+J）模型拟合
- AIC, AICc, 权重的模型比较
- 似然比测试
- 参数可视化
- 祖先范围绘制
- 可通过YAML参数定制

references/

biogeobears_details.md - 综合参考包括：
- 详细的模型描述
- 输入文件格式规范
- 参数解释指南
- 绘图选项和定制
- 引用和进一步阅读
- 计算考虑

加载此参考时：

当用户需要模型、参数或解释的详细说明时
用于故障排除复杂问题
帮助用户为出版物编写方法部分

最佳实践

始终在分析前验证输入文件 - 节省后续调试时间
在专用目录中组织分析 - 保持一切在一起并可重现
默认运行全部6个模型 - 模型比较对生物地理推断至关重要
记录参数和决策 - 分析README有助于可重现性
生成两种可视化类型 - 饼图用于不确定性，文本标签用于清晰度
保存中间结果 - RMarkdown模板自动执行此操作
检查参数估计 - 不切实际的值表明数据或模型问题
用可视化提供上下文 - 解释扩散/灭绝率对用户系统的意义

输出解释

向用户呈现结果时，解释：

模型选择

AIC权重代表每个模型是最好的概率
ΔAIC < 2：模型基本等价
ΔAIC 2-7：支持度较低
ΔAIC > 10：几乎没有支持

参数估计

d（扩散率）：较高=更多范围扩展
e（灭绝率）：较高=更多局部灭绝
j（创始事件率）：较高=更多跳跃扩散在物种形成时
比率d/e：> 1支持扩展，< 1支持收缩

祖先范围

饼图：较大的部分=更高概率
颜色：代表区域（单一区域=亮色，多个区域=混合色）
节点标签：最可能的祖先范围
分裂事件（在角落）：物种形成时范围变化

统计测试

LRT p < 0.05：+J参数显著提高拟合度
高AIC权重（>0.7）：一个模型的强有力证据
相似的AIC权重：模型不确定性 - 报告多个模型的结果

示例用法

用户：“我有30种鸟类的系统发育树和它们在5个岛屿上的分布。你能帮我弄清楚它们的祖先生活在哪里吗？”

Claude（使用此技能）：
1. 询问树和分布文件路径
2. 验证树文件（检查30个树尖，有根）
3. 验证/重新格式化地理文件（5个区域）
4. 询问max_range_size（建议4个区域）
5. 询问模型（建议全部6个）
6. 设置biogeobears_analysis/目录结构
7. 复制带有参数的模板RMarkdown脚本
8. 生成README.md和run_analysis.sh
9. 提供清晰的运行分析说明
10. 解释预期输出及如何解释它们

结果：用户拥有完整、准备就绪的分析和文档

归属

此技能是基于以下创建的：

BioGeoBEARS包由Nicholas Matzke创建
教程资源来自http://phylo.wikidot.com/biogeobears
来自BioGeoBEARS GitHub仓库的示例工作流程

附加说明

技能执行时间估计：

文件验证：1-2分钟
目录设置：< 1分钟
总设置时间：5-10分钟

分析运行时间（与技能执行分开）：

取决于树大小和区域数量
小数据集（<50树尖，≤5区域）：10-30分钟
大数据集（>100树尖，>5区域）：1-6小时

安装要求（用户必须有）：

R（≥4.0）
BioGeoBEARS R包
支持包：ape, rmarkdown, knitr, kableExtra
Python 3（用于验证脚本）

何时参考：

当用户需要模型、参数或解释的详细说明时
用于故障排除复杂问题
帮助用户为出版物编写方法部分