Displaying 1 result from an estimated 1 matches for "logitll".
Did you mean:
logitel
2007 Jun 01
0
Metropolis code help
...t(y)
if (log(runif(1)) < newlogalpha - oldlogalpha) {
sims[i,] <- y
oldlogalpha <- newlogalpha
accepts <- accepts + 1
} else {
sims[i,] <- sims[i-1,]
}
}
cat('Accepted ',100*accepts/(R-1),'%\n')
sims
}
# Use the binomial likelihood
logitll=function(beta,y,X)
{
X<- cbind(1,DF$nsaid,DF$diuretic,DF$diuretic*DF$nsaid)
y<- DF$Var3
lF1=plogis(X%*%as.vector(beta),log.p=TRUE)
lF2=plogis(-X%*%as.vector(beta),log.p=TRUE)
sum(y*lF1+(1-y)*lF2)
}
# Use a uniform prior for p
logprior <- function(beta,y,X) 0
# The log posterior is the su...