eeg_motivated_contagion_analysis_eegdata.Rmd 11.5 KB
Newer Older
1 2 3 4 5 6 7
---
title: "EEG_analysis"
author: "Amit"
date: "October 21, 2016"
output: html_document
---

Amit Goldenberg's avatar
Amit Goldenberg committed
8
```{r packages and functions, include=FALSE}
9

Amit Goldenberg's avatar
Amit Goldenberg committed
10 11
if (!require("pacman")) install.packages("pacman")
pacman::p_load("ggplot2","lme4",'lmerTest','plyr','dplyr','psych','car',"tidyr",'Rmisc','ez',"visreg",'GGally','coefplot')
12

13
pl1way <- function(dat,y,x) {
14 15 16 17
	pld <- summarySEwithin(dat,y,withinvars = x, idvar = "subject")
	tmp <- dat
	colnames(tmp) <- gsub(as.name(y),"y",colnames(tmp))
	colnames(tmp) <- gsub(as.name(x),"x",colnames(tmp))
18
	res <- ezANOVA(data=tmp,dv=y, wid=subject, within=.(x), type=3)$ANOVA                # ANOVA
19 20 21 22 23
	result <- paste0("F(",res$DFn,",",res$DFd,")=",round(res$F,2), ", GES=", round(res$ges,2), ", p=", round(res$p,3))
	ggplot(pld,aes(x=pld[,x],y=pld[,y], ymin=pld[,y]-ci,ymax=pld[,y]+ci)) + geom_pointrange(fatten = 8) +
		labs(x=x,y=y, title = result)
}

24 25
```

26 27 28 29 30
###some dataset definision
###social - second phase
###individual - first phase. 


Amit Goldenberg's avatar
Amit Goldenberg committed
31
```{r get and convert the file, include=FALSE}
Amit Goldenberg's avatar
Amit Goldenberg committed
32 33
 setwd("C:/Users/Amit/Dropbox/Research/emotional conformity-non conformity/motivated_contation_EEG_2017/Results")
#setwd("C:/Users/ro/OneDrive - Tartu likool/Teadus/REGULATION/16 motivated contagion")
34
d <- read.csv("eeg_data_motivated_contagion.csv")
Amit Goldenberg's avatar
Amit Goldenberg committed
35
d=d[,c(1,3:5,9,11,14,17:18,21:22,24:26)]
36
d$grate = ifelse(d$grate == 0, NA, d$grate)
37

38 39
##remove participants 
d = subset (d, subject != "5787127" & subject != "5858527")
40

41 42 43 44 45
# separate phases
a <- d[d$phase=="indiv",]
b <- d[d$phase=="social",]

# make a wide format with each picture in a separate row and the phases in separate columns
Amit Goldenberg's avatar
Amit Goldenberg committed
46 47
d1 <- merge(a[,c("subject", "photo", "order","rate", "rt", "P3", "LPP", "SLW")],b,by=c("subject", "photo"))
rm(a,b)
48

Amit Goldenberg's avatar
Amit Goldenberg committed
49 50
colnames(d1) <- gsub(".x", ".indiv", colnames(d1), fixed =T)
colnames(d1) <- gsub(".y", ".social", colnames(d1), fixed =T)
51

Amit Goldenberg's avatar
Amit Goldenberg committed
52
# remove "none"" and fix comparison levels  
53
#d1 <- d1[d1$condition!="none",]
Amit Goldenberg's avatar
Amit Goldenberg committed
54 55

#I just want to make sure that this just organizes the condition 
Amit Goldenberg's avatar
Amit Goldenberg committed
56

57
d1$condition <- factor(as.character(d1$condition), levels = c( "low", "high","same", "none"))
Amit Goldenberg's avatar
Amit Goldenberg committed
58 59 60


# add phase differences
61 62 63 64
d1$P3.dif = d1$P3.social-d1$P3.indiv; hist (d1$P3.dif)
d1$LPP.dif = d1$LPP.social-d1$LPP.indiv
d1$SLW.dif = d1$SLW.social-d1$SLW.indiv
d1$rate.dif = d1$rate.social-d1$rate.indiv
Amit Goldenberg's avatar
Amit Goldenberg committed
65
# add difference between group and individual ratings
66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96
d1$rateDifCon <- d1$rate.social-d1$grate 


d1e = subset (d1,photoType =="emo" &(rate.indiv<8 &rate.indiv ))
d1n = subset (d1,photoType =="ntr" )


d1$photoType <- relevel(d1$photoType,"ntr")

```

