---
title: "007_VerbAnalysis_3_3_20"
author: "Katharine Aveni"
date: "October 1, 2021"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

```{r Load Libraries & Data}

sentence_window_clean <- readRDS("007_AVT_data_clean_with_anticipated.rds", refhook = NULL)

#load required packages

library("Matrix")
library("plyr")
library("tidyr")
library("lme4")
library(lmerTest)
library("ggplot2")
library("eyetrackingR")
library("dplyr")
```


#Subset into verb window (offset by 200ms)
```{r}
verb_window <- subset_by_window(sentence_window_clean,
                                window_start_time = 658,
                                window_end_time = 1158,
                                rezero = TRUE, remove = TRUE)

#names(verb_window)[names(verb_window) == "file"] <- "sentence"
```
Extended verb window 37 ms into 2nd article in order to have even-size bins

##Verb window
```{r Create verb window time series}

verb_time_seq <-make_time_sequence_data(verb_window, time_bin_size = 50,
                                   predictor_columns = c("Group", "trial_name"),
                                   aois = c("Target", "ActionDis", "AgentDis", "Unrelated")
)


```


###Interpolate values for set 9 trials for subjects 2-6.
007_02: list 9; 007_03 & 007_04: list 10; 007_05: list 1; 007_06: list 2

Each chunk stores vectors for: Proportions for each bin in surrounding 2 trials (unless surrounding trial was a question trial or pre-anticipated trial, in which case the next closest trial was stored).
Stored vectors were averaged and written to new trial.
```{r Create vectors for interpolated trials}

verb_time_seq <- as.data.frame(verb_time_seq)

#SUBJECT 007_05
#interpolate trial 9_3
`vi_05_3_2` <- subset(verb_time_seq, Subject == "007_05" & trial_name == "ExpA_Set3_Item2", 
                    select = c("Prop"))
`vi_05_7_3` <- subset(verb_time_seq, Subject == "007_05" & trial_name == "ExpA_Set7_Item3", 
                    select = c("Prop"))
verb_time_seq[which(verb_time_seq$Subject == "007_05" & verb_time_seq$trial_name == "ExpA_Set9_Item3"), "Prop"] <- ((`vi_05_3_2` + `vi_05_7_3`)/2)


#interpolate trial 9_2
"vi_05_7_2" <- subset(verb_time_seq, Subject == "007_05" & trial_name == "ExpA_Set7_Item2", 
                    select = c("Prop"))
"vi_05_12_1" <- subset(verb_time_seq, Subject == "007_05" & trial_name == "ExpA_Set12_Item1", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_05" & verb_time_seq$trial_name == "ExpA_Set9_Item2"), "Prop"] <- ((vi_05_7_2 + vi_05_12_1) / 2)

#interpolate trial 9_4
"vi_05_11_4" <- subset(verb_time_seq, Subject == "007_05" & trial_name == "ExpA_Set11_Item4", 
                    select = c("Prop"))
"vi_05_12_2" <- subset(verb_time_seq, Subject == "007_05" & trial_name == "ExpA_Set12_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_05" & verb_time_seq$trial_name == "ExpA_Set9_Item4"), "Prop"] <- ((vi_05_11_4 + vi_05_12_2) / 2)

#SUBJECT 007_02
#interpolate trial 9_4
"vi_02_10_2" <- subset(verb_time_seq, Subject == "007_02" & trial_name == "ExpA_Set10_Item2", 
                    select = c("Prop"))
"vi_02_3_2" <- subset(verb_time_seq, Subject == "007_02" & trial_name == "ExpA_Set3_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_02" & verb_time_seq$trial_name == "ExpA_Set9_Item4"), "Prop"] <- ((vi_02_10_2 + vi_02_3_2) / 2)

#interpolate trial 9_1
"vi_02_3_2" <- subset(verb_time_seq, Subject == "007_02" & trial_name == "ExpA_Set3_Item2", 
                    select = c("Prop"))
"vi_02_5_3" <- subset(verb_time_seq, Subject == "007_02" & trial_name == "ExpA_Set5_Item3", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_02" & verb_time_seq$trial_name == "ExpA_Set9_Item1"), "Prop"] <- ((vi_02_3_2 + vi_02_5_3) / 2)

#interpolate trial 9_2
"vi_02_3_4" <- subset(verb_time_seq, Subject == "007_02" & trial_name == "ExpA_Set3_Item4", 
                    select = c("Prop"))
"vi_02_5_1" <- subset(verb_time_seq, Subject == "007_02" & trial_name == "ExpA_Set5_Item1", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_02" & verb_time_seq$trial_name == "ExpA_Set9_Item2"), "Prop"] <- ((vi_02_3_4 + vi_02_5_1) / 2)

#SUBJECT 007_06
#interpolate trial 9_4
"vi_06_2_1" <- subset(verb_time_seq, Subject == "007_06" & trial_name == "ExpA_Set2_Item1", 
                    select = c("Prop"))
"vi_06_7_2" <- subset(verb_time_seq, Subject == "007_06" & trial_name == "ExpA_Set7_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_06" & verb_time_seq$trial_name == "ExpA_Set9_Item4"), "Prop"] <- ((vi_06_2_1 + vi_06_7_2) / 2)

#interpolate trial 9_2
"vi_06_7_4" <- subset(verb_time_seq, Subject == "007_06" & trial_name == "ExpA_Set7_Item4", 
                    select = c("Prop"))
"vi_06_4_2" <- subset(verb_time_seq, Subject == "007_06" & trial_name == "ExpA_Set4_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_06" & verb_time_seq$trial_name == "ExpA_Set9_Item2"), "Prop"] <- ((vi_06_7_4 + vi_06_4_2) / 2)

#SUBJECT 007_03
#interpolate trial 9_4
"vi_03_11_2" <- subset(verb_time_seq, Subject == "007_03" & trial_name == "ExpA_Set11_Item2", 
                    select = c("Prop"))
"vi_03_12_1" <- subset(verb_time_seq, Subject == "007_03" & trial_name == "ExpA_Set12_Item1", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_03" & verb_time_seq$trial_name == "ExpA_Set9_Item4"), "Prop"] <- ((vi_03_11_2 + vi_03_12_1) / 2)

#interpolate trial 9_2
"vi_03_6_2" <- subset(verb_time_seq, Subject == "007_03" & trial_name == "ExpA_Set6_Item2", 
                    select = c("Prop"))
"vi_03_5_1" <- subset(verb_time_seq, Subject == "007_03" & trial_name == "ExpA_Set5_Item1", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_03" & verb_time_seq$trial_name == "ExpA_Set9_Item2"), "Prop"] <- ((vi_03_6_2 + vi_03_5_1) / 2)

#interpolate trial 9_3
"vi_03_7_4" <- subset(verb_time_seq, Subject == "007_03" & trial_name == "ExpA_Set7_Item4", 
                    select = c("Prop"))
"vi_03_4_2" <- subset(verb_time_seq, Subject == "007_03" & trial_name == "ExpA_Set4_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_03" & verb_time_seq$trial_name == "ExpA_Set9_Item3"), "Prop"] <- ((vi_03_7_4 + vi_03_4_2) / 2)

#SUBJECT 007_04
#interpolate trial 9_4
"vi_04_11_2" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set11_Item2", 
                    select = c("Prop"))
"vi_04_7_3" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set7_Item3", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_04" & verb_time_seq$trial_name == "ExpA_Set9_Item4"), "Prop"] <- ((vi_04_11_2 + vi_04_7_3) / 2)

#interpolate trial 9_1
"vi_04_10_4" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set10_Item4", 
                    select = c("Prop"))
"vi_04_1_2" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set1_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_04" & verb_time_seq$trial_name == "ExpA_Set9_Item1"), "Prop"] <- ((vi_04_10_4 + vi_04_1_2) / 2)

#interpolate trial 9_2
"vi_04_6_2" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set6_Item2", 
                    select = c("Prop"))
"vi_04_5_1" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set5_Item1", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_04" & verb_time_seq$trial_name == "ExpA_Set9_Item2"), "Prop"] <- ((vi_04_6_2 + vi_04_5_1) / 2)

#interpolate trial 9_3
"vi_04_1_3" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set1_Item3", 
                    select = c("Prop"))
"vi_04_4_2" <- subset(verb_time_seq, Subject == "007_04" & trial_name == "ExpA_Set4_Item2", 
                    select = c("Prop"))
verb_time_seq[ which(verb_time_seq$Subject == "007_04" & verb_time_seq$trial_name == "ExpA_Set9_Item3"), "Prop"] <- ((vi_04_1_3 + vi_04_4_2) / 2)
```


