R语言 plantlist 0.6.5 + ape 包:14个物种系统发育树构建与性状数据整合实战 R语言 plantlist 0.6.5 ape 包14个物种系统发育树构建与性状数据整合实战在生态学和进化生物学研究中系统发育树Phylogenetic Tree是理解物种间进化关系的重要工具。传统构建系统发育树需要基因序列数据但对于许多研究场景特别是植物学研究我们往往只需要基于分类学关系构建一个合理的系统发育框架。本文将详细介绍如何利用R语言中的plantlist、ape和phylosignal等包实现从物种名录到系统发育树构建再到性状数据整合可视化的完整流程。1. 环境准备与数据获取1.1 安装必要的R包在开始之前确保已安装以下R包。如果尚未安装可以使用以下命令安装install.packages(c(plantlist, ape, phylosignal, phylobase, readxl, dplyr))这些包的主要功能如下plantlist用于植物学名查询和验证ape系统发育树分析和可视化phylosignal系统发育信号分析phylobase系统发育数据整合readxl读取Excel文件dplyr数据整理和操作1.2 准备物种名录数据假设我们需要为以下14个植物物种构建系统发育树species_list - c( Chimonanthus praecox, Michelia figo, Cycas revoluta, Ginkgo biloba, Litsea rotundifolia, Bombax ceiba, Delonix regia, Ficus microcarpa, Mangifera indica, Dimocarpus longan, Triadica sebifera, Liquidambar formosana, Alstonia scholaris, Camphora officinarum )提示实际研究中物种名录通常存储在Excel或CSV文件中。可以使用readxl::read_excel()或read.csv()导入。2. 物种分类学信息查询2.1 使用plantlist包验证学名plantlist包的TPL()函数可以验证植物学名并获取分类学信息library(plantlist) taxa_check - TPL(species_list)查看前几行结果head(taxa_check)输出将包含每个物种的科、属、种信息以及APG III分类系统中的科编号等。2.2 转换分类学格式为了后续构建系统发育树我们需要将查询结果转换为科/属/种格式taxa_table - taxa.table(taxa_check) print(taxa_table)输出示例[1] Calycanthaceae/Chimonanthus/Chimonanthus_praecox [2] Magnoliaceae/Michelia/Michelia_figo [3] Cycadaceae/Cycas/Cycas_revoluta ...3. 系统发育树构建3.1 使用Phylomatic生成Newick格式树虽然可以直接在R中完成所有步骤但使用Phylomatic工具可以获得更准确的系统发育关系。以下是完整流程将分类学信息保存为文本文件writeLines(taxa_table, taxa_list.txt)使用Phylomatic的在线服务或本地awk版本生成Newick格式树文件保存为tree.newick3.2 在R中读取和可视化系统发育树使用ape包读取和绘制系统发育树library(ape) phy_tree - read.tree(tree.newick) # 或直接使用字符串 plot(phy_tree, cex 0.8)基础绘图效果可能不够理想我们可以进行美化plot(phy_tree, type fan, edge.width 2, tip.color blue, cex 0.8, main 14种植物系统发育树)4. 性状数据整合与可视化4.1 准备性状数据假设我们有以下性状数据通常来自实验测量或文献speciesbiomassleaf_sizegrowth_rateChimonanthus praecox12.515.20.45Michelia figo8.712.10.38............将数据导入R并整理library(dplyr) trait_data - read.csv(trait_data.csv) %% mutate(species gsub( , _, species)) %% group_by(species) %% summarise(across(everything(), mean, na.rm TRUE)) %% column_to_rownames(species)4.2 性状数据与系统发育树整合使用phylobase包将性状数据映射到系统发育树上library(phylobase) library(phylosignal) p4d - phylo4d(phy_tree, trait_data)可视化整合结果barplot.phylo4d(p4d, tree.type phylo, tree.ladderize TRUE, trait.bg.col lightblue, scale FALSE)4.3 系统发育信号分析phylosignal包可以量化性状的系统发育信号phylo_signal - phyloSignal(p4d) print(phylo_signal)结果将显示每个性状的Blombergs K和Pagels λ等指标帮助我们理解性状在系统发育树上的分布模式。5. 高级分析与可视化技巧5.1 系统发育主成分分析Phylogenetic PCAlibrary(phytools) phylo_pca - phyl.pca(phy_tree, trait_data) biplot(phylo_pca)5.2 性状演化模型拟合比较不同演化模型对性状数据的拟合优度models - c(BM, OU, EB) aic_scores - sapply(models, function(m) { fitContinuous(phy_tree, trait_data$biomass, model m)$opt$aic })5.3 使用ggtree增强可视化library(ggtree) library(ggplot2) ggtree(phy_tree) %% trait_data geom_tiplab(aes(color biomass), size 3) scale_color_gradient(low blue, high red) theme_tree2()6. 完整可复现代码示例以下是一个端到端的可执行脚本整合了上述所有步骤# 加载所需包 library(plantlist) library(ape) library(phylobase) library(phylosignal) library(dplyr) library(ggtree) # 1. 物种名录 species - c(Chimonanthus praecox, Michelia figo, Cycas revoluta, Ginkgo biloba, Litsea rotundifolia, Bombax ceiba, Delonix regia, Ficus microcarpa, Mangifera indica, Dimocarpus longan, Triadica sebifera, Liquidambar formosana, Alstonia scholaris, Camphora officinarum) # 2. 验证学名并获取分类信息 taxa_check - TPL(species) taxa_table - taxa.table(taxa_check) # 3. 构建系统发育树 writeLines(taxa_table, taxa_list.txt) # 此处需运行Phylomatic生成tree.newick phy_tree - read.tree(tree.newick) # 4. 性状数据整合 trait_data - data.frame( species gsub( , _, species), biomass runif(14, 5, 15), leaf_size runif(14, 10, 20) ) %% column_to_rownames(species) # 5. 可视化 p4d - phylo4d(phy_tree, trait_data) barplot.phylo4d(p4d) # 6. 系统发育信号分析 phyloSignal(p4d) # 7. ggtree高级可视化 ggtree(phy_tree) %% as.data.frame(trait_data) geom_tiplab(aes(color biomass), size 3, offset 0.1) geom_tippoint(aes(size leaf_size), shape 16) scale_color_gradient(low blue, high red) theme_tree2()在实际项目中这套流程已经帮助我快速分析了多个植物功能性状的系统发育模式。特别是在处理50-100个物种的中等规模数据集时这种方法既保证了科学性又大大节省了时间。