### perliminary analysis, get raw means
```{r}


means = d1e%>%
  group_by(photoType,condition) %>%
  summarise(
    P3s =mean(P3.social, na.rm = T),
    LPPs = mean(LPP.social, na.rm = T),
    SLWs = mean(SLW.social , na.rm =T),
    P3i = mean(P3.indiv, na.rm = T),
    LPPi = mean(LPP.indiv, na.rm = T),
    SLWi = mean(SLW.indiv, na.rm = T),
    P3.dif = mean(P3.dif, na.rm = T),
    LPP.dif = mean(LPP.dif, na.rm = T),
    SLW.dif = mean(SLW.dif, na.rm = T)
    )


data.frame(means)
97

98
rm (means)
99

Amit Goldenberg's avatar
Amit Goldenberg committed
100 101
```

102 103


Amit Goldenberg's avatar
Amit Goldenberg committed
104 105
### A couple of descriptive pictures
```{r descriptives, warning=F}
106

Amit Goldenberg's avatar
Amit Goldenberg committed
107
# frequency table of trials
108 109 110 111 112
a=table(d1$subject, d1$condition)

length(table(mean(data.frame(apply(a,1,sum))[,1])))
table(d1e$subject, d1e$condition)

113

Amit Goldenberg's avatar
Amit Goldenberg committed
114
# independent variables
Amit Goldenberg's avatar
Amit Goldenberg committed
115
ggpairs(d1e[,c("rate.social","grate", "rateDifCon")]) 
116

Amit Goldenberg's avatar
Amit Goldenberg committed
117 118 119 120
# dependent variables
ggpairs(d[,c("P3","LPP", "SLW")])

# differences 
Amit Goldenberg's avatar
Amit Goldenberg committed
121 122 123 124
ggpairs(d1[,c("P3.indiv","P3.social", "P3.dif")])
ggpairs(d1[,c("LPP.indiv","LPP.social", "LPP.dif")])
ggpairs(d1[,c("SLW.indiv","SLW.social", "SLW.dif")])

125 126
ggpairs(d1e[,c("P3.indiv","P3.social")])
ggpairs(d1e[,c("rate.indiv","LPP.social")])
Amit Goldenberg's avatar
Amit Goldenberg committed
127
ggpairs(d1[,c("SLW.indiv","SLW.social")])
128
ggpairs(d1e[,c("rate.indiv","rate.social")])
Amit Goldenberg's avatar
Amit Goldenberg committed
129 130 131

# Do we need to control for order? No
cor(d[,c("P3","LPP", "SLW", "rate", "grate")],d$order, use="pairwise.complete.obs")
132 133 134

```

Amit Goldenberg's avatar
Amit Goldenberg committed
135
I've been thinking about some of the researcher degrees of freedom we have in constructing the models. Here's a summary:
136

Amit Goldenberg's avatar
Amit Goldenberg committed
137 138 139
### 1) We have a choice between taking a "brain-centric" and "task-centric" approach
Brain-centric would look somethink like this: P3.social ~ P3.indiv + condition + (1|subject)
Task-centric something like this: P3.social ~ photoType + condition + (1|subject) + (1|photo)
140

Amit Goldenberg's avatar
Amit Goldenberg committed
141
Both approaches aim to isolate brain responses to the social manipulation.
142

Amit Goldenberg's avatar
Amit Goldenberg committed
143
The brain-centric model assumes that the first presentation response (P3.indiv term) captures all brain responses to both the affective and nonaffective aspects of the stimuli. In this model, we shouldn't include other stuff like individual ratings, photoType, or photo random factor, because processes triggered by all these things are already contained in the P3.indiv variance. 
144

