---
title: "Main analysis: naming task"
author: "Marlou Rasenberg"
date: '(this version: `r format(Sys.Date())`)'
output:
  github_document:
    toc: true
  html_document:
    toc: true
    toc_float: true
editor_options:
  chunk_output_type: console
---

```{r global_options, include=FALSE}
knitr::opts_chunk$set(echo=TRUE, warning=FALSE, message=FALSE, fig.path='3_out/')
```

```{r packages, results='hide',include=F}
# Packages
list.of.packages <- c("dplyr", "tidyr","lsa", "ggplot2","cowplot","ggpubr", "stringr", "RColorBrewer")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)>0) install.packages(new.packages)
lapply(list.of.packages, require, character.only=T)
rm(list.of.packages,new.packages)

```

This script supports the paper 'The primacy of multimodal alignment in converging on shared symbols for novel referents' by Marlou Rasenberg, Asli Özyürek, Sara Bögels and Mark Dingemanse (for the CABB team). To increase readability, not all code chunks present in the .Rmd source are shown in the output.


## Aim
With this script we analyse the similarity of the labels that dyad members used to name the Fribbles in the naming task.  
We compare:  
* pre to post: are the names more similar after compared to before the interaction?  
* real to random dyads: is the naming similarity (post) contingent on the interaction, or simply a result of spending time with the stimuli/referring to them over multiple rounds?  

## Data
The input data contains the names that participants provided for the Fribble in the naming task. They are cleaned, i.e. lemmatized and compounds split.

```{r read in data}
df <- read.csv2("1_data/naming_task/naming_task_clean.csv")
```

```{r preliminaries, include=F}
#reorder pre and post (pre first)
df$Session <- as.factor(df$Session)
df$Session <- relevel(df$Session, "pre")
```

## Processing

### Shuffle the names (non-dyad comparisons)

```{r shuffle names}
# simple shifter function
shiftthis <- function(x, n = 1) { if (n == 0) x else c(tail(x, -n), head(x, n)) }

# shift names
df <- df %>%
  group_by(Session,Fribble_nr) %>%
  mutate(Name_b_shifted = shiftthis(Name_b,1))

#change order columns
df <- df[,c(1:5,8,6:7)]

# sanity check: does shifting work?
df %>% filter(Fribble_nr==9,Session=="post") %>%
  dplyr::select(Pair,Session,Fribble_nr,Name_a,Name_b,Name_b_shifted)
```

### Calculate cosine similarity
Here we calculate the cosine similarity of the names that participants provided

First for real dyads:
```{r cosine similarity real dyads}
#create column for cosine similarity
df$name_sim_realdyads = NA
df$name_sim_realdyads <- as.numeric(df$name_sim_realdyads)

#loop over the files to calculate cosine similarity
for (x in 1:nrow(df)){
  
#make a list of words from both A and B
a_list = as.list(unlist(strsplit(df[x,]$Name_a, ' ')))
b_list = as.list(unlist(strsplit(df[x,]$Name_b, ' ')))

#make a list of all unique words from both A and B
all_words = unique(c(a_list,b_list))

#make two new lists for a and b
a_numlist = c()
b_numlist = c()
  #add numbers to the lists for a and b (1 if in the list, 0 if not)
  for (y in all_words){
    if(y %in% a_list){
      a_numlist = c(a_numlist,1)
    } else {a_numlist = c(a_numlist, 0)}
    if (y %in% b_list){
      b_numlist = c(b_numlist,1)
    } else {b_numlist = c(b_numlist, 0)}
  }

#calculate cosine similarity of a and b
df[x,]$name_sim_realdyads <- cosine(a_numlist, b_numlist)[1,1]
}

```

Then for non-dyads:
```{r cosine similarity non-dyads}
#create column for cosine similarity
df$name_sim_nondyads = NA
df$name_sim_nondyads <- as.numeric(df$name_sim_nondyads)

#loop over the files to calculate cosine similarity

for (x in 1:nrow(df)){
  
#make a list of words from both A and B
a_list = as.list(unlist(strsplit(df[x,]$Name_a, ' ')))
b_s_list = as.list(unlist(strsplit(df[x,]$Name_b_shifted, ' ')))

#make a list of all unique words from both A and B
all_words = unique(c(a_list,b_s_list))

#make two new lists for a and b
a_numlist = c()
b_s_numlist = c()
  #add numbers to the lists for a and b (1 if in the list, 0 if not)
  for (y in all_words){
    if(y %in% a_list){
      a_numlist = c(a_numlist,1)
    } else {a_numlist = c(a_numlist, 0)}
    if (y %in% b_s_list){
      b_s_numlist = c(b_s_numlist,1)
    } else {b_s_numlist = c(b_s_numlist, 0)}
  }

#calculate cosine similarity of a and b
df[x,]$name_sim_nondyads <- cosine(a_numlist, b_s_numlist)[1,1]
}
```

## Analysis (§3.1 in ms)

