1 Set up R environment

library(tidyverse)
library(ggplot2)
library(ggpubr)
library(plyr)
library(magick)
library(png)
library(EBImage)
library(lme4)
library(lmerTest)

2 Set the R working drectory to the main experiment directory.

setwd("/Users/adambarnas/Box/MetaAwareness/data/")

3 Read in the individual subject files.

tbl_all <- list.files(path = "./Ma_Recontact", pattern = "*.csv", full.names = T, ignore.case = F) %>%
  map_df(~read.csv(., colClasses=c("gender..m.f."="character")))
tbl_all = subset(tbl_all, select = c(user_resp.keys,user_resp.rt,workerId,image_a))
col_idx <- grep("workerId", names(tbl_all))
tbl_all <- tbl_all[, c(col_idx, (1:ncol(tbl_all))[-col_idx])]
tbl_all <- data.frame(na.omit(tbl_all))
tbl_all <- tbl_all %>%
separate(image_a,into=c('database', 'image', NA), sep = "([\\_\\-])")

4 Compute average likelihood rating.

tbl_all_subj_avg <- tbl_all %>%
  group_by(workerId,image) %>%
  dplyr::summarize(average = mean(user_resp.keys)) %>%
  spread(image,average) %>% 
  mutate(subj_avg = rowMeans(.[-1], na.rm = TRUE))
mean(tbl_all_subj_avg$subj_avg)
## [1] 3.030042
tbl_all_img_avg <- data.frame(img_avg = colMeans(tbl_all_subj_avg[,2:70], na.rm = TRUE))
tbl_all_img_avg <- tibble::rownames_to_column(tbl_all_img_avg, "image")

5 Merge Mudsplash and Meta-Awareness data files.

ma_RTs_raw <- read_csv("Ma_RTs_raw.csv")
ma_RTs_raw <- ma_RTs_raw[, -c(2,3,6:16,18,19)]
ma_RTs_raw <- ma_RTs_raw[(ma_RTs_raw$workerId %in% tbl_all_subj_avg$workerId),]

tbl_all <- tbl_all[order(tbl_all$workerId, tbl_all$image), , drop = FALSE]
ma_RTs_raw <- ma_RTs_raw[order(ma_RTs_raw$workerId, ma_RTs_raw$image), , drop = FALSE]
ma_RTs_raw <- ma_RTs_raw %>%  
    mutate(image = as.character(image))

ma_RTs_likelihood <- left_join(tbl_all, ma_RTs_raw, by = c("workerId", "image"))
colnames(ma_RTs_likelihood)[2] <- "likelihood_rating"
colnames(ma_RTs_likelihood)[3] <- "likelihood_rt"
colnames(ma_RTs_likelihood)[7] <- "detection_rt"
ma_RTs_likelihood <- ma_RTs_likelihood[, c(-4)]
colnames(ma_RTs_likelihood)[5] <- "stim_set"
ma_RTs_likelihood <- ma_RTs_likelihood[,c(1,4,5,6,2,3)]

6 Add change_type to dataframe.

ma_change_type<- read_csv("Ma_change_type.csv")
ma_change_type <- ma_change_type %>%  
    mutate(image = as.character(image))
ma_RTs_likelihood <- left_join(ma_RTs_likelihood, ma_change_type, by = "image")

7 Compute likelihood rating for each image.

ma_RTs_likelihood %>%
  ggbarplot(x = "image", y = "likelihood_rating", ylab = "Likelihood of Detecting Change", title = "All images (30 per subject)", fill = "#f7a800", add = "mean_se", font.xtickslab = 8, sort.val = c("asc")) + rotate_x_text() + theme(legend.position = "none")

ma_RTs_likelihood_no_NA <- ma_RTs_likelihood %>%
  drop_na()

ma_RTs_likelihood_no_NA %>% 
  ggbarplot(x = "image", y = "likelihood_rating", ylab = "Likelihood of Detecting Change", title = "'Correct' images", fill = "#f7a800", add = "mean_se", font.xtickslab = 8, sort.val = c("asc")) + rotate_x_text() + theme(legend.position = "none")

write.csv(ma_RTs_likelihood,'Ma_RTs_likelihood.csv', row.names=FALSE)
write.csv(ma_RTs_likelihood_no_NA,'Ma_RTs_likelihood_no_NA.csv', row.names=FALSE)

