Loading the dataset

data.test4 <- read.csv("~/Final_Adult_Study_R_Docs/adult_study011615.csv")
# Load the psych package
library(psych)
items <- grep("GRIT[0-9]*", names(data.test4), value=TRUE)
scaleKey <- c(-1,-1,-1,-1,1,1,1,1)
data.test4[,items] <- apply(data.test4[,items], 2, as.numeric)
data.test4$meanGRIT <- scoreItems(scaleKey, items = data.test4[, items], delete = FALSE)$score

library(reshape2); library(car); library(Amelia);library(mitools);library(nlme);library(predictmeans)
## 
## Attaching package: 'car'
## 
## The following object is masked from 'package:psych':
## 
##     logit
## 
## Loading required package: Rcpp
## ## 
## ## Amelia II: Multiple Imputation
## ## (Version 1.7.3, built: 2014-11-14)
## ## Copyright (C) 2005-2015 James Honaker, Gary King and Matthew Blackwell
## ## Refer to http://gking.harvard.edu/amelia/ for more information
## ## 
## Loading required package: lme4
## Loading required package: Matrix
## 
## Attaching package: 'lme4'
## 
## The following object is masked from 'package:nlme':
## 
##     lmList
#Remove the meanGRIT and ID Group and wave from data.test4 and create a new #dataset with only those variables.
data <- data.test4[,c("ID", "GROUP", "wave", "meanGRIT")]
#Use dcast to cnage from long-format data to wide format data
data <- dcast(data, ID + GROUP ~ wave, mean, value.var = "meanGRIT")
# Changing all NaNs to NA
data[,3:5] <- apply(data[,3:5],2,function(x) recode(x, "NaN = NA") )

Unsing the mapply function we create a new data set with ID Group baseline meanGRIT and wave 2 and 3 of meanGRIT. So we have a Baseline, which is Time 1 (placed in column 3 one on top of the other) to compare to both Time 2 and 3 (placed in column 4 one on top of the other). In the next line of code we then create a separate column called “wave” which calls the first 89 (which compares Time 2 to Baseline) “wave 1” and then the second 89 we call “wave 2” which compares Time 3 to Baseline. In the third line of code we add names to the new columns of the dataset.

data2 <- as.data.frame(mapply(c,data[,1:4], data[,c(1:3,5)]))
data2$wave <- rep(1:2, each=89)
names(data2) <- c("ID", "GROUP", "BASELINE", "meanGRIT", "WAVE")

