
解码布兰特鹅的饮食密码用SIMMR实现稳定同位素混合模型全流程解析生态学家们常常面临一个看似简单却极其复杂的问题野生动物究竟吃什么传统的研究方法如粪便分析或直接观察往往存在局限性而稳定同位素分析技术为我们打开了一扇新窗。本文将带您走进稳定同位素混合模型的世界通过R语言中的SIMMR包以布兰特鹅的饮食研究为例展示从原始数据到可视化结果的完整分析流程。1. 稳定同位素分析的基础原理当布兰特鹅在湿地觅食时它们摄入的食物中的稳定同位素会在其组织中留下独特的化学指纹。碳(δ13C)和氮(δ15N)稳定同位素比值能够反映生物体的饮食来源因为δ13C值主要指示食物的基础来源如陆生vs水生植物δ15N值反映生物在食物链中的营养级位置这些同位素比值在食物链传递过程中会发生可预测的变化使我们能够追溯消费者的饮食组成。SIMMRStable Isotope Mixing Model in R正是利用这一原理通过贝叶斯统计方法解决以下核心问题给定消费者的同位素组成和潜在食物源的同位素特征各种食物源在消费者饮食中的比例是多少关键概念速览表术语解释典型值范围δ13C碳同位素比值-30‰ 到 -10‰δ15N氮同位素比值0‰ 到 20‰TEF营养富集因子δ13C: ~0.5-1‰; δ15N: ~2-4‰凸包源在同位素空间形成的多边形-# 查看SIMMR包内置数据集 library(simmr) data(packagesimmr) # 输出示例 Data sets in package simmr: geese_data Geese stable isotope mixing data set geese_data_day1 A smaller version of the geese dataset2. 数据准备与SIMMR环境搭建2.1 系统环境配置在开始分析前需要确保系统已安装以下组件R版本≥3.5.0JAGS(Just Another Gibbs Sampler)SIMMR依赖的MCMC引擎R必要包R2jags、ggplot2、simmr# 安装SIMMR依赖包 install.packages(c(R2jags, ggplot2, simmr)) # 验证JAGS安装 library(R2jags) test_jags() # 应返回TRUE2.2 数据加载与初步探索我们将使用SIMMR内置的geese_data_day1数据集包含9个鹅组织样本的稳定同位素测量值4种潜在食物源Zostera大叶藻、Grass草、U.lactuca石莼、Enteromorpha浒苔2种同位素δ13C和δ15N校正因子营养富集因子(TEFs)的均值和标准差浓度依赖各食物源中元素的相对比例# 加载数据并创建simmr输入对象 data(geese_data_day1) simmr_in - with( geese_data_day1, simmr_load( mixtures mixtures, source_names source_names, source_means source_means, source_sds source_sds, correction_means correction_means, correction_sds correction_sds, concentration_means concentration_means ) )2.3 同位素空间可视化数据质量检查是分析的关键第一步。plot.simmr_input函数生成同位素空间图直观展示混合物鹅组织与食物源的关系# 绘制同位素空间图 plot(simmr_in, xlab expression(paste(delta^13, C (\u2030))), ylab expression(paste(delta^15, N (\u2030))))理想的建模数据应满足大多数混合物点位于食物源构成的凸包内各食物源在同位素空间中有明显区分数据点分布无明显异常值注意如果混合物点大量落在凸包外可能需要重新考虑食物源选择或检查数据质量3. 模型运行与诊断3.1 MCMC采样过程SIMMR使用马尔可夫链蒙特卡洛(MCMC)方法估计膳食比例的后验分布。核心函数simmr_mcmc的关键参数包括iter总迭代次数默认10000burn老化期迭代数默认1000thin抽样间隔默认10n.chain并行链数默认4# 运行MCMC采样 simmr_out - simmr_mcmc(simmr_in) # 查看模型概要 print(simmr_out)3.2 收敛性诊断MCMC结果的可靠性取决于采样是否充分收敛。SIMMR提供多种诊断工具# 收敛诊断Gelman-Rubin统计量 summary(simmr_out, type diagnostics) # 各参数应满足 # - 点估计≈1 # - 97.5%分位数1.1若诊断结果不理想可尝试增加迭代次数iter50000延长老化期burn5000增大细化间隔thin203.3 后验统计量提取获取膳食比例的定量估计# 获取分位数统计 summary(simmr_out, type quantiles) # 获取均值与标准差 summary(simmr_out, type statistics) # 获取参数间相关性 summary(simmr_out, type correlations)典型输出示例Summary for quantiles: 2.5% 25% 50% 75% 97.5% Zostera 0.421 0.549 0.617 0.687 0.817 Grass 0.035 0.060 0.073 0.088 0.116 U.lactuca 0.019 0.072 0.124 0.189 0.322 Enteromorpha 0.022 0.082 0.151 0.239 0.4274. 结果可视化与解读4.1 膳食比例分布可视化SIMMR提供多种图形展示后验分布# 箱线图展示中位数和四分位距 plot(simmr_out, type boxplot) # 密度图展示概率分布形态 plot(simmr_out, type density) # 矩阵图揭示源间的相关性 plot(simmr_out, type matrix)4.2 食物源组合分析当两个食物源在同位素空间位置相近时可考虑合并以减少不确定性# 合并U.lactuca和Enteromorpha simmr_combined - combine_sources( simmr_out, to_combine c(U.lactuca, Enteromorpha), new_source_name Algae ) # 绘制合并后的结果 plot(simmr_combined, type boxplot)4.3 食物源贡献比较统计比较不同食物源的相对贡献# 比较Zostera和Grass的贡献 compare_sources(simmr_out, source_names c(Zostera, Grass)) # 输出示例 Prob ( proportion of Zostera proportion of Grass ) 1.005. 高级应用与技巧5.1 多组比较分析对于包含多个组别的数据如不同季节或种群SIMMR可进行组间比较# 加载多组数据 data(geese_data) simmr_multi - with( geese_data, simmr_load( mixtures mixtures, source_names source_names, source_means source_means, source_sds source_sds, correction_means correction_means, correction_sds correction_sds, concentration_means concentration_means, group groups ) ) # 组间比较例如比较Zostera在组1和组2中的贡献 compare_groups(simmr_multi_out, source Zostera, groups 1:2)5.2 先验信息整合当有先验知识时如饲养实验数据可通过simmr_elicit整合# 设置先验分布假设Zostera占50%±8% prior - simmr_elicit( n_sources 4, proportion_means c(0.5, 0.2, 0.2, 0.1), proportion_sds c(0.08, 0.02, 0.01, 0.02) ) # 运行含先验的模型 simmr_prior_out - simmr_mcmc( simmr_in, prior_control list( means prior$mean, sd prior$sd ) )5.3 营养富集因子估计对于缺乏可靠TEF数据的研究可使用饲喂实验数据估计# 假设已知膳食比例p p_matrix - matrix( c(0.4, 0.3, 0.2, 0.1), ncol 4, nrow nrow(simmr_in$mixtures), byrow TRUE ) # 估计TEFs simmr_tef_out - simmr_mcmc_tdf( simmr_in, p p_matrix ) # 获取TEF估计 summary(simmr_tef_out, type quantiles)6. 实际应用中的注意事项食物源选择应包含所有潜在重要食物源遗漏关键源会导致估计偏差同位素基线确保食物源数据与消费者处于相同时空尺度模型假设SIMMR假设线性混合和恒定TEF实际情况可能更复杂结果解释膳食比例是统计学估计需结合生态学知识综合判断# 常见问题检查清单 checklist - data.frame( Issue c(混合物点超出凸包, 高后验相关性, 收敛不佳), Solution c(检查源是否遗漏, 考虑合并相似源, 增加MCMC迭代) ) knitr::kable(checklist)通过本案例我们不仅解答了布兰特鹅吃什么的具体问题更掌握了一套可推广的稳定同位素分析方法。这种技术框架同样适用于鱼类、土壤有机质来源解析等多样化的生态学研究场景。