如何修复glm的截距值答案

【问题标题】：How to fix intercept value of glm如何修复glm的截距值
【发布时间】：2014-12-29 11:32:21
【问题描述】：

这是我工作的例子：

data2 = data.frame( X = c(0,2,4,6,8,10),
                    Y = c(300,220,210,90,80,10))
attach(data2)
model <- glm(log(Y)~X)
model

Call:  glm(formula = log(Y) ~ X)

Coefficients:
(Intercept)            X  
     6.0968      -0.2984  

Degrees of Freedom: 5 Total (i.e. Null);  4 Residual
Null Deviance:      7.742 
Residual Deviance: 1.509    AIC: 14.74

我的问题是：

glm 函数上有一个选项可以让我用我想要的值修复截距系数？并预测 x 值？

例如：我希望我的曲线以较高的“Y”值开始 ==> 我想用 log(300) 更改截距

【问题讨论】：

model <- lm(log(Y)~X-1+offset(rep(log(300),nrow(data2))),data2) ; ?predict
我可以在这里假设我有： Y = AX+B , B = log(300), A = model$Coefficients ？ ,我有生物学背景，这一切对我来说都是新的#谢谢

标签： r glm intercept coefficients

【解决方案1】：

您错误地使用了glm(...)，这在 IMO 中是一个比偏移量更大的问题。

最小二乘回归的主要基本假设是响应中的误差呈正态分布，方差不变。如果Y 中的错误是正态分布的，那么log(Y) 肯定不是。因此，虽然您可以在log(Y)~X 上“运行数字”，但结果将没有意义。广义线性建模理论就是为了解决这个问题而发展起来的。所以使用glm，而不是适合log(Y) ~X，你应该适合Y~X和family=poisson。前者适合

log(Y) = b₀ + b₁x

后者适合

Y = exp(b₀ + b₁x)

在后一种情况下，如果Y 中的误差呈正态分布，并且模型有效，则残差将按要求呈正态分布。请注意，这两种方法对于 b₀ 和 b₁ 给出了非常不同的结果。

fit.incorrect <- glm(log(Y)~X,data=data2)
fit.correct   <- glm(Y~X,data=data2,family=poisson)
coef(summary(fit.incorrect))
#               Estimate Std. Error  t value     Pr(>|t|)
# (Intercept)  6.0968294 0.44450740 13.71592 0.0001636875
# X           -0.2984013 0.07340798 -4.06497 0.0152860490
coef(summary(fit.correct))
#               Estimate Std. Error   z value     Pr(>|z|)
# (Intercept)  5.8170223 0.04577816 127.06982 0.000000e+00
# X           -0.2063744 0.01122240 -18.38951 1.594013e-75

尤其是，X 的系数在使用正确的方法时几乎减小了 30%。

注意模型有何不同：

plot(Y~X,data2)
curve(exp(coef(fit.incorrect)[1]+x*coef(fit.incorrect)[2]),
      add=T,col="red")
curve(predict(fit.correct,  type="response",newdata=data.frame(X=x)),
      add=T,col="blue")

正确拟合的结果（蓝色曲线）或多或少随机地通过数据，而错误拟合的结果严重高估了较小的X 的数据，而低估了较大的X 的数据。我想知道这是否就是您要“修复”拦截的原因。查看另一个答案，您可以看到当您修复 Y₀ = 300 时，整个拟合都被低估了。

相比之下，让我们看看当我们正确使用 glm 修复 Y₀ 时会发生什么。

data2$b0 <- log(300)   # add the offset as a separate column
# b0 not fixed
fit <- glm(Y~X,data2,family=poisson)
plot(Y~X,data2)
curve(predict(fit,type="response",newdata=data.frame(X=x)),
      add=TRUE,col="blue")
# b0 fixed so that Y0 = 300
fit.fixed <-glm(Y~X-1+offset(b0), data2,family=poisson)
curve(predict(fit.fixed,type="response",newdata=data.frame(X=x,b0=log(300))),
      add=TRUE,col="green")

这里，蓝色曲线是无约束拟合（正确完成），绿色曲线是约束 Y₀ = 300 的拟合。您可以看到它们差别不大，因为正确（无约束）的拟合已经相当不错了。

【讨论】：

我认为对响应进行对数转换，即假设对数正态分布的残差，不一定一定是错误。
是的，当然。如果您知道错误是对数正态分布的，那么正确的方法是使用lm(log(Y)~X)，而像我一样使用glm(...) 将是不正确的方法。您需要查看 QQ 图，看看哪种方法是正确的。但是看看 OP 的数据 - 使用 poisson glm 的拟合要好得多。
我不这么认为。比较fit.correct <- glm(Y~X,data=data2, family=quasipoisson);coefplot2::coefplot2(list(fit.correct,fit.incorrect),intercept=TRUE)（对不起，coefplot2 不再出现在 CRAN 上）表明系数并没有那么不同。 “小 30%”，但 -0.3 (SE 0.07) 与 -0.2 (SE 0.04) 实际上并没有显着差异......（注意我在这里使用quasipoisson）
这对于评论线程来说有点密集，但我认为 X 的系数是 NSD，因为从拟合到 log(Y)~X 的估计值相对较差（p=0.015），而使用泊松 glm 从拟合到 Y~X 的估计要好得多（plog(Y)~X 拟合给出 SSE=29081。
确实如此，但对于它的价值而言，p<2e-16 完全是由于过度分散——这确实需要一个允许过度分散的模型。

【解决方案2】：

data2 <- data.frame( X = c(0,2,4,6,8,10),
                Y = c(300,220,210,90,80,10)
m1 <- lm(log(Y)~X-1+offset(rep(log(300),nrow(data2))),data2)

有一个predict() 函数，但在这里它可能更容易手动预测。

par(las=1,bty="l")
plot(Y~X,data=data2)
curve(300*exp(coef(m1)*x),add=TRUE)

对于它的价值，如果你想比较 log-Normal 和 Poisson 模型，你可以通过

library("ggplot2")
 theme_set(theme_bw())
 ggplot(data2,aes(X,Y))+geom_point()+
    geom_smooth(method="glm",family=quasipoisson)+
    geom_smooth(method="glm",family=quasi(link="log",var="mu^2"),
      colour="red",fill="red")

【讨论】：