### Plots
```{r violin plots - Fig. 5 in ms}
#set colours
nb.cols <- 10
mycolors <- colorRampPalette(brewer.pal(8, "Set2"))(nb.cols)


#add 'datapoint' index to be able to plot lines from pre to post
df$PairFribble <- as.character("")
for (x in 1:nrow(df)) {
  df[x,]$PairFribble <- paste0(as.character(df[x,]$Pair), "+",as.character(df[x,]$Fribble_nr))
}


violin_realdyads_lines <- df %>%
  mutate(Fribble_nr = as.factor(Fribble_nr)) %>%
  ggplot(aes(Session,name_sim_realdyads,color=Session)) +
  ggtitle("Real dyads") + 
  ylab("naming similarity score")+
  xlab("session")+
  scale_color_manual(values=mycolors)+
  geom_violin(show.legend = F, colour="darkgrey", fill="lightgrey", alpha=.2, lwd=0.2) +
  geom_jitter(aes(colour=as.factor(Pair)),show.legend=F, height=0.015, width=.015, size=0.8)+
  geom_line(aes(group=as.factor(PairFribble), colour=as.factor(Pair), y = jitter(name_sim_realdyads, .3), x = jitter(as.numeric(Session), .2)), alpha=.4)+
  theme_light()+
        theme(panel.grid.major = element_blank(),
              panel.grid.minor = element_blank(),
              panel.border = element_blank(),
              axis.line.x = element_line(colour="grey", size=0.2), 
              axis.line.y=element_line(colour="grey", size=0.2))+
  theme(plot.title = element_text(hjust = 0.5, size=7))+
  theme(axis.title.y=element_text(size=6), axis.text.y=element_text(size=6, colour="black"),axis.title.x=element_text(size=6), axis.text.x=element_text(size=6, colour="black"), title=element_text(size=7))+
  theme(legend.position="none")



violin_nondyads_lines <- df %>%
  mutate(Fribble_nr = as.factor(Fribble_nr)) %>%
  ggplot(aes(Session,name_sim_nondyads,color=Session)) +
  ggtitle("Non-dyads") + 
  ylab("naming similarity score")+
  xlab("session")+
  scale_color_manual(values=c("darkgrey","darkgrey"))+
  geom_violin(show.legend = F, colour="darkgrey", fill="lightgrey", alpha=.2, lwd=0.2) +
  geom_line(aes(group=as.factor(PairFribble), y = jitter(name_sim_nondyads, .3), x = jitter(as.numeric(Session), .2)), alpha=.4, colour="darkgrey")+
  geom_jitter(show.legend=F, height=0.015, width=.015, size=0.8)+
  theme_light()+
        theme(panel.grid.major = element_blank(),
              panel.grid.minor = element_blank(),
              panel.border = element_blank(),
              axis.line.x = element_line(colour="grey", size=0.2), 
              axis.line.y=element_line(colour="grey", size=0.2))+
  theme(plot.title = element_text(hjust = 0.5, size=7))+
  theme(axis.title.y=element_text(size=6, colour="black"), axis.text.y=element_text(size=6),axis.title.x=element_text(size=6), axis.text.x=element_text(size=6, colour="black"), title=element_text(size=7))+
  theme(legend.position="none")


p_violins_withlines <- plot_grid(violin_realdyads_lines, violin_nondyads_lines, labels=c("A","B"), label_size=7)
p_violins_withlines

```

### Hypothesis test
We test whether name similarity increased from pre to post (for real and non-real dyads).

#### Real dyads
First, we inspect the distribution of the data
```{r distribution pre-post difference scores - real dyads}
#check if difference pre vs. post is normally distributed

d_diff_real <- density((df[df$Session=="post",]$name_sim_realdyads)-(df[df$Session=="pre",]$name_sim_realdyads), na.rm=T)
plot(d_diff_real, main="")
ggqqplot((df[df$Session=="post",]$name_sim_realdyads)-(df[df$Session=="pre",]$name_sim_realdyads), na.rm=T)

shapiro.test(df[df$Session=="post",]$name_sim_realdyads-df[df$Session=="pre",]$name_sim_realdyads)

```
The density and quantile-quantile plot indicate that the difference scores (naming similarity post - pre) are not normally distributed.  
This is confirmed by the Shapiro-Wilk normality test.  
Hence, we use the Wilcoxon signed-rank test.

```{r Wilcoxon signed-rank test - real dyads}
wtest_real <- wilcox.test(df[df$Session=="post",]$name_sim_realdyads, df[df$Session=="pre",]$name_sim_realdyads, paired = TRUE, alternative = "two.sided") 
wtest_real

Zstat_real<-qnorm(wtest_real$p.value/2)
Zeffectsize_real <- abs(Zstat_real)/sqrt(80) #0.7

print(paste0("Effect size: ", round(Zeffectsize_real, digits=2)))

#https://www.sheffield.ac.uk/polopoly_fs/1.714576!/file/stcp-marquier-WilcoxonR.pdf

summary(df[df$Session=="pre",]$name_sim_realdyads) #mean=0.07, median=0
summary(df[df$Session=="post",]$name_sim_realdyads) #mean=0.46, meidan=0.41

```