```{r Save clean, interpolated verb time sequence}
#save the cleaned file to run in verb/object window analysis scripts. Or plot/analyse overall sentence below.
#saveRDS(verb_time_seq, file = "007_verb_data_interpolated.rds")
saveRDS(verb_time_seq, file = "007_verb_data_interpolated_with_anticipated.rds")

```
Saved as of 10/1/21 (just with pre-anticipated trials included)

```{r load clean, interpolated verb time sequence}
#verb_time_seq <- readRDS("007_verb_data_clean_interpolated.rds", refhook = NULL)
verb_time_seq <- readRDS("007_verb_data_interpolated_with_anticipated.rds", refhook = NULL)

verb_time_seq$Group <- as.factor(verb_time_seq$Group)
verb_time_seq$AOI <- as.factor(verb_time_seq$AOI)


#load required packages

library("Matrix")
library("plyr")
library("tidyr")
library("lme4")
library(lmerTest)
library("ggplot2")
library("eyetrackingR")
library("dplyr")
#library(sjstats)
```

##Final prep for modeling
```{r final prep for modeling}
#label high vs. low action verbs
verb_time_seq %>% 
  mutate(mt_content = ifelse(file == "boxer_hits_punchingbag" | file == "driver_hits_brakes" |
                               file == "cat_catches_mouse" | file == "	child_catches_football" |
                               file == "dolphin_jumps_waves" | file == "horse_jumps_fence" |
                               file == "lifeguard_saves_swimmer" | file == "detective_hunts_fugitive" |
                               file == "lion_hunts_zebra" | file == "lion_crosses_jungle" |
                               file == "detective_crosses_courtroom" | file == "mechanic_organizes_tools" |
                               file == "waitress_organizes_dishes", "Hi", "Lo")) -> verb_time_seq

verb_time_seq$mt_content <- as.factor(verb_time_seq$mt_content)

#create column for natural quadratic time polynomial, if using
verb_time_seq %>% 
  mutate(Time2 = Time^2) %>% 
#create column of proportions rounded to 0 or 1
  mutate(BiProp = ifelse(Prop > .5, 1, 0)) -> verb_time_seq

#center factors
contrasts(verb_time_seq$mt_content) <- c(.5, -.5)
contrasts(verb_time_seq$Group) <- c(-.5, .5)

#create dataframe with average proportions from each subject, splits AOIs into different columns
verb_time_seq %>%
  filter(!is.na(Prop)) %>%
  group_by(Subject, AOI) %>% 
  summarise(MeanProp = mean(Prop), y=sum(Prop), N=length(Prop), Time=min(Time)) %>% 
  spread(AOI, MeanProp)-> verb_time_seq_avg

#create dataframe that excludes set 9 trials for early participants, for use in motion content analysis
verb_time_seq_9rm <- verb_time_seq[!((verb_time_seq$Subject=="007_02" | verb_time_seq$Subject=="007_03" | verb_time_seq$Subject=="007_04" | verb_time_seq$Subject=="007_05" | verb_time_seq$Subject=="007_06") & (verb_time_seq$trial_name=="ExpA_Set9_Item1" | verb_time_seq$trial_name=="ExpA_Set9_Item2" | verb_time_seq$trial_name=="ExpA_Set9_Item3" | verb_time_seq$trial_name=="ExpA_Set9_Item4")),]

#subset the data into target AOIs only
`target_verb_time_seq` <- subset(verb_time_seq, AOI == "Target", 
                    select = c("Prop", "Subject", "Group", "ot1", "ot2", "trial_name", "Time", "Time2", "BiProp", "mt_content"))

#target AOIs only, minus early Ss set 9
`target_verb_time_seq_9rm` <- subset(verb_time_seq_9rm, AOI == "Target", 
                    select = c("Prop", "Subject", "Group", "ot1", "ot2", "trial_name", "Time", "Time2", "BiProp", "mt_content"))

#PD only, target AOIs only
`PD_target_verb_time_seq` <- subset(verb_time_seq, (AOI == "Target" & Group == "PD"), 
                    select = c("Prop", "Subject", "Group", "ot1", "ot2", "trial_name", "Time", "Time2", "BiProp", "mt_content"))

#subset the data into agent-distractor AOIs only
`agentdis_verb_time_seq` <- subset(verb_time_seq, AOI == "AgentDis", 
                    select = c("Prop", "Subject", "Group", "ot1", "ot2", "trial_name", "Time", "Time2", "BiProp", "mt_content"))

#subset the data into action-distractor AOIs only
`actiondis_verb_time_seq` <- subset(verb_time_seq, AOI == "ActionDis", 
                    select = c("Prop", "Subject", "Group", "ot1", "ot2", "trial_name", "Time", "Time2", "BiProp", "mt_content"))

#create dataframe with average proportions from each subject, splits AOIs into different columns
target_verb_time_seq %>%
  filter(!is.na(Prop)) %>%
  group_by(Subject, Time, Group) %>% 
  summarise(MeanProp = mean(Prop), y=sum(Prop), N=length(Prop), Time=min(Time)) -> target_verb_time_seq_avg

#add new column that converts probabilities to logits
target_verb_time_seq_avg$Logit <- qlogis(target_verb_time_seq_avg$MeanProp)
```


