我正在尝试对比例数据(0到1之间的数字(运行GAM。但我收到了的警告
In eval(family$initialize) : non-integer #successes in a binomial glm!
基本上,我正在根据海面温度对温暖适应物种的出现次数与温暖和寒冷适应物种的总出现次数进行建模,并使用另一个天气系统(NAO(的数据作为随机效应,以及其他三个分类、参数变量。
m5 <- gam(prop ~ s(SST_mean) + s(NAO, bs="re") + WarmCold + Cycle6 + Region,
family=binomial, data=DAT_WC, method = "REML")
prop=发生的比例,WarmCold=物种是温暖适应还是寒冷适应,周期6=6年时间段,区域=4个区域之一。我的数据集样本位于下方
structure(list(WarmCold = structure(c(1L, 1L, 1L, 1L, 2L, 2L), .Label = c("Cold",
"Warm"), class = "factor"), Season = structure(c(2L, 2L, 2L,
2L, 2L, 2L), .Label = c("Autumn", "Spring", "Summer", "Winter"
), class = "factor"), Region = structure(c(1L, 2L, 3L, 4L, 1L,
2L), .Label = c("OSPARII_N", "OSPARII_S", "OSPARIII_N", "OSPARIII_S"
), class = "factor"), Cycle6 = structure(c(1L, 1L, 1L, 1L, 1L,
1L), .Label = c("1990-1995", "1996-2001", "2002-2007", "2008-2013",
"2014-2019"), class = "factor"), WC.Strandings = c(18L, 10L,
0L, 3L, 5L, 25L), SST_mean = c(7.4066298185553, 7.49153086390094,
9.28247524767124, 10.8654859624361, 7.4066298185553, 7.49153086390094
), NAO = c(0.542222222222222, 0.542222222222222, 0.542222222222222,
0.542222222222222, 0.542222222222222, 0.542222222222222), AMO = c(-0.119444444444444,
-0.119444444444444, -0.119444444444444, -0.119444444444444, -0.119444444444444,
-0.119444444444444), Total.Strandings = c(23, 35, 5, 49, 23,
35), prop = c(0.782608695652174, 0.285714285714286, 0, 0.0612244897959184,
0.217391304347826, 0.714285714285714)), row.names = c(NA, 6L), class = "data.frame")
从文献(Zuur,2009(来看,二项式分布似乎是用于比例数据的最佳方法。但它似乎不起作用。它正在运行,但发出了上述警告,并输出了没有意义的内容。我在这里做错了什么?
这是一个警告,而不是错误,但它确实表明有些地方不正确;二项式分布支持非负整数值,因此在没有形成比例的样本总数的情况下传递非整数值是没有意义的。
您可以使用weights
参数来完成此操作,在本例中,该参数应采用一个整数向量,该向量包含计算比例的每个观测值的计数总和。
或者,如果平均方差关系适合您的数据,可以考虑使用family = quasibinomial
;气候变暖会消失,但届时你将无法使用AIC和相关工具来预测真正的可能性。
如果比例为真比例,则考虑family = betar
以拟合贝塔回归模型,其中响应的条件分布支持单位区间(0,1(上的实数值(但从技术上讲不是0或1;如果响应中有0或1个值,mgcv将加或减一个小数字来调整数据(。
我还发现,不是计算总数,而是使用cbind((和2列感兴趣的列删除了警告,例如
m8 <- gam(cbind(WC.Strandings, Total.Strandings) ~ s(x1) + x2,
family=binomial(link="logit"), data=DAT, method = "REML")