1. Worry and Rumination Survey
A. Calculate the average ratings across thought type per topic per
subject
Worry
#create empty matrix to fill in with NAs - 36 obesrvations and 15 cols ( 7 freq 7 intensity and 1 subid col)
#may need to update # of obs if using the updated doc where we removed fake respones/999s
AvgWorry <- as.data.frame(matrix(NA, nrow = 87, ncol = 13))
#assign the columns variable names
colnames(AvgWorry) <- c("ID", "Worry1avgFreq", "Worry1avgInt", "Worry2avgFreq","Worry2avgInt", "Worry3avgFreq", "Worry3avgInt", "Worry4avgFreq", "Worry4avgInt", "Worry5avgFreq","Worry5avgInt", "Worry6avgFreq","Worry6avgInt")
#pull subids from original data set
AvgWorry$ID <- surv$Q1
#create within-persons avgs
AvgWorry$Worry1avgFreq <-rowSums(surv[,c(47:50)])/4 # Q62_1-4
AvgWorry$Worry1avgInt <- rowSums(surv[,c(53:56)])/4 # Q63_1-4
AvgWorry$Worry2avgFreq <- rowSums(surv[,c(67:70)])/4
AvgWorry$Worry2avgInt <- rowSums(surv[,c(73:76)])/4
AvgWorry$Worry3avgFreq <- rowSums(surv[,c(87:90)])/4
AvgWorry$Worry3avgInt <- rowSums(surv[,c(93:96)])/4
AvgWorry$Worry4avgFreq <-rowSums(surv[,c(107:110)])/4
AvgWorry$Worry4avgInt <- rowSums(surv[,c(113:116)])/4
AvgWorry$Worry5avgFreq <-rowSums(surv[,c(127:130)])/4
AvgWorry$Worry5avgInt <- rowSums(surv[,c(133:136)])/4
AvgWorry$Worry6avgFreq <- rowSums(surv[,c(147:150)])/4
AvgWorry$Worry6avgInt <- rowSums(surv[,c(153:156)])/4
Rumination
#create empty matrix to fill in with NAs - 25 obesrvations and 15 cols ( 7 freq 7 intensity and 1 subid col)
AvgRum <- as.data.frame(matrix(NA, nrow=87, ncol=13))
#assign the columns variable names
colnames(AvgRum) <- c("ID", "Rum1avgFreq", "Rum1avgInt", "Rum2avgFreq", "Rum2avgInt", "Rum3avgFreq", "Rum3avgInt" , "Rum4avgFreq", "Rum4avgInt", "Rum5avgFreq", "Rum5avgInt", "Rum6avgFreq", "Rum6avgInt")
#pull subids from original data set
AvgRum$ID <- surv$Q1
#Fill in new matrix directly with mean calcs using rowSums and division
AvgRum$Rum1avgFreq <- rowSums(surv[,c(196:199)])/4 #Q40_1-4
AvgRum$Rum1avgInt <- rowSums(surv[,c(202:205)])/4
AvgRum$Rum2avgFreq <- rowSums(surv[,c(216:219)])/4
AvgRum$Rum2avgInt <- rowSums(surv[,c(222:225)])/4
AvgRum$Rum3avgFreq <- rowSums(surv[,c(236:239)])/4
AvgRum$Rum3avgInt <- rowSums(surv[,c(242:245)])/4
AvgRum$Rum4avgInt <- rowSums(surv[,c(256:259)])/4
AvgRum$Rum4avgFreq <- rowSums(surv[,c(262:265)])/4
AvgRum$Rum5avgInt <- rowSums(surv[,c(276:279)])/4
AvgRum$Rum5avgFreq <- rowSums(surv[,c(282:285)])/4
AvgRum$Rum6avgInt <- rowSums(surv[,c(296:299)])/4
AvgRum$Rum6avgFreq <- rowSums(surv[,c(302:305)])/4
B. Combine Worry and Rum Averages into single dataset (Supplemental
Figure 1)
AvgWorryRum <- cbind(AvgRum,AvgWorry)
AvgWorryRum <- AvgWorryRum[,-c(14)] #remove duplicate ID variable
### For scan subs only:
scansubs <- unique(wide$Subject)
AvgRum_int <- AvgRum[c(which(AvgRum$ID %in% scansubs)),c(1,seq(3,13,by=2))]
AvgRum_freq <- AvgRum[c(which(AvgRum$ID %in% scansubs)),c(1,seq(2,13,by=2))]
AvgWorry_int <- AvgWorry[c(which(AvgRum$ID %in% scansubs)),c(1,seq(3,13,by=2))]
AvgWorry_freq <- AvgWorry[c(which(AvgRum$ID %in% scansubs)),c(1,seq(2,13,by=2))]
AvgRum_int <- AvgRum_int[order(AvgRum_int$Rum1avgInt),]
AvgRum_freq <- AvgRum_freq[order(AvgRum_freq$Rum1avgFreq),]
AvgWorry_int <- AvgWorry_int[order(AvgWorry_int$Worry1avgInt),]
AvgWorry_freq <- AvgWorry_freq[order(AvgWorry_freq$Worry1avgFreq),]
AvgRum_int$ID <- factor(AvgRum_int$ID, levels=(AvgRum_int$ID)[order(AvgRum_int$Rum1avgInt)])
AvgRum_freq$ID <- factor(AvgRum_freq$ID, levels=(AvgRum_freq$ID)[order(AvgRum_freq$Rum1avgFreq)])
AvgWorry_int$ID <- factor(AvgWorry_int$ID, levels=(AvgWorry_int$ID)[order(AvgWorry_int$Worry1avgInt)])
AvgWorry_freq$ID <- factor(AvgWorry_freq$ID, levels=(AvgWorry_freq$ID)[order(AvgWorry_freq$Worry1avgFreq)])
AvgRum_int_mltd <- melt(data=AvgRum_int)
AvgRum_freq_mltd <- melt(data=AvgRum_freq)
AvgWorry_int_mltd <- melt(data=AvgWorry_int)
AvgWorry_freq_mltd <- melt(data=AvgWorry_freq)
ari <- ggplot(AvgRum_int_mltd,
aes(x = variable, y = ID, fill = value))+
geom_tile()+scale_fill_gradient(high = "coral", low = "slateblue3") +
scale_x_discrete(labels=c("1", "2", "3","4","5","6"))+
labs(x="Topic Rank")
arf <- ggplot(AvgRum_freq_mltd,
aes(x = variable, y = ID, fill = value))+
geom_tile()+scale_fill_gradient(high = "coral", low = "slateblue3") +
scale_x_discrete(labels=c("1", "2", "3","4","5","6"))+
labs(x="Topic Rank")
awi <- ggplot(AvgWorry_int_mltd,
aes(x = variable, y = ID, fill = value))+
geom_tile()+scale_fill_gradient(high = "coral", low = "slateblue3") +
scale_x_discrete(labels=c("1", "2", "3","4","5","6"))+
labs(x="Topic Rank")
awf <- ggplot(AvgWorry_freq_mltd,
aes(x = variable, y = ID, fill = value))+
geom_tile()+scale_fill_gradient(high = "coral", low = "slateblue3") +
scale_x_discrete(labels=c("1", "2", "3","4","5","6"))+
labs(x="Topic Rank")
ari
arf
awi
awf
#write.csv(AvgWorryRum,"/Users/nikki/Desktop/Research/HealthyU_Scanning_Local/Emily_HUresources/WRGsurvey_Avg_31422.csv")
AvgWorryRum <- read.csv("/Users/nikki/Desktop/Research/HealthyU_Scanning_Local/Emily_HUresources/WRGsurvey_Avg_31422.csv")
AvgWorryRum <- AvgWorryRum[,-c(1)]
C. Examine and save out the topic rankings (Figure not included in
Manuscript)
#First pull necessary columns from big dataset
worryrank<-surv[,c(18,25:37)]
#rum rank should be different columns from worry rank:
rumrank<- surv[,c(18,174:186)]
colnames(worryrank) <- c("ID",
"School",
"Friend",
"Job",
"Health",
"Finances",
"Romance",
"WorldEvent",
"parents",
"Family",
"Safety",
"Roommates",
"Achievements",
"Other")
colnames(rumrank) <- c("ID",
"School",
"Friend",
"Job",
"Health",
"Finances",
"Romance",
"WorldEvent",
"parents",
"Family",
"Safety",
"Roommates",
"Achievements",
"Other")
# worrysummary <- psych::describe(worryrank)
# rumsummary <- psych::describe(rumrank)
worryrank_red <- subset(worryrank, worryrank$ID %in% unique(wide$Subject))
worrysummary_dis <- as.data.frame(psych::describe(worryrank_red))
worrysummary_dis<- worrysummary_dis[c(3,4)]
colnames(worrysummary_dis) <- c("w_mean","w_sd")
rumrank_red <- subset(rumrank, rumrank$ID %in% unique(wide$Subject))
rumsummary_dis <- as.data.frame(psych::describe(rumrank_red))
rumsummary_dis<- rumsummary_dis[c(3,4)]
colnames(rumsummary_dis) <- c("r_mean","r_sd")
wr_sum <- cbind(worrysummary_dis,rumsummary_dis)
wr_sum <- wr_sum[-c(1),-c(2,4)]
wr_sum$topic <- row.names(wr_sum)
wr_sum_ord<- wr_sum[order(wr_sum$w_mean),]
# create a long form data frame for plotting
wr_sum_long <- gather(wr_sum_ord, RNTtype, AvgRank, w_mean:r_mean, factor_key=TRUE)
wr_sum_long$topic <- as.factor(wr_sum_long$topic)
wr_sum_long <- wr_sum_long[order(wr_sum_long$AvgRank),]
# cisualize the topic rankings
ggplot(data=wr_sum_long, aes(x=factor(topic, level=c("Other","Roommates","WorldEvent","Safety","Job","Romance","Family","Finances","Health","Friend","parents","Achievements","School")), y=AvgRank, fill=RNTtype)) +
geom_bar(stat="identity", position=position_dodge())+
coord_flip()
#fill in rankings per subject with a for loop
# rows = topic
# in the cell = rank
for (i in 1:length(worryrank$Q19_1)){
worryrank$Schoollrank[i] <- which(worryrank[i,]==1)
worryrank$Friendrank[i] <- which(worryrank[i,]==2)
worryrank$Jobrank[i] <- which(worryrank[i,]==3)
worryrank$Healthrank[i] <- which(worryrank[i,]==4)
worryrank$Financesrank[i] <- which(worryrank[i,]==5)
worryrank$Romancerank[i] <- which(worryrank[i,]==6)
worryrank$Worldeventrank[i] <- which(worryrank[i,]==7)
worryrank$Parentsrank[i] <- which(worryrank[i,]==8)
worryrank$Familyrank[i] <-which(worryrank[i,]==9)
worryrank$Safetyrank[i] <- which(worryrank[i,]==10)
worryrank$Roommatesrank[i] <- which(worryrank[i,]==11)
worryrank$Achievementsrank[i] <- which(worryrank[i,]==12)
worryrank$Otherrank[i] <- which(worryrank[i,]==13)
}
for (i in 1:length(rumrank$Q54_1)){
rumrank$Schoolrank[i] <- which(rumrank[i,]==1)
rumrank$Friendrank[i] <- which(rumrank[i,]==2)
rumrank$Jobrank[i]<- which(rumrank[i,]==3)
rumrank$Healthrank[i] <- which(rumrank[i,]==4)
rumrank$Financesrank[i] <- which(rumrank[i,]==5)
rumrank$Romancerank[i] <- which(rumrank[i,]==6)
rumrank$Worldeventrank[i] <- which(rumrank[i,]==7)
rumrank$Parentsrank[i] <- which(rumrank[i,]==8)
rumrank$Familyrank[i] <- which(rumrank [i,]==9)
rumrank$Safetyrank[i] <- which(rumrank[i,]==10)
rumrank$Roommatesrank[i] <- which(rumrank[i,]==11)
rumrank$Achievementsrank[i] <- which(rumrank[i,]==12)
rumrank$Otherrank[i] <- which(rumrank[i,]==13)
}
ranksummary <- as.data.frame(matrix(NA, nrow=13, ncol=6))
colnames(ranksummary) <- c("topic","avg_worry_rank", "worry_rankrank", "avg_rum_rank", "rum_rankrank")
#manually enter in topics
ranksummary$topic <- c("School",
"Friend",
"Job",
"Health",
"Finances",
"Romance",
"WorldEvent",
"parents",
"Family",
"Safety",
"Roommates",
"Achievements",
"Other")
ranksummary$avg_worry_rank[1] <- mean(worryrank$Schoollrank, na.rm=TRUE)
ranksummary$avg_worry_rank[2] <- mean(worryrank$Friendrank, na.rm=TRUE)
ranksummary$avg_worry_rank[3] <- mean(worryrank$Jobrank, na.rm=TRUE)
ranksummary$avg_worry_rank[4] <- mean(worryrank$Healthrank, na.rm=TRUE)
ranksummary$avg_worry_rank[5] <- mean(worryrank$Financesrank, na.rm=TRUE)
ranksummary$avg_worry_rank[6] <- mean(worryrank$Romancerank, na.rm=TRUE)
ranksummary$avg_worry_rank[7] <- mean(worryrank$Worldeventrank, na.rm=TRUE)
ranksummary$avg_worry_rank[8] <- mean(worryrank$Parentsrank, na.rm=TRUE)
ranksummary$avg_worry_rank[9] <- mean(worryrank$Familyrank, na.rm=TRUE)
ranksummary$avg_worry_rank[10] <- mean(worryrank$Safetyrank, na.rm=TRUE)
ranksummary$avg_worry_rank[11] <- mean(worryrank$Roommatesrank, na.rm=TRUE)
ranksummary$avg_worry_rank[12] <- mean(worryrank$Achievementsrank, na.rm=TRUE)
ranksummary$avg_worry_rank[13] <- mean(worryrank$Otherrank, na.rm=TRUE)
#order the data by avg worry rank in order to assign rankrank
order <- order(ranksummary$avg_worry_rank, decreasing = FALSE)
ranksummary$worry_rankrank[c(order)] <- c(1:13)
##Repeat for rumination
#fill in avg. values
ranksummary$avg_rum_rank[1] <- mean(rumrank$Schoolrank, na.rm=TRUE)
ranksummary$avg_rum_rank[2] <- mean(rumrank$Friendrank, na.rm=TRUE)
ranksummary$avg_rum_rank[3] <- mean(rumrank$Jobrank, na.rm=TRUE)
ranksummary$avg_rum_rank[4] <- mean(rumrank$Healthrank, na.rm=TRUE)
ranksummary$avg_rum_rank[5] <- mean(rumrank$Financesrank, na.rm=TRUE)
ranksummary$avg_rum_rank[6] <- mean(rumrank$Romancerank, na.rm=TRUE)
ranksummary$avg_rum_rank[7] <- mean(rumrank$Worldeventrank, na.rm=TRUE)
ranksummary$avg_rum_rank[8] <- mean(rumrank$Parentsrank, na.rm=TRUE)
ranksummary$avg_rum_rank[9] <- mean(rumrank$Familyrank, na.rm=TRUE)
ranksummary$avg_rum_rank[10] <- mean(rumrank$Safetyrank, na.rm=TRUE)
ranksummary$avg_rum_rank[11] <- mean(rumrank$Roommatesrank, na.rm=TRUE)
ranksummary$avg_rum_rank[12] <- mean(rumrank$Achievementsrank, na.rm=TRUE)
ranksummary$avg_rum_rank[13] <- mean(rumrank$Otherrank, na.rm=TRUE)
#order the data by avg rum rank in order to assign rankrank
order <- order(ranksummary$avg_rum_rank, decreasing = FALSE)
ranksummary$rum_rankrank[c(order)] <- c(1:13)
#writeout
write.csv(ranksummary, "/Users/nikki/Desktop/Research/HealthyU_Scanning_Local/Emily_HUresources/ranksummary_31422.csv")
E. SURVEY Multilevel models: INTENSITY~FREQ and BIGGER
RATING~TOPIC+RNTtype+RATEtype
Num_statements<-aggregate(long_surv_reduc$value,list(long_surv_reduc$ID),length) #confirm everyone has 24 statements
# multilevel regression:
# Does Surv Freq predict Int?
long_surv_reduc <- long_surv_reduc[c(1,3,4,6,7)]
surv_reduc <- spread(long_surv_reduc, "Ratingtype", "value")
avgF<-aggregate(surv_reduc$Freq,list(surv_reduc$ID),mean,na.rm=T)
names(avgF)<-c("ID","avgF") #rename variables
avgI<-aggregate(surv_reduc$Int,list(surv_reduc$ID),mean,na.rm=T)
names(avgI)<-c("ID","avgI") #rename variables
avgs <- merge(avgF,avgI,by = "ID")
surv_reduc <- merge(surv_reduc, avgs,by = "ID")
surv_reduc$Freq_c <- surv_reduc$Freq - surv_reduc$avgF
surv_reduc$Int_c <- surv_reduc$Int - surv_reduc$avgI
#Surv Freq significantly predicts Int rating
summary(lmer(Int~Freq+(1+Freq|ID), data=surv_reduc))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Int ~ Freq + (1 + Freq | ID)
## Data: surv_reduc
##
## REML criterion at convergence: 808.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.9717 -0.5649 -0.0358 0.6372 2.9081
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## ID (Intercept) 0.23724 0.4871
## Freq 0.01711 0.1308 -0.87
## Residual 0.29774 0.5457
## Number of obs: 456, groups: ID, 39
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.93136 0.12268 34.71467 7.592 7.06e-09 ***
## Freq 0.60341 0.04364 36.64969 13.828 4.07e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## Freq -0.916
# Visualized (not included in the manuscipt)
surv_reduc$ID <- as.numeric(as.character(surv_reduc$ID))
ss_agree<- ggplot(data = surv_reduc, aes(x = Freq_c, y = Int_c, group = ID, colour = ID)) +
geom_smooth(method = "lm", se = F, size = .35, alpha = .65) +
scale_color_gradient(low = "mediumpurple4", high = "plum3") +
geom_smooth(aes(group = 1), size = 1, colour = "black", method = "lm", se = F) +
labs(title="", x = "Within-person Survey Frequency Ratings", y = "Within-person Survey Intensity Ratings") +
xlim(-2, 2) +
ylim(-2, 2)
ss_agree + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"))
G. Assess whether survey responses correlate with trait responses
(Figure 1 and Supplemental Figure 3)
# pull in the trait measures to this long form dataframe with survey information
l_srv_trt <- merge(long_surv_reduc,meas, by.x = "ID", by.y ="id",all.x = T)
l_srv_trt$RNT_1_sc <- scale(l_srv_trt$PTQ_tot_T1)
l_srv_trt$RNT_3_sc <- scale(l_srv_trt$PTQ_tot_T3)
l_srv_trt$TrWorry_1_sc <- scale(l_srv_trt$PSWQ_tot_T1)
l_srv_trt$TrWorry_3_sc <- scale(l_srv_trt$PSWQ_tot_T3)
l_srv_trt$TrRum_1_sc <- scale(l_srv_trt$RRS_tot_T1)
l_srv_trt$TrRum_3_sc <- scale(l_srv_trt$RRS_tot_T3)
# model whether survey ratings are significanly related to various trait measure scores
###only End of semester RNT is related to both W and R survey scores
###
#trending (p=0.072) worry ~ PTQ relationship
srv_trt_1 <- lmer(value~RNT_1_sc*RNTtype+(1|ID), data=l_srv_trt)
summary(srv_trt_1)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ RNT_1_sc * RNTtype + (1 | ID)
## Data: l_srv_trt
##
## REML criterion at convergence: 1799.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3427 -0.6989 0.0004 0.6784 3.2332
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 0.3027 0.5501
## Residual 0.4040 0.6356
## Number of obs: 870, groups: ID, 37
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.20113 0.09552 38.95969 23.043 < 2e-16 ***
## RNT_1_sc 0.07886 0.09553 38.88480 0.826 0.4141
## RNTtypeWorry 0.27543 0.04323 831.77261 6.371 3.11e-10 ***
## RNT_1_sc:RNTtypeWorry 0.09734 0.04301 831.47376 2.263 0.0239 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) RNT_1_s RNTtyW
## RNT_1_sc -0.001
## RNTtypeWrry -0.229 0.003
## RNT_1_:RNTW 0.003 -0.228 -0.003
simple_slopes(srv_trt_1)
## RNT_1_sc RNTtype Test Estimate Std. Error df t value Pr(>|t|) Sig.
## 1 -1.00393 sstest 0.1777 0.0612 831.9524 2.9042 0.00378 **
## 2 0.00068 sstest 0.2755 0.0432 831.7722 6.3726 3.076e-10 ***
## 3 1.005291 sstest 0.3733 0.0610 831.2911 6.1146 1.488e-09 ***
## 4 sstest Rum 0.0789 0.0955 38.8848 0.8255 0.41411
## 5 sstest Worry 0.1762 0.0954 38.6925 1.8468 0.07244 .
#sig worry ~ PTQ AND rum ~ PTQ relationships
srv_trt_1.5 <- lmer(value~RNT_3_sc*RNTtype+(1|ID), data=l_srv_trt)
summary(srv_trt_1.5)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ RNT_3_sc * RNTtype + (1 | ID)
## Data: l_srv_trt
##
## REML criterion at convergence: 1694.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.4461 -0.6422 0.0110 0.6787 3.2669
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 0.2498 0.4998
## Residual 0.3837 0.6194
## Number of obs: 842, groups: ID, 36
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.17999 0.08876 38.55422 24.561 < 2e-16 ***
## RNT_3_sc 0.21267 0.08864 38.27131 2.399 0.0214 *
## RNTtypeWorry 0.30132 0.04287 805.06729 7.028 4.47e-12 ***
## RNT_3_sc:RNTtypeWorry 0.07118 0.04344 805.13135 1.639 0.1017
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) RNT_3_s RNTtyW
## RNT_3_sc -0.002
## RNTtypeWrry -0.247 0.004
## RNT_3_:RNTW 0.004 -0.237 0.008
simple_slopes(srv_trt_1.5)
## RNT_3_sc RNTtype Test Estimate Std. Error df t value Pr(>|t|) Sig.
## 1 -0.996298 sstest 0.2304 0.0607 805.1777 3.7971 0.0001574 ***
## 2 -0.007732 sstest 0.3008 0.0429 805.0685 7.0155 4.862e-12 ***
## 3 0.980835 sstest 0.3711 0.0607 805.0222 6.1162 1.494e-09 ***
## 4 sstest Rum 0.2127 0.0886 38.2713 2.3993 0.0214019 *
## 5 sstest Worry 0.2839 0.0890 38.8879 3.1894 0.0028159 **
###no main effects, only End of semester PSWQ is related to Worry survey scores
###
#Same as PTQ -> trending (p=0.083) worry ~ PSWQ relationship
srv_trt_2 <- lmer(value~TrWorry_1_sc*RNTtype+(1|ID), data=l_srv_trt)
summary(srv_trt_2)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ TrWorry_1_sc * RNTtype + (1 | ID)
## Data: l_srv_trt
##
## REML criterion at convergence: 1803.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3678 -0.7031 -0.0028 0.6942 3.1705
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 0.2947 0.5428
## Residual 0.4063 0.6374
## Number of obs: 870, groups: ID, 37
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.20073 0.09441 39.08875 23.309 < 2e-16 ***
## TrWorry_1_sc 0.14489 0.09449 39.11749 1.533 0.133
## RNTtypeWorry 0.27544 0.04336 831.80603 6.352 3.49e-10 ***
## TrWorry_1_sc:RNTtypeWorry 0.02263 0.04327 831.87967 0.523 0.601
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrW_1_ RNTtyW
## TrWrry_1_sc -0.003
## RNTtypeWrry -0.232 0.007
## TrW_1_:RNTW 0.007 -0.234 -0.016
simple_slopes(srv_trt_2)
## TrWorry_1_sc RNTtype Test Estimate Std. Error df t value Pr(>|t|) Sig.
## 1 -0.991705 sstest 0.2530 0.0615 832.5081 4.1140 4.276e-05 ***
## 2 0.011079 sstest 0.2757 0.0434 831.7910 6.3587 3.352e-10 ***
## 3 1.013864 sstest 0.2984 0.0612 831.1498 4.8770 1.291e-06 ***
## 4 sstest Rum 0.1449 0.0945 39.1175 1.5335 0.13321
## 5 sstest Worry 0.1675 0.0943 38.8021 1.7766 0.08347 .
#sig worry ~ PSWQ relationship
srv_trt_2.5 <- lmer(value~TrWorry_3_sc*RNTtype+(1|ID), data=l_srv_trt)
summary(srv_trt_2.5)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ TrWorry_3_sc * RNTtype + (1 | ID)
## Data: l_srv_trt
##
## REML criterion at convergence: 1675
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.4011 -0.6746 0.0224 0.6995 3.1209
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 0.2449 0.4948
## Residual 0.3750 0.6124
## Number of obs: 842, groups: ID, 36
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.18617 0.08788 38.67642 24.877 < 2e-16 ***
## TrWorry_3_sc 0.14954 0.08845 39.57091 1.691 0.0988 .
## RNTtypeWorry 0.29509 0.04241 805.11251 6.958 7.15e-12 ***
## TrWorry_3_sc:RNTtypeWorry 0.20127 0.04358 807.29852 4.618 4.51e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrW_3_ RNTtyW
## TrWrry_3_sc -0.009
## RNTtypeWrry -0.247 0.019
## TrW_3_:RNTW 0.019 -0.263 -0.032
simple_slopes(srv_trt_2.5)
## TrWorry_3_sc RNTtype Test Estimate Std. Error df t value Pr(>|t|) Sig.
## 1 -0.961432 sstest 0.1016 0.0606 807.4051 1.6769 0.0939490 .
## 2 0.020652 sstest 0.2992 0.0424 805.0579 7.0593 3.616e-12 ***
## 3 1.002736 sstest 0.4969 0.0599 804.9419 8.2950 4.548e-16 ***
## 4 sstest Rum 0.1495 0.0885 39.5709 1.6906 0.0987693 .
## 5 sstest Worry 0.3508 0.0877 38.3056 3.9995 0.0002802 ***
### no main effects, end of semester RRS is related to Worry survey scores
###
# RRS sig predicts worry ratings
srv_trt_3 <- lmer(value~TrRum_1_sc*RNTtype+(1|ID), data=l_srv_trt)
summary(srv_trt_3)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ TrRum_1_sc * RNTtype + (1 | ID)
## Data: l_srv_trt
##
## REML criterion at convergence: 1801.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3830 -0.7051 -0.0079 0.6844 3.1704
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 0.2875 0.5362
## Residual 0.4060 0.6371
## Number of obs: 870, groups: ID, 37
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.20037 0.09337 39.11017 23.565 < 2e-16 ***
## TrRum_1_sc 0.15702 0.09333 38.94643 1.683 0.10
## RNTtypeWorry 0.27652 0.04335 831.78464 6.379 2.96e-10 ***
## TrRum_1_sc:RNTtypeWorry 0.04226 0.04331 831.55831 0.976 0.33
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrR_1_ RNTtyW
## TrRum_1_sc 0.003
## RNTtypeWrry -0.235 -0.006
## TrR_1_:RNTW -0.006 -0.231 0.024
simple_slopes(srv_trt_3)
## TrRum_1_sc RNTtype Test Estimate Std. Error df t value Pr(>|t|) Sig.
## 1 -1.019066 sstest 0.2335 0.0611 830.9688 3.8193 0.0001438 ***
## 2 -0.018591 sstest 0.2757 0.0433 831.7594 6.3626 3.274e-10 ***
## 3 0.981884 sstest 0.3180 0.0614 832.3360 5.1757 2.847e-07 ***
## 4 sstest Rum 0.1570 0.0933 38.9464 1.6825 0.1004717
## 5 sstest Worry 0.1993 0.0934 39.0544 2.1338 0.0391920 *
#
srv_trt_3.5 <- lmer(value~TrRum_3_sc*RNTtype+(1|ID), data=l_srv_trt)
summary(srv_trt_3.5)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: value ~ TrRum_3_sc * RNTtype + (1 | ID)
## Data: l_srv_trt
##
## REML criterion at convergence: 1685.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6069 -0.6537 0.0065 0.6871 3.2779
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 0.2966 0.5446
## Residual 0.3770 0.6140
## Number of obs: 842, groups: ID, 36
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.17495 0.09573 37.63412 22.720 < 2e-16 ***
## TrRum_3_sc 0.05648 0.09622 38.30374 0.587 0.561
## RNTtypeWorry 0.30519 0.04252 804.80862 7.177 1.62e-12 ***
## TrRum_3_sc:RNTtypeWorry 0.17525 0.04341 806.35535 4.037 5.94e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrR_3_ RNTtyW
## TrRum_3_sc 0.006
## RNTtypeWrry -0.227 -0.013
## TrR_3_:RNTW -0.013 -0.242 0.031
simple_slopes(srv_trt_3.5)
## TrRum_3_sc RNTtype Test Estimate Std. Error df t value Pr(>|t|) Sig.
## 1 -1.008096 sstest 0.1285 0.0601 804.5621 2.1396 0.03269 *
## 2 -0.020808 sstest 0.3015 0.0425 804.7650 7.0947 2.844e-12 ***
## 3 0.966479 sstest 0.4746 0.0607 806.5476 7.8241 1.603e-14 ***
## 4 sstest Rum 0.0565 0.0962 38.3037 0.5870 0.56066
## 5 sstest Worry 0.2317 0.0955 37.2243 2.4261 0.02023 *
##End of Semester RNT
Surv_tr <- ggplot(data = l_srv_trt, aes(x = RNT_3_sc, y = value)) +
stat_summary(geom = "point", fun.x = "mean", col = "palegreen3",size = 2, shape = 24, fill = "palegreen3")+
#geom_point(aes(x = RNT_3_sc, group=ID), color = "darkolivegreen", size = 1, alpha = .5) +
geom_smooth(aes(x = RNT_3_sc), color = "palegreen4", method = "lm", size = 1,se = F) +
#geom_smooth(aes(x = TrWorry_3_sc), color = "dodgerblue3", method = "lm", size = 1,se = F) +
#geom_smooth(aes(x = TrRum_3_sc), color = "mediumpurple3", method = "lm", size = 1,se = F) +
ylim(1.25, 3.75) +
xlim(-3.75, 3.75) +
facet_wrap(~RNTtype) +
labs(title="", x = "End of Semester Trait RNT Scores", y = "Average Survey Ratings") +
theme(strip.text.x = element_text(size = 15)) +
theme(axis.title = element_text(size = 15))
Surv_tr + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"))
#End of Semester worry rum
ggplot(data = l_srv_trt, aes(y = value)) +
geom_smooth(aes(x = RNT_3_sc), color = "darkolivegreen3", method = "lm", size = 1,se = F) +
geom_smooth(aes(x = TrWorry_3_sc), color = "dodgerblue3", method = "lm", size = 1,se = F) +
geom_smooth(aes(x = TrRum_3_sc), color = "mediumpurple3", method = "lm", size = 1,se = F) +
ylim(1.75, 3.25) +
xlim(-3, 3) +
facet_wrap(~RNTtype) +
labs(title="", x = "End of Semester: Normalized Trait Questionnaire Scores", y = "Survey Rating") +
theme(strip.text.x = element_text(size = 15)) +
theme(axis.title = element_text(size = 15))
#Beginning of Semester worry rum
ggplot(data = l_srv_trt, aes(y = value)) +
geom_smooth(aes(x = RNT_1_sc), color = "darkolivegreen3", method = "lm", size = .75,linetype = 5, se = F) +
geom_smooth(aes(x = TrWorry_1_sc), color = "dodgerblue3", method = "lm", size = .75,linetype = 5, se = F) +
geom_smooth(aes(x = TrRum_1_sc), color = "mediumpurple3", method = "lm",size = .75,linetype = 5, se = F) +
ylim(1.75, 3.25) +
xlim(-3, 3) +
facet_wrap(~RNTtype) +
labs(title="", x = "Start of Semester: Normalized Trait Questionnaire Scores", y = "Survey Rating") +
theme(strip.text.x = element_text(size = 15)) +
theme(axis.title = element_text(size = 15))
2. Worry and Rumination fMRI paradigm
H. Assess whether scanner and survey intensity ratings agree (Figure
1)
# subset the long fmri by intensity rating (ignore freq) and then create both within and between person effects for the model
# the "long" data frame contains both survey and scanner ratings for each statement
long[long==99]<-NA
int <- subset(long, long$DataType == "IntensityRate")
surv.int.bw <- aggregate(qualInt ~ Subject, long, mean, na.action = na.omit )
colnames(surv.int.bw) <- c("Subject","surv.int.bw")
int<- merge(int,surv.int.bw,by="Subject",all.x=T)
int$surv.int.win <- calc.mcent(qualInt, Subject, data=int)
int$surv.int.bw <- calc.mean(qualInt, Subject, data=int, expand=TRUE)
# without random slope
agreement_mod <- lmer(Response ~ int$surv.int.win + int$surv.int.bw + (1 | Subject), data = int)
summary(agreement_mod)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Response ~ int$surv.int.win + int$surv.int.bw + (1 | Subject)
## Data: int
##
## REML criterion at convergence: 3701.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.99320 -0.67290 0.00191 0.70977 2.69684
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subject (Intercept) 0.1676 0.4094
## Residual 0.6037 0.7770
## Number of obs: 1541, groups: Subject, 39
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.798e+00 2.699e-01 3.701e+01 6.663 7.99e-08 ***
## int$surv.int.win 2.632e-01 2.375e-02 1.501e+03 11.082 < 2e-16 ***
## int$surv.int.bw 2.624e-01 9.612e-02 3.700e+01 2.730 0.00964 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) int$srv.nt.w
## int$srv.nt.w 0.000
## int$srv.nt.b -0.967 0.000
# with random slope
agreement_mod2 <- lmer(Response ~ surv.int.win + surv.int.bw + (1 + surv.int.win| Subject), data = int)
summary(agreement_mod2)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Response ~ surv.int.win + surv.int.bw + (1 + surv.int.win | Subject)
## Data: int
##
## REML criterion at convergence: 3680.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.71044 -0.63280 -0.03213 0.69443 2.69198
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subject (Intercept) 0.16799 0.4099
## surv.int.win 0.03048 0.1746 0.31
## Residual 0.58365 0.7640
## Number of obs: 1541, groups: Subject, 39
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.83559 0.26455 37.70222 6.939 3.10e-08 ***
## surv.int.win 0.27345 0.03804 32.52298 7.189 3.33e-08 ***
## surv.int.bw 0.24864 0.09408 37.57513 2.643 0.0119 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) srv.nt.w
## surv.int.wn 0.085
## surv.int.bw -0.966 -0.029
# including the random slope fits better
anova(agreement_mod2,agreement_mod)
## Data: int
## Models:
## agreement_mod: Response ~ int$surv.int.win + int$surv.int.bw + (1 | Subject)
## agreement_mod2: Response ~ surv.int.win + surv.int.bw + (1 + surv.int.win | Subject)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## agreement_mod 5 3699.1 3725.8 -1844.6 3689.1
## agreement_mod2 7 3682.9 3720.3 -1834.5 3668.9 20.175 2 4.16e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# visualize
ggplot(data = int, aes(y = Response, group = Subject, color = Subject)) +
geom_smooth(aes(x= qualInt), method="lm",se=F)
int$norm.surv.int.win <- scale(int$surv.int.win)
int$norm.surv.int.bw <- scale(int$surv.int.bw)
win <- ggplot(data = int, aes(y = Response)) +
geom_smooth(aes(x = norm.surv.int.win, group = Subject, color = Subject), alpha = .4, size = .75, method = "lm", se = F) +
geom_smooth(aes(x = norm.surv.int.win, group = 1), color = "black",size=1, method = "lm", se = F) +
scale_color_gradient(low = "indianred4", high = "indianred2") +
xlim(-3,3) +
ylim(.75,4) +
xlab("Normalized Within-person Survey Intensity Ratings") +
ylab("Scanner Intensity Ratings")
win <- win + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"))
bw <- ggplot(data = int, aes(y = Response)) +
geom_smooth(aes(x = norm.surv.int.bw), color = "black", size = 1.5,method = "lm", se = F) +
stat_summary(aes(x = norm.surv.int.bw,group = Subject, color = Subject,fill=Subject), geom = "point", fun.y = "mean",
size = 3, shape = 17, alpha = .75) +
scale_color_gradient(low = "indianred4", high = "indianred2") +
scale_fill_gradient(low = "indianred4", high = "indianred2") +
xlim(-3,3) +
ylim(.75,4) +
xlab("Normalized Person Mean Survey Intensity Ratings") +
ylab("Scanner Intensity Ratings")
bw <- bw + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"))
grid.arrange(win,bw,nrow=2)
I. Confirm the worry and rumination are more negative and intense
than neutral stim in the scanner (Figure 1)
options(scipen = 100)
fr <- subset(long, DataType=="FeelingRate")
fr$Response <- as.numeric(fr$Response)
fr <- fr[-which(is.na(fr$Response)),]
fr1 <- subset(long, DataType=="IntensityRate")
#this makes r the intercept
fr$EventType <- factor(fr$EventType,
levels=c("rum","neutral","worry"),
labels=c("R","N","W"))
scansanity_mod1 <- lmer(Response ~ EventType + (1 | Subject), data = fr)
summary(scansanity_mod1)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Response ~ EventType + (1 | Subject)
## Data: fr
##
## REML criterion at convergence: 2986.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.0323 -0.5348 -0.3012 -0.0217 3.8104
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subject (Intercept) 0.02365 0.1538
## Residual 0.23377 0.4835
## Number of obs: 2095, groups: Subject, 39
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.29822 0.03026 61.55969 42.903 < 0.0000000000000002 ***
## EventTypeN 0.91834 0.02678 2054.23337 34.295 < 0.0000000000000002 ***
## EventTypeW -0.14512 0.02479 2054.34392 -5.853 0.0000000056 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) EvntTN
## EventTypeN -0.381
## EventTypeW -0.412 0.465
#this makes r the intercept
fr1$EventType <- factor(fr1$EventType,
levels=c("rum","neutral","worry"),
labels=c("R","N","W"))
scansanity_mod2 <- lmer(Response ~ EventType + (1 | Subject), data = fr1)
summary(scansanity_mod2)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Response ~ EventType + (1 | Subject)
## Data: fr1
##
## REML criterion at convergence: 4887
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.8324 -0.6143 -0.0884 0.6442 3.2386
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subject (Intercept) 0.1225 0.3500
## Residual 0.5567 0.7461
## Number of obs: 2121, groups: Subject, 39
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.34050 0.06215 48.97226 37.656 <0.0000000000000002 ***
## EventTypeN -1.08126 0.04103 2080.06172 -26.352 <0.0000000000000002 ***
## EventTypeW 0.34197 0.03800 2080.06604 8.999 <0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) EvntTN
## EventTypeN -0.283
## EventTypeW -0.306 0.463
##FREQUENCY
Feels<-ggplot(fr, aes(EventType, Response))+
geom_violin(scale = "area", fill = "palegreen3", color="palegreen3")+
stat_summary(fun=mean, color="palegreen4", geom="point", size = 1) +
#ggtitle("Valence Rating by Statement Type")+
theme(plot.title = element_text(hjust = 0.5)) + labs(y="Valence Rating", x="Statement Type")
Feels <- Feels + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"))
##INTESITY
Ints<-ggplot(fr1, aes(EventType, Response))+
geom_violin(scale = "area", fill = "indianred2",color= NA)+
#ggtitle("Intensity Rating by Statement Type")+
stat_summary(fun=mean, color="indianred4", geom="point", size = 1) +
theme(plot.title = element_text(hjust = 0.5)) + labs(y="Intensity Rating", x="Statement Type")
Ints <- Ints+ theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"))
grid.arrange(Feels,Ints, nrow = 2)