```{r Visualize verb window results}

#change order of levels in AOI factor in order to change plot default
levels(as.factor(verb_time_seq$AOI))
#verb_time_seq$AOI <- factor(verb_time_seq$AOI, levels = c("Target", "AgentDis", "ActionDis", "Unrelated"))

#tiff("Verb_plot_7_21.tiff", units="in", width=7.5, height=5, res=300)
ggplot(verb_time_seq_9rm, aes(x = Time, y = BiProp, col= AOI)) +
  stat_summary(fun.data = mean_se, geom = "errorbar", aes(color=paste("mean_se", AOI))) +
  stat_summary(geom="line", size = 2) + 
  #geom_smooth() +
  theme_minimal() +
  theme(legend.text=element_text(size=rel(1.2))) +
  theme(legend.position="bottom", legend.box = "horizontal") +
  facet_grid(. ~ Group) +
  #scale_color_manual(labels = c("Target", "AgentDistractor", "ActionDistractor", "Unrelated"), values=c('darkorchid2','turquoise3', 'red2', 'chartreuse3')) +
  labs(title="Predictive Sentences-- Verb Window",
        x ="Time (ms)", y = "Proportion of Looks")
```

#model glmer

## 1) run null model to find variance components
```{r ICC models}
M0S.verb <- glmer(BiProp ~ (1|Subject), data=target_verb_time_seq, family = binomial, control = glmerControl(optimizer = "bobyqa"))
summary(M0S.verb)
performance::icc(M0S.verb)

M0T.verb <- glmer(BiProp ~ (1|trial_name), data=target_verb_time_seq, family = binomial, control = glmerControl(optimizer = "bobyqa"))
summary(M0T.verb)
performance::icc(M0T.verb)

M0 <- glmer(BiProp ~ (1|Subject) + (1|trial_name), data=target_verb_time_seq, family = binomial, control = glmerControl(optimizer = "bobyqa"))
summary(M0)
performance::icc(M0)
```
ICC (SUBJECT): 0.062. Low correlation within subjects, so variability occurs within, not between, subjects; but still a nonindependence issue, so proceed with MLM. .071 w/all trials in
ICC (TRIAL): 0.085. (.071 all trials.)  Keep both variance components in model.
ICC (total): 0.147 (.151 all trials)

