## READ DATA from file!

library(lme4)

fruitdat1 <- read.csv("../DATA/wildflower_pods_complete.csv")

sessionInfo()

fruitdat1$State = factor(fruitdat1$State, levels = c("F", "L", "M", "S", "D"))
t1 <- system.time(nfit <- glmer(toF ~ State + Pods + (1|PlantID) + (-1 + Pods|PlantID) + (1|Year), family = binomial, data = fruitdat1))["elapsed"]

if (class(nfit)=="try-error")  {  # or if (inherits(nfit,"try-error"))
  results <- list(obj=NA,coef=rep(NA,4),
                  sd=rep(NA,4))
} else {
  ## manipulate VarCorr to get properly named results
  dropvc <- unlist(lapply(VarCorr(nfit),
                function(x) {
                  s <- drop(x)
                  names(s) <- rownames(x)
                  s
                }))
  coef.vc <- dropvc[c("PlantID.(Intercept)","PlantID.Pods","Year.(Intercept)")]
  coefvec <- c(fixef(nfit),coef.vc)
  sdvec <- c(coef(summary(nfit))[,2],rep(NA,length(coef.vc)))
  results <- list(obj=-logLik(nfit),  ## objective function/Neg log likelihood (NULL for MCMC)
                  coef=coefvec,
                sd=sdvec, ## no 
                vcov=vcov(nfit),
                  confint.quad=cbind(coefvec-1.96*sdvec,
                    coefvec+1.96*sdvec),
 #               confint.profile=confint(nfit),
               time=c(fit=unname(t1)))
}

save("results",file="fit.RData")
