用户定义的“负指数”链接glm-Java 学习之路

我试着按照这个例子modify glm... user specificed link function in r但是收到错误 . 我有二进制数据，并希望将链接功能从"logit"更改为负指数链接 . 我想预测一下

成功概率（p）= 1-exp（线性预测因子）

我需要此链接而不是其中一个内置链接的原因是p以0和0.5之间的凸起方式增加，但"logit"，"cloglog"，"probit"和"cauchy"仅允许凹形 . 参见附图以供参考：predicted p vs binned observations

模拟数据

location<-as.character(LETTERS[rep(seq(from=1,to=23),30)])
success<-rbinom(n=690, size=1, prob=0.15)
df<-data.frame(location,success)
df$random_var<-rnorm(690,5,3)
df$seedling_size<-abs((0.1+df$success)^(1/df$random_var))
df<-df[order(df$location)]

创建自定义链接功能 . 注意：eta =线性预测值，mu =概率

negex<-function(){
##link
linkfun<-function(mu) log(-mu+1)
linkinv<-function(eta) 1-exp(eta)
## derivative of inverse link with respect to eta
mu.eta<-function(eta)-exp(eta)
valideta<-function(eta) TRUE
link<-"log(-mu+1)"
structure(list(linkfun=linkfun,linkinv=linkinv,
             mu.eta=mu.eta,valideta=valideta,
             name=link),
        class="link-glm")
}

模型成功作为幼苗大小的函数

negexp<-negex()

model1<-glm(success~seedling_size,family=binomial(link=negexp),data=df)

错误：未找到有效的系数集：请提供起始值

使用glmer的模型（我的终极目标）

model2<-glmer(success~seedling_size+ (1|location),family=binomial(link=negexp),data=df)

错误（函数（fr，X，reTrms，family，nAGQ = 1L，verbose = 0L，maxit = 100L，:( maxstephalfit）PIRLS步长减半无法减少pwrssUpdate中的偏差

我得到不同的错误消息，但我认为无论是使用glmer还是glm，问题都是一样的，那就是我的链接功能在某种程度上是错误的 .

1 回答

我找到了答案 . 最有帮助的是R thread from 2016 . 有2个问题 . 首先，我的链接功能是错误的 . 我把它修改为：

negex <- function() 
 { 
 linkfun <- function(mu) -log(1-mu) 
 linkinv <- function(eta) 1-exp(-eta) 
 mu.eta <- function(eta) exp(-eta) 
 valideta <- function(eta) all(is.finite(eta)&eta>0) 
 link <- paste0("negexp") 
 structure(list(linkfun = linkfun, linkinv = linkinv, 
             mu.eta = mu.eta, valideta = valideta, name = link), 
        class = "link-glm") 
}

其次，该模型需要特定的起始值 . 这些将是您的数据所特有的 . 以下是我实际找到解决方案的数据的前几行：

site plot sub_plot oak_success oak_o1_gt05ft..1
  0001   10        3           1                0
  0001   12        2           0                0
  0001   12        3           0                0
  0001   12        4           0                0
  0001   13        4           0                0

我不知道如何将完整数据发布到这个网站，但如果有人想要它运行这个例子，请给我发一封电子邮件：lake.graboski@gmail.com

negexp<-negex()

希望这有助于未来的某些人，因为我发现没有其他的例子可以在堆栈溢出或在线上解决 . 使用新的起始值，我能够运行模型：

starting_values<-c(1,0) #1 for the intercept and 0 for the slope
h_gt05_solo_negex2<-glm(oak_success~ oak_o1_gt05ft..1 , 
                    family=binomial(link=negexp),start=starting_values,data=rocdf)
summary(h_gt05_solo_negex2)
Call:
glm(formula = oak_success ~ oak_o1_gt05ft..1, family = binomial(link = negexp), 
data = lt40, start = starting_values)

Deviance Residuals: 
Min       1Q   Median       3Q      Max  
-1.3808  -0.4174  -0.2637  -0.2637   2.5985  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)      0.034774   0.005484   6.341 2.28e-10 ***
oak_o1_gt05ft..1 0.023253   0.002187  10.635  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 1416.9  on 2078  degrees of freedom
Residual deviance: 1213.5  on 2077  degrees of freedom
AIC: 1217.5

Number of Fisher Scoring iterations: 6

收敛存在一些问题 . 由于幼苗高度（oak_o1_gt05ft..1）高于40英尺，参数估计值变得不可靠convergence issues . 我在这个范围内的观察很少，所以我将数据限制在预测值<40英尺的观察值并重新运行模型 . 我还包括"site"（与模拟数据中的"location"相同）） . 您在this figure中看到的是橡树成功的预测，关于每个地点/位置的橡树幼苗高度（黑色圆圈），成功/样本（大绿点）的分类观察以及没有场地因素的成功概率的预测（蓝线） . 当考虑遗址时，看起来幼苗大小变量的斜率更准确 .

不幸的是，我无法让这个模型在glmer中运行，因此必须将站点作为固定效果包含在内，因此，橡树幼苗高度的标准误差和坡度估计可能有点保守 .

回复于 2024-05-06T22:54:06+08:00

用户定义的“负指数”链接glm

1 回答

相关问题