In this study, we compared the eye-movements of participants as they completed a joint attention game with a virtual partner or avatar. Half of the participants believed that the avatar was controlled by another human. The other half were correctly informed that the avatar was controlled by a computer algorithm.
The code imports trial and interest area reports and combines them into a single dataframe describing each trial. Accuracy and eyetracking data are presented via boxplots. Eyetracking data are analysed via ANOVAs and t-tests. We use non-parmetric stats to explore subject ratings of the task.
Caruana, N., Spirou, D., & Brock, J. (2017). Human agency beliefs influence joint attention behaviour. Manuscript submitted for publication.
library(reshape2)
library(tables)
library(ez)
library(knitr)
library(plyr)
library(gridExtra)
library(ReporteRs)
library(stringr)
sessionInfo()
## R version 3.3.2 (2016-10-31)
## Platform: x86_64-apple-darwin13.4.0 (64-bit)
## Running under: OS X El Capitan 10.11.6
##
## locale:
## [1] en_AU.UTF-8/en_AU.UTF-8/en_AU.UTF-8/C/en_AU.UTF-8/en_AU.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] stringr_1.2.0 ReporteRs_0.8.8 ReporteRsjars_0.0.2
## [4] gridExtra_2.2.1 plyr_1.8.4 knitr_1.15.1
## [7] ez_4.4-0 tables_0.8 Hmisc_4.0-2
## [10] ggplot2_2.2.1 Formula_1.2-1 survival_2.40-1
## [13] lattice_0.20-34 reshape2_1.4.2
##
## loaded via a namespace (and not attached):
## [1] Rcpp_0.12.9 png_0.1-7 assertthat_0.1
## [4] rprojroot_1.2 digest_0.6.12 mime_0.5
## [7] R6_2.2.0 backports_1.0.5 acepack_1.4.1
## [10] MatrixModels_0.4-1 evaluate_0.10 gdtools_0.1.3
## [13] lazyeval_0.2.0 minqa_1.2.4 data.table_1.10.4
## [16] SparseM_1.74 car_2.1-4 nloptr_1.0.4
## [19] R.utils_2.5.0 R.oo_1.21.0 rpart_4.1-10
## [22] Matrix_1.2-7.1 checkmate_1.8.2 rmarkdown_1.3
## [25] splines_3.3.2 lme4_1.1-12 foreign_0.8-67
## [28] htmlwidgets_0.8 munsell_0.4.3 shiny_1.0.0
## [31] httpuv_1.3.3 base64enc_0.1-3 mgcv_1.8-15
## [34] rvg_0.1.2 htmltools_0.3.5 nnet_7.3-12
## [37] tibble_1.2 htmlTable_1.9 MASS_7.3-45
## [40] R.methodsS3_1.7.1 grid_3.3.2 xtable_1.8-2
## [43] nlme_3.1-128 gtable_0.2.0 magrittr_1.5
## [46] scales_0.4.1 stringi_1.1.2 latticeExtra_0.6-28
## [49] xml2_1.1.1 RColorBrewer_1.1-2 tools_3.3.2
## [52] parallel_3.3.2 pbkrtest_0.4-6 yaml_2.1.14
## [55] colorspace_1.3-2 cluster_2.0.5 rJava_0.9-8
## [58] quantreg_5.29
Draws boxplot with overlaid datapoints
get_Boxplot <- function (Data, DV, yMin=0, yMax=1, yBreaks=0.2, yLabel="Test") {
ggplot(Data, aes_string(x="factor(Group)", y=DV, fill="Condition")) +
geom_boxplot(outlier.colour=NA) +
scale_fill_brewer(palette="BuGn") +
theme_minimal(base_size = 16, base_family = "") +
scale_x_discrete("") +
geom_point(position=position_jitterdodge(dodge.width=0.75, jitter.width=0.3), alpha=0.2, size=3) +
scale_y_continuous(yLabel,limits=c(yMin,yMax), breaks=seq(yMin,yMax,yBreaks))
}
Grabs the legend from one of the figures so it can be placed in a multi-panel figure Source: http://www.sthda.com/english/wiki/print.php?id=177
get_legend<-function(myggplot){
tmp <- ggplot_gtable(ggplot_build(myggplot))
leg <- which(sapply(tmp$grobs, function(x) x$name) == "guide-box")
legend <- tmp$grobs[[leg]]
return(legend)
}
We used Eyelink DataViewer software to export trial and interest area reports as text files.
Extract accuracy data from trial report.
AccData <- read.delim("data/Trial_Report_Accuracy.txt", stringsAsFactors=FALSE)
Data <- AccData[,c("RECORDING_SESSION_LABEL", "trial_id_num", "INDEX", "ERROR_TYPE")]
Data$SubNum <- as.numeric(lapply(Data$RECORDING_SESSION_LABEL, function(x) substr(x, 1, 3)))
Data$UTI <- Data$SubNum*1000 + Data$trial_id_num
Extract saccadic reaction times from interest area report and add to trial data.
SRTData <- read.delim("data/IA_Report_RJA_SRT.txt", stringsAsFactors=FALSE)
SRTData <- SRTData [(SRTData$IA_LABEL=="IA_BURGLAR "),]
SRTData$IPStartTime <- SRTData$IP_START_TIME - SRTData$TRIAL_START_TIME
SRTData$IA_FIRST_SACCADE_START_TIME[SRTData$IA_FIRST_SACCADE_START_TIME=="."] <- NA
SRTData$IA_FIRST_SACCADE_START_TIME <- as.numeric(SRTData$IA_FIRST_SACCADE_START_TIME)
SRTData$SRT <- SRTData$IA_FIRST_SACCADE_START_TIME - (SRTData$IP_START_TIME - SRTData$TRIAL_START_TIME)
SRTData$SubNum <- as.numeric(lapply(SRTData$RECORDING_SESSION_LABEL, function(x) substr(x, 1, 3)))
SRTData$UTI <- SRTData$SubNum*1000 + SRTData$trial_id_num
Data$FirstSaccadeSRT <- SRTData$SRT[match(Data$UTI, SRTData$UTI)]
Extract dwell times from interest area report and add to trial data.
DwellData <- read.delim("data/IA_Report_IJA_Dwell.txt", stringsAsFactors=FALSE)
DwellData <- DwellData [(DwellData$IA_LABEL=="IA_BURGLAR "),]
DwellData$SubNum <- as.numeric(lapply(DwellData$RECORDING_SESSION_LABEL, function(x) substr(x, 1, 3)))
DwellData$UTI <- DwellData$SubNum*1000 + DwellData$trial_id_num
Data$BurglarDwellTime <- DwellData$IA_DWELL_TIME[match(Data$UTI, DwellData$UTI)]
Extract premature saccade data from interest area report and add to trial data.
PremiData <- read.delim("data/IA_Report_IJA_Premi.txt", stringsAsFactors=FALSE)
PremiData <- PremiData [(PremiData$IA_LABEL=="IA_BURGLAR "),]
PremiData$PrematureSaccade <- lapply(PremiData$IA_FIXATION_COUNT, function(x) if(x==0){0}else{1})
PremiData$SubNum <- as.numeric(lapply(PremiData$RECORDING_SESSION_LABEL, function(x) substr(x, 1, 3)))
PremiData$UTI <- PremiData$SubNum*1000 + PremiData$trial_id_num
Data$PrematureSaccade <- PremiData$PrematureSaccade[match(Data$UTI, PremiData$UTI)]
Data$PrematureSaccade <- vapply(Data$PrematureSaccade, paste, collapse = ", ", character(1L)) # For some reason PrematureSaccade is identified as a list. This command "flattens" the list, allowing it to be written to csv.
Decode correct response / error type
Data$CorrectResponse <- ifelse(Data$ERROR_TYPE==0,1,0)
Data$SearchError <- ifelse(Data$ERROR_TYPE==1,1,0)
Data$TimeOut <- ifelse(Data$ERROR_TYPE==2,1,0)
Data$LocationError <- ifelse(Data$ERROR_TYPE==3,1,0)
Data$RecalibrationError <- ifelse(Data$ERROR_TYPE==4,1,0)
Add information about trials.
TrialInfo <- read.csv("data/TrialInfo.csv", stringsAsFactors=FALSE)
Data$Condition <- TrialInfo$Condition[match(Data$trial_id_num, TrialInfo$trial_id_num)]
Data$Interface <- TrialInfo$Interface[match(Data$trial_id_num, TrialInfo$trial_id_num)]
Data$SubjectRole <- TrialInfo$SubjectRole[match(Data$trial_id_num, TrialInfo$trial_id_num)]
Add information about subjects.
SubjectData <- read.csv("data/SubjectData.csv", stringsAsFactors=FALSE)
Data$Group <- SubjectData$Group[match(Data$SubNum, SubjectData$SubNum)]
Exclude two participants who weren’t deceived.
Data$Exclude <- SubjectData$Exclude[match(Data$SubNum, SubjectData$SubNum)]
Data <- Data [(Data$Exclude==0),]
Rename dataframe columns for consistency with other datasets.
names(Data)[names(Data)=="trial_id_num"] <- "TrialID"
names(Data)[names(Data)=="INDEX"] <- "TrialNum"
Remove columns that are now redundant.
Data <- subset(Data, select = -c(Exclude, RECORDING_SESSION_LABEL, ERROR_TYPE))
Write trial data to CSV file.
write.csv(Data, file="output/data/Data.csv", row.names=FALSE)
Identify factors
FactorNames <- c("SubNum", "TrialID", "Interface", "SubjectRole", "Group", "Condition")
Data[,FactorNames] <- colwise(as.factor)(Data[,FactorNames])
str(Data)
## 'data.frame': 10368 obs. of 16 variables:
## $ TrialID : Factor w/ 432 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ TrialNum : int 1 2 3 4 5 6 7 8 9 10 ...
## $ SubNum : Factor w/ 48 levels "101","102","103",..: 3 3 3 3 3 3 3 3 3 3 ...
## $ UTI : num 103001 103002 103003 103004 103005 ...
## $ FirstSaccadeSRT : num NA NA 486 NA NA 504 627 NA NA NA ...
## $ BurglarDwellTime : int 0 2956 0 1256 698 0 0 636 624 1938 ...
## $ PrematureSaccade : chr "0" "1" "0" "1" ...
## $ CorrectResponse : num 1 1 1 1 1 1 1 0 1 1 ...
## $ SearchError : num 0 0 0 0 0 0 0 0 0 0 ...
## $ TimeOut : num 0 0 0 0 0 0 0 1 0 0 ...
## $ LocationError : num 0 0 0 0 0 0 0 0 0 0 ...
## $ RecalibrationError: num 0 0 0 0 0 0 0 0 0 0 ...
## $ Condition : Factor w/ 4 levels "IJA","IJAc","RJA",..: 4 2 4 2 2 4 3 1 1 1 ...
## $ Interface : Factor w/ 2 levels "Control","Test": 1 1 1 1 1 1 2 2 2 2 ...
## $ SubjectRole : Factor w/ 2 levels "Initiate","Respond": 2 1 2 1 1 2 2 1 1 1 ...
## $ Group : Factor w/ 2 levels "Computer","Human": 1 1 1 1 1 1 1 1 1 1 ...
CurrentData <- Data [(Data$SubjectRole=="Respond"), ]
CurrentData$Condition <- factor(CurrentData$Condition)
CurrentData$ScreenIn <-ifelse((CurrentData$RecalibrationError==0 & CurrentData$SearchError==0 & !is.na(CurrentData$CorrectResponse)), 1, 0)
mean(CurrentData$ScreenIn)
## [1] 0.9876543
ScreenedData <- CurrentData [(CurrentData$ScreenIn==1), ]
BySubjectsData <- ddply(ScreenedData, .(Group, SubNum, Condition), summarise, ProportionCorrect = mean(CorrectResponse))
RJA_Acc_Boxplot <- get_Boxplot(BySubjectsData, "ProportionCorrect", yMin=0, yMax=1, yBreaks=0.2, yLabel="Proportion Correct")
RJA_Acc_Boxplot
CurrentData$ScreenIn <-ifelse((CurrentData$CorrectResponse==1 & CurrentData$FirstSaccadeSRT>=150 & CurrentData$FirstSaccadeSRT<=3000 & !is.na(CurrentData$FirstSaccadeSRT)), 1, 0)
mean(CurrentData$ScreenIn)
## [1] 0.8661265
ScreenedData <- CurrentData [(CurrentData$ScreenIn==1), ]
BySubjectsData <- ddply(ScreenedData, .(Group, SubNum, Condition), summarise, FirstSaccadeSRT = mean(FirstSaccadeSRT))
RJA_SRT_Boxplot <- get_Boxplot(BySubjectsData, "FirstSaccadeSRT", yMin=0, yMax=1200, yBreaks=200, yLabel="Reaction Time")
RJA_SRT_Boxplot
kable(ezANOVA(BySubjectsData, dv = FirstSaccadeSRT, wid = SubNum, within = .(Condition), between = Group, type = 3), row.names=FALSE, caption="First saccade reaction time", digits=3)
|
CurrentData <- Data [(Data$SubjectRole=="Initiate"), ]
CurrentData$Condition <- factor(CurrentData$Condition)
CurrentData$ScreenIn <-ifelse((CurrentData$RecalibrationError==0 & CurrentData$SearchError==0 & !is.na(CurrentData$CorrectResponse)), 1, 0)
mean(CurrentData$ScreenIn)
## [1] 0.9907407
ScreenedData <- CurrentData [(CurrentData$ScreenIn==1), ]
BySubjectsData <- ddply(ScreenedData, .(Group, SubNum, Condition), summarise, ProportionCorrect = mean(CorrectResponse))
IJA_Acc_Boxplot <- get_Boxplot(BySubjectsData, "ProportionCorrect", yMin=0, yMax=1, yBreaks=0.2, yLabel="Proportion Correct")
IJA_Acc_Boxplot
CurrentData$ScreenIn <-ifelse((CurrentData$CorrectResponse==1 & CurrentData$BurglarDwellTime>=150 & CurrentData$BurglarDwellTime<=3000 & !is.na(CurrentData$BurglarDwellTime)), 1, 0)
mean(CurrentData$ScreenIn)
## [1] 0.9371142
ScreenedData <- CurrentData [(CurrentData$ScreenIn==1), ]
BySubjectsData <- ddply(ScreenedData, .(Group, SubNum, Condition), summarise, BurglarDwellTime = mean(BurglarDwellTime))
IJA_Dwell_Boxplot <- get_Boxplot (BySubjectsData, "BurglarDwellTime", yMin=0, yMax=1800, yBreaks=200, yLabel="Dwell Time / ms")
IJA_Dwell_Boxplot
kable(ezANOVA(BySubjectsData, dv = BurglarDwellTime, wid = SubNum, within = .(Condition), between = Group, type = 3), row.names=FALSE, caption="ANOVA: Burglar Dwell Time", digits=3)
|
TTest.Dwell <- data.frame(Measure=c("Dwell Time", "Dwell Time"),
Group=c("Computer", "Human"),
df=numeric(2),
t=numeric(2),
p=numeric(2))
for (i in 1:2) {
Group <- TTest.Dwell$Group[i]
TTest <- t.test(BurglarDwellTime ~ Condition, BySubjectsData[(BySubjectsData$Group==Group),], paired=TRUE)
TTest.Dwell$df[i] <- TTest$parameter
TTest.Dwell$t[i] <- TTest$statistic
TTest.Dwell$p[i] <- TTest$p.value
}
kable(TTest.Dwell, row.names=FALSE, caption="T-tests: Effect of Condition on Dwell Times", digits=3)
Measure | Group | df | t | p |
---|---|---|---|---|
Dwell Time | Computer | 23 | 0.888 | 0.383 |
Dwell Time | Human | 23 | 5.581 | 0.000 |
CurrentData$ScreenIn <-ifelse((CurrentData$RecalibrationError==0 & CurrentData$SearchError==0 & !is.na(CurrentData$PrematureSaccade)), 1, 0)
mean(CurrentData$ScreenIn)
## [1] 0.9907407
ScreenedData <- CurrentData [(CurrentData$ScreenIn==1), ]
ScreenedData$PrematureSaccade <- as.numeric(ScreenedData$PrematureSaccade)
BySubjectsData <- ddply(ScreenedData, .(Group, SubNum, Condition), summarise, PrematureSaccades = mean(PrematureSaccade))
IJA_Premi_Boxplot <- get_Boxplot (BySubjectsData, "PrematureSaccades", yMin=0, yMax=1, yBreaks=0.2, yLabel="Proportion of Premature Saccades")
IJA_Premi_Boxplot
kable(ezANOVA(BySubjectsData, dv = PrematureSaccades, wid = SubNum, within = .(Condition), between = Group, type = 3), row.names=FALSE, caption="PrematureSaccade", digits=3)
|
TTest.Premi <- data.frame(Measure=c("Premature Saccades", "Premature Saccades"),
Group=c("Computer", "Human"),
df=numeric(2),
t=numeric(2),
p=numeric(2))
for (i in 1:2) {
Group <- TTest.Premi$Group[i]
TTest <- t.test(PrematureSaccades ~ Condition, BySubjectsData[(BySubjectsData$Group==Group),], paired=TRUE)
TTest.Premi$df[i] <- TTest$parameter
TTest.Premi$t[i] <- TTest$statistic
TTest.Premi$p[i] <- TTest$p.value
}
kable(TTest.Premi, row.names=FALSE, caption="T-tests: Effect of Condition on Premature Saccades", digits=3)
Measure | Group | df | t | p |
---|---|---|---|---|
Premature Saccades | Computer | 23 | 2.109 | 0.046 |
Premature Saccades | Human | 23 | 6.145 | 0.000 |
RatingsData <- subset(SubjectData[(SubjectData$Exclude==0),], select = -c(Exclude, Sex, Age, EdHAnd, Block.1, Block.2, ET_Accuracy) )
RatLong <- melt(RatingsData, id=c("Group", "SubNum"), variable.name="Measure", value.name="Rating")
RatLong$Question <- as.factor(str_split_fixed(RatLong$Measure, "_", 2)[,1])
RatLong$Condition <- as.factor(str_split_fixed(RatLong$Measure, "_", 2)[,2])
RatLong$Question <- mapvalues(RatLong$Question, from = c("PreferVirtual", "PreferAlone", "Feel", "Behave", "Appear"), to = c("Prefer Virtual", "Prefer Alone", "Feel Humanlike", "Behave Humanlike", "Appear Humanlike"))
QuestionList <- c("Prefer Virtual", "Prefer Alone", "Behave Humanlike", "Appear Humanlike", "Feel Humanlike", "Cooperative", "Pleasant","Intuitive","Natural","Difficult")
RatingsBoxplot <- ggplot(RatLong[(RatLong$Condition!="Arrow"),], aes(x=factor(Question), y=Rating, fill=Group)) +
geom_boxplot(outlier.colour=NA) +
theme_minimal(base_size = 16, base_family = "") +
xlim(QuestionList) +
xlab("Question") +
scale_y_continuous("Rating",limits=c(1,10), breaks=seq(1,10,1)) +
scale_fill_brewer(palette="Blues")+
coord_flip()
RatingsBoxplot
Wilcoxon.DF <- data.frame(Question=rev(QuestionList),
W=numeric(length(QuestionList)),
p=numeric(length(QuestionList)))
for (i in 1:length(QuestionList)) {
QuestionData <- RatLong[(RatLong$Question==Wilcoxon.DF$Question[i]) & (RatLong$Condition=="Alan" | RatLong$Condition=="Both"),]
Wilcox <- wilcox.test(Rating ~ Group, data = QuestionData, paired=FALSE, exact=FALSE)
Wilcoxon.DF$W[i] <- Wilcox$statistic
Wilcoxon.DF$p[i] <- Wilcox$p.value
}
kable(Wilcoxon.DF, row.names=FALSE, caption="Wilcoxon rank sum test with continuity correction", digits=3)
Question | W | p |
---|---|---|
Difficult | 346.5 | 0.209 |
Natural | 216.5 | 0.139 |
Intuitive | 384.5 | 0.045 |
Pleasant | 188.5 | 0.038 |
Cooperative | 193.0 | 0.039 |
Feel Humanlike | 206.5 | 0.090 |
Appear Humanlike | 284.5 | 0.950 |
Behave Humanlike | 214.0 | 0.124 |
Prefer Alone | 358.5 | 0.142 |
Prefer Virtual | 299.5 | 0.817 |
TTest.DF <- rbind(TTest.Dwell, TTest.Premi)
doc <- docx()
doc <- addTitle(doc, "Table 2: ", level=2)
doc <- addFlexTable( doc, vanilla.table(TTest.DF))
writeDoc(doc, file = "output/doc/Tables.docx")
legend <- get_legend(RJA_Acc_Boxplot)
png(filename="output/figures/Respond.png", width=650, height=400)
grid.arrange (RJA_Acc_Boxplot + theme(legend.position="none"), RJA_SRT_Boxplot + theme(legend.position="none"), legend,
ncol=3, nrow =1, widths=c(2.3, 2.3, 0.8), heights=c(2.3))
dev.off()
## quartz_off_screen
## 2
legend <- get_legend(IJA_Acc_Boxplot)
png(filename="output/figures/Initiate.png", width=950, height=400)
grid.arrange (IJA_Acc_Boxplot + theme(legend.position="none"), IJA_Dwell_Boxplot + theme(legend.position="none"), IJA_Premi_Boxplot + theme(legend.position="none"), legend,
ncol=4, nrow =1, widths=c(2.3, 2.3, 2.3, 0.8), heights=c(2.3))
dev.off()
## quartz_off_screen
## 2
png(filename="output/figures/Ratings.png", width=600, height=600)
RatingsBoxplot
dev.off()
## quartz_off_screen
## 2
```