daily pastebin goal
69%
SHARE
TWEET

Untitled

simondp May 30th, 2018 139 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. ### This script calculates the probability of a chimpanzee pulling the left lever
  2. ### Using rstanarm and leave one out cv (on actor level)
  3.  
  4. # Import libraries and data
  5. library(rethinking)
  6. library(rstanarm)
  7. library(loo)
  8. options(mc.cores = parallel::detectCores())
  9. data(chimpanzees)
  10. d <- chimpanzees
  11.  
  12. # How many actors?
  13. nractors <- length(unique(d$actor))
  14.  
  15. ### Define null model
  16. m0 <- stan_glmer(pulled_left~+(1|actor),
  17.                  data=d, family=binomial)
  18.  
  19. # Define model with one predictor
  20. m1 <- update(m0,.~.+condition)
  21.  
  22. ### Compare using "classic loo"
  23. loo_m0 <- loo(m0)
  24. loo_m1 <- loo(m1)
  25. compare(loo_m0,loo_m1)
  26.  
  27. ### Compare using grouped k-fold cv
  28. folds <- kfold_split_stratified(K = nractors, x = d$actor)
  29. kf_m0 <- kfold(m0, K = nractors, folds = folds)
  30. kf_m1 <- kfold(m1, K = nractors, folds = folds)
  31. compare(kf_m0,kf_m1)
  32.  
  33. ### Compare respecting grouping structure
  34. m0par<- matrix(NA,nractors,1)
  35. m1par<- matrix(NA,nractors,1)
  36. for(i in 1:nractors){
  37.   m0par[i] <- sum(kf_m0$pointwise[(folds==i),])
  38.   m1par[i] <- sum(kf_m1$pointwise[(folds==i),])
  39. }
  40.  
  41. round(sum(m1par-m0par),2) # Diff
  42. round(sd(m1par-m0par),2) # se
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top