R语言优雅的进行广义可加模型泊松回归分析

R语言优雅的进行广义可加模型泊松回归分析

引用

CSDN

1.

https://blog.csdn.net/dege857/article/details/140666143

泊松回归是一种用于分析计数数据的回归方法，在许多领域都有广泛的应用。本文将介绍如何使用R语言中的tidygam包进行广义可加模型（GAM）的泊松回归分析。通过一个具体的婴儿手势数据集，详细展示了模型建立、预测值生成以及可视化图表的绘制过程。

数据准备

首先，我们需要导入所需的R包和数据。数据使用的是tidygam自带的gest数据集，这是一个关于10、11和12个月婴儿手势数量的数据集。

library(tidygam)
library(mgcv)
library(dplyr)
library(ggplot2)
data("gest")
theme_set(theme_light())

这个数据表包含了来自孟加拉、中国和英国背景的60名婴儿所做的三种手势的计数。变量说明如下：

• dyad：ID标识
• background：表明来自哪个国家
• months：月份，只有10，11，12共3个月
• gesture：手势的类型
• count：手势的计数，这是结局变量

基本模型

我们先建立一个基本的GAM模型，分析月份与手势计数的关系。

gs <- gam(
  count ~ s(months, k = 3),
  data = gest,
  family = poisson
)

生成预测值：

gs_pred <- predict_gam(gs)

绘制图形：

gs_pred %>%
  plot(series = "months")

由于生成的Y轴是count的对数值，我们需要转换一下：

predict_gam(gs, tran_fun = exp) %>%
  plot(series = "months")

按国家分类的模型

如果想了解不同国家的差异，可以在模型中加入background变量：

gs_by <- gam(
  count ~ s(months, by = background, k = 3),
  data = gest,
  family = poisson
)

绘制分类图形：

gs_by %>%
  predict_gam(length_out = 20, series = "months", tran_fun = exp) %>%
  plot(comparison = "background")

可以使用ggplot2的方法进一步美化图形：

gs_by %>%
  predict_gam(length_out = 20, series = "months", tran_fun = exp) %>%
  plot(comparison = "background") +
  scale_color_brewer(type = "qual") + scale_fill_brewer(type = "qual")

包含两个非连续变量的模型

如果考虑两个非连续变量的影响，可以在模型中定义：

gs_by_2 <- gam(
  count ~ s(months, by = background, k = 3) +
    s(months, by = gesture, k = 3),
  data = gest,
  family = poisson
)

绘制图形：

gs_by_2 %>%
  predict_gam(length_out = 20, series = "months", tran_fun = exp) %>%
  plot(comparison = "gesture") +
  scale_color_brewer(type = "qual") + scale_fill_brewer(type = "qual") +
  facet_grid(~ background)

交互效应模型

如果想了解background和gesture的交互关系，需要先生成一个交互变量：

gest <- gest %>%
  mutate(back_gest = interaction(background, gesture))

然后使用这个交互变量进行模型分类：

gs_i <- gam(
  count ~ s(months, by = back_gest, k = 3),
  data = gest,
  family = poisson
)
summary(gs_i)

绘制交互效应图形：

predict_gam(
  gs_i, tran_fun = exp,
  separate = list(back_gest = c("background", "gesture"))
) %>%
  plot(series = "months", comparison = "gesture") +
  facet_grid(~ background)

这样一个亚组交互图形就生成好了。