title: “EG_Happiness” output: html_document

Load Packages

library(arm); library(lmerTest); library(psych)
## Loading required package: MASS
## Loading required package: Matrix
## Loading required package: lme4
## Loading required package: Rcpp
## 
## arm (Version 1.7-03, built: 2014-4-27)
## 
## Working directory is /private/var/folders/9c/88_yf73x0ys02zvn5hkvvq6r0000gn/T/RtmpmxWTU5
## 
## KernSmooth 2.23 loaded
## Copyright M. P. Wand 1997-2009
## 
## Attaching package: 'lmerTest'
## 
## The following object is masked from 'package:lme4':
## 
##     lmer
## 
## The following object is masked from 'package:stats':
## 
##     step
## 
## 
## Attaching package: 'psych'
## 
## The following objects are masked from 'package:arm':
## 
##     logit, rescale, sim

Loading Data

setwd("/Volumes/TOSHIBA EXT/Dropbox/Schools Study Data/Emily Griffith Data and R Scripts")

data <- read.csv("EmilyGriffith_all.csv")

Creating School ID as ID

data$ID <- data$Q1

#Create scale scores
data$meanHAPPI <- apply(data[, c("HAPPI.1_1", "HAPPI.2_1", "HAPPI.3_1", "HAPPI.4_1")], 1, mean, na.rm = TRUE)
#Means or plotting
data$baseline <- ifelse(data$Time<4,0,1) 
pdata <- tapply(data[,"meanHAPPI"], data[,3], mean, na.rm=TRUE)
plot(pdata, type="l")

plot of chunk unnamed-chunk-4 Baseline Model

M0 <- lmer(meanHAPPI ~ 1 + (1|ID), data=data)
fixef(M0)
## (Intercept) 
##       4.904

Confidence Intervals

confint(M0)
## Computing profile confidence intervals ...
##              2.5 % 97.5 %
## .sig01      0.0000 0.5517
## .sigma      0.8642 1.1717
## (Intercept) 4.7073 5.1018

Model 1

M1 <- update(M0, .~. + Time, REML=FALSE)
fixef(M1)
## (Intercept)        Time 
##     5.09514    -0.08583

Confidence Intervals for Model 1

confint(M1)
## Computing profile confidence intervals ...
## Warning: convergence code 3 from bobyqa: bobyqa -- a trust region step
## failed to reduce q
##               2.5 % 97.5 %
## .sig01       0.0000 0.5508
## .sigma       0.8612 1.1678
## (Intercept)  4.6210 5.5694
## Time        -0.2801 0.1083

Model 2

M2 <- update(M1, .~. + baseline)
fixef(M2)
## (Intercept)        Time    baseline 
##      5.4187     -0.2730      0.8164
confint(M2)
## Computing profile confidence intervals ...
##                2.5 %   97.5 %
## .sig01       0.00000  0.54475
## .sigma       0.84421  1.14613
## (Intercept)  4.86356  5.97426
## Time        -0.53250 -0.01427
## baseline     0.05271  1.58209

Model 3

M3 <- update(M2, .~. + I(Time^2))
fixef(M3)
## (Intercept)        Time    baseline   I(Time^2) 
##      7.8677     -3.1113     -1.6384      0.7103

Confidence Intervals for Model 3

confint(M3)
## Computing profile confidence intervals ...
##               2.5 %  97.5 %
## .sig01       0.0000  0.6077
## .sigma       0.7620  1.0497
## (Intercept)  6.5266  9.1947
## Time        -4.5463 -1.6586
## baseline    -3.0464 -0.2130
## I(Time^2)    0.3519  1.0641
Pdata <- tapply(data[,"meanHAPPI"], data[,3], mean, na.rm=TRUE)
# Add random noise to time to better see the points of interest
data$TimeJIT <- data$Time+runif(126, min=-.1, max=.1)
with(data, plot(TimeJIT, meanHAPPI, col="grey", pch="*"))
lines(pdata, col="red", lwd=2)

plot of chunk unnamed-chunk-11

Intervention was given between 3 and 4 (Time 1 is not relevant because the data is not clean)

boxplot(data$meanHAPPI~Time, data=data, notch=F, col=(c("red","blue", "green")), main="HAPPINESS", xlab="HAPPINESS", xlim = c(), ylim = c(2, 9 ), yaxs = "i")                                                                                                 

plot of chunk unnamed-chunk-12 ```

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.