如何在R中构建含交互项的GLM模型，验证活动期相关假设？

一、你的现有代码问题 #

你给出的模型代码有两处明显问题：

1. 重复写入了percip和time变量，冗余变量会干扰模型估计，必须删除重复项；
2. 交互项的设置没有匹配核心假设——你认为高活动期环境变量无影响，低活动期才受环境变量作用，但activity * temp是让所有活动期都与温度产生交互，没有针对性区分不同活动期的效应差异。

二、匹配假设的GLM模型构建 #

你的因变量是计数数据，优先选择泊松GLM；如果数据存在过度离散（残差方差远大于均值），改用负二项GLM（需加载MASS包）。

1. 先设置参考水平 #

把activity的参考组设为"高"，这样模型输出中，中、低活动期的系数都是相对于高活动期的差异，方便直接检验假设：

# 将activity转为因子并设置参考水平
data$activity <- relevel(factor(data$activity), ref = "高")

2. 核心模型代码 #

模型需要包含activity主效应、环境变量（temp、percip、time）主效应，以及activity与每个环境变量的交互项——这样能分别检验不同活动期内环境变量的效应：

# 泊松GLM（基础版本）
glm_pois <- glm(
  formula = counts ~ activity + temp + percip + time + activity:temp + activity:percip + activity:time,
  data = your_data,  # 替换成你的数据集名称
  family = poisson(link = "log")
)

# 若存在过度离散，改用负二项GLM
library(MASS)
glm_nb <- glm.nb(
  formula = counts ~ activity + temp + percip + time + activity:temp + activity:percip + activity:time,
  data = your_data
)

注：activity:temp表示仅保留交互项，若想简化写法，也可以用activity*temp（等价于activity + temp + activity:temp），但手动列出所有项更清晰。

3. 对应假设的检验方法 #

• 假设1：月份间差异无法用环境变量解释：
  拟合仅含activity的空模型，和上面的全模型做似然比检验：

  glm_null <- glm(counts ~ activity, data = your_data, family = poisson)
  anova(glm_null, glm_pois, test = "Chisq")
  

  如果检验的p值不显著，说明环境变量无法解释月份间的差异，支持假设1。

• 假设2：高活动期内环境变量不影响计数：
  查看summary(glm_pois)的输出，temp、percip、time的主效应系数对应的p值如果都不显著，说明高活动期（参考组）内这些环境变量对计数无影响，支持假设2。

• 假设3：低活动期内计数受环境变量影响：
  查看activity低:temp、activity低:percip、activity低:time的系数p值，若其中任意一个显著不为0，说明低活动期内对应的环境变量对计数有影响，支持假设3。

三、R中运行含交互项与非交互项GLM的通用步骤 #

1. 数据预处理：确保分类变量转为因子（factor()），连续变量为数值型，无缺失值（或用na.omit()处理）；
2. 选择模型族：根据因变量类型选对应的family参数：
   • 计数数据：poisson()或negbin()（负二项）；
   • 二分类数据：binomial(link = "logit")；
   • 连续数据：gaussian()（默认）；
3. 构建模型公式：
   • 非交互项用+连接，如y ~ x1 + x2；
   • 交互项两种写法：
     • x1*x2：自动包含x1、x2的主效应和交互项x1:x2；
     • x1 + x2 + x1:x2：手动列出主效应和交互项，和x1*x2效果一致；
4. 拟合与诊断：用glm()拟合模型，summary()看系数显著性，plot()做残差诊断，anova()做模型比较。

示例（用模拟数据）：

# 生成模拟数据
set.seed(123)
sim_data <- data.frame(
  counts = rpois(150, lambda = c(20, 10, 5)[sample(1:3, 150, replace = TRUE)]),
  activity = sample(c("高", "中", "低"), 150, replace = TRUE),
  temp = rnorm(150, 25, 4),
  percip = rnorm(150, 60, 15)
)

# 设置参考水平
sim_data$activity <- relevel(factor(sim_data$activity), ref = "高")

# 拟合含交互项的GLM
model <- glm(counts ~ activity + temp + percip + activity:temp + activity:percip,
             data = sim_data,
             family = poisson)

# 查看结果
summary(model)

内容的提问来源于stack exchange，提问作者Benjamin Colbert