#### Non-dyads
Again, we inspect the distribution of the data
```{r distribution pre-post difference scores - non-dyads}
#check if difference pre vs. post is normally distributed

d_diff_non <- density((df[df$Session=="post",]$name_sim_nondyads)-(df[df$Session=="pre",]$name_sim_nondyads), na.rm=T)
plot(d_diff_non, main="")
ggqqplot((df[df$Session=="post",]$name_sim_nondyads)-(df[df$Session=="pre",]$name_sim_nondyads), na.rm=T)

shapiro.test(df[df$Session=="post",]$name_sim_nondyads-df[df$Session=="pre",]$name_sim_nondyads)
```
The density and quantile-quantile plot indicate that the difference scores (naming similarity post - pre) are not normally distributed.  
This is confirmed by the Shapiro-Wilk normality test.  
Hence, we use the Wilcoxon signed-rank test.

```{r Wilcoxon signed-rank test - non-dyads}
wtest_non <- wilcox.test(df[df$Session=="post",]$name_sim_nondyads, df[df$Session=="pre",]$name_sim_nondyads, paired = TRUE, alternative = "two.sided") 
wtest_non

summary(df[df$Session=="pre",]$name_sim_nondyads) #mean=0.11, median=0
summary(df[df$Session=="post",]$name_sim_nondyads) #mean=0.05, median=0

```



```{r time on task, results='hide',include=F}
df_ToT <- df %>%
  group_by(Pair, Session) %>%
  summarize(ToT_a = sum(RT_a), ToT_b = sum(RT_b))

#wide to long
df_ToT <- gather(df_ToT, participant, ToT, ToT_a:ToT_b)

df_ToT$participant <- gsub("ToT_","",df_ToT$participant)

mean_ToT <- df_ToT %>%
  group_by(Session) %>%
  summarize(ToT_mean_msec = mean(ToT), ToT_min_msec = min(ToT), ToT_max_msec = max(ToT))

mean_ToT$ToT_mean_min = mean_ToT$ToT_mean_msec/60000
mean_ToT$ToT_min_min = mean_ToT$ToT_min_msec/60000
mean_ToT$ToT_max_min = mean_ToT$ToT_max_msec/60000
  
```

```{r save plot, results='hide',include=T}
ggsave(p_violins_withlines, file="2_figures/Fig5.png",  bg="white", width=143, height=60, units="mm", dpi=600)
ggsave(p_violins_withlines, file="2_figures/Fig5.pdf",  bg="white", width=143, height=60, units="mm")
```

```{r export df, results='hide', include=F}
#export df with similarity scores ==> used in script supplementary_analysis_alignment-namingtask
write.csv2(df, file="1_data/naming_task/naming_task_Routput.csv", row.names=F)
```

```{r session info, include=F}
sessionInfo()

# R version 4.0.2 (2020-06-22)
# Platform: x86_64-w64-mingw32/x64 (64-bit)
# Running under: Windows 10 x64 (build 19043)
# 
# Matrix products: default
# 
# locale:
# [1] LC_COLLATE=Dutch_Netherlands.1252  LC_CTYPE=Dutch_Netherlands.1252    LC_MONETARY=Dutch_Netherlands.1252 LC_NUMERIC=C                      
# [5] LC_TIME=Dutch_Netherlands.1252    
# 
# attached base packages:
# [1] stats     graphics  grDevices utils     datasets  methods   base     
# 
# other attached packages:
# [1] RColorBrewer_1.1-2 stringr_1.4.0      ggpubr_0.4.0       cowplot_1.0.0      ggplot2_3.3.1      lsa_0.73.2         SnowballC_0.7.0    tidyr_1.1.0       
# [9] dplyr_1.0.0       
# 
# loaded via a namespace (and not attached):
#  [1] tidyselect_1.1.0  xfun_0.14         purrr_0.3.4       haven_2.3.1       lattice_0.20-41   carData_3.0-4     colorspace_1.4-1  vctrs_0.3.0      
#  [9] generics_0.0.2    htmltools_0.4.0   yaml_2.2.1        utf8_1.1.4        rlang_0.4.11      pillar_1.4.4      foreign_0.8-80    glue_1.4.1       
# [17] withr_2.2.0       readxl_1.3.1      lifecycle_0.2.0   munsell_0.5.0     ggsignif_0.6.0    gtable_0.3.0      cellranger_1.1.0  zip_2.0.4        
# [25] evaluate_0.14     labeling_0.3      knitr_1.28        rio_0.5.16        forcats_0.5.0     curl_4.3          fansi_0.4.1       broom_0.5.6      
# [33] Rcpp_1.0.7        scales_1.1.1      backports_1.1.7   abind_1.4-5       farver_2.0.3      hms_0.5.3         digest_0.6.25     stringi_1.4.6    
# [41] openxlsx_4.1.5    rstatix_0.6.0     grid_4.0.2        cli_2.0.2         tools_4.0.2       magrittr_1.5      tibble_3.0.1      crayon_1.3.4     
# [49] car_3.0-8         pkgconfig_2.0.3   ellipsis_0.3.1    data.table_1.12.8 assertthat_0.2.1  rmarkdown_2.2     rstudioapi_0.11   R6_2.4.1         
# [57] nlme_3.1-148      compiler_4.0.2 
```