Amit Goldenberg's avatar
Amit Goldenberg committed
145 146 147
The alternative is of course to use what we manipulated to regress out brain variance not related to the social manipulation. If we do that, we shouldn't any more include the first phase EEG, because we then run the risk of regressing the same thing out twice. 

### 2) For brain-centric approach, should we use a difference score or residualization?
148
In principle using the difference between P3.social and P.indiv as a depent variable is also a way to do the brain-centric approach. There is some reason to think residualization works better, because variance from the first phase will "leak" into difference scores but not into residualized scores:
Amit Goldenberg's avatar
Amit Goldenberg committed
149

Amit Goldenberg's avatar
Amit Goldenberg committed
150 151 152 153 154
```{r individual and social , warning=F}
head(d1)

rsd <- lmer(P3.social ~ P3.indiv + (1|subject), data = d1) ; summary (rsd)

155 156 157 158 159 160 161 162
r = lmer (P3.social ~ rate.indiv +(1|subject) , d1e); summary(r)
r = lmer (P3.social ~ rate.social+grate+(1|subject) , d1e); summary(r)


r = lmer (LPP.social ~ grate+rate.social +(1|subject) +(1|photo), d1e); summary(r)

r = lmer (SLW.social ~ rate.indiv +(1|subject) , d1e); summary(r)
r = lmer (SLW.social ~ rate.social+grate  +(1|subject) , d1e); summary(r)
Amit Goldenberg's avatar
Amit Goldenberg committed
163
#rsd <- lmer(P3.social ~ P3.indiv + (1|subject), data = d1e) ; summary (rsd)
Amit Goldenberg's avatar
Amit Goldenberg committed
164

Amit Goldenberg's avatar
Amit Goldenberg committed
165
ggpairs(cbind(d1[,c("P3.indiv", "P3.social","P3.dif")],P3.residuals=resid(rsd)))
Amit Goldenberg's avatar
Amit Goldenberg committed
166

Amit Goldenberg's avatar
Amit Goldenberg committed
167
ggpairs(cbind(d1[,c("P3.indiv", "P3.social","P3.dif")],P3.residuals=resid(rsd)))
Amit Goldenberg's avatar
Amit Goldenberg committed
168 169


Amit Goldenberg's avatar
Amit Goldenberg committed
170 171
rsd <- lmer(rate.social ~ rate.indiv + (1|subject), data = d1) ; summary (rsd)

Amit Goldenberg's avatar
Amit Goldenberg committed
172

Amit Goldenberg's avatar
Amit Goldenberg committed
173 174 175 176 177 178 179
```

### 3) Should we use condition as a factor or as a continuous difference between group and individual ratings?
One reaosn to prefer the factor is that this allows us to detect nonlinear effects (e.g. low differing from both same and high, which is what happens for P3 amplitudes)


Based on various reasons, including some saturday-night-p-hacking, I stuck to the brain-centric models using group ratings as factor: 
180

181 182
### P3 analysis

Amit Goldenberg's avatar
Amit Goldenberg committed
183
```{r P3 dif, warning=F}
Amit Goldenberg's avatar
Amit Goldenberg committed
184 185

# simplest model
186 187 188
summary(fit <- lmer(P3.social ~  condition +grate+P3.indiv+order.indiv+order.social+ (1|subject)+(1|photo), data = d1e))


Amit Goldenberg's avatar
Amit Goldenberg committed
189 190 191 192 193

# get p-values for the random effects
rand(fit)
# get a significance test for all terms
anova(fit)
194
coefplot(fit)
Amit Goldenberg's avatar
Amit Goldenberg committed
195
#visreg(fit, xvar = "condition")
196

Amit Goldenberg's avatar
Amit Goldenberg committed
197 198
# plot means using raw data
# pl1way(D3,"P3.dif", "condition")
199

Amit Goldenberg's avatar
Amit Goldenberg committed
200 201 202 203 204 205 206 207 208 209
# plot condition effect using residuals from the model above
rsd <- lmer(P3.social ~ P3.indiv + (1|subject), data = D3)
pl1way(cbind(D3,dep=resid(rsd)),"dep", "condition")

```

