Examining the data for 5 ADHD participants (AD01, AD04, AD08, AD10, and AD15) before and after medication
This section examines whether taking ADHD medicine before or after —- has an effect on the response variables while controlling for other factors.
# data can be easily uploaded via the 'import dataset' tab in R studio...
library(readxl)
adhd.med <- read_excel(file.choose(), sheet = "ADHDmed")
adhd.med2 <- subset(adhd.med, select=c(
# independent variables
`ID entered in RETeval device`,age,gender,VERT,
Ethnicity_Simple,`Tested eye`,iris,
Strength,Before_After,
# dependent variables
a_time, a_amp, b_time, b_amp, bdiva_amp, p72, Tmin, BT,
p_ratio, w_ratio))
# checking the data's structure...
str(adhd.med2)
## tibble [355 x 19] (S3: tbl_df/tbl/data.frame)
## $ ID entered in RETeval device: chr [1:355] "AD01" "AD01" "AD01" "AD01" ...
## $ age : num [1:355] 15.5 15.5 15.5 15.5 15.5 15.5 15.5 15.5 15.5 15.5 ...
## $ gender : num [1:355] 0 0 0 0 0 0 0 0 0 0 ...
## $ VERT : num [1:355] 2 2 2 2 2 2 2 2 2 2 ...
## $ Ethnicity_Simple : num [1:355] 4 4 4 4 4 4 4 4 4 4 ...
## $ Tested eye : num [1:355] 1 1 1 1 1 1 1 1 1 1 ...
## $ iris : num [1:355] 1.38 1.38 1.38 1.38 1.38 ...
## $ Strength : num [1:355] -0.367 -0.367 -0.119 0.114 0.398 0.398 0.477 0.477 0.602 0.602 ...
## $ Before_After : num [1:355] 0 0 0 0 0 0 0 0 0 0 ...
## $ a_time : num [1:355] 14.9 15.1 13.3 11.3 12.6 ...
## $ a_amp : num [1:355] -4.22 -3.01 -3.3 -5.79 -7.2 ...
## $ b_time : num [1:355] 20.2 21.8 21.7 24.2 24.5 ...
## $ b_amp : num [1:355] 14.2 10.7 17.2 25.9 38.1 ...
## $ bdiva_amp : num [1:355] 3.37 3.56 5.2 4.48 5.29 ...
## $ p72 : num [1:355] 1.74 -1.02 -2.67 -7.69 -4.42 ...
## $ Tmin : num [1:355] 99.9 99.4 100.1 100.1 85.9 ...
## $ BT : num [1:355] -2.63 -6.23 -4.41 -11.78 -8.71 ...
## $ p_ratio : num [1:355] -0.175 0.133 0.192 0.382 0.143 ...
## $ w_ratio : num [1:355] 0.888 1.301 1.064 1.231 1.04 ...
dim(adhd.med2)
## [1] 355 19
# rendering variables into their right type
adhd.med2$`ID entered in RETeval device` <- as.factor(adhd.med2$`ID entered in RETeval device`)
adhd.med2$gender <- as.factor(adhd.med2$gender)
adhd.med2$VERT <- as.factor(adhd.med2$VERT)
adhd.med2$VERT<-relevel(adhd.med2$VERT, ref='2') # ref level in this variable
adhd.med2$Ethnicity_Simple <- as.factor(adhd.med2$Ethnicity_Simple)
adhd.med2$`Tested eye` <- as.factor(adhd.med2$`Tested eye`)
adhd.med2$Strength <- as.factor(adhd.med2$Strength)
adhd.med2$Before_After <- as.factor(adhd.med2$Before_After)
# changing the name of the columns
colnames(adhd.med2) <- c('Participant','Age','Gender','Vert',
'Ethnicity','Eye','Iris','Strength','med before-after',
'a_time', 'a_amp', 'b_time', 'b_amp','bdiva_amp',
'p72', 'Tmin', 'BT', 'p_ratio', 'w_ratio')
# relabeling some of the variables' levels
library(plyr)
adhd.med2$Gender <-
mapvalues(adhd.med2$Gender,
from = c('0','1'),
to = c('Male','Female'))
adhd.med2$Ethnicity <-
mapvalues(adhd.med2$Ethnicity,
from = c('1','2','3','4','5'),
to = c('Caucasian','Asian','Afro-Caribbean','Latino','Mixed'))
adhd.med2$Eye <-
mapvalues(adhd.med2$Eye,
from = c('0','1'),
to = c('Right','Left'))
adhd.med2$`med before-after` <-
mapvalues(adhd.med2$`med before-after`,
from = c('0','1'),
to = c('Before','After'))
# checking if there are missing values...
# see https://cran.r-project.org/web/packages/naniar/vignettes/naniar-visualisation.html
# for more visualisations
#library(naniar)
#gg_miss_upset(adhd.med2)
Yep definitely a good thing to not have NAs!.
Descriptive statistics
This section shows some descriptive statistics along with robust estimations of location (median) and scale (MAD = median absolute deviations) for all dependent variables regardless of Strength and Eye.
library(tidyverse)
# Note 'Eye' is left out in order to avoid duplicates
# VERT and Strength are also not included
some_stats_by_participants.adhd <- adhd.med2 %>%
group_by(Participant, Age, Gender, Ethnicity, `med before-after`) %>%
summarize(median.A.time = median(a_time),
MAD.A.time=mad(a_time),
CV.A.time = mad(a_time)/median(a_time),
median.A.amp = median(a_amp),
MAD.A.amp=mad(a_amp),
CV.A.amp = mad(a_amp)/median(a_amp),
median.B.time = median(b_time),
MAD.B.time=mad(b_time),
CV.B.time = mad(b_time)/median(b_time),
median.B.amp = median(b_amp),
MAD.B.amp=mad(b_amp),
CV.B.amp = mad(b_amp)/median(b_amp),
median.bdiva_amp = median(bdiva_amp),
MAD.bdiva_amp=mad(bdiva_amp),
CV.bdiva_amp = mad(bdiva_amp)/median(bdiva_amp),
median.p72 = median(p72),
MAD.p72 =mad(p72),
CV.p72 = mad(p72)/median(p72),
median.Tmin = median(Tmin),
MAD.Tmin =mad(Tmin),
CV.Tmin = mad(Tmin)/median(Tmin),
median.BT = median(BT),
MAD.BT=mad(BT),
CV.BT = mad(BT)/median(BT),
median.p_ratio = median(p_ratio),
MAD.p_ratio =mad(p_ratio),
CV.p_ratio = mad(p_ratio)/median(p_ratio),
median.w_ratio = median(w_ratio),
MAD.w_ratio =mad(w_ratio),
CV.w_ratio = mad(w_ratio)/median(w_ratio))
# filtering the data for demographics only
demographics.adhd <- some_stats_by_participants.adhd[ , c(1:5)]
# number of observations per gender, ethnicity and group
# note: observations can come from the same participant
library(ggplot2)
library(ggpubr)
theme_set(theme_pubr())
count.data.adhd<-with(demographics.adhd,
ftable(Gender, `med before-after`, Ethnicity))
count.data.adhd1<-as.data.frame(count.data.adhd)
names(count.data.adhd1)[4] <- c('Count')
As in the above analyses, an important aspect here is that there should be one iris measure per eye per participant. That is, each participant should have two iris measurements; one for the left eye and another for the right eye.
# here b-time is simply an excuse to obtain unique values for the other variables
# the way to do is by setting the function to 'class'
iris.aggregate.adhd<-aggregate(data=adhd.med2,
b_time ~ Participant + Gender + Age + Ethnicity + `med before-after` + Eye + Iris, FUN='class')
iris.aggregate.sorted.adhd<-
data.frame(sort(table(iris.aggregate.adhd$Participant)))
# those with more than four observations
iris.aggregate.sorted.adhd$Var1[iris.aggregate.sorted.adhd$Freq>4]
## factor(0)
## Levels: AD01 AD04 AD08 AD10 AD15
# those with less than 2 observations
iris.aggregate.sorted.adhd$Var1[iris.aggregate.sorted.adhd$Freq<4]
## [1] AD01
## Levels: AD01 AD04 AD08 AD10 AD15
Modelling the effect of ‘before’ vs ‘after’ medicine
The effect of ‘before’ vs ‘after’ medicine in five ADHD participants will be assessed on each dependent variable while controlling for other variables.
a-time
A GEE model gives error. According to some sites it’s due to the design matrix not being invertible and this happens due to high multicolinearity or from the \(n < p\) issue; i.e. less observations than variables . Recall, the data is from just 5 participants.
A robust linear fixed effects models will be used .
library(robustbase)
library(knitr)
a.time.model.adhd <- lmrob(a_time ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
library(car)
anova.a.time.model.adhd <- Anova(a.time.model.adhd,
test.statistic='Chisq', type=3)
anova.a.time.model.adhd %>% kable()
| (Intercept) |
1 |
4.4418294 |
0.0350687 |
med before-after |
1 |
0.0501923 |
0.8227290 |
| Strength |
9 |
136.0113643 |
0.0000000 |
| Vert |
2 |
12.5818053 |
0.0018531 |
| Ethnicity |
2 |
4.4761873 |
0.1066616 |
| Iris |
1 |
0.9164068 |
0.3384202 |
| Gender |
1 |
0.0385071 |
0.8444283 |
| Age |
1 |
0.1029108 |
0.7483641 |
| Eye |
1 |
3.5289429 |
0.0603062 |
med before-after:Strength |
9 |
2.8845030 |
0.9686984 |
a-amp
a.amp.model.adhd <- lmrob(a_amp ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.a.amp.model.adhd <- Anova(a.amp.model.adhd,
test.statistic='Chisq', type=3)
anova.a.amp.model.adhd %>% kable()
| (Intercept) |
1 |
13.2633800 |
0.0002706 |
med before-after |
1 |
0.1431994 |
0.7051209 |
| Strength |
9 |
170.9715197 |
0.0000000 |
| Vert |
2 |
10.0895148 |
0.0064430 |
| Ethnicity |
2 |
12.8083152 |
0.0016547 |
| Iris |
1 |
0.8250441 |
0.3637095 |
| Gender |
1 |
23.0495922 |
0.0000016 |
| Age |
1 |
7.8366010 |
0.0051199 |
| Eye |
1 |
1.4759473 |
0.2244097 |
med before-after:Strength |
9 |
5.8626030 |
0.7535870 |
b-time
b.time.model.adhd <- lmrob(b_time ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.b.time.model.adhd <- Anova(b.time.model.adhd,
test.statistic='Chisq', type=3)
anova.b.time.model.adhd %>% kable()
| (Intercept) |
1 |
46.4339618 |
0.0000000 |
med before-after |
1 |
1.7290445 |
0.1885330 |
| Strength |
9 |
4785.2874645 |
0.0000000 |
| Vert |
2 |
1.6592069 |
0.4362222 |
| Ethnicity |
2 |
91.9502639 |
0.0000000 |
| Iris |
1 |
0.0108317 |
0.9171096 |
| Gender |
1 |
6.0001896 |
0.0143043 |
| Age |
1 |
5.3688885 |
0.0204990 |
| Eye |
1 |
2.4566369 |
0.1170291 |
med before-after:Strength |
9 |
5.6545118 |
0.7739360 |
b-amp
b.amp.model.adhd <- lmrob(b_amp ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.b.amp.model.adhd <- Anova(b.amp.model.adhd,
test.statistic='Chisq', type=3)
anova.b.amp.model.adhd %>% kable()
| (Intercept) |
1 |
89.0104283 |
0.0000000 |
med before-after |
1 |
0.7861413 |
0.3752695 |
| Strength |
9 |
672.2421771 |
0.0000000 |
| Vert |
2 |
14.6578668 |
0.0006563 |
| Ethnicity |
2 |
113.2517539 |
0.0000000 |
| Iris |
1 |
0.0001358 |
0.9907008 |
| Gender |
1 |
155.6425186 |
0.0000000 |
| Age |
1 |
87.6112895 |
0.0000000 |
| Eye |
1 |
7.3469220 |
0.0067178 |
med before-after:Strength |
9 |
12.9368234 |
0.1654858 |
bdiva_amp
bdiva_amp.model.adhd <- lmrob(bdiva_amp ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.bdiva_amp.model.adhd <- Anova(bdiva_amp.model.adhd,
test.statistic='Chisq', type=3)
anova.bdiva_amp.model.adhd %>% kable()
| (Intercept) |
1 |
0.2576293 |
0.6117536 |
med before-after |
1 |
1.9772521 |
0.1596802 |
| Strength |
9 |
93.7781270 |
0.0000000 |
| Vert |
2 |
2.5760579 |
0.2758139 |
| Ethnicity |
2 |
0.8432307 |
0.6559863 |
| Iris |
1 |
15.2468690 |
0.0000943 |
| Gender |
1 |
0.8547777 |
0.3552043 |
| Age |
1 |
2.0547066 |
0.1517365 |
| Eye |
1 |
0.5243572 |
0.4689889 |
med before-after:Strength |
9 |
7.7402809 |
0.5605167 |
p72
p72.model.adhd <- lmrob(p72 ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.p72.model.adhd <- Anova(p72.model.adhd,
test.statistic='Chisq', type=3)
anova.p72.model.adhd %>% kable()
| (Intercept) |
1 |
7.992669 |
0.0046967 |
med before-after |
1 |
1.735135 |
0.1877564 |
| Strength |
9 |
77.338688 |
0.0000000 |
| Vert |
2 |
2.491695 |
0.2876970 |
| Ethnicity |
2 |
6.040445 |
0.0487904 |
| Iris |
1 |
5.524869 |
0.0187480 |
| Gender |
1 |
5.426438 |
0.0198341 |
| Age |
1 |
3.874992 |
0.0490105 |
| Eye |
1 |
1.655538 |
0.1982069 |
med before-after:Strength |
9 |
9.859407 |
0.3619866 |
Tmin
Tmin.model.adhd <- lmrob(Tmin ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.Tmin.model.adhd <- Anova(Tmin.model.adhd,
test.statistic='Chisq', type=3)
anova.Tmin.model.adhd %>% kable()
| (Intercept) |
1 |
9.2165863 |
0.0023983 |
med before-after |
1 |
8.4136186 |
0.0037242 |
| Strength |
9 |
15.2563402 |
0.0841312 |
| Vert |
2 |
0.3297030 |
0.8480196 |
| Ethnicity |
2 |
116.1197838 |
0.0000000 |
| Iris |
1 |
0.7209381 |
0.3958364 |
| Gender |
1 |
8.2661757 |
0.0040390 |
| Age |
1 |
4.5640509 |
0.0326498 |
| Eye |
1 |
0.0003914 |
0.9842168 |
med before-after:Strength |
9 |
9.6684205 |
0.3779891 |
BT
BT.model.adhd <- lmrob(BT ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.BT.model.adhd <- Anova(BT.model.adhd,
test.statistic='Chisq', type=3)
anova.BT.model.adhd %>% kable()
| (Intercept) |
1 |
6.927651 |
0.0084873 |
med before-after |
1 |
1.803102 |
0.1793380 |
| Strength |
9 |
52.678194 |
0.0000000 |
| Vert |
2 |
3.168554 |
0.2050960 |
| Ethnicity |
2 |
10.077884 |
0.0064806 |
| Iris |
1 |
6.754120 |
0.0093531 |
| Gender |
1 |
3.352360 |
0.0671087 |
| Age |
1 |
2.098274 |
0.1474655 |
| Eye |
1 |
3.746515 |
0.0529177 |
med before-after:Strength |
9 |
12.310007 |
0.1963914 |
p_ratio
p_ratio.model.adhd <- lmrob(p_ratio ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.p_ratio.model.adhd <- Anova(p_ratio.model.adhd,
test.statistic='Chisq', type=3)
anova.p_ratio.model.adhd %>% kable()
| (Intercept) |
1 |
0.2532536 |
0.6147934 |
med before-after |
1 |
1.3677893 |
0.2421921 |
| Strength |
9 |
44.3435452 |
0.0000012 |
| Vert |
2 |
0.8132433 |
0.6658961 |
| Ethnicity |
2 |
0.8700534 |
0.6472474 |
| Iris |
1 |
1.6817016 |
0.1946985 |
| Gender |
1 |
0.6581762 |
0.4172045 |
| Age |
1 |
0.1000497 |
0.7517700 |
| Eye |
1 |
0.2108639 |
0.6460911 |
med before-after:Strength |
9 |
10.8657120 |
0.2850262 |
w_ratio
w_ratio.model.adhd <- lmrob(w_ratio ~
`med before-after` +
Strength + Vert + Ethnicity + Iris + Gender + Age + Eye
+ `med before-after`:Strength,
data = adhd.med2)
anova.w_ratio.model.adhd <- Anova(w_ratio.model.adhd,
test.statistic='Chisq', type=3)
anova.w_ratio.model.adhd %>% kable()
| (Intercept) |
1 |
7.330515e+01 |
0.0000000 |
med before-after |
1 |
3.492242e+03 |
0.0000000 |
| Strength |
9 |
4.715867e+04 |
0.0000000 |
| Vert |
2 |
8.735018e+00 |
0.0126828 |
| Ethnicity |
2 |
4.182045e+00 |
0.1235607 |
| Iris |
1 |
2.532768e+00 |
0.1115044 |
| Gender |
1 |
9.956382e-01 |
0.3183682 |
| Age |
1 |
7.701000e-04 |
0.9778614 |
| Eye |
1 |
5.023976e+00 |
0.0249987 |
med before-after:Strength |
9 |
3.718770e+03 |
0.0000000 |