An experiment on the effects of anti-anxiety medicine on memory recall when being primed with happy or sad memories. The participants were done on novel Islanders whom mimic real-life humans in response to external factors.
Drugs of interest (known-as) [Dosage 1, 2, 3]:
A - Alprazolam (Xanax, Long-term) [1mg/3mg/5mg]
T - Triazolam (Halcion, Short-term) [0.25mg/0.5mg/0.75mg]
S- Sugar Tablet (Placebo) [1 tab/2tabs/3tabs]
Participants - all genders above 25+ years old to ensure a fully developed pre-frontal cortex, a region responsible for higher level cognition and memory recall.
pacman::p_load(tidyverse, GGally, ggExtra, reshape2, AICcmodavg, RColorBrewer)
memtest <- read.csv('Islander_data.csv', header = TRUE)
summary(memtest)
## first_name last_name age Happy_Sad_group
## Length:198 Length:198 Min. :24.00 Length:198
## Class :character Class :character 1st Qu.:30.00 Class :character
## Mode :character Mode :character Median :37.00 Mode :character
## Mean :39.53
## 3rd Qu.:48.00
## Max. :83.00
## Dosage Drug Mem_Score_Before Mem_Score_After
## Min. :1.00 Length:198 Min. : 27.20 Min. : 27.10
## 1st Qu.:1.00 Class :character 1st Qu.: 46.52 1st Qu.: 47.17
## Median :2.00 Mode :character Median : 54.80 Median : 56.75
## Mean :1.99 Mean : 57.97 Mean : 60.92
## 3rd Qu.:3.00 3rd Qu.: 68.40 3rd Qu.: 73.25
## Max. :3.00 Max. :110.00 Max. :120.00
## Diff
## Min. :-40.400
## 1st Qu.: -3.175
## Median : 1.700
## Mean : 2.955
## 3rd Qu.: 5.925
## Max. : 49.000
str(memtest)
## 'data.frame': 198 obs. of 9 variables:
## $ first_name : chr "Bastian" "Evan" "Florencia" "Holly" ...
## $ last_name : chr "Carrasco" "Carrasco" "Carrasco" "Carrasco" ...
## $ age : int 25 52 29 50 52 37 35 38 29 36 ...
## $ Happy_Sad_group : chr "H" "S" "H" "S" ...
## $ Dosage : int 1 1 1 1 1 1 1 1 1 1 ...
## $ Drug : chr "A" "A" "A" "A" ...
## $ Mem_Score_Before: num 63.5 41.6 59.7 51.7 47 66.4 44.1 76.3 56.2 54.8 ...
## $ Mem_Score_After : num 61.2 40.7 55.1 51.2 47.1 58.1 56 74.8 45 75.9 ...
## $ Diff : num -2.3 -0.9 -4.6 -0.5 0.1 -8.3 11.9 -1.5 -11.2 21.1 ...
sum(is.na(memtest))
## [1] 0
There is no missing value from our data set but we need to convert variables into categorical format for further analysis.
memtest$Happy_Sad_group <- as.factor(memtest$Happy_Sad_group)
memtest$Drug <- as.factor(memtest$Drug)
table(memtest$Dosage, memtest$Drug)
##
## A S T
## 1 23 22 22
## 2 22 22 22
## 3 22 22 21
The table shows the number of participants for each drug and dosage. It can be seen that the participants in the clinical trial were grouped into equal groups.
table(memtest$Dosage, memtest$Drug, memtest$Happy_Sad_group)
## , , = H
##
##
## A S T
## 1 11 11 11
## 2 11 11 11
## 3 11 11 11
##
## , , = S
##
##
## A S T
## 1 12 11 11
## 2 11 11 11
## 3 11 11 10
The participants in each group were allocated evenly among the the happy and sad memory groups.
memtest %>% ggplot(aes(x = age, color = 4)) +
geom_histogram(bins = 45, fill = '#2196f3', position = "dodge", show.legend = FALSE) +
labs(x = 'Age', y = 'Number of participants') +
ggtitle('Distribution of participants by age') +
theme(plot.title = element_text(hjust = 0.5))
ggpairs(memtest[,c('age','Dosage','Mem_Score_Before','Mem_Score_After','Diff')],
title = 'Correlation between the variables',
aes(color = memtest$Happy_Sad_group, alpha = 0.25),
columnLabels = c('age','Dosage','Mem_Score_Before','Mem_Score_After','Diff'),
lower = list(continuous = "smooth")) +
theme(plot.title = element_text(hjust = 0.5))
The age distribution shows a non normal distribution among the participants. The correlation plot shows positive correlation between the memory score before and after the test also we can observe a slight positive correlation between the memory score after and the difference between the scores.
hist <- ggplot(memtest, aes(x = Mem_Score_Before, y = Mem_Score_After, color = Happy_Sad_group)) +
geom_point() +
scale_y_continuous(breaks = seq(25, 125, by = 10)) +
scale_x_continuous(breaks = seq(25, 125, by = 10)) +
labs(x = 'Memory score before', y = 'Memory score after') +
ggtitle('Memory score by block comarison')
ggMarginal(hist, type = 'histogram', xparams = list(fill = 4), yparams = list(fill = 3))
happy_grp <- select(memtest, first_name, last_name, Mem_Score_Before, Mem_Score_After, Happy_Sad_group) %>% filter(Happy_Sad_group == 'H') %>% melt(id.vars = c('first_name', 'last_name', 'Happy_Sad_group'))
sad_grp <- select(memtest, first_name, last_name, Mem_Score_Before, Mem_Score_After, Happy_Sad_group) %>% filter(Happy_Sad_group == 'S') %>% melt(id.vars = c('first_name', 'last_name', 'Happy_Sad_group'))
ggplot(happy_grp, aes(x = variable, y = value, color = variable)) +
geom_point(cex = 1.5, pch = 1.0, position = position_jitter(w = 0.1, h = 0)) +
geom_boxplot(alpha = 0.5) +
geom_point(data = sad_grp ,cex = 1.5, pch = 1.0, position = position_jitter(w = 0.1, h = 0)) +
geom_boxplot(data = sad_grp, alpha = 0.5) +
facet_wrap(~Happy_Sad_group) +
labs(x = 'Memory score before vs after', y = 'Score') +
ggtitle('Memory score by block before and after the test comarison') +
theme(plot.title = element_text(hjust = 0.5))
The correlation between the scores before and after the test among the happy and sad blocks. On the second plot can be seen that the means of the groups are negligible between the two blocks indicating that the happy and sad memories that were primed 10 minutes prior to testing has no significant contribution.
memtest %>% ggplot(aes(x = Mem_Score_Before, y = Mem_Score_After, color = Dosage)) +
geom_point() +
facet_wrap(~Drug) +
labs(x = 'Memory score before', y = 'Memory score after') +
ggtitle('Memory score by drugs and dosage comparison') +
theme(plot.title = element_text(hjust = 0.5))
The correlation between the two different drugs and the control group among the two blocks.
memtest %>% ggplot(aes(x = Dosage, y = Mem_Score_Before, group = Dosage, color = Dosage)) +
geom_point(cex = 1.5, pch = 1.0, position = position_jitter(w = 0.1, h = 0)) +
geom_boxplot(alpha = 0.5) +
facet_wrap(~Drug) +
labs(x = 'Dosage', y = 'Memory score before') +
ggtitle('Memory score by drug groups and dosage before the test') +
theme(plot.title = element_text(hjust = 0.5))
memtest %>% ggplot(aes(x = Dosage, y = Mem_Score_After, group = Dosage, color = Dosage)) +
geom_point(cex = 1.5, pch = 1, position = position_jitter(w = 0.1, h = 0)) +
geom_boxplot(alpha = 0.5) +
facet_wrap(~Drug) +
labs(x = 'Dosage', y = 'Memory score after') +
ggtitle('Memory score by drug groups and dosage after the test') +
theme(plot.title = element_text(hjust = 0.5))
The Alprazolam has the biggest effect on memory score with the increase of dosage there is a significant increase on scores. The sugar and Triazolam didn’t show any significant increase in scores only a slight improvement with the Triazolam.
A_ordered <- memtest[order(memtest$Mem_Score_Before, decreasing = TRUE),] %>% filter(Drug == 'A') %>% select(first_name, last_name, Mem_Score_Before, Mem_Score_After) %>%
head(n = 10)
S_ordered <- memtest[order(memtest$Mem_Score_Before, decreasing = TRUE),] %>% filter(Drug == 'S') %>% select(first_name, last_name, Mem_Score_Before, Mem_Score_After) %>%
head(n = 10)
T_ordered <- memtest[order(memtest$Mem_Score_Before, decreasing = TRUE),] %>% filter(Drug == 'T') %>% select(first_name, last_name, Mem_Score_Before, Mem_Score_After) %>%
head(n = 10)
melt(A_ordered, id.vars = c('first_name', 'last_name')) %>%
ggplot(aes(x = value, y = first_name)) +
geom_line() +
geom_point(aes(color = variable), size = 3) +
theme(legend.position = 'bottom', plot.title = element_text(hjust = 0.5)) +
ggtitle('Top 10 performer before vs after - Drug A') +
labs(x = 'Performance score', y = 'Top 10 performer')
melt(S_ordered, id.vars = c('first_name', 'last_name')) %>%
ggplot(aes(x = value, y = first_name)) +
geom_line() +
geom_point(aes(color = variable), size = 3) +
theme(legend.position = 'bottom', plot.title = element_text(hjust = 0.5)) +
ggtitle('Top 10 performer before vs after - Drug S') +
labs(x = 'Performance score', y = 'Top 10 performer')
melt(T_ordered, id.vars = c('first_name', 'last_name')) %>%
ggplot(aes(x = value, y = first_name)) +
geom_line() +
geom_point(aes(color = variable), size = 3) +
theme(legend.position = 'bottom', plot.title = element_text(hjust = 0.5)) +
ggtitle('Top 10 performer before vs after - Drug T') +
labs(x = 'Performance score', y = 'Top 10 performer')
The orange dots show the top 10 performers score before they received their dosage compared with the score after the dosage shown as blue for each Drug.
Performing an ANNOVA test to estimate how the quantitative dependent variable changes according to the independent variables. Our null hypothesis (\(H_0\)) of the ANOVA is no difference in means and the alternative hypothesis (\(H_a\)) is that the means are differ from each other. We convert the dosage into a factorial variable before performing the test.
memtest$Dosage <- as.factor(memtest$Dosage)
anova_two_way <- aov(formula = Diff~Drug+Dosage, data = memtest)
summary(anova_two_way)
## Df Sum Sq Mean Sq F value Pr(>F)
## Drug 2 4305 2152.4 24.082 4.6e-10 ***
## Dosage 2 1231 615.4 6.885 0.00129 **
## Residuals 193 17250 89.4
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova_interaction <- aov(formula = Diff~Drug*Dosage, data = memtest)
summary(anova_interaction)
## Df Sum Sq Mean Sq F value Pr(>F)
## Drug 2 4305 2152.4 33.596 3.28e-13 ***
## Dosage 2 1231 615.4 9.606 0.000106 ***
## Drug:Dosage 4 5141 1285.3 20.062 8.74e-14 ***
## Residuals 189 12109 64.1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova_inter_confounding <- aov(formula = Diff~Drug*Dosage+Happy_Sad_group, data = memtest)
summary(anova_inter_confounding)
## Df Sum Sq Mean Sq F value Pr(>F)
## Drug 2 4305 2152.4 33.443 3.75e-13 ***
## Dosage 2 1231 615.4 9.562 0.000111 ***
## Happy_Sad_group 1 8 8.0 0.124 0.725429
## Drug:Dosage 4 5142 1285.5 19.973 1.01e-13 ***
## Residuals 188 12100 64.4
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova_three_inter <- aov(formula = Diff~Drug*Dosage*Happy_Sad_group, data = memtest)
summary(anova_three_inter)
## Df Sum Sq Mean Sq F value Pr(>F)
## Drug 2 4305 2152.4 33.957 3.07e-13 ***
## Dosage 2 1231 615.4 9.709 9.90e-05 ***
## Happy_Sad_group 1 8 8.0 0.126 0.7234
## Drug:Dosage 4 5142 1285.5 20.280 8.33e-14 ***
## Drug:Happy_Sad_group 2 106 53.2 0.839 0.4338
## Dosage:Happy_Sad_group 2 9 4.7 0.075 0.9280
## Drug:Dosage:Happy_Sad_group 4 574 143.6 2.266 0.0639 .
## Residuals 180 11409 63.4
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
From the above summary we can decide whether any of the group means differ from the overall mean. We conducted four different ANOVA test, first a two-way ANOVA, secondly we included the interaction effect, third time the Happy and Sad grouping was included and lastly the interaction between the three independent variables were tested.
From the summary we can conclude that the Happy_Sad_group variable and it’s interaction with the other variables is not significant. It has a low sum of squares and a high p-value which means that it’s not effecting the variation in the dependent variable. The only significant variables are the Drug, Dosage and the interaction between them.
Finding the best fit model that explains the variation in the model, we use the Akaike Information Criterion that calculates the information value on each model. We choose the model with the lowest value as it explains more variation in the data.
model.set <- list(anova_two_way, anova_interaction, anova_inter_confounding, anova_three_inter)
model.names <- c('anova_two_way', 'anova_interaction', 'anova_inter_confounding', 'anova_three_inter')
aictab(model.set, modnames = model.names)
##
## Model selection based on AICc:
##
## K AICc Delta_AICc AICcWt Cum.Wt LL
## anova_interaction 10 1397.53 0.00 0.74 0.74 -688.18
## anova_inter_confounding 11 1399.63 2.10 0.26 0.99 -688.10
## anova_three_inter 19 1406.85 9.32 0.01 1.00 -682.29
## anova_two_way 6 1458.86 61.33 0.00 1.00 -723.21
The ANOVA model with the interaction between the independent variables Drug and Dosage has the lowest AIC score and AIC weight meaning that it explains 74% of the total variation in the dependent variable.
We plot the model diagnostic
par(mfrow = c(2,2))
plot(anova_interaction)
par(mfrow = c(1,1))
The residuals plot shows the unexplained variance in the data. The Residuals vs Fitted looks good, it’s centered around zero. The Q-Q plot shows the regression between theoretical and the actual residuals of the model. It shows some issues at both tails but seems sufficient. The Scale-Location is centered around at 1 meaning no large outliers.
From the ANOVA test we concluded that there are differences between the group means but we don’t know which groups are statistically different from one other. We perform the Tukey’s HSD post-hoc test to do a pairwise comparison between each group and determine which groups are statistically significant.
TukeyHSD(anova_interaction, ordered = TRUE, conf.level = 0.95)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
## factor levels have been ordered
##
## Fit: aov(formula = Diff ~ Drug * Dosage, data = memtest)
##
## $Drug
## diff lwr upr p adj
## S-T 0.4164802 -2.887694 3.720654 0.9523138
## A-T 10.0578416 6.765925 13.349758 0.0000000
## A-S 9.6413614 6.362128 12.920595 0.0000000
##
## $Dosage
## diff lwr upr p adj
## 2-1 1.709833 -1.5694003 4.989067 0.4359693
## 3-1 5.943694 2.6517782 9.235611 0.0000932
## 3-2 4.233861 0.9296874 7.538035 0.0078887
##
## $`Drug:Dosage`
## diff lwr upr p adj
## T:3-S:3 0.4216450 -7.24131446 8.084605 1.0000000
## T:1-S:3 0.9045455 -6.66878568 8.477877 0.9999887
## S:2-S:3 1.3909091 -6.18242205 8.964240 0.9996978
## A:1-S:3 2.4498024 -5.04075764 9.940362 0.9829030
## T:2-S:3 3.2954545 -4.27787659 10.868786 0.9092533
## S:1-S:3 4.5318182 -3.04151296 12.105149 0.6298208
## A:2-S:3 8.0272727 0.45394159 15.600604 0.0286613
## A:3-S:3 24.7863636 17.21303250 32.359695 0.0000000
## T:1-T:3 0.4829004 -7.18005905 8.145860 0.9999999
## S:2-T:3 0.9692641 -6.69369541 8.632224 0.9999823
## A:1-T:3 2.0281573 -5.55300962 9.609324 0.9954813
## T:2-T:3 2.8738095 -4.78914995 10.536769 0.9604468
## S:1-T:3 4.1101732 -3.55278632 11.773133 0.7560868
## A:2-T:3 7.6056277 -0.05733177 15.268587 0.0534433
## A:3-T:3 24.3647186 16.70175914 32.027678 0.0000000
## S:2-T:1 0.4863636 -7.08696750 8.059695 0.9999999
## A:1-T:1 1.5452569 -5.94530309 9.035817 0.9992873
## T:2-T:1 2.3909091 -5.18242205 9.964240 0.9863640
## S:1-T:1 3.6272727 -3.94605841 11.200604 0.8530021
## A:2-T:1 7.1227273 -0.45060387 14.696058 0.0833773
## A:3-T:1 23.8818182 16.30848704 31.455149 0.0000000
## A:1-S:2 1.0588933 -6.43166673 8.549453 0.9999583
## T:2-S:2 1.9045455 -5.66878568 9.477877 0.9970566
## S:1-S:2 3.1409091 -4.43242205 10.714240 0.9298886
## A:2-S:2 6.6363636 -0.93696750 14.209695 0.1379000
## A:3-S:2 23.3954545 15.82212341 30.968786 0.0000000
## T:2-A:1 0.8456522 -6.64490783 8.336212 0.9999927
## S:1-A:1 2.0820158 -5.40854420 9.572576 0.9941342
## A:2-A:1 5.5774704 -1.91308965 13.068030 0.3257899
## A:3-A:1 22.3365613 14.84600126 29.827121 0.0000000
## S:1-T:2 1.2363636 -6.33696750 8.809695 0.9998752
## A:2-T:2 4.7318182 -2.84151296 12.305149 0.5726311
## A:3-T:2 21.4909091 13.91757795 29.064240 0.0000000
## A:2-S:1 3.4954545 -4.07787659 11.068786 0.8772990
## A:3-S:1 20.2545455 12.68121432 27.827877 0.0000000
## A:3-A:2 16.7590909 9.18575977 24.332422 0.0000000
The Tukey’s HSD test shows that the interaction of Alprazolam and Triazolam, Alprazolam and Sugar are significant however the Triazolam and Sugar comparison found to be not significant with a p-value of 0.95. The significant groups are:
plot(TukeyHSD(anova_interaction, conf.level = 0.95),las = 2.5)
The significant group-wise differences are where the confidence interval doesn’t include zero, meaning that the p-value is <0.05.
ggplot(memtest, aes(x = Dosage, y = Diff, color = Drug)) +
geom_point(cex = 1.5, pch = 1.0, position = position_jitter(w = 0.1, h = 0)) +
geom_boxplot(alpha = 0.2) +
stat_summary(fun.data = 'mean_se', geom = 'errorbar', width = 0.2, color = 'red') +
stat_summary(fun.data = 'mean_se', geom = 'pointrange', color = 'red') +
facet_wrap(~Drug) +
labs(x = 'Dosage from 1-3', y = 'Difference in memory score') +
ggtitle('Mean difference of Drugs by dosage with standard error') +
theme(plot.title = element_text(hjust = 0.5))
The mean difference between the Drugs and Dosages shown above with standard errors for each group.
m <- matrix(data = 0, nrow = 3, ncol = 3, dimnames = list(c('1','2','3'),c('A','T','S')))
for (i in 1:3){
for (j in c('A', 'T', 'S')){
x <- memtest %>% group_by(Dosage) %>% filter(Dosage == i, Drug == j) %>% summarise(mean = mean(Diff))
m[i,j] <- x[[2]]
}
}
m <- as.data.frame(m)
ggplot(m, aes(x = seq(1,3, by = 1))) +
geom_line(aes(y = A), color = 1, linewidth = 0.5) +
geom_line(aes(y = T), color = 2, linewidth = 0.5) +
geom_line(aes(y = S), color = 3, linewidth = 0.5) +
ggtitle('Interaction plot of A,T,S by dosage') +
labs(x = "Dosage", y = "Mean Difference") +
theme(plot.title = element_text(hjust = 0.5))
m
## A T S
## 1 0.3043478 -1.240909 2.3863636
## 2 5.8818182 1.150000 -0.7545455
## 3 22.6409091 -1.723810 -2.1454545
The interaction plot between the two drugs and the placebo shown above. The Alprazolam has the largest range of variance but it didn’t interact with the Triazolam while the control group (Sugar) has interactions with both group.