##Then the full model, below. 
Use either random intercepts only or random slopes only depending on ANOVA results from prior step.
in this case, AIM.verb sig. differed, so use random slopes.
```{r full verb model}
#analyze, with all trials in
verb.target.binom.model.alltrials <- glmer(BiProp ~ (ot1+ot2)*Group + (1+ot1 | Subject) + (1+ot1 | trial_name), data=target_verb_time_seq, family = binomial, glmerControl(optimizer = c("bobyqa")))
summary(verb.target.binom.model.alltrials)
saveRDS(verb.target.binom.model.alltrials, "./verb.target.binom.model.alltrials.rds")

#reload results
verb.target.binom.model.alltrials <- readRDS("./verb.target.binom.model.alltrials.rds")
print(summary(verb.target.binom.model.alltrials))


#OR without anticipated trials included:
# verb.target.binom.model <- glmer(BiProp ~ (ot1+ot2)*Group + (1+ot1 | Subject) + (1+ot1 | trial_name), data=target_verb_time_seq, family = binomial, glmerControl(optimizer = c("bobyqa")))
# summary(verb.target.binom.model)
# 
# saveRDS(verb.target.binom.model, "./verb.target.binom.model.rds")

```


```{r target plot data vs model}
#reload model if needed
verb.target.binom.model.alltrials <- readRDS("./verb.target.binom.model.alltrials.rds")

new.verb <- target_verb_time_seq
new.verb$predicted.verb <- predict(verb.target.binom.model.alltrials, newdata=target_verb_time_seq)

tiff("verb_plot_modelfit_12_2.tiff", units="in", width=3, height=5, res=600)
ggplot() +
  geom_smooth(data = target_verb_time_seq_avg, fun.data = mean_se, geom="line", alpha = 0.15, size = 1, aes(x = Time, y = Logit, color='darkorchid'), fill= 'darkorchid') +
  stat_summary(data = new.verb, fun.data = mean_se, geom="line", size = 1, aes(x = Time, y =
                predicted.verb, color='black')) +
  facet_grid(Group ~ .) +
  theme_minimal() +
  theme(legend.text=element_text(size=rel(1.2))) +
  theme(legend.position = "bottom") +
  scale_color_manual(
    name=NULL,
    values=c('black', 'darkorchid2'),
    labels = c('Model prediction', 'Data')) +
  scale_fill_manual(
    name = NULL,
    values=c('white', 'darkorchid'),
    labels = c('Model prediction', 'Data')) +
  expand_limits(y = c(-1.2, 1.2)) +
  labs(title="Verb window: Model fit",
        x ="Time (ms)", y = "Log odds of fixations on target AOI") +
  guides(color=guide_legend(override.aes=list(fill=c('white', 'darkorchid'))))

dev.off()
```