Intention to treat model (ITT) where we keep the cases who dropped out and did not complete the study (http://en.wikipedia.org/wiki/Intention-to-treat_analysis).

data2[which(data2$GROUP ==2), "GROUP"] <- 1

Make GROUP and ID a factor

data2$GROUP <-as.factor(data2$GROUP)
data2$ID <-as.factor(data2$ID)

Imputing missing data

MI <- amelia(data2, 50, idvars = c("ID"), ords = "GROUP")
## -- Imputation 1 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 2 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 3 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 4 --
## 
##   1  2  3  4  5  6  7  8  9 10
## 
## -- Imputation 5 --
## 
##   1  2  3  4  5  6
## 
## -- Imputation 6 --
## 
##   1  2  3  4  5  6  7  8  9 10 11
## 
## -- Imputation 7 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 8 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 9 --
## 
##   1  2  3  4  5  6  7  8  9 10 11
## 
## -- Imputation 10 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 11 --
## 
##   1  2  3  4  5  6  7  8  9
## 
## -- Imputation 12 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 13 --
## 
##   1  2  3  4  5  6
## 
## -- Imputation 14 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 15 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 16 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 17 --
## 
##   1  2  3  4  5  6  7  8  9 10
## 
## -- Imputation 18 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 19 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 20 --
## 
##   1  2  3  4  5  6  7  8  9
## 
## -- Imputation 21 --
## 
##   1  2  3  4  5  6  7  8  9
## 
## -- Imputation 22 --
## 
##   1  2  3  4  5  6  7  8  9
## 
## -- Imputation 23 --
## 
##   1  2  3  4  5  6
## 
## -- Imputation 24 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 25 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 26 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 27 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 28 --
## 
##   1  2  3  4  5  6
## 
## -- Imputation 29 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 30 --
## 
##   1  2  3  4  5  6  7  8  9 10
## 
## -- Imputation 31 --
## 
##   1  2  3  4  5  6  7  8  9
## 
## -- Imputation 32 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 33 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 34 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 35 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 36 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 37 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 38 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 39 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 40 --
## 
##   1  2  3  4  5  6  7  8  9 10 11
## 
## -- Imputation 41 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 42 --
## 
##   1  2  3  4  5
## 
## -- Imputation 43 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 44 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 45 --
## 
##   1  2  3  4  5  6  7
## 
## -- Imputation 46 --
## 
##   1  2  3  4  5  6
## 
## -- Imputation 47 --
## 
##   1  2  3  4  5  6  7  8
## 
## -- Imputation 48 --
## 
##   1  2  3  4  5  6
## 
## -- Imputation 49 --
## 
##   1  2  3  4  5  6  7  8  9
## 
## -- Imputation 50 --
## 
##   1  2  3  4  5  6  7  8  9

Creating new dataset with missing data imputed

data(MI$imputations)
## Warning in data(MI$imputations): data set 'MI$imputations' not found
allimplogreg<-lapply(MI$imputations,function(X) {lme(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = X, method = "ML", na.action = "na.omit")})
betas<-MIextract(allimplogreg, fun=fixef)
vars<-MIextract(allimplogreg, fun=vcov)
summary(MIcombine(betas,vars))
## Multiple imputation results:
##       MIcombine.default(betas, vars)
##                 results         se     (lower    upper) missInfo
## (Intercept)  1.22826513 0.24706865  0.7422800 1.7142503     38 %
## GROUP1       0.01813835 0.21455642 -0.4032384 0.4395151     29 %
## WAVE        -0.09674956 0.10372741 -0.3008356 0.1073365     40 %
## BASELINE     0.70836868 0.05131513  0.6074234 0.8093139     39 %
## GROUP1:WAVE  0.08140390 0.13438879 -0.1825876 0.3453954     30 %

Check results with Imputations using Zelig

library("Zelig")
## Loading required package: boot
## 
## Attaching package: 'boot'
## 
## The following object is masked from 'package:car':
## 
##     logit
## 
## The following object is masked from 'package:psych':
## 
##     logit
## 
## Loading required package: MASS
## Loading required package: sandwich
## ZELIG (Versions 4.2-1, built: 2013-09-12)
## 
## +----------------------------------------------------------------+
## |  Please refer to http://gking.harvard.edu/zelig for full       |
## |  documentation or help.zelig() for help with commands and      |
## |  models support by Zelig.                                      |
## |                                                                |
## |  Zelig project citations:                                      |
## |    Kosuke Imai, Gary King, and Olivia Lau.  (2009).            |
## |    ``Zelig: Everyone's Statistical Software,''                 |
## |    http://gking.harvard.edu/zelig                              |
## |   and                                                          |
## |    Kosuke Imai, Gary King, and Olivia Lau. (2008).             |
## |    ``Toward A Common Framework for Statistical Analysis        |
## |    and Development,'' Journal of Computational and             |
## |    Graphical Statistics, Vol. 17, No. 4 (December)             |
## |    pp. 892-913.                                                |
## |                                                                |
## |   To cite individual Zelig models, please use the citation     |
## |   format printed with each model run and in the documentation. |
## +----------------------------------------------------------------+
## 
## 
## 
## Attaching package: 'Zelig'
## 
## The following objects are masked from 'package:psych':
## 
##     alpha, describe, sim
## 
## The following object is masked from 'package:utils':
## 
##     cite
zelig.fit <- zelig(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = MI$imputations,  model = "ls", cite = FALSE)
summary(zelig.fit)
## 
##   Model: ls
##   Number of multiply imputed data sets: 50 
## 
## Combined results:
## 
## Call:
## lm(formula = formula, weights = weights, model = F, data = data)
## 
## Coefficients:
##                   Value Std. Error      t-stat      p-value
## (Intercept)  1.22485148 0.24827870  4.93337323 1.253360e-06
## GROUP1       0.01835897 0.22820040  0.08045108 9.358996e-01
## WAVE        -0.09674956 0.11075504 -0.87354541 3.828779e-01
## BASELINE     0.70931953 0.04944241 14.34637887 1.474250e-35
## GROUP1:WAVE  0.08139392 0.14483680  0.56196988 5.743102e-01
## 
## For combined results from datasets i to j, use summary(x, subset = i:j).
## For separate results, use print(summary(x), subset = i:j).

Check assumptions with Random Computations

data1=MI$imputations[[1]]
library("Zelig")
zelig.fitdata1 <- zelig(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data1,  model = "ls", cite = FALSE)
summary(zelig.fitdata1)
## 
## Call:
## lm(formula = formula, weights = weights, model = F, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.98709 -0.27480  0.02221  0.24586  1.61164 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.99700    0.18953   5.260 4.21e-07 ***
## GROUP1      -0.03539    0.19040  -0.186    0.853    
## WAVE        -0.12463    0.08651  -1.441    0.151    
## BASELINE     0.77842    0.03654  21.301  < 2e-16 ***
## GROUP1:WAVE  0.16146    0.12033   1.342    0.181    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4011 on 173 degrees of freedom
## Multiple R-squared:  0.7251, Adjusted R-squared:  0.7187 
## F-statistic: 114.1 on 4 and 173 DF,  p-value: < 2.2e-16

Describe the meanGRIT variable by the GROUP variable

describeBy(data1[,3:4], group = data1$GROUP)
## group: 0
##          vars  n mean   sd median trimmed  mad  min max range  skew
## BASELINE    1 86 3.59 0.75   3.62    3.61 0.93 2.00   5  3.00 -0.24
## meanGRIT    2 86 3.60 0.73   3.56    3.62 0.77 1.88   5  3.12 -0.11
##          kurtosis   se
## BASELINE    -0.71 0.08
## meanGRIT    -0.63 0.08
## -------------------------------------------------------- 
## group: 1
##          vars  n mean   sd median trimmed  mad  min  max range  skew
## BASELINE    1 92 3.37 0.89   3.25    3.42 0.93 1.00 4.88  3.88 -0.47
## meanGRIT    2 92 3.64 0.78   3.75    3.68 0.76 1.25 5.11  3.86 -0.57
##          kurtosis   se
## BASELINE    -0.11 0.09
## meanGRIT     0.48 0.08

Create a plot that visualizes meanGRIT variable by the GROUP variable

library(ggplot2)
## 
## Attaching package: 'ggplot2'
## 
## The following object is masked from 'package:psych':
## 
##     %+%
library(influence.ME)
## 
## Attaching package: 'influence.ME'
## 
## The following object is masked from 'package:stats':
## 
##     influence

Take a look at the residuals

residual <- lm(meanGRIT ~ BASELINE, data=data1)$residual

Plot the residuals to see that they are random

plot(density(residual))# A density plot

qqnorm(residual) # A quantile normal plot to checking normality
qqline(residual)

Checking the different between intervention and control groups residuals. This allows us to control for individual unsystematic differences.

data2$residual <- NA
sel1 <- which(!is.na(data1$meanGRIT)) 
sel2 <- which(!is.na(data1$BASELINE))
data1$residual[intersect(sel1,sel2)] <- residual
qplot(GROUP, meanGRIT, data=data1, geom="boxplot")

Plot of the difference between intervention and control groups.

qplot(GROUP, residual, data=data1, geom="boxplot")

Two way repeated measures ======================================================== Graphing the Two-Way Interaction. Both meanGRIT and the Residuals

# Load the nlme package
library(nlme)
with(data1, boxplot(meanGRIT ~ WAVE + GROUP))

with(data1, boxplot(residual ~ WAVE + GROUP))
Linear Mixed-Effects Model

Comparing Basline to Wave 2 and 3 by Group.

fullModeldata1 <- lme(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data1, method = "ML", na.action = "na.omit")

Cooks Distence

CookD(fullModeldata1)

Plot Cook’s distance:

plot(fullModeldata1, which="cook")
Check results on this random Imputation model
Results

Explanation of significance:

We asses the significance of our models by comparing them from the baseline model using the anova() function.
(Intercept): Where everything is 0
GROUP1: Is there a difference between group. If it is significant than there is a difference and the treatment had an effect.
WAVE: Asseses whether the effects gets bigger beteen time 2 and 3 (does not have to be significant)
BASELINE: Should not be significant. If it is then it shows that there is a difference between groups before the treatment.
GROUP1:WAVE: If this is significant then it means that the effect was either fleeting or it happened after the treatment i.e. between time 2 and 3.

summary(fullModeldata1)
## Linear mixed-effects model fit by maximum likelihood
##  Data: data1 
##        AIC      BIC    logLik
##   188.3381 210.6106 -87.16906
## 
## Random effects:
##  Formula: ~1 | ID
##         (Intercept) Residual
## StdDev:   0.1100524 0.379832
## 
## Fixed effects: meanGRIT ~ GROUP * WAVE + BASELINE 
##                  Value  Std.Error DF   t-value p-value
## (Intercept)  1.0016777 0.18989769 87  5.274828  0.0000
## GROUP1      -0.0356515 0.18442659 87 -0.193310  0.8472
## WAVE        -0.1246345 0.08309202 86 -1.499957  0.1373
## BASELINE     0.7771130 0.03789675 86 20.506060  0.0000
## GROUP1:WAVE  0.1614363 0.11557952 86  1.396755  0.1661
##  Correlation: 
##             (Intr) GROUP1 WAVE   BASELI
## GROUP1      -0.530                     
## WAVE        -0.656  0.676              
## BASELINE    -0.716  0.041  0.000       
## GROUP1:WAVE  0.468 -0.940 -0.719  0.005
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -2.47134813 -0.64321135  0.02777506  0.59012091  3.72950435 
## 
## Number of Observations: 178
## Number of Groups: 89

Another random imputation

data10=MI$imputations[[10]]
library("Zelig")
zelig.fitdata10 <- zelig(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data10,  model = "ls", cite = FALSE)
summary(zelig.fitdata10)
## 
## Call:
## lm(formula = formula, weights = weights, model = F, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.07600 -0.24692 -0.01129  0.26417  1.44324 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.08547    0.19513   5.563 9.93e-08 ***
## GROUP1       0.28950    0.19614   1.476    0.142    
## WAVE        -0.11260    0.08905  -1.264    0.208    
## BASELINE     0.74494    0.03763  19.797  < 2e-16 ***
## GROUP1:WAVE -0.04115    0.12388  -0.332    0.740    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4129 on 173 degrees of freedom
## Multiple R-squared:  0.6963, Adjusted R-squared:  0.6893 
## F-statistic: 99.17 on 4 and 173 DF,  p-value: < 2.2e-16

Describe the meanGRIT variable by the GROUP variable

describeBy(data10[,3:4], group = data10$GROUP)
## group: 0
##          vars  n mean   sd median trimmed  mad  min max range  skew
## BASELINE    1 86 3.59 0.75   3.62    3.61 0.93 2.00   5  3.00 -0.24
## meanGRIT    2 86 3.59 0.72   3.51    3.59 0.82 1.88   5  3.12 -0.01
##          kurtosis   se
## BASELINE    -0.71 0.08
## meanGRIT    -0.71 0.08
## -------------------------------------------------------- 
## group: 1
##          vars  n mean   sd median trimmed  mad  min  max range  skew
## BASELINE    1 92 3.37 0.89   3.25    3.43 0.93 1.00 4.88  3.88 -0.47
## meanGRIT    2 92 3.66 0.76   3.75    3.67 0.78 1.25 5.35  4.10 -0.34
##          kurtosis   se
## BASELINE    -0.10 0.09
## meanGRIT     0.33 0.08

Create a plot that visualizes meanGRIT variable by the GROUP variable

library(ggplot2)
library(influence.ME)

Take a look at the residuals

residual <- lm(meanGRIT ~ BASELINE, data=data10)$residual

Plot the residuals to see that they are random

plot(density(residual))# A density plot

qqnorm(residual) # A quantile normal plot to checking normality
qqline(residual)

Checking the different between intervention and control groups residuals. This allows us to control for individual unsystematic differences.

data2$residual <- NA
sel1 <- which(!is.na(data10$meanGRIT)) 
sel2 <- which(!is.na(data10$BASELINE))
data10$residual[intersect(sel1,sel2)] <- residual
qplot(GROUP, meanGRIT, data=data10, geom="boxplot")

Plot of the difference between intervention and control groups.

qplot(GROUP, residual, data=data10, geom="boxplot")

Two way repeated measures ======================================================== Graphing the Two-Way Interaction. Both meanGRIT and the Residuals

# Load the nlme package
library(nlme)
with(data10, boxplot(meanGRIT ~ WAVE + GROUP))

with(data10, boxplot(residual ~ WAVE + GROUP))
Linear Mixed-Effects Model

Comparing Basline to Wave 2 and 3 by Group.

fullModeldata10 <- lme(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data10, method = "ML", na.action = "na.omit")

Cooks Distence

CookD(fullModeldata10)

Plot Cook’s distance:

plot(fullModeldata10, which="cook")
Check results on this random Imputation model
Results

Explanation of significance:

We asses the significance of our models by comparing them from the baseline model using the anova() function.
(Intercept): Where everything is 0
GROUP1: Is there a difference between group. If it is significant than there is a difference and the treatment had an effect.
WAVE: Asseses whether the effects gets bigger beteen time 2 and 3 (does not have to be significant)
BASELINE: Should not be significant. If it is then it shows that there is a difference between groups before the treatment.
GROUP1:WAVE: If this is significant then it means that the effect was either fleeting or it happened after the treatment i.e. between time 2 and 3.

summary(fullModeldata10)
## Linear mixed-effects model fit by maximum likelihood
##  Data: data10 
##        AIC     BIC    logLik
##   193.6405 215.913 -89.82024
## 
## Random effects:
##  Formula: ~1 | ID
##         (Intercept)  Residual
## StdDev:   0.2018203 0.3535264
## 
## Fixed effects: meanGRIT ~ GROUP * WAVE + BASELINE 
##                  Value  Std.Error DF   t-value p-value
## (Intercept)  1.0909092 0.19634933 87  5.555961  0.0000
## GROUP1       0.2890766 0.17593802 87  1.643059  0.1040
## WAVE        -0.1125967 0.07733739 86 -1.455916  0.1491
## BASELINE     0.7434288 0.04189798 86 17.743786  0.0000
## GROUP1:WAVE -0.0410892 0.10758757 86 -0.381914  0.7035
##  Correlation: 
##             (Intr) GROUP1 WAVE   BASELI
## GROUP1      -0.512                     
## WAVE        -0.591  0.659              
## BASELINE    -0.766  0.067  0.000       
## GROUP1:WAVE  0.437 -0.918 -0.719 -0.016
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -2.33162837 -0.59014635  0.02886088  0.58783729  2.66029176 
## 
## Number of Observations: 178
## Number of Groups: 89

Another random imputation

data15=MI$imputations[[15]]
library("Zelig")
zelig.fitdata15 <- zelig(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data15,  model = "ls", cite = FALSE)
summary(zelig.fitdata15)
## 
## Call:
## lm(formula = formula, weights = weights, model = F, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.96523 -0.27878  0.02434  0.29788  1.45914 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.36435    0.19672   6.936 7.74e-11 ***
## GROUP1      -0.09872    0.19821  -0.498    0.619    
## WAVE        -0.01278    0.09004  -0.142    0.887    
## BASELINE     0.65498    0.03781  17.322  < 2e-16 ***
## GROUP1:WAVE  0.09867    0.12525   0.788    0.432    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4175 on 173 degrees of freedom
## Multiple R-squared:  0.6365, Adjusted R-squared:  0.6281 
## F-statistic: 75.74 on 4 and 173 DF,  p-value: < 2.2e-16

Describe the meanGRIT variable by the GROUP variable

describeBy(data15[,3:4], group = data15$GROUP)
## group: 0
##          vars  n mean   sd median trimmed  mad  min max range  skew
## BASELINE    1 86 3.59 0.75   3.62    3.61 0.93 2.00   5  3.00 -0.24
## meanGRIT    2 86 3.70 0.68   3.62    3.71 0.74 1.88   5  3.12 -0.14
##          kurtosis   se
## BASELINE    -0.71 0.08
## meanGRIT    -0.48 0.07
## -------------------------------------------------------- 
## group: 1
##          vars  n mean   sd median trimmed  mad  min  max range  skew
## BASELINE    1 92 3.39 0.90   3.25    3.45 0.93 1.00 4.88  3.88 -0.49
## meanGRIT    2 92 3.62 0.69   3.75    3.64 0.74 1.25 5.05  3.80 -0.42
##          kurtosis   se
## BASELINE    -0.13 0.09
## meanGRIT     0.29 0.07

Create a plot that visualizes meanGRIT variable by the GROUP variable

library(ggplot2)
library(influence.ME)

Take a look at the residuals

residual <- lm(meanGRIT ~ BASELINE, data=data15)$residual

Plot the residuals to see that they are random

plot(density(residual))# A density plot

qqnorm(residual) # A quantile normal plot to checking normality
qqline(residual)

Checking the different between intervention and control groups residuals. This allows us to control for individual unsystematic differences.

data2$residual <- NA
sel1 <- which(!is.na(data15$meanGRIT)) 
sel2 <- which(!is.na(data15$BASELINE))
data15$residual[intersect(sel1,sel2)] <- residual
qplot(GROUP, meanGRIT, data=data15, geom="boxplot")

Plot of the difference between intervention and control groups.

qplot(GROUP, residual, data=data15, geom="boxplot")

Two way repeated measures ======================================================== Graphing the Two-Way Interaction. Both meanGRIT and the Residuals

# Load the nlme package
library(nlme)
with(data15, boxplot(meanGRIT ~ WAVE + GROUP))

with(data15, boxplot(residual ~ WAVE + GROUP))
Linear Mixed-Effects Model

Comparing Basline to Wave 2 and 3 by Group.

fullModeldata15 <- lme(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data15, method = "ML", na.action = "na.omit")

Cooks Distence

CookD(fullModeldata15)

Plot Cook’s distance:

plot(fullModeldata15, which="cook")
Check results on this random Imputation model
Results

Explanation of significance:

We asses the significance of our models by comparing them from the baseline model using the anova() function.
(Intercept): Where everything is 0
GROUP1: Is there a difference between group. If it is significant than there is a difference and the treatment had an effect.
WAVE: Asseses whether the effects gets bigger beteen time 2 and 3 (does not have to be significant)
BASELINE: Should not be significant. If it is then it shows that there is a difference between groups before the treatment.
GROUP1:WAVE: If this is significant then it means that the effect was either fleeting or it happened after the treatment i.e. between time 2 and 3.

summary(fullModeldata15)
## Linear mixed-effects model fit by maximum likelihood
##  Data: data15 
##        AIC      BIC    logLik
##   201.3456 223.6181 -93.67279
## 
## Random effects:
##  Formula: ~1 | ID
##         (Intercept)  Residual
## StdDev:   0.1543421 0.3815704
## 
## Fixed effects: meanGRIT ~ GROUP * WAVE + BASELINE 
##                  Value  Std.Error DF   t-value p-value
## (Intercept)  1.3655590 0.19743727 87  6.916420  0.0000
## GROUP1      -0.0987940 0.18678130 87 -0.528929  0.5982
## WAVE        -0.0127837 0.08347230 86 -0.153149  0.8786
## BASELINE     0.6546425 0.04035735 86 16.221146  0.0000
## GROUP1:WAVE  0.0986778 0.11610934 86  0.849870  0.3978
##  Correlation: 
##             (Intr) GROUP1 WAVE   BASELI
## GROUP1      -0.523                     
## WAVE        -0.634  0.670              
## BASELINE    -0.734  0.049  0.000       
## GROUP1:WAVE  0.461 -0.933 -0.719 -0.006
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -2.22439626 -0.60147121  0.06578369  0.66220700  3.06933848 
## 
## Number of Observations: 178
## Number of Groups: 89

Another random imputation

data25=MI$imputations[[25]]

library("Zelig")
zelig.fitdata25 <- zelig(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data25,  model = "ls", cite = FALSE)
summary(zelig.fitdata25)
## 
## Call:
## lm(formula = formula, weights = weights, model = F, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.94971 -0.24064  0.01018  0.24552  1.53378 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.10490    0.19246   5.741 4.14e-08 ***
## GROUP1       0.07309    0.19317   0.378    0.706    
## WAVE        -0.07478    0.08774  -0.852    0.395    
## BASELINE     0.72989    0.03716  19.642  < 2e-16 ***
## GROUP1:WAVE  0.05193    0.12204   0.426    0.671    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4068 on 173 degrees of freedom
## Multiple R-squared:  0.6908, Adjusted R-squared:  0.6837 
## F-statistic: 96.64 on 4 and 173 DF,  p-value: < 2.2e-16

Describe the meanGRIT variable by the GROUP variable

describeBy(data25[,3:4], group = data25$GROUP)
## group: 0
##          vars  n mean   sd median trimmed  mad  min max range  skew
## BASELINE    1 86 3.59 0.75   3.62    3.61 0.93 2.00   5  3.00 -0.24
## meanGRIT    2 86 3.61 0.73   3.62    3.62 0.74 1.88   5  3.12 -0.12
##          kurtosis   se
## BASELINE    -0.71 0.08
## meanGRIT    -0.71 0.08
## -------------------------------------------------------- 
## group: 1
##          vars  n mean   sd median trimmed  mad  min  max range  skew
## BASELINE    1 92 3.38 0.89   3.25    3.43 0.93 1.00 4.88  3.88 -0.49
## meanGRIT    2 92 3.61 0.72   3.75    3.66 0.74 1.05 4.88  3.83 -0.92
##          kurtosis   se
## BASELINE    -0.05 0.09
## meanGRIT     1.30 0.08

Create a plot that visualizes meanGRIT variable by the GROUP variable

library(ggplot2)
library(influence.ME)

Take a look at the residuals

residual <- lm(meanGRIT ~ BASELINE, data=data25)$residual

Plot the residuals to see that they are random

plot(density(residual))# A density plot

qqnorm(residual) # A quantile normal plot to checking normality
qqline(residual)

Checking the different between intervention and control groups residuals. This allows us to control for individual unsystematic differences.

data2$residual <- NA
sel1 <- which(!is.na(data25$meanGRIT)) 
sel2 <- which(!is.na(data25$BASELINE))
data25$residual[intersect(sel1,sel2)] <- residual
qplot(GROUP, meanGRIT, data=data25, geom="boxplot")

Plot of the difference between intervention and control groups.

qplot(GROUP, residual, data=data25, geom="boxplot")

Two way repeated measures ======================================================== Graphing the Two-Way Interaction. Both meanGRIT and the Residuals

# Load the nlme package
library(nlme)
with(data25, boxplot(meanGRIT ~ WAVE + GROUP))

with(data25, boxplot(residual ~ WAVE + GROUP))
Linear Mixed-Effects Model

Comparing Basline to Wave 2 and 3 by Group.

fullModeldata25 <- lme(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data25, method = "ML", na.action = "na.omit")

Cooks Distence

CookD(fullModeldata25)

Plot Cook’s distance:

plot(fullModeldata25, which="cook")
Check results on this random Imputation model
Results

Explanation of significance:

We asses the significance of our models by comparing them from the baseline model using the anova() function.
(Intercept): Where everything is 0
GROUP1: Is there a difference between group. If it is significant than there is a difference and the treatment had an effect.
WAVE: Asseses whether the effects gets bigger beteen time 2 and 3 (does not have to be significant)
BASELINE: Should not be significant. If it is then it shows that there is a difference between groups before the treatment.
GROUP1:WAVE: If this is significant then it means that the effect was either fleeting or it happened after the treatment i.e. between time 2 and 3.

summary(fullModeldata25)
## Linear mixed-effects model fit by maximum likelihood
##  Data: data25 
##        AIC      BIC    logLik
##   192.2751 214.5476 -89.13754
## 
## Random effects:
##  Formula: ~1 | ID
##         (Intercept)  Residual
## StdDev:   0.1467658 0.3732422
## 
## Fixed effects: meanGRIT ~ GROUP * WAVE + BASELINE 
##                  Value  Std.Error DF   t-value p-value
## (Intercept)  1.1055585 0.19326470 87  5.720437  0.0000
## GROUP1       0.0730443 0.18258380 87  0.400059  0.6901
## WAVE        -0.0747832 0.08165042 86 -0.915895  0.3623
## BASELINE     0.7297109 0.03955780 86 18.446702  0.0000
## GROUP1:WAVE  0.0519338 0.11357540 86  0.457263  0.6486
##  Correlation: 
##             (Intr) GROUP1 WAVE   BASELI
## GROUP1      -0.526                     
## WAVE        -0.634  0.671              
## BASELINE    -0.735  0.053  0.000       
## GROUP1:WAVE  0.461 -0.933 -0.719 -0.007
## 
## Standardized Within-Group Residuals:
##          Min           Q1          Med           Q3          Max 
## -2.272648874 -0.633121213 -0.007359014  0.602637748  3.289055657 
## 
## Number of Observations: 178
## Number of Groups: 89

Check assumptions on model without any imputations

Describe the meanGRIT variable by the GROUP variable

describeBy(data2[,3:4], group = data2$GROUP)
## group: 0
##          vars  n mean   sd median trimmed  mad  min max range  skew
## BASELINE    1 86 3.59 0.75   3.62    3.61 0.93 2.00   5  3.00 -0.24
## meanGRIT    2 59 3.72 0.73   3.62    3.74 0.74 1.88   5  3.12 -0.25
##          kurtosis   se
## BASELINE    -0.71 0.08
## meanGRIT    -0.52 0.09
## -------------------------------------------------------- 
## group: 1
##          vars  n mean   sd median trimmed  mad  min  max range  skew
## BASELINE    1 88 3.37 0.90   3.25    3.42 0.93 1.00 4.88  3.88 -0.46
## meanGRIT    2 54 3.62 0.66   3.75    3.65 0.74 1.25 4.88  3.62 -0.72
##          kurtosis   se
## BASELINE    -0.13 0.10
## meanGRIT     1.35 0.09

Create a plot that visualizes meanGRIT variable by the GROUP variable

library(ggplot2)
library(influence.ME)

Take a look at the residuals

residual <- lm(meanGRIT ~ BASELINE, data=data2)$residual

Plot the residuals to see that they are random

plot(density(residual))# A density plot

qqnorm(residual) # A quantile normal plot to checking normality
qqline(residual)

Checking the different between intervention and control groups residuals. This allows us to control for individual unsystematic differences.

data2$residual <- NA
sel1 <- which(!is.na(data2$meanGRIT)) 
sel2 <- which(!is.na(data2$BASELINE))
data2$residual[intersect(sel1,sel2)] <- residual
qplot(GROUP, meanGRIT, data=data2, geom="boxplot")
## Warning: Removed 65 rows containing non-finite values (stat_boxplot).

Plot of the difference between intervention and control groups.

qplot(GROUP, residual, data=data2, geom="boxplot")
## Warning: Removed 69 rows containing non-finite values (stat_boxplot).

Two way repeated measures ======================================================== Graphing the Two-Way Interaction. Both meanGRIT and the Residuals

# Load the nlme package
library(nlme)
with(data2, boxplot(meanGRIT ~ WAVE + GROUP))

with(data2, boxplot(residual ~ WAVE + GROUP))
Linear Mixed-Effects Model

Comparing Basline to Wave 2 and 3 by Group.

fullModel <- lme(meanGRIT ~ GROUP * WAVE + BASELINE, random = ~1 | ID, data = data2, method = "ML", na.action = "na.omit")

Cooks Distence

CookD(fullModel)

Plot Cook’s distance:

plot(fullModel, which="cook")
Results on Model with data that contains no imputations
Results

Explanation of significance:

We asses the significance of our models by comparing them from the baseline model using the anova() function.
(Intercept): Where everything is 0
GROUP1: Is there a difference between group. If it is significant than there is a difference and the treatment had an effect.
WAVE: Asseses whether the effects gets bigger beteen time 2 and 3 (does not have to be significant)
BASELINE: Should not be significant. If it is then it shows that there is a difference between groups before the treatment.
GROUP1:WAVE: If this is significant then it means that the effect was either fleeting or it happened after the treatment i.e. between time 2 and 3.

summary(fullModel)
## Linear mixed-effects model fit by maximum likelihood
##  Data: data2 
##       AIC      BIC    logLik
##   124.824 143.6634 -55.41201
## 
## Random effects:
##  Formula: ~1 | ID
##         (Intercept)  Residual
## StdDev:   0.2562764 0.3253572
## 
## Fixed effects: meanGRIT ~ GROUP * WAVE + BASELINE 
##                  Value  Std.Error DF   t-value p-value
## (Intercept)  1.1440052 0.24862785 66  4.601275  0.0000
## GROUP1       0.0533591 0.21189136 66  0.251823  0.8020
## WAVE        -0.0648842 0.09340478 38 -0.694656  0.4915
## BASELINE     0.7199886 0.05651616 66 12.739516  0.0000
## GROUP1:WAVE  0.0438335 0.13678769 38  0.320449  0.7504
##  Correlation: 
##             (Intr) GROUP1 WAVE   BASELI
## GROUP1      -0.456                     
## WAVE        -0.494  0.607              
## BASELINE    -0.818  0.082 -0.031       
## GROUP1:WAVE  0.350 -0.898 -0.682  0.005
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -2.04700338 -0.49141196  0.02432851  0.50910118  2.50432546 
## 
## Number of Observations: 109
## Number of Groups: 69
Table with P-value

|             |       Value|  Std.Error|  DF|     t-value|    p-value|
|:------------|-----------:|----------:|---:|-----------:|----------:|
|(Intercept)  |   1.1440052|  0.2486279|  66|   4.6012753|  0.0000196|
|GROUP1       |   0.0533591|  0.2118914|  66|   0.2518231|  0.8019603|
|WAVE         |  -0.0648842|  0.0934048|  38|  -0.6946563|  0.4914965|
|BASELINE     |   0.7199886|  0.0565162|  66|  12.7395156|  0.0000000|
|GROUP1:WAVE  |   0.0438335|  0.1367877|  38|   0.3204491|  0.7503835|

Table with confidence intervals

est. lower upper
(Intercept) 1.1440052 0.6591227 1.6288877
GROUP1 0.0533591 -0.3598786 0.4665969
WAVE -0.0648842 -0.2495845 0.1198161
BASELINE 0.7199886 0.6097688 0.8302083
GROUP1:WAVE 0.0438335 -0.2266530 0.3143199