8 Count number of ratings.

ma_RTs_likelihood_count <- ma_RTs_likelihood_no_NA %>% 
  group_by(workerId,image) %>% 
  dplyr::summarize(counts = n()) %>%
  spread(image,counts) %>%
  mutate(sum = rowSums(.[-1], na.rm = TRUE))
#head(tbl_all_counts,10)

ma_RTs_likelihood_count <- data.frame(count = colSums(ma_RTs_likelihood_count[,2:70], na.rm = TRUE))
ma_RTs_likelihood_count <- tibble::rownames_to_column(ma_RTs_likelihood_count, "image")
ma_RTs_likelihood_count
##       image count
## 1  10504629     8
## 2  10810329     8
## 3   1191801     9
## 4  12115280     8
## 5  12178414     8
## 6  13141692     8
## 7   1383096     8
## 8  13873251    10
## 9  16527526     5
## 10 18169626     6
## 11 18345691     9
## 12 22020472     7
## 13 23024660     9
## 14 23199105     8
## 15 24383097     7
## 16 25107991     8
## 17 27857618     9
## 18  3099758    10
## 19 31236119     9
## 20 32289063     6
## 21 38466626     6
## 22 38546029     4
## 23 42429798     3
## 24  4247084     8
## 25 44993860     7
## 26 45525109     7
## 27 46475259     8
## 28 46635293     5
## 29 48384711     5
## 30 48486405    10
## 31 51537628     7
## 32 51856108     7
## 33 55174490     8
## 34 56835136     7
## 35 57861456     7
## 36 61118260    12
## 37 62096551     9
## 38 62224663     8
## 39 67862299     7
## 40 69128765     4
## 41 70687495     7
## 42 72488522     8
## 43 73637203     9
## 44 74173745     8
## 45 75081153     5
## 46 75958241    10
## 47 77345858     9
## 48 77574131     6
## 49 77793328     9
## 50 79191795     5
## 51 79222679     7
## 52 79241011     8
## 53 79573638     7
## 54  8197559     4
## 55 81993755     1
## 56 83536470     5
## 57 83691215     5
## 58 83785171     8
## 59 85741618     6
## 60 86520382     4
## 61 87983207     3
## 62 88767165     7
## 63 89354846     7
## 64  8974554     9
## 65 90405028     8
## 66 95091295     4
## 67 97475929     3
## 68 98156944     9
## 69 98265889     3

9 Mixed effects model and correlation.

fit0 <- lmer(detection_rt ~ likelihood_rating + (1 | workerId) + (1 | image) + (1 | change_type), data=ma_RTs_likelihood_no_NA)
summary(fit0)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: detection_rt ~ likelihood_rating + (1 | workerId) + (1 | image) +  
##     (1 | change_type)
##    Data: ma_RTs_likelihood_no_NA
## 
## REML criterion at convergence: 2852.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5387 -0.3390 -0.1115  0.0754  9.3376 
## 
## Random effects:
##  Groups      Name        Variance Std.Dev.
##  image       (Intercept)  0.08619 0.2936  
##  workerId    (Intercept)  9.05703 3.0095  
##  change_type (Intercept)  0.23028 0.4799  
##  Residual                19.28909 4.3919  
## Number of obs: 483, groups:  image, 69; workerId, 20; change_type, 5
## 
## Fixed effects:
##                   Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)         9.6279     0.9151  39.5192  10.521 5.07e-13 ***
## likelihood_rating  -0.1631     0.1672 244.3875  -0.975     0.33    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## liklhd_rtng -0.575
corr <- ma_RTs_likelihood_no_NA %>% 
  group_by(image) %>% 
  dplyr::summarize(detection_rt = mean(detection_rt), likelihood_rating = mean(likelihood_rating))
corr %>%
  ggscatter(y = "detection_rt", x = "likelihood_rating", ylab = "Change Detection RT (sec)", xlab = "Likelihood of Detecting Change", title = "N = 20", fill = "#f7a800", color = "#f7a800", add = "reg.line", cor.coef = TRUE, cor.coeff.args = list(method = "pearson", label.x = 1, label.y = 5, label.sep = "\n"), xlim = c(1, 5), ylim = c(0, 25), label = "image", font.label = c(5, "plain", "black"))