### LPP
```{r, warning=F}

# simplest model
210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230
summary(fit <- lmer(LPP.social ~  condition +LPP.indiv+ (1|subject), data = d1e))

270-470
470-800

200/530
330/530

d1e$a = d1e$LPP.social*0.3773585
d1e$b = d1e$LPP.social*0.6226415

d1e$c = d1e$LPP.indiv*0.3773585
d1e$d = d1e$LPP.indiv*0.6226415


d1e$early.social= (rowMeans(d1e[,c('a', 'b')], na.rm=T))
d1e$early.indiv= (rowMeans(d1e[,c('c', 'd')], na.rm=T))


summary(fit <- lmer(early.social ~condition +early.indiv+(1|photo)+(1|subject), data = d1e))

Amit Goldenberg's avatar
Amit Goldenberg committed
231 232 233 234 235

# get p-values for the random effects
rand(fit)
# get a significance test for all terms
anova(fit)
236
coefplot(fit)
Amit Goldenberg's avatar
Amit Goldenberg committed
237 238 239 240
#visreg(fit, xvar = "condition")

# plot means using raw data
# pl1way(D3,"P3.dif", "condition")
241

Amit Goldenberg's avatar
Amit Goldenberg committed
242
# plot condition effect using residuals from the model above
Amit Goldenberg's avatar
Amit Goldenberg committed
243
rsd <- lmer(LPP.social ~ LPP.indiv + (1|subject), data = d1)
Amit Goldenberg's avatar
Amit Goldenberg committed
244
pl1way(cbind(D3,dep=resid(rsd)),"dep", "condition")
245

Amit Goldenberg's avatar
Amit Goldenberg committed
246
```
247

Amit Goldenberg's avatar
Amit Goldenberg committed
248 249
### SLW
```{r, warning=F}
250

Amit Goldenberg's avatar
Amit Goldenberg committed
251
# simplest model
252
summary(fit <- lmer(SLW.social ~  condition +SLW.indiv+ order.indiv+order.social+ (1|subject)+(1|photo), data = d1e))
253

254
 # get p-values for the random effects
Amit Goldenberg's avatar
Amit Goldenberg committed
255 256 257 258 259
rand(fit)
# get a significance test for all terms
anova(fit)
coefplot(fit)
#visreg(fit, xvar = "condition")
260

Amit Goldenberg's avatar
Amit Goldenberg committed
261 262
# plot means using raw data
# pl1way(D3,"P3.dif", "condition")
263

Amit Goldenberg's avatar
Amit Goldenberg committed
264
# plot condition effect using residuals from the model above
Amit Goldenberg's avatar
Amit Goldenberg committed
265
rsd <- lmer(SLW.social ~ SLW.indiv + (1|subject), data = d1)
Amit Goldenberg's avatar
Amit Goldenberg committed
266
pl1way(cbind(D3,dep=resid(rsd)),"dep", "condition")
267 268 269

```

Amit Goldenberg's avatar
Amit Goldenberg committed
270 271 272 273 274 275

Based on these analyses, we have a small but significant effect on P3 and no effects in later stages. I'm pretty sure we can build a story out of this alone. But I also am getting increasing eager to take a more explorative approach to ask where in the brain the social cues actually have an impact. Given our design, LIMO is probably still the best way to do that, we can of course discuss this further.



```{r dif, include=FALSE}
276 277


278
summary(fit <- lmer(scale(P3.dif, scale = F) ~ condition+ (1|subject) + (1|photo), data = d1e))
279

280
summary(fit <- lmer(scale(LPP.dif, scale = F) ~  condition +order.social+order.indiv+ (1|subject) + (1|photo), d1e))
281 282


283 284
summary(fit <- lmer(scale(SLW.dif, scale = F)  ~  condition + (1|subject) + (1|photo), data = d1[d1$photoType=="emo",]))
```
285 286


