如何按站点分类运行litterfitter模型并提取参数？

你遇到的错误核心是嵌套后未正确引用每个站点的子集数据：

• nest(-site_code)后生成的是包含site_code和data两列的数据框，每个data单元格对应一个站点的子集数据
• 你的代码里用map(decomp_test, ...)，相当于遍历原数据集的3列，而非6个站点的子集，导致fit列长度与数据框行数不匹配
• 同时fit_litter里直接调用decomp_test$days_between，用的是整个数据集的数据，没有按站点拆分

library(litterfitter)
library(tidyverse)

decomp_test <- structure(list(site_code = c("CCPp1a", "CCPp1a", "CCPp1a", "CCPp1a", 
"CCPp1a", "CCPp1b", "CCPp1b", "CCPp1b", "CCPp1b", "CCPp1b", "CCPp1c", 
"CCPp1c", "CCPp1c", "CCPp1c", "CCPp1c", "CCPp1d", "CCPp1d", "CCPp1d", 
"CCPp1d", "CCPp1d", "CCPp1e", "CCPp1e", "CCPp1e", "CCPp1e", "CCPp1e", 
"CCPp1f", "CCPp1f", "CCPp1f", "CCPp1f", "CCPp1f"), days_between = c(0L, 
118L, 229L, 380L, 572L, 0L, 118L, 229L, 380L, 572L, 0L, 118L, 
229L, 380L, 572L, 0L, 118L, 229L, 380L, 572L, 0L, 118L, 229L, 
380L, 572L, 0L, 118L, 229L, 380L, 572L), mass_remaining = c(1, 
0.7587478816, 0.7366473295, 0.6038150404, 0.6339366063, 1, 0.7609346914, 
0.7487194938, 0.7336179508, 0.6595702348, 1, 0.777213425, 0.734006734, 
0.6963752241, 0.5827854154, 1, 0.7716566866, 0.7002094345, 0.6913555798, 
0.7519095328, 1, 0.7403565314, 0.6751289171, 0.6572164948, 0.620339994, 
1, 0.8126440236, 0.7272999401, 0.7223268259, 0.6805293006)), row.names = c(NA, 
-30L), class = "data.frame")

# 按站点拟合离散平行模型，并提取参数、AIC等信息
discrete_parallel <- decomp_test %>%
  # 按site_code嵌套，每个站点对应独立数据集
  nest(data = -site_code) %>%
  mutate(
    # 用possibly包裹，避免单个站点拟合失败中断全流程
    fit = map(data, possibly(~ fit_litter(
      time = .x$days_between,
      mass.remaining = .x$mass_remaining,
      model = 'discrete.parallel',
      iters = 1000
    ), otherwise = NA)),
    # 提取模型参数（离散平行模型含p、k1、k2三个参数）
    coefs = map(fit, ~ if(!is.na(.x)) coef(.x) else NA),
    # 提取AIC值
    aic = map_dbl(fit, ~ if(!is.na(.x)) .x$fitAIC else NA),
    # 提取模型名称
    model_type = map_chr(fit, ~ if(!is.na(.x)) .x$model else NA)
  ) %>%
  # 将参数列表拆分为单独列
  unnest_wider(coefs, names_sep = "_")

# 查看结果
print(discrete_parallel)

• 引用子集数据：用map(data, ...)遍历每个站点的子集，lambda函数中用.x指代当前子集，调用.x$days_between和.x$mass_remaining传入模型
• 容错处理：possibly()包裹拟合函数，避免单个站点拟合失败导致整个流程终止
• 参数提取：离散平行模型返回的coef是包含p（快速分解比例）、k1（快速分解速率）、k2（慢速分解速率）的列表，用unnest_wider()拆分为单独列
• 指标提取：直接从拟合对象中提取fitAIC和model属性，用map_dbl处理数值型AIC，map_chr处理字符型模型名称

• 若要对比多个模型，可使用crossing()生成站点+模型的组合，再批量拟合
• 若需查看拟合预测值，可在mutate中添加preds = map(fit, ~ if(!is.na(.x)) .x$predicted else NA)

内容的提问来源于stack exchange，提问作者Daniel Fishburn