#Verb action model
in this case, use random slopes only for ot1. ot2 p = 0.4618, and resulted in singular boundary fit
```{r full MT model}

#high vs low motion content verbs with anticipated trials
MT.verb.target.binom.model.alltrials.9rm <- glmer(BiProp ~ (ot1+ot2) * Group * mt_content + (1+ot1| Subject) + (1+ot1 | trial_name), data=target_verb_time_seq_9rm, family = binomial, control=glmerControl(optimizer="bobyqa",optCtrl=list(maxfun=100000)))
summary(MT.verb.target.binom.model.alltrials.9rm)
saveRDS(MT.verb.target.binom.model.alltrials.9rm, "./MT.verb.target.binom.model.alltrials.9rm.rds")

#reload
MT.verb.target.binom.model.alltrials.9rm <- readRDS("./MT.verb.target.binom.model.alltrials.9rm.rds")
print(summary(MT.verb.target.binom.model.alltrials.9rm))

#high vs low motion content verbs no anticipated trials
# MT.verb.target.binom.model.9rm <- glmer(BiProp ~ (ot1+ot2) * Group * mt_content + (1+ot1| Subject) + (1+ot1 | trial_name), data=target_verb_time_seq_9rm, family = binomial, control=glmerControl(optimizer="bobyqa",
#             optCtrl=list(maxfun=100000)))
# summary(MT.verb.target.binom.model.9rm)
# saveRDS(MT.verb.target.binom.model.9rm, "./MT.verb.target.binom.model.9rm.rds")

#did inclusion of MT content sig. improve model fit? 
#anova(verb.target.binom.model, NP.verb.target.binom.model)


```
MT model: anova results, p = 0.117

