# If you haven't already done it, you'll need to install the ggplot2 library and all its dependencies before # running the next command. library(ggplot2) ###################################################################################################### # A thorough understanding of ggplot is beyond what I can show here, but here are a couple of examples # that may prove useful. The first is a plot of myr1 from the simple predictors example. ###################################################################################################### # Note, EVERY variable in the model being plotted must be included in the newdat grid as well as referenced # in the ggplot commands below. newdat <- expand.grid( logMagRatio = c(0, .5, 1, 1.5, 2, 2.5) , logDelayRatio =c(.8873) , Waited = 0 ) modelToPlot<-myr1 # The following code was obtained on-line to extract the correct error bars (plo/phi) and put them in newdat. mm <- model.matrix(terms(modelToPlot),newdat) newdat$Waited <- mm %*% fixef(modelToPlot) pvar1 <- diag(mm %*% tcrossprod(vcov(modelToPlot),mm)) newdat <- data.frame( newdat , ProbWaited = inv.logit(newdat$Waited) , plo = inv.logit(newdat$Waited-sqrt(pvar1)) , phi = inv.logit(newdat$Waited+sqrt(pvar1)) ) newdat # if you want to look at what we're going to plot. # Plotting behavior as a function of logMagRatio. Note, the LLdelay is three times longer than the SSdelay # so it isn't surprising that when the LLmagn = SSmagn (log(MagRatio) = 0) that people were only waiting 25% of the time g1 <- ggplot(newdat, aes(x=logMagRatio, y= ProbWaited))+geom_line()+geom_point()+ ylim(0,1) g1 + geom_errorbar(aes(ymin = plo, ymax = phi), linetype=1, width=.02)+ xlab("Log(LLMagnitude/SSMagnitude)")+ylab("P(Waited)\n")+theme_bw()+theme(text = element_text(size=18)) ###################################################################################################### # Now, we're going to plot myr2.1 from the simple predictions example. ###################################################################################################### newdat <- expand.grid( logMagRatio = c(0, .5, 1, 1.5, 2, 2.5) , logDelayRatio =c(.8873) , c.Sex = c(-1, 1) , Waited = 0 ) modelToPlot<-myr2.1 mm <- model.matrix(terms(modelToPlot),newdat) newdat$Waited <- mm %*% fixef(modelToPlot) pvar1 <- diag(mm %*% tcrossprod(vcov(modelToPlot),mm)) newdat <- data.frame( newdat , ProbWaited = inv.logit(newdat$Waited) , plo = inv.logit(newdat$Waited-sqrt(pvar1)) , phi = inv.logit(newdat$Waited+sqrt(pvar1)) ) newdat$Sex <- ifelse(newdat$c.Sex==1, "Female", "Male") newdat$logDelay <- as.factor(newdat$logDelay) # Plotting behavior as a function of MagRatio g1 <- ggplot(newdat, aes(x=logMagRatio, y= ProbWaited, group=Sex, linetype=Sex))+geom_line()+geom_point()+ ylim(0,1) g1 + geom_errorbar(aes(ymin = plo, ymax = phi), linetype=1, width=.02)+ xlab("Log(LLMagnitude/SSMagnitude)")+ylab("P(Waited)\n")+theme_bw()+theme(text = element_text(size=18)) + facet_grid(. ~ logDelay) # Repeat of above but using position_dodge to avoid error bar overlap. pd <- position_dodge(width = 0.06) g1 <- ggplot(newdat, aes(x=logMagRatio, y= ProbWaited, group=Sex, linetype=Sex))+geom_line(position = pd)+geom_point(position = pd)+ ylim(0,1) g1 + geom_errorbar(aes(ymin = plo, ymax = phi), linetype=1, width=.02, position = pd)+ xlab("Log(LLMagnitude/SSMagnitude)")+ylab("P(Waited)\n")+theme_bw()+theme(text = element_text(size=18)) + facet_grid(. ~ logDelay) ###################################################################################################### # Here I'm plotting the results from myr1 from the Example one group.R file. ###################################################################################################### # Problem arises in creating newdat with the complex predictor names, so redoing the analysis. myd<-read.csv("Example one group data glmer.dat") myd$logMagRatio <- log(myd$LLmagn/myd$SSmagn) myd$logDelay <- log(myd$LLdelay+1) myr1<-glmer(Waited~logMagRatio+logDelay-1+(logMagRatio+logDelay-1|Subject), data=myd, family=binomial) newdat <- expand.grid( logMagRatio = c(0, .1, .2, .3, .4, .5, .6, .7, .8, .9, 1) , logDelay =c(2, 3.5, 5) , Waited = 0 ) modelToPlot<-myr1 mm <- model.matrix(terms(modelToPlot),newdat) newdat$Waited <- mm %*% fixef(modelToPlot) pvar1 <- diag(mm %*% tcrossprod(vcov(modelToPlot),mm)) newdat <- data.frame( newdat , ProbWaited = inv.logit(newdat$Waited) , plo = inv.logit(newdat$Waited-sqrt(pvar1)) , phi = inv.logit(newdat$Waited+sqrt(pvar1)) ) # Plotting behavior as a function of logMagRatio. Note that the beta for logMagRatio is much greater than # that for logDelay pd <- position_dodge(width = 0.06) g1 <- ggplot(newdat, aes(x=logMagRatio, y= ProbWaited, group=logDelay, colour=logDelay))+geom_line(position = pd)+geom_point(position = pd)+ ylim(0,1) g1 + geom_errorbar(aes(ymin = plo, ymax = phi), linetype=1, width=.02, position = pd)+ xlab("Log(LLMagnitude/SSMagnitude)")+ylab("P(Waited)\n")+theme_bw()+theme(text = element_text(size=18)) # Plotting same results in a different way g2 <- ggplot(newdat, aes(x=logMagRatio, y= ProbWaited))+geom_line(position = pd)+geom_point(position = pd)+ ylim(0,1) # Want to break it out by a factor, so need to change logDelay into a factor newdat$logDelay <- factor(newdat$logDelay) g2 + geom_errorbar(aes(ymin = plo, ymax = phi), linetype=1, width=.02, position = pd)+ xlab("Log(LLMagnitude/SSMagnitude)")+ylab("P(Waited)\n")+theme_bw()+theme(text = element_text(size=18))+facet_grid(logDelay~.) # I don't think in terms of logs, so let's create new columns with the backtransformed values newdat$MagRatio <- exp(newdat$logMagRatio) newdat$Delay <- exp(as.numeric(as.character(newdat$logDelay))) - 1 # complicated because I changed it to a factor earlier newdat$Delay <- as.factor(round(newdat$Delay)) # removes ugly decimal places g3 <- ggplot(newdat, aes(x=MagRatio, y= ProbWaited))+geom_line(position = pd)+geom_point(position = pd)+ ylim(0,1) g3 + geom_errorbar(aes(ymin = plo, ymax = phi), linetype=1, width=.02, position = pd)+ xlab("LLMagnitude/SSMagnitude")+ylab("P(Waited)\n")+theme_bw()+theme(text = element_text(size=18))+facet_grid(Delay~.)