287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355
```{r order}
summary(fit <- lmer(LPP.social ~  order.social+ (1|subject) + (1|photo)-1, d1e))


summary(fit <- lmer(LPP.indiv ~  order.indiv+ (1|subject) + (1|photo)-1, d1e))
```


###actual ratings
```{r actual ratings}

r = lmer(rate.indiv~ condition+(1|subject) +(1|photo),d1e);summary (r)



r = lmer(scale(rate.dif, scale =F)~ condition+order.indiv+order.social -1+(1|subject) +(1|photo),d1e);summary (r)


head(d1e)
r = lmer (rate.social ~ condition +order.indiv+order.social + (1|subject) + (1|photo), d1e) ; summary (r)



```

####reaction time
```{r}


r = lmer(rt.social~ condition+order.indiv+order.social +(1|subject) +(1|photo),d1e);summary (r)


##correlate reaction time with lpp

r = lmer(LPP.social~ rt.social+order.indiv+order.social +(1|subject) +(1|photo),d1e);summary (r)

```




### analysis of individual participants
```{r}


addmargins(table(d1e$subject, d1e$condition))


models = d1e%>%
  group_by(subject)  %>% 
  do(r = print(summary(lm(LPP.social ~ condition +P3.indiv, data = .))))


coef = models %>%
  do(data.frame (low = (.$r)$coefficients[[1]],
               high = (.$r)$coefficients[[2]],
               same = (.$r)$coefficients[[1]]
               )) 
coef =cbind(models[,1], coef)

coef$test= ifelse(coef$high>1,1,2)

table(coef$test)

coef=arrange (coef,coef[,1])



rm(models)
356

Amit Goldenberg's avatar
Amit Goldenberg committed
357
```
358 359 360 361




362
###compute differences from neutral 
Amit Goldenberg's avatar
Amit Goldenberg committed
363 364 365 366 367 368 369 370 371 372 373 374 375
```{r compute difference from neutral}
head(d)
##creeat dataset of the averages. 
dd = d1%>%
  group_by(subject,photoType,condition) %>%
  summarise(
    P3s =mean(P3.social, na.rm = T),
    LPPs = mean(LPP.social, na.rm = T),
    SLWs = mean(SLW.social , na.rm =T),
    P3i = mean(P3.indiv, na.rm = T),
    LPPi = mean(LPP.indiv, na.rm = T),
    SLWi = mean(SLW.indiv, na.rm = T)
    )
376

Amit Goldenberg's avatar
Amit Goldenberg committed
377 378
dde = subset (dd, photoType == "emo")
ddn = subset (dd, photoType == "ntr")
379

Amit Goldenberg's avatar
Amit Goldenberg committed
380 381
dd= left_join(dde, ddn, by = c("subject", "condition"))
head(dd)
382 383


Amit Goldenberg's avatar
Amit Goldenberg committed
384 385
colnames(dd) <- gsub(".x", ".e", colnames(dd), fixed =T)
colnames(dd) <- gsub(".y", ".n", colnames(dd), fixed =T)
386

Amit Goldenberg's avatar
Amit Goldenberg committed
387
rm (dde, ddn)
388

Amit Goldenberg's avatar
Amit Goldenberg committed
389
dd = dd[,-2]
390

Amit Goldenberg's avatar
Amit Goldenberg committed
391 392
dd$P3 = dd$P3s.e- dd$P3s.n
hist (dd$P3)
393

Amit Goldenberg's avatar
Amit Goldenberg committed
394 395 396 397
dd$LPP = dd$LPPs.e- dd$LPPs.n
hist (dd$LPP)
dd$SLW = dd$SLWs.e- dd$SLWs.n
hist (dd$SLW)
398 399


Amit Goldenberg's avatar
Amit Goldenberg committed
400 401 402
r= lm (P3~ condition, dd); summary (r)
r= lm (LPP~ condition+LPPi.e+LPPi.n, dd); summary (r)
r= lm (SLW~ condition+SLWi.e+SLWi.n, dd); summary (r)
403

404 405


Amit Goldenberg's avatar
Amit Goldenberg committed
406
```