123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169 |
- # _____ _ _ _ ______
- # / ____| | | (_) | | | ____|
- # | (___ _ __ ___ ___ | | ___ _ __ __ _ __ _ _ __ __| | | |__ ___ __ _ _ __
- # \___ \| '_ ` _ \ / _ \| |/ / | '_ \ / _` | / _` | '_ \ / _` | | __/ _ \/ _` | '__|
- # ____) | | | | | | (_) | <| | | | | (_| | | (_| | | | | (_| | | | | __/ (_| | |
- # |_____/|_| |_| |_|\___/|_|\_\_|_| |_|\__, | \__,_|_| |_|\__,_| |_| \___|\__,_|_|
- # __/ |
- # |___/
- #Analysis script Fear ratings and US expectancy ratings (day1 and day 2)
- #written by Madeleine Mueller, February 2024
- library("lme4")
- library("car")
- library("effects")
- library("emmeans")
- # __ _ _
- # / _| | | (_)
- # | |_ ___ __ _ _ __ _ __ __ _| |_ _ _ __ __ _ ___
- # | _/ _ \/ _` | '__| | '__/ _` | __| | '_ \ / _` / __|
- # | || __/ (_| | | | | | (_| | |_| | | | | (_| \__ \
- # |_| \___|\__,_|_| |_| \__,_|\__|_|_| |_|\__, |___/
- # _ _ _ _ __/ |
- # | | | | | | | | |___/
- # | |__ ___ | |_| |__ __| | __ _ _ _ ___
- # | '_ \ / _ \| __| '_ \ / _` |/ _` | | | / __|
- # | |_) | (_) | |_| | | | | (_| | (_| | |_| \__ \
- # |_.__/ \___/ \__|_| |_| \__,_|\__,_|\__, |___/
- # __/ |
- # |___/
- #read in data of fear rating for both days
- onstud_fear<-read.csv(file='/Madeleine/Online Studie/analysis/R/FearRatingsOnstudforR.csv',head=TRUE,sep=",")
- #all four groups (1=non-smoker, 2=smoker with 6h break after ACQ, 3=smoker with cigarette after ACQ, 4=smoker with cigarette before ACQ)
- onstud_fear$group<-as.factor(onstud_fear$group)
- #all 4 groups, see above, but smokers from group 2 that did smoke during the 6h smoke break were regrouped into group 3
- onstud_fear$regroup<-as.factor(onstud_fear$regroup)
- #stimulus: CSP=CS+, CSM=CS-
- onstud_fear$stim<-as.factor(onstud_fear$stim)
- #time: pre ACQ/pre GEN (depending on day)=1; post ACQ/ post GEN = 2
- onstud_fear$time<-as.factor(onstud_fear$prepost)
- #smoker(group2,3 and 4)= 1, non-smoker(group1) = 0
- onstud_fear$smoker<-as.factor(onstud_fear$smoker)
- #subject number
- onstud_fear$sub<-onstud_fear$ID
- #fear ratings day 1 models with different groups
- both.f1s<-(lmer(ratingd1~(1|sub)+stim*time*smoker, data=onstud_fear,control = lmerControl(optimizer = "bobyqa")))
- both.f1g<-(lmer(ratingd1~(1|sub)+stim*time*regroup, data=onstud_fear,control = lmerControl(optimizer = "bobyqa")))
- #fear ratings day 1 models with different groups
- both.f2s<-(lmer(ratingd2~(1|sub)+stim*time*smoker, data=onstud_fear,control = lmerControl(optimizer = "bobyqa")))
- both.f2g<-(lmer(ratingd2~(1|sub)+stim*time*regroup, data=onstud_fear,control = lmerControl(optimizer = "bobyqa")))
- #anova based on models, include different model name for all results
- Anova(both.f1s, type="3", test="F" )
- Anova(both.f2s, type="3", test="F" )
- #post-hoc tests, include effects in which you are interested
- emmeans(both.f2s, list(pairwise ~ stim:time:smoker), adjust = "none")
- #plot effects, change effect name and model name as needed
- plot(effect("stim*smoker",both.f1s))
- plot(effect("stim*time*smoker",both.f2s))
- # _ _ _____ _
- # | | | |/ ____| | |
- # | | | | (___ _____ ___ __ ___ ___| |_ __ _ _ __ ___ _ _
- # | | | |\___ \ / _ \ \/ / '_ \ / _ \/ __| __/ _` | '_ \ / __| | | |
- # | |__| |____) | | __/> <| |_) | __/ (__| || (_| | | | | (__| |_| |
- # \____/|_____/ \___/_/\_\ .__/ \___|\___|\__\__,_|_| |_|\___|\__, |
- # | | /_ | | | __/ |
- # __| | __ _ _ _ | | |_| |___/
- # / _` |/ _` | | | | | |
- # | (_| | (_| | |_| | | |
- # \__,_|\__,_|\__, | |_|
- # __/ |
- # |___/
- #read in data
- data_ex<-read.csv(file='/Madeleine/Online Studie/analysis/R/D1USexponstudR.csv',head=TRUE,sep=",")
- #all four groups (1=non-smoker, 2=smoker with 6h break after ACQ, 3=smoker with cigarette after ACQ, 4=smoker with cigarette before ACQ)
- data_ex$group<-as.factor(data_ex$group)
- #all 4 groups, see above, but smokers from group 2 that did smoke during the 6h smoke break were regrouped into group 3
- data_ex$regroup<-as.factor(data_ex$regroup)
- #stimulus: CSP=CS+, CSM=CS-
- data_ex$stim<-as.factor(data_ex$stim)
- # trial 1-12 per stimulus, so in total 24 trials
- data_ex$trial<-as.factor(data_ex$trial)
- #block 1-3, one block consists of 4 trials per stimulus, so 8 trials in total
- data_ex$block<-as.factor(data_ex$block)
- #difference rated between pleasantness US and nUS (mean US - mean nUS)
- data_ex$usdiff<-as.factor(data_ex$usdiff)
- # avoidance of screen when US is presented, higher rating means more avoidance
- data_ex$avoid<-as.factor(data_ex$avoid)
- # CS-US contingency, 1=not aware of correct CS-US contingency, 2=aware of correct CS-US contingency
- data_ex$aware<-as.factor(data_ex$aware)
- #smoker(group2,3 and 4)= 1, non-smoker(group1) = 0
- data_ex$smoker<-as.factor(data_ex$smoker)
- #subject number
- data_ex$sub<-data_ex$sub
- #########################################
- #models
- #########################################
- USex1a<-lmer(rating~(1|sub)+stim*block*regroup,data = data_ex, control = lmerControl(optimizer = "bobyqa"))
- Anova(USex1a, type="3", test="F" )
- USex1sm<-lmer(rating~(1|sub)+stim*block*smoker,data = data_ex,control = lmerControl(optimizer = "bobyqa"))
- Anova(USex1sm, type="3", test="F" )
- emmeans(USex1a, list(pairwise ~ stim), adjust = "none")
- plot(effect("block*regroup",USex1sa))
- plot(effect("stim*block*smoker",USex1sm))
- plot(effect("regroup",USex1a))
- # _ _ _____ _
- # | | | |/ ____| | |
- # | | | | (___ _____ ___ __ ___ ___| |_ __ _ _ __ ___ _ _
- # | | | |\___ \ / _ \ \/ / '_ \ / _ \/ __| __/ _` | '_ \ / __| | | |
- # | |__| |____) | | __/> <| |_) | __/ (__| || (_| | | | | (__| |_| |
- # \____/|_____/ \___/_/\_\ .__/ \___|\___|\__\__,_|_| |_|\___|\__, |
- # | | |__ \ | | __/ |
- # __| | __ _ _ _ ) ||_| |___/
- # / _` |/ _` | | | | / /
- # | (_| | (_| | |_| | / /_
- # \__,_|\__,_|\__, | |____|
- # __/ |
- # |___/
- d_ex<-read.csv(file='/Madeleine/Online Studie/analysis/R/D2USexponstudR.csv',head=TRUE, sep=",")
- #all four groups (1=non-smoker, 2=smoker with 6h break after ACQ, 3=smoker with cigarette after ACQ, 4=smoker with cigarette before ACQ)
- d_ex$group<-as.factor(d_ex$group)
- #all 4 groups, see above, but smokers from group 2 that did smoke during the 6h smoke break were regrouped into group 3
- d_ex$regroup<-as.factor(d_ex$regroup)
- #stimulus: CSP=CS+, GS2-GS9=generalized stimuli, CSM=CS-
- d_ex$stim<-as.factor(d_ex$stim)
- #block 1 and 2, one block consists of 1 trial per stimulus, so 10 trials in total
- d_ex$block<-as.factor(d_ex$block)
- # CS-US contingency, 1=not aware of correct CS-US contingency, 2=aware of correct CS-US contingency
- d_ex$aware<-as.factor(d_ex$aware)
- #smoker(group2,3 and 4)= 1, non-smoker(group1) = 0
- d_ex$smoker<-as.factor(d_ex$smoker)
- #subject number
- d_ex$sub<-d_ex$sub
- #########################################
- #models
- #########################################
- #models for either smoker vs non-smoker or regroup
- D2onstud<-lmer(rating~(1|sub)+stim*block*smoker,data = d_ex,control = lmerControl(optimizer = "bobyqa"))
- D2bonstud<-lmer(rating~(1|sub)+stim*block*regroup,data = d_ex,control = lmerControl(optimizer = "bobyqa"))
- Anova(D2onstud, type="3", test="F" )
- plot(effect("stim",D2onstud))
- plot(effect("regroup",D2bonstud))
- plot(effect("block",D2onstud))
- plot(effect("stim*regroup",D2bonstud))
- plot(effect("stim*block*smoker",D2onstud))
- plot(effect("stim*smoker",D2onstud))
- #posthoc
- emmeans(D2onstud, list(pairwise ~ stim:block:smoker), adjust = "none")
|