Compare the impact of different smooths on overall fit

hello

I want to estimate the relative impact of different covariates (or groups of covariates) on the overall fit of a brms model. Does the approach below using the change in “log-likelihood” make sense, is it well-behaved in general, is there a better approach?

thanks

Greg

library(brms)
library(ggplot2)
library(gridExtra)

set.seed(1)
x1<-runif(200)
x2<-runif(200)
x3<-runif(200)

f1<-function(x1) sin(8*x1)*exp(-4*x1)
f2<-function(x2) 2*sin(8*x2)*exp(-2*x2)
f3<-function(x3) 0.5*x3

par(mfrow=c(1,3))
plot(f1,ylim=c(-1.5,1.5))
plot(f2,ylim=c(-1.5,1.5))
plot(f3,ylim=c(-1.5,1.5))
par(mfrow=c(1,1))

lp<--0.5+f1(x1)+f2(x2)+f3(x3)
h<-function(t) (1+exp(-t))^-1
p<-h(lp)
T<-rpois(200,50)

y<-rbinom(200,size=T,prob=p)
dat1<-data.frame(x1,x2,x3,T,y)
bf1<-bf(y|trials(T)~s(x1)+s(x2)+s(x3),family="binomial")
mod1<-brm(bf1,family="binomial",data=dat1,file="test1",
cores = getOption("mc.cores", 2),control=list(adapt_delta=0.95))
summary(mod1)

cs1<-conditional_smooths(mod1,smooths="s(x1)")
cs2<-conditional_smooths(mod1,smooths="s(x2)")
cs3<-conditional_smooths(mod1,smooths="s(x3)")
g1<-plot(cs1)[[1]]+ylim(-1.2,1.2)
g2<-plot(cs2)[[1]]+ylim(-1.2,1.2)
g3<-plot(cs3)[[1]]+ylim(-1.2,1.2)
grid.arrange(g1,g2,g3,nrow=1)

mat<-as.matrix(mod1)
dim(mat)
colnames(mat)
sdata<-standata(mod1)
names(sdata)
sdata$K

dfx<-dim(sdata$X)[2]
for (j in 1:dfx)
{
assign("musamj",outer(mat[,j],sdata$X[,j]))
assign(paste0("con",j),apply(musamj,2,mean))
rm(musamj)
}

fir=-1
dsm<-dim(sdata$Xs)[2]
for (j in 1:dsm)
{
assign("musamj",outer(mat[,dfx+j],sdata$Xs[,j])+mat[,fir+8*j+1:8]%*%t(eval(parse(text=paste0("sdata$Zs_",j,"_1")))))
assign(paste0("con",dfx+j),apply(musamj,2,mean))
rm(musamj)
}

conmat<-con1

for (j in 2:(dfx+dsm))
{
conmat<-cbind(conmat,eval(parse(text=paste0("con",j))))
}

dim(conmat)
pfit<-h(apply(conmat,1,sum))

LL<-sum(y*log(pfit)+(T-y)*log(1-pfit))
pfx1<-h(apply(conmat[,-2],1,sum))
pfx2<-h(apply(conmat[,-3],1,sum))
pfx3<-h(apply(conmat[,-4],1,sum))
ddx1<--2*sum(y*log(pfx1/pfit)+(T-y)*log((1-pfx1)/(1-pfit)))
ddx2<--2*sum(y*log(pfx2/pfit)+(T-y)*log((1-pfx2)/(1-pfit)))
ddx3<--2*sum(y*log(pfx3/pfit)+(T-y)*log((1-pfx1)/(1-pfit)))

S<-ddx1+ddx2+ddx3
round(c(ddx1/S,ddx2/S,ddx3/S),3)
#0.087 0.840 0.073
#f1 and f3 are roughly comparable, and f2 has around 10x the impact.

  • Operating System: Linux Mint 19.3 Cinnamon
  • brms Version: 2.17.0