Last updated: 2017-03-04

Code version: 5d0fa13282db4a97dc7d62e2d704e88a5afdb824


include the most complex concepts required to understand the material.


Suppose we have a logistic regression \(Y_i | X_i \sim Bern(p_i)\) where \[log(p_i/(1-p_i)) = \mu + \theta X_i.\]

We will assume that \(X_i \in {-1,+1}\), and assume priors for \(\mu\) and \(\theta\): \[\mu \sim N(0,100)\] \[\theta \sim N(0,1)\]

For illustration we simulate data where \(\mu=\theta=0\):

x = sample(c(-1,1),1000,replace=TRUE)
y = rbinom(1000,1,0.5)

#b is a vector b=(mu,theta)
#loglikelihood for logistic regression
loglik = function(b){
  eta = b[1]+b[2]*x
  p = exp(eta)/(1+exp(eta))

#b is a vector b=(mu,theta)
log_prior = function(b){
  return(dnorm(b[1],0,10, log=TRUE)+dnorm(b[2],0,1,log=TRUE))

#b is a vector b=(mu,theta)
log_post = function(b){

Let’s compute a 95% CI for \(\theta\). First try a discrete grid

Note: This is still a work in progress.

m = seq(-10,10,length=100)
t = seq(-2,2,length=100)
df = expand.grid(m=m,t=t)
#df = c(df,dplyr::ddply(df,log_post))



Session information

R version 3.3.0 (2016-05-03)
Platform: x86_64-apple-darwin13.4.0 (64-bit)
Running under: OS X 10.10.5 (Yosemite)

[1] en_NZ.UTF-8/en_NZ.UTF-8/en_NZ.UTF-8/C/en_NZ.UTF-8/en_NZ.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] tidyr_0.4.1     dplyr_0.5.0     ggplot2_2.1.0   knitr_1.15.1   
 [5] MASS_7.3-45     expm_0.999-0    Matrix_1.2-6    viridis_0.3.4  
 [9] workflowr_0.3.0 rmarkdown_1.3  

loaded via a namespace (and not attached):
 [1] Rcpp_0.12.5      git2r_0.18.0     plyr_1.8.4       tools_3.3.0     
 [5] digest_0.6.9     evaluate_0.10    tibble_1.1       gtable_0.2.0    
 [9] lattice_0.20-33  shiny_0.13.2     DBI_0.4-1        yaml_2.1.14     
[13] gridExtra_2.2.1  stringr_1.2.0    gtools_3.5.0     rprojroot_1.2   
[17] grid_3.3.0       R6_2.1.2         reshape2_1.4.1   magrittr_1.5    
[21] backports_1.0.5  scales_0.4.0     htmltools_0.3.5  assertthat_0.1  
[25] mime_0.5         colorspace_1.2-6 xtable_1.8-2     httpuv_1.3.3    
[29] labeling_0.3     stringi_1.1.2    munsell_0.4.3   

This site was created with R Markdown