```{r reload motion content models}
#re-load model(s)
MT.verb.target.binom.model <- readRDS("./MT.verb.target.binom.model.rds")
print(summary(MT.verb.target.binom.model))
#OR
MT.verb.target.binom.model.9rm <- readRDS("./MT.verb.target.binom.model.9rm.rds")
print(summary(MT.verb.target.binom.model.9rm))
#OR
MT.verb.target.binom.model.alltrials <- readRDS("./MT.verb.target.binom.model.alltrials.rds")
print(summary(MT.verb.target.binom.model.alltrials))
#OR
MT.verb.target.binom.model.alltrials.9rm <- readRDS("./MT.verb.target.binom.model.alltrials.9rm.rds")
print(summary(MT.verb.target.binom.model.alltrials.9rm))
```

#modeling looks to agent-related distractor

##full agent-related model
```{r}

verb.agentdis.binom.model <- glmer(BiProp ~ (ot1+ot2)*Group + (1+ot1 | Subject) + (1+ot1 | trial_name), data=agentdis_verb_time_seq, family = binomial, glmerControl(optimizer = c("bobyqa")))
summary(verb.agentdis.binom.model)
saveRDS(verb.agentdis.binom.model, "./verb.agentdis.binom.model.rds")

#reload
verb.agentdis.binom.model.alltrials <- readRDS("./verb.agentdis.binom.model.rds")
print(summary(verb.agentdis.binom.model))


#or without anticipated trials:
verb.agentdis.binom.model.reducedtrials <- glmer(BiProp ~ (ot1+ot2)*Group + (1+ot1 | Subject) + (1+ot1 | trial_name), data=agentdis_verb_time_seq, family = binomial, glmerControl(optimizer = c("bobyqa")))
summary(verb.agentdis.binom.model.reducedtrials)
saveRDS(verb.agentdis.binom.model.reducedtrials, "./verb.agentdis.binom.model.reducedtrials.rds")

#reload
verb.agentdis.binom.model.alltrials <- readRDS("./verb.agentdis.binom.model.rds")
print(summary(verb.agentdis.binom.model))


#test difference in mean proportions to agent/target
var.test(verb_time_seq_avg$Target,verb_time_seq_avg$AgentDis) #variances unequal, use Welch
t.test(verb_time_seq_avg$Target,verb_time_seq_avg$AgentDis)
sd(verb_time_seq_avg$Target, na.rm=TRUE)
sd(verb_time_seq_avg$AgentDis, na.rm=TRUE)

```


```{r values of interest from raw data}

PD.initial.target<-target_verb_time_seq[which(target_verb_time_seq$Group.x=="PD" & target_verb_time_seq$Time==0), ]
mean(PD.initial.target$BiProp, na.rm=TRUE)

PD.end.target<-target_verb_time_seq[which(target_verb_time_seq$Group.x=="PD" & target_verb_time_seq$Time==450), ]
mean(PD.end.target$BiProp, na.rm=TRUE)

Control.initial.target<-target_verb_time_seq[which(target_verb_time_seq$Group.x=="CONTROL" & target_verb_time_seq$Time==0), ]
mean(Control.initial.target$BiProp, na.rm=TRUE)

Control.end.target<-target_verb_time_seq[which(target_verb_time_seq$Group.x=="CONTROL" & target_verb_time_seq$Time==450), ]
mean(Control.end.target$BiProp, na.rm=TRUE)

```
0 - 500ms
probability of looks to target for PD       @ timebin 0 = 34.7% @ time 450 = 44.1%. linear increase = 9.4%
probability of looks to target for controls @ timebin 0 = 36.4% @ time 450 = 48.7%. linear increase = 12.3%

```{r view model results for writeup}
summary(MT.agent.target.binom.model.9rm)

summary(MT.verb.target.binom.model.9rm)

summary(MT.obj.target.binom.model.9rm)

summary(end.target.binom.model)
```