OSCR

The neurocomputational mechanisms underlying the impact of social comparison on effort investment.

Code ↔ Paper

The paper beside its authors' code: matches between them have not been computed for this paper yet.

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

R Markdown · 1,925 lines · 71 KB · no license

  1. ---
  2. title: "Exp1 and Exp2"
  3. author: "JiaRui"
  4. date: "2026-2-11"
  5. output: html_document
  6. ---
  8. ```{r setup, include=FALSE}
  9. library(Rmisc)
  10. library(tidyverse)
  11. library(bruceR)
  12. library(ggpubr)
  13. library(afex)
  14. library(emmeans)
  15. library(sjPlot)
  16. library(readxl)
  17. library(dplyr)
  18. library(tidyr)
  19. library(xlsx)
  20. library(ggplot2)
  21. library(ggtext)
  22. library(psych)
  23. library(effects)
  24. library(brms)
  25. library(parameters)
  26. library(simr)
  28. gg.rr <-theme(axis.line = element_line(colour = "black"),text = element_text(family='Arial',size=11),
  29. axis.title.y = element_markdown(margin = margin(t = 20, r = 5, b = 20, l = 0)),
  30. axis.title.x = element_markdown(),
  31. panel.grid.major = element_blank(),
  32. panel.grid.minor = element_blank(),
  33. panel.border = element_blank(),
  34. panel.background = element_blank(),
  35. strip.background = element_blank(),
  36. plot.title = element_markdown(hjust = 0.5,size=12))
  38. gg.side <-theme(axis.line = element_line(colour = "black"),text =element_text(family='serif',size=16),
  39. axis.title.y = element_markdown(margin = margin(t = 0, r = 8, b = 0, l = 0)), ##up right down left
  40. axis.title.x = element_markdown(margin = margin(t = 8, r = 0, b = 0, l = 0)),axis.text = element_text(colour = "black"),
  41. panel.grid.major = element_blank(),
  42. panel.grid.minor = element_blank(),
  43. panel.border = element_blank(),
  44. panel.background = element_blank(),
  45. strip.background = element_blank(),
  46. plot.title = element_markdown(hjust = 0.5,size=12))
  47. ```
  50. # Exp1 and Exp2 behavior analysis
  51. ```{r import exp1 and exp2 efficacy data}
  52. ###Set working directory###
  53. set.wd("")
  55. ###exp1###
  56. exp1_interval = read_csv("/beha/exp1_interval.csv")
  57. exp1_interval_learn = read_csv("/beha/exp1_interval_learnHB.csv")
  58. exp1_check = read_excel("/beha/exp1_check.xlsx")
  60. exp1_interval_learn = exp1_interval_learn %>%
  61. group_by(SubID) %>%
  62. mutate(Cong_prev = lag(meanCongruency, 1)) %>%
  63. dplyr::mutate(un_PE_prev = lag(unsigned_PE_efficacy, 1)) %>%
  64. dplyr::mutate(PE_prev = lag(signed_PE_efficacy, 1)) %>%
  65. dplyr::mutate(zCRPS = scale(CRPS)) %>%
  66. dplyr::mutate(zmeanRT = scale(meanRT)) %>%
  67. dplyr::mutate(zmbased_efficacy_prev = scale(mbased_efficacy_prev)) %>%
  68. dplyr::mutate(zun_PE_prev = scale(un_PE_prev)) %>%
  69. dplyr::mutate(zPE_prev = scale(PE_prev))
  71. exp1_interval_learn$duration = as.numeric(exp1_interval_learn$duration)
  72. exp1_interval$gender = as.factor(exp1_interval$gender)
  73. exp1_interval_learn$gender = as.factor(exp1_interval_learn$gender)
  74. exp1_interval_learn$FeedCode = as.factor(exp1_interval_learn$FeedCode)
  75. exp1_interval_learn$FeedPrev = as.factor(exp1_interval_learn$FeedPrev)
  76. exp1_interval_learn$mbased_efficacy = as.numeric(exp1_interval_learn$mbased_efficacy)
  77. exp1_interval_learn$mbased_efficacy_prev = as.numeric(exp1_interval_learn$mbased_efficacy_prev)
  78. exp1_interval_learn$unsigned_PE_efficacy = as.numeric(exp1_interval_learn$unsigned_PE_efficacy)
  79. exp1_interval_learn$un_PE_prev = as.numeric(exp1_interval_learn$un_PE_prev)
  80. exp1_interval_learn$signed_PE_efficacy = as.numeric(exp1_interval_learn$signed_PE_efficacy)
  81. exp1_interval_learn$PE_prev = as.numeric(exp1_interval_learn$PE_prev)
  82. exp1_interval$FeedPrev = as.factor(exp1_interval$FeedPrev)
  84. ###exp2###
  85. exp2_interval = read_csv("/beha/exp2_interval.csv")
  86. exp2_interval_learn = read_csv("/beha/exp2_interval_learnHB.csv")
  88. exp2_check = read_excel("/beha/exp2_check.xlsx")
  89. exp2_interval_learn = exp2_interval_learn %>%
  90. group_by(SubID) %>%
  91. mutate(Cong_prev = lag(meanCongruency, 1)) %>%
  92. mutate(un_PE_prev = lag(unsigned_PE_efficacy, 1)) %>%
  93. mutate(PE_prev = lag(signed_PE_efficacy, 1))
  95. exp2_interval_learn$duration = as.numeric(exp2_interval_learn$duration)
  96. exp2_interval$gender = as.factor(exp2_interval$gender)
  97. exp2_interval_learn$gender = as.factor(exp2_interval_learn$gender)
  98. exp2_interval_learn$FeedCode = as.factor(exp2_interval_learn$FeedCode)
  99. exp2_interval_learn$FeedPrev = as.factor(exp2_interval_learn$FeedPrev)
  100. exp2_interval_learn$mbased_efficacy = as.numeric(exp2_interval_learn$mbased_efficacy)
  101. exp2_interval_learn$mbased_efficacy_prev = as.numeric(exp2_interval_learn$mbased_efficacy_prev)
  102. exp1_interval_learn$unsigned_PE_efficacy = as.numeric(exp1_interval_learn$unsigned_PE_efficacy)
  103. exp2_interval_learn$un_PE_prev = as.numeric(exp2_interval_learn$un_PE_prev)
  104. exp1_interval_learn$signed_PE_efficacy = as.numeric(exp1_interval_learn$signed_PE_efficacy)
  105. exp2_interval_learn$PE_prev = as.numeric(exp2_interval_learn$PE_prev)
  106. exp2_interval$FeedPrev = as.factor(exp2_interval$FeedPrev)
  107. ```
  109. ### Exp1 and Exp2 post rating
  110. ```{r exp1 and exp2 post rating des}
  111. exp1_check_emotion_sub = exp1_check %>%
  112. dplyr::summarise(emotion_downward = mean(emotion1),
  113. emotion_upward = mean(emotion2),
  114. emotion_lateral = mean(emotion3))
  115. exp1_check_emotion_sub
  117. exp2_check_emotion_sub = exp2_check %>%
  118. dplyr::summarise(emotion_downward = mean(emotion1),
  119. emotion_upward = mean(emotion2),
  120. emotion_lateral = mean(emotion3))
  121. exp2_check_emotion_sub
  122. ```
  124. ```{r exp1 check post emotion}
  125. exp1_check_emotion = exp1_check %>%
  126. pivot_longer(cols = starts_with("emotion"),
  127. names_to = "cond",
  128. values_to = "emotion") %>%
  129. mutate(cond = factor(cond,
  130. levels = c("emotion1","emotion2","emotion3"),
  131. labels = c("emotion_downward","emotion_upward","emotion_lateral")))
  133. exp1_check_emotion_avo = aov_ez(id = "SubID", dv = "emotion",
  134. within = "cond", data = exp1_check_emotion)
  135. exp1_check_emotion_avo
  137. exp1_check_emotion_post = pairs(emmeans(exp1_check_emotion_avo, ~ cond),
  138. adjust = "bonferroni")
  139. exp1_check_emotion_post
  141. ```
  144. ```{r exp2 check post emotion}
  145. exp2_check_emotion = exp2_check %>%
  146. pivot_longer(cols = starts_with("emotion"),
  147. names_to = "cond",
  148. values_to = "emotion") %>%
  149. mutate(cond = factor(cond,
  150. levels = c("emotion1","emotion2","emotion3"),
  151. labels = c("emotion_downward","emotion_upward","emotion_lateral")))
  153. exp2_check_emotion_avo = aov_ez(id = "SubID", dv = "emotion",
  154. within = "cond", data = exp2_check_emotion)
  155. exp2_check_emotion_avo
  157. exp2_check_emotion_post = pairs(emmeans(exp2_check_emotion_avo, ~ cond),
  158. adjust = "bonferroni")
  159. exp2_check_emotion_post
  160. ```
  165. ### Exp1 descriptive statistics
  166. ```{r exp1 rate_score and model_based efficacy descriptive}
  167. exp1_interval_efficacy = exp1_interval_learn %>%
  168. group_by(SubID) %>%
  169. mutate(FeedCode = factor(FeedCode, levels = c("1", "2"),
  170. labels=c("Downward","Upward")))
  172. contrasts(exp1_interval_efficacy$gender) = contr.sum
  173. contrasts(exp1_interval_efficacy$FeedCode) = contr.sum
  175. ###exp1 rate_score###
  176. exp1_rate_feedback = exp1_interval_learn %>%
  177. dplyr::group_by(SubID) %>%
  178. fill(rate_score, .direction = "up")%>%
  179. drop_na(rate_score)%>%
  180. dplyr::group_by(FeedCode) %>%
  181. dplyr::summarise(rate_score_mean = round(mean(rate_score),3),
  182. rate_score_sd = round(sd(rate_score),3),
  183. rate_score_se=rate_score_sd/sqrt(32))
  184. exp1_rate_feedback
  186. ###exp1 model_based efficacy###
  187. exp1_efficacy_feedback = exp1_interval_learn %>%
  188. dplyr::group_by(SubID) %>%
  189. fill(mbased_efficacy, .direction = "up")%>%
  190. drop_na(mbased_efficacy)%>%
  191. dplyr::group_by(FeedCode) %>%
  192. dplyr::summarise(mbased_efficacy_mean = round(mean(mbased_efficacy),3),
  193. mbased_efficacy_sd = round(sd(mbased_efficacy),3),
  194. mbased_efficacy_se=mbased_efficacy_sd/sqrt(32))
  196. exp1_efficacy_feedback
  198. ```
  200. ```{r exp1 behavior descriptive}
  201. ###exp1 CRPS ###
  202. exp1_interval_CRPS = exp1_interval %>%
  203. dplyr::group_by(SubID, FeedPrev) %>%
  204. filter(CRPS > mean(CRPS) - 3 * sd(CRPS),
  205. CRPS < mean(CRPS) + 3 * sd(CRPS))
  207. exp1_CRPS_feedback = exp1_interval_CRPS %>%
  208. dplyr::group_by(FeedPrev) %>%
  209. dplyr::summarise(CRPS_mean = round(mean(CRPS),3),
  210. CRPS_sd = round(sd(CRPS),3),
  211. CRPS_se=CRPS_sd/sqrt(32))
  213. exp1_CRPS_feedback
  215. ###exp1 RT ###
  216. exp1_interval_RT = exp1_interval %>%
  217. dplyr::group_by(SubID, FeedPrev) %>%
  218. filter(meanRT > mean(meanRT) - 3 * sd(meanRT),
  219. meanRT < mean(meanRT) + 3 * sd(meanRT))
  221. exp1_RT_feedback = exp1_interval_RT %>%
  222. dplyr::group_by(FeedPrev) %>%
  223. dplyr::summarise(RT_mean = round(mean(meanRT)),
  224. RT_sd = round(sd(meanRT)),
  225. RT_se=RT_sd/sqrt(32))
  227. exp1_RT_feedback
  229. ###exp1 ACC ###
  230. exp1_ACC_feedback = exp1_interval %>%
  231. drop_na(FeedPrev) %>%
  232. dplyr::group_by(FeedPrev) %>%
  233. dplyr::summarise(ACC_mean = round(mean(ACC),3),
  234. ACC_sd = round(sd(ACC),3),
  235. ACC_se=ACC_sd/sqrt(32))
  236. exp1_ACC_feedback
  238. ###exp1 feedback count###
  239. exp1_feedback_sub = exp1_interval %>%
  240. group_by(SubID,FeedCode) %>%
  241. dplyr::summarise(counts = n()) %>%
  242. ungroup()
  244. exp1_feedback_count = exp1_feedback_sub %>%
  245. dplyr::group_by(FeedCode) %>%
  246. dplyr::summarise(feed_mean = round(mean(counts),3),
  247. feed_sd = round(sd(counts),3),
  248. feed_se=feed_sd/sqrt(32))
  249. exp1_feedback_count
  250. ```
  253. ### Exp1 efficacy analysis
  254. ```{r exp1 rate_score and feedback}
  255. exp1_rate.feedback_lmm = lmerTest::lmer(rate_score ~ age+gender+FeedCode+
  256. (1+FeedCode|SubID),
  257. exp1_interval_efficacy, REML=FALSE)
  258. exp1_rate.feedback_anova = anova(exp1_rate.feedback_lmm)
  259. exp1_rate.feedback_anova
  260. summary(exp1_rate.feedback_lmm)
  263. exp1_rate.feedback_std_ci = model_parameters(
  264. exp1_rate.feedback_lmm,
  265. standardize = "refit",
  266. ci_method = "wald",
  267. ci = 0.95,
  268. effects = "all",
  269. iterations = 1000,
  270. summary = getOption("parameters_mixed_summary", FALSE),
  271. digits = 3)
  273. exp1_rate.feedback_std_ci
  274. ```
  276. ```{r exp1 rate_score and feedback plot}
  277. exp1_sub_rate_feedback = exp1_interval_efficacy %>%
  278. dplyr::group_by(SubID) %>%
  279. fill(rate_score, .direction = "up")%>%
  280. drop_na(rate_score)%>%
  281. dplyr::group_by(SubID, FeedCode) %>%
  282. dplyr::summarise(rate_score_mean = round(mean(rate_score),3),
  283. rate_score_sd = round(sd(rate_score),3),
  284. rate_score_se=rate_score_sd/sqrt(32))
  286. exp1_rate.feedback_plot = ggbarplot(exp1_sub_rate_feedback, x = "FeedCode",
  287. y = "rate_score_mean",
  288. alpha = 0.6, ylab= "Rate",
  289. xlab = "Social Comparison",
  290. color = "black",fill = "FeedCode",
  291. add = c("mean_se", "jitter"),
  292. add.params = list(color = "FeedCode"),
  293. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  294. legend = 'right', position = position_dodge(width =0.4))+
  295. scale_y_continuous(breaks = seq(0, 100, by = 20), expand = c(0, 0)) +
  296. coord_cartesian(ylim = c(0, 105))+
  297. theme(axis.text = element_text(size = 16, family = "serif"),
  298. axis.title = element_text(size = 20, family = "serif"))
  300. exp1_rate.feedback_plot = ggpar(exp1_rate.feedback_plot, legend = 'right') +
  301. theme(legend.text = element_text(size = 12, family = "serif"),
  302. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  304. exp1_rate.feedback_plot
  306. ```
  310. ```{r exp1 rate_score and efficacy}
  311. exp1_efficacy.rate_lmm = lmerTest::lmer(rate_score ~ age+gender+ mbased_efficacy+
  312. (1+mbased_efficacy|SubID),
  313. exp1_interval_efficacy, REML=FALSE)
  314. exp1_efficacy.rate_anova = anova(exp1_efficacy.rate_lmm)
  315. exp1_efficacy.rate_anova
  316. summary(exp1_efficacy.rate_lmm)
  318. exp1_efficacy.rate_std_ci = model_parameters(
  319. exp1_efficacy.rate_lmm,
  320. standardize = "refit",
  321. ci_method = "wald",
  322. ci = 0.95,
  323. effects = "all",
  324. iterations = 1000,
  325. summary = getOption("parameters_mixed_summary", FALSE),
  326. digits = 3)
  328. exp1_efficacy.rate_std_ci
  330. ```
  333. ```{r exp1 rate_score and efficacy plot}
  334. eff_df0 <- Effect(c("mbased_efficacy"), exp1_efficacy.rate_lmm,
  335. xlevels = list(mbased_efficacy = seq(min(exp1_interval_efficacy$mbased_efficacy, na.rm = TRUE),
  336. max(exp1_interval_efficacy$mbased_efficacy, na.rm = TRUE),
  337. by = 0.1)))
  339. exp1.rate.efficacy = as.data.frame(eff_df0)
  340. head(exp1.rate.efficacy)
  342. len1 = length(exp1_interval_efficacy$mbased_efficacy)
  344. exp1.plot.rate.efficacy <- ggplot(exp1.rate.efficacy,
  345. aes(x = mbased_efficacy, y = fit)) +
  346. geom_line(color = "#F8766D", linewidth = 1, linetype = 1) +
  347. geom_ribbon(aes(ymin = fit - se, ymax = fit + se),
  348. fill = "#F8766D", alpha = 0.1) +
  349. gg.side +
  350. scale_x_continuous(limits = c(0, 1))
  352. exp1.plot.rate.efficacy
  354. ```
  356. ```{r exp1 rate_score and efficacy plot}
  357. exp1_sub_rate_efficacy = exp1_interval_learn %>%
  358. group_by(RoundID) %>%
  359. dplyr::summarise(
  360. efficacy_mean = round(mean(mbased_efficacy, na.rm = TRUE), 3),
  361. efficacy_sd = round(sd(mbased_efficacy, na.rm = TRUE), 3),
  362. efficacy_se = efficacy_sd / sqrt(32),
  363. Rate_mean = round(mean(EfficacyProbeRespLin, na.rm = TRUE),3),
  364. Rate_sd = round(sd(EfficacyProbeRespLin, na.rm = TRUE), 3),
  365. Rate_se = Rate_sd / sqrt(32))%>%
  366. drop_na(Rate_mean)
  368. exp1_rate.efficacy_plot = ggplot(data = exp1_sub_rate_efficacy, aes(x = RoundID)) +
  369. # Plot the "rate" line and confidence interval
  370. geom_line(aes(y = Rate_mean, color = "rate"), size = 1) +
  371. geom_ribbon(aes(ymin =Rate_mean - Rate_se,
  372. ymax = Rate_mean + Rate_se, fill = "rate"),
  373. alpha = 0.2) +
  374. # Plot the "CRPS" line and confidence interval
  375. geom_line(aes(y = efficacy_mean, color = "efficacy"), size = 1) +
  376. geom_ribbon(aes(ymin = efficacy_mean - efficacy_se,
  377. ymax = efficacy_mean + efficacy_se, fill = "efficacy"),
  378. alpha = 0.2) +
  379. # Labels and scales
  380. labs(x = "RoundID", y = "Value", color = "Legend", fill = "Legend") +
  381. scale_x_continuous(breaks = seq(0, max(exp1_sub_rate_efficacy$RoundID), by = 24)) +
  382. scale_color_manual(values = c("rate" = "#808080", "efficacy" = "#749857")) + # Specify colors
  383. scale_fill_manual(values = c("rate" = "#A9A9A9", "efficacy" = "#749857")) + # Specify fill colors
  384. ylim(0.25, 0.75) +
  385. theme_minimal() +
  386. theme(legend.position = "top")+gg.side # Place the legend at the top
  388. exp1_rate.efficacy_plot
  390. ```
  395. ```{r exp1 efficacy and social comparison}
  396. exp1_efficacy.feedback_lmm = lmerTest::lmer(mbased_efficacy ~age+ gender+FeedCode+
  397. (1+FeedCode|SubID),exp1_interval_efficacy, REML=FALSE)
  398. exp1_efficacy.feedback_anova = anova(exp1_efficacy.feedback_lmm)
  399. exp1_efficacy.feedback_anova
  400. summary(exp1_efficacy.feedback_lmm)
  402. exp1_efficacy.feedback_std_ci = model_parameters(
  403. exp1_efficacy.feedback_lmm,
  404. standardize = "refit",
  405. ci_method = "wald",
  406. ci = 0.95,
  407. effects = "all",
  408. iterations = 1000,
  409. summary = getOption("parameters_mixed_summary", FALSE),
  410. digits = 3)
  412. exp1_efficacy.feedback_std_ci
  413. ```
  415. ```{r exp1 efficacy and social comparison plot}
  416. exp1_sub_efficacy_feedback = exp1_interval_efficacy %>%
  417. dplyr::filter(FeedCode == "Downward"|FeedCode == "Upward") %>%
  418. dplyr::group_by(SubID, FeedCode) %>%
  419. dplyr::summarise(mbased_efficacy_mean = round(mean(mbased_efficacy),3),
  420. mbased_efficacy_sd =round(sd(mbased_efficacy),3),
  421. mbased_efficacy_se=mbased_efficacy_sd/sqrt(32))
  423. exp1_efficacy.feedback_plot= ggbarplot(exp1_sub_efficacy_feedback, x = "FeedCode",
  424. y = "mbased_efficacy_mean",
  425. alpha = 0.6,ylab= "Efficacy",
  426. xlab = "Social Comparison",
  427. color = "black",fill = "FeedCode",
  428. add = c("mean_se", "jitter"),
  429. add.params = list(color = "FeedCode"),
  430. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  431. legend = 'right',
  432. position = position_dodge(width =0.4))+
  433. scale_y_continuous(breaks = seq(0, 1, by = 0.2), expand = c(0, 0)) +
  434. coord_cartesian(ylim = c(0, 1.05))+
  435. theme(axis.text = element_text(size = 16, family = "serif"),
  436. axis.title = element_text(size = 20, family = "serif"))
  438. exp1_efficacy.feedback_plot = ggpar(exp1_efficacy.feedback_plot, legend = 'right') +
  439. theme(legend.text = element_text(size = 12, family = "serif"),
  440. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  441. exp1_efficacy.feedback_plot
  443. ```
  445. ## Exp1 behavior analysis
  446. ### Exp1 CRPS analysis
  447. ```{r exp1 CRPS and social comparison}
  448. exp1_feedback.CRPS = exp1_interval_CRPS %>%
  449. dplyr::filter(FeedPrev == "1"|FeedPrev == "2") %>%
  450. dplyr::mutate(FeedPrev = factor(recode(FeedPrev,"1" = "Downward", "2" = "Upward"),
  451. levels = c("Downward","Upward")))
  453. contrasts(exp1_feedback.CRPS$FeedPrev) = contr.sum
  454. contrasts(exp1_feedback.CRPS$gender) = contr.sum
  456. exp1_feedback.CRPS_lmm = lmerTest::lmer(CRPS ~ age+gender+FeedPrev+meanCongruency+
  457. (1+FeedPrev+meanCongruency|SubID),exp1_feedback.CRPS, REML=FALSE)
  458. exp1_feedback.CRPS_anova = anova(exp1_feedback.CRPS_lmm)
  459. exp1_feedback.CRPS_anova
  460. summary(exp1_feedback.CRPS_lmm)
  462. exp1_feedback.CRPS_std_ci = model_parameters(
  463. exp1_feedback.CRPS_lmm,
  464. standardize = "refit",
  465. ci_method = "wald",
  466. ci = 0.95,
  467. effects = "all",
  468. iterations = 1000,
  469. summary = getOption("parameters_mixed_summary", FALSE),
  470. digits = 3)
  472. exp1_feedback.CRPS_std_ci
  473. ```
  475. ```{r exp1 CRPS and social comparison plot}
  476. exp1_sub_CRPS_feedback = exp1_feedback.CRPS %>%
  477. dplyr::group_by(SubID,FeedPrev) %>%
  478. dplyr::summarise(CRPS_mean = round(mean(CRPS),3),
  479. CRPS_sd = round(sd(CRPS),3),
  480. CRPS_se=CRPS_sd/sqrt(32))
  482. exp1_CRPS.feedback_plot = ggbarplot(exp1_sub_CRPS_feedback, x = "FeedPrev",
  483. y = "CRPS_mean", alpha = 0.6,
  484. ylab= "CRPS", xlab = "Social Comparison",
  485. color = "black",fill = "FeedPrev",
  486. add = c("mean_se", "jitter"),
  487. add.params = list(color = "FeedPrev"),
  488. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  489. legend = 'right', position = position_dodge(width =0.4)) +
  490. scale_y_continuous(breaks = seq(0, 1.25, by = 0.2), expand = c(0, 0)) +
  491. coord_cartesian(ylim = c(0, 1.25))+
  492. annotate("text",x=1.35,y=0.97,label="",size = 6) +
  493. theme(axis.text = element_text(size = 16, family = "serif"),
  494. axis.title = element_text(size = 20, family = "serif"))
  496. exp1_CRPS.feedback_plot = ggpar(exp1_CRPS.feedback_plot, legend = 'right') +
  497. theme(legend.text = element_text(size = 12, family = "serif"),
  498. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  500. exp1_CRPS.feedback_plot
  501. ```
  503. ```{r exp1 CRPS and efficacy}
  504. ###exp1 CRPS and efficacy LMM### un_PE_prev PE_prev
  505. exp1_efficacy.CRPS = exp1_interval_learn %>%
  506. group_by(SubID,FeedPrev) %>%
  507. filter(CRPS > mean(CRPS) - 3 * sd(CRPS),
  508. CRPS < mean(CRPS) + 3 * sd(CRPS))
  511. contrasts(exp1_efficacy.CRPS$gender) = contr.sum
  513. exp1_efficacy.CRPS_lmm =lmerTest::lmer(CRPS ~age+gender+mbased_efficacy_prev+meanCongruency+
  514. (1+mbased_efficacy_prev+meanCongruency|SubID),
  515. exp1_efficacy.CRPS, REML=FALSE)
  516. exp1_efficacy.CRPS_anova = anova(exp1_efficacy.CRPS_lmm)
  517. summary(exp1_efficacy.CRPS_lmm)
  518. exp1_efficacy.CRPS_anova
  521. exp1_efficacy.CRPS_std_ci = model_parameters(
  522. exp1_efficacy.CRPS_lmm,
  523. standardize = "refit",
  524. ci_method = "wald",
  525. ci = 0.95,
  526. effects = "all",
  527. iterations = 1000,
  528. summary = getOption("parameters_mixed_summary", FALSE),
  529. digits = 3)
  531. exp1_efficacy.CRPS_std_ci
  532. ```
  537. ```{r CRPS and efficacy plot}
  538. eff_df1 = Effect(c("mbased_efficacy_prev"), exp1_efficacy.CRPS_lmm,
  539. xlevels = list(mbased_efficacy_prev = seq(min(exp1_efficacy.CRPS$mbased_efficacy_prev, na.rm = TRUE),
  540. max(exp1_efficacy.CRPS$mbased_efficacy_prev, na.rm = TRUE),
  541. by = 0.1)))
  543. exp1.CRPS.efficacy = as.data.frame(eff_df1)
  544. head(exp1.CRPS.efficacy)
  546. len1 = length(exp1.CRPS.efficacy$mbased_efficacy_prev)
  548. exp1.plot.CRPS.efficacy = ggplot(exp1.CRPS.efficacy,
  549. aes(x = mbased_efficacy_prev, y = fit)) +
  550. geom_line(color = "#F8766D", linewidth = 1, linetype = 1) +
  551. geom_ribbon(aes(ymin = fit - se, ymax = fit + se),
  552. fill = "#F8766D", alpha = 0.1) +
  553. gg.side +
  554. scale_x_continuous(limits = c(0, 1)) +
  555. scale_y_continuous(
  556. limits = c(0.92, 1.01),
  557. breaks = seq(0, 2, by = 0.02) # y轴每0.2一个刻度
  558. )
  560. exp1.plot.CRPS.efficacy
  562. ```
  565. ### Exp1 RT analysis
  566. ```{r exp1 RT and social comparison}
  567. ###exp1 RT LMM###
  568. exp1_feedback.RT = exp1_interval_RT %>%
  569. filter(FeedPrev == "1"|FeedPrev == "2") %>%
  570. dplyr::mutate(FeedPrev = factor(recode(FeedPrev,"1" = "Downward", "2" = "Upward"),
  571. levels = c("Downward","Upward")))
  573. contrasts(exp1_feedback.RT$FeedPrev)<-contr.sum
  574. contrasts(exp1_feedback.RT$gender)<-contr.sum
  576. exp1_feedback.RT_lmm = lmerTest::lmer(meanRT~age+gender+FeedPrev+meanCongruency+
  577. (1+FeedPrev+meanCongruency|SubID),exp1_feedback.RT, REML=FALSE)
  578. exp1_feedback.RT_anova = anova(exp1_feedback.RT_lmm)
  579. exp1_feedback.RT_anova
  580. summary(exp1_feedback.RT_lmm)
  582. exp1_feedback.RT_std_ci = model_parameters(
  583. exp1_feedback.RT_lmm,
  584. standardize = "refit",
  585. ci_method = "wald",
  586. ci = 0.95,
  587. effects = "all",
  588. iterations = 1000,
  589. summary = getOption("parameters_mixed_summary", FALSE),
  590. digits = 3)
  592. exp1_feedback.RT_std_ci
  593. ```
  598. ```{r exp1 RT and social comparison plot}
  599. exp1_sub_RT_feedback = exp1_interval_RT %>%
  600. filter(FeedPrev == "1"|FeedPrev == "2") %>%
  601. mutate(FeedPrev = factor(recode(FeedPrev, "1" = "Downward", "2" = "Upward"),
  602. levels = c("Downward", "Upward"))) %>%
  603. dplyr::group_by(SubID, FeedPrev) %>%
  604. dplyr::summarise(RT_mean = round(mean(meanRT),3),
  605. RT_sd = round(sd(meanRT),3),
  606. RT_se=RT_sd/sqrt(32))
  608. exp1_RT.feedback_plot= ggbarplot(exp1_sub_RT_feedback, x = "FeedPrev",
  609. y = "RT_mean", alpha = 0.6,
  610. ylab= "Reaction time (ms)",
  611. xlab = "Social Comparison",
  612. color = "black",fill = "FeedPrev",
  613. add = c("mean_se", "jitter"),
  614. add.params = list(color = "FeedPrev"),
  615. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  616. legend = 'right', position = position_dodge(width =0.4)) +
  617. scale_y_continuous(breaks = c(0, seq(0, 1150, by = 200)), expand = c(0, 0))+
  618. coord_cartesian(ylim = c(0, 1150)) +
  619. annotate("text",x=1.35,y=710,label="",size = 6) +
  620. theme(axis.text = element_text(size = 16, family = "serif"),
  621. axis.title = element_text(size = 20, family = "serif"))
  623. exp1_RT.feedback_plot = ggpar(exp1_RT.feedback_plot, legend = 'right') +
  624. theme(legend.text = element_text(size = 12, family = "serif"),
  625. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  627. exp1_RT.feedback_plot
  629. ```
  631. ```{r exp1 RT and efficacy}
  632. ###exp1 CRPS and efficacy LMM### un_PE_prev PE_prev
  633. exp1_efficacy.RT = exp1_interval_learn %>%
  634. dplyr::group_by(SubID,FeedPrev) %>%
  635. filter(meanRT > mean(meanRT) - 3 * sd(meanRT),
  636. meanRT < mean(meanRT) + 3 * sd(meanRT))
  638. contrasts(exp1_efficacy.RT$gender) = contr.sum
  640. exp1_efficacy.RT_lmm = lmerTest::lmer(meanRT ~ age+gender+mbased_efficacy_prev+meanCongruency+
  641. (1+mbased_efficacy_prev+meanCongruency|SubID),
  642. data = exp1_efficacy.RT, REML=FALSE)
  643. exp1_efficacy.RT_anova = anova(exp1_efficacy.RT_lmm)
  644. exp1_efficacy.RT_anova
  645. summary(exp1_efficacy.RT_lmm)
  647. exp1_efficacy.RT_std_ci = model_parameters(
  648. exp1_efficacy.RT_lmm,
  649. standardize = "refit",
  650. df_method = "satterthwaite",
  651. ci_method = "wald",
  652. ci = 0.95,
  653. effects = "all",
  654. iterations = 1000,
  655. summary = getOption("parameters_mixed_summary", FALSE),
  656. digits = 3)
  658. exp1_efficacy.RT_std_ci
  659. ```
  661. ```{r exp1 RT and efficacy plot}
  662. eff_df2 = Effect(c("mbased_efficacy_prev"),exp1_efficacy.RT_lmm,
  663. xlevels = list(mbased_efficacy_prev = seq(min(exp1_efficacy.RT$mbased_efficacy_prev, na.rm = TRUE),
  664. max(exp1_efficacy.RT$mbased_efficacy_prev, na.rm = TRUE),
  665. by = 0.1)))
  667. exp1.RT.efficacy = as.data.frame(eff_df2)
  668. head(exp1.RT.efficacy)
  669. len2<-length(exp1.RT.efficacy$mbased_efficacy_prev)
  671. exp1.plot.RT.efficacy <- ggplot()+
  672. geom_line(data=exp1.RT.efficacy,
  673. aes(x=mbased_efficacy_prev, y=fit),
  674. color="#F8766D",size=1,linetype=1)+
  675. geom_ribbon(data=exp1.RT.efficacy,
  676. aes(x=mbased_efficacy_prev, max = fit + se, min = fit- se),
  677. fill = "#F8766D",
  678. alpha=0.1,
  679. inherit.aes = FALSE)+ylim(c(600, 740))+ scale_x_continuous(limits = c(0, 1)) +gg.side
  681. exp1.plot.RT.efficacy
  683. ```
  687. ###Exp1 LMM power
  688. ```{r exp1 efficacy and behavior power}
  689. ###rate & feedack###
  690. exp1_rate.feedback_lmm = lmerTest::lmer(rate_score ~ age+gender+FeedCode+
  691. (1+FeedCode|SubID),
  692. exp1_interval_efficacy, REML=FALSE)
  693. exp1_rate.feedback_anova = anova(exp1_rate.feedback_lmm)
  694. exp1_rate.feedback_anova
  695. summary(exp1_rate.feedback_lmm)
  697. model_rate_fb_exp1=powerSim(exp1_rate.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  700. ###efficacy & rate###
  701. exp1_efficacy.rate_lmm = lmerTest::lmer(rate_score ~ age+gender+ mbased_efficacy+
  702. (1+mbased_efficacy|SubID),
  703. exp1_interval_efficacy, REML=FALSE)
  704. exp1_efficacy.rate_anova = anova(exp1_efficacy.rate_lmm)
  705. exp1_efficacy.rate_anova
  706. summary(exp1_efficacy.rate_lmm)
  709. model_rate_eff_exp1=powerSim(exp1_efficacy.rate_lmm,fixed("mbased_efficacy", "t"),nsim=1000)
  711. ###efficacy & feedback###
  712. exp1_efficacy.feedback_lmm = lmerTest::lmer(mbased_efficacy ~ age+gender+FeedCode+
  713. (1+FeedCode|SubID),exp1_interval_efficacy, REML=FALSE)
  714. exp1_efficacy.feedback_anova = anova(exp1_efficacy.feedback_lmm)
  715. exp1_efficacy.feedback_anova
  716. summary(exp1_efficacy.feedback_lmm)
  718. model_eff_fb_exp1=powerSim(exp1_efficacy.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  720. # ###CRPS & feedback###
  721. exp1_feedback.CRPS_lmm = lmerTest::lmer(CRPS ~ age+gender+FeedPrev+meanCongruency+
  722. (1+FeedPrev+meanCongruency|SubID),exp1_feedback.CRPS, REML=FALSE)
  723. exp1_feedback.CRPS_anova = anova(exp1_feedback.CRPS_lmm)
  724. exp1_feedback.CRPS_anova
  725. summary(exp1_feedback.CRPS_lmm)
  727. model_CRPS_fb_exp1=powerSim(exp1_feedback.CRPS_lmm,fixed("FeedPrev", "f"),nsim=1000)
  729. ###CRPS & efficacy###
  730. exp1_efficacy.CRPS_lmm =lmerTest::lmer(CRPS ~ age+gender+mbased_efficacy_prev+meanCongruency+
  731. (1+mbased_efficacy_prev+meanCongruency|SubID),
  732. exp1_efficacy.CRPS, REML=FALSE)
  733. exp1_efficacy.CRPS_anova = anova(exp1_efficacy.CRPS_lmm)
  734. exp1_efficacy.CRPS_anova
  735. summary(exp1_efficacy.CRPS_lmm)
  737. model_CRPS_eff_exp1=powerSim(exp1_efficacy.CRPS_lmm,
  738. fixed("mbased_efficacy_prev", "t"),nsim=1000)
  740. # ###RT & feedback###
  741. exp1_feedback.RT_lmm = lmerTest::lmer(meanRT~age+gender+FeedPrev+meanCongruency+
  742. (1+FeedPrev+meanCongruency|SubID),
  743. exp1_feedback.RT, REML=FALSE)
  744. exp1_feedback.RT_anova = anova(exp1_feedback.RT_lmm)
  745. exp1_feedback.RT_anova
  746. summary(exp1_feedback.RT_lmm)
  748. model_RT_fb_exp1=powerSim(exp1_feedback.RT_lmm,fixed("FeedPrev", "f"),nsim=1000)
  750. ###RT & efficacy###
  751. exp1_efficacy.RT_lmm = lmerTest::lmer(meanRT ~ age+gender+mbased_efficacy_prev+meanCongruency+
  752. (1+mbased_efficacy_prev+meanCongruency|SubID),
  753. data = exp1_efficacy.RT, REML=FALSE)
  754. exp1_efficacy.RT_anova = anova(exp1_efficacy.RT_lmm)
  755. exp1_efficacy.RT_anova
  756. summary(exp1_efficacy.RT_lmm)
  758. model_RT_eff_exp1=powerSim(exp1_efficacy.RT_lmm,
  759. fixed("mbased_efficacy_prev", "t"),nsim=1000)
  760. ```
  762. ## Exp2 model analysis
  763. ### Exp2 descriptive statistics
  764. ```{r exp2 rate_score and social comparison descriptive}
  765. exp2_interval_efficacy = exp2_interval_learn %>%
  766. group_by(SubID) %>%
  767. mutate(FeedCode = factor(FeedCode, levels = c("1", "2"),
  768. labels=c("Downward","Upward")))
  770. contrasts(exp2_interval_efficacy$FeedCode) = contr.sum
  771. contrasts(exp2_interval_efficacy$gender) = contr.sum
  773. ###exp2 rate_score###
  774. exp2_rate_feedback = exp2_interval_learn %>%
  775. dplyr::group_by(SubID) %>%
  776. fill(rate_score, .direction = "up")%>%
  777. drop_na(rate_score)%>%
  778. dplyr::group_by(FeedCode) %>%
  779. dplyr::summarise(rate_score_mean = round(mean(rate_score),3),
  780. rate_score_sd = round(sd(rate_score),3),
  781. rate_score_se=rate_score_sd/sqrt(34))
  783. exp2_rate_feedback
  785. ###exp2 model_based efficacy###
  786. exp2_efficacy_feedback = exp2_interval_learn %>%
  787. dplyr::group_by(SubID) %>%
  788. fill(mbased_efficacy, .direction = "up")%>%
  789. drop_na(mbased_efficacy)%>%
  790. dplyr::group_by(FeedCode) %>%
  791. dplyr::summarise(mbased_efficacy_mean = round(mean(mbased_efficacy),3),
  792. mbased_efficacy_sd = round(sd(mbased_efficacy),3),
  793. mbased_efficacy_se=mbased_efficacy_sd/sqrt(34))
  795. exp2_efficacy_feedback
  797. ```
  800. ```{r exp2 CRPS and RT descriptive}
  801. ###exp2 CRPS ###
  802. exp2_interval_CRPS = exp2_interval_learn %>%
  803. dplyr::group_by(SubID, FeedPrev) %>%
  804. filter(CRPS > mean(CRPS) - 3 * sd(CRPS),
  805. CRPS < mean(CRPS) + 3 * sd(CRPS))
  807. exp2_sub_CRPS_feedback = exp2_interval %>%
  808. dplyr::group_by(SubID, FeedPrev) %>%
  809. filter(CRPS > mean(CRPS) - 3 * sd(CRPS),
  810. CRPS < mean(CRPS) + 3 * sd(CRPS)) %>%
  811. dplyr::group_by(SubID, FeedPrev) %>%
  812. dplyr::summarise(CRPS = round(mean(CRPS),3))
  814. exp2_CRPS_feedback = exp2_sub_CRPS_feedback %>%
  815. dplyr::group_by(FeedPrev) %>%
  816. dplyr::summarise(CRPS_mean = round(mean(CRPS),3),
  817. CRPS_sd = round(sd(CRPS),3),
  818. CRPS_se=CRPS_sd/sqrt(34))
  819. exp2_CRPS_feedback
  821. ###exp2 RT ###
  822. exp2_interval_RT = exp2_interval_learn %>%
  823. dplyr::group_by(SubID, FeedPrev) %>%
  824. filter(meanRT > mean(meanRT) - 3 * sd(meanRT),
  825. meanRT < mean(meanRT) + 3 * sd(meanRT))
  827. exp2_sub_RT_feedback=exp2_interval_learn %>%
  828. dplyr::group_by(SubID, FeedPrev) %>%
  829. filter(meanRT > mean(meanRT) - 3 * sd(meanRT),
  830. meanRT < mean(meanRT) + 3 * sd(meanRT))%>%
  831. dplyr::group_by(SubID, FeedPrev) %>%
  832. dplyr::summarise(RT = round(mean(meanRT)))
  834. exp2_RT_feedback = exp2_sub_RT_feedback %>%
  835. dplyr::group_by(FeedPrev) %>%
  836. dplyr::summarise(RT_mean = round(mean(RT)),
  837. RT_sd = round(sd(RT)),
  838. RT_se=RT_sd/sqrt(34))
  839. exp2_RT_feedback
  841. ###exp2 ACC ###
  842. exp2_ACC_feedback = exp2_interval %>%
  843. drop_na(FeedPrev) %>%
  844. dplyr::group_by(FeedPrev) %>%
  845. dplyr::summarise(ACC_mean = round(mean(ACC),3),
  846. ACC_sd = round(sd(ACC),3),
  847. ACC_se=ACC_sd/sqrt(32))
  848. exp2_ACC_feedback
  850. ###exp2 feedback count###
  851. exp2_feedback_sub = exp2_interval %>%
  852. group_by(SubID,FeedCode) %>%
  853. dplyr::summarise(counts = n()) %>%
  854. ungroup()
  856. exp2_feedback_count = exp2_feedback_sub %>%
  857. dplyr::group_by(FeedCode) %>%
  858. dplyr::summarise(feed_mean = round(mean(counts),3),
  859. feed_sd = round(sd(counts),3),
  860. feed_se=feed_sd/sqrt(32))
  861. exp2_feedback_count
  862. ```
  864. ###Exp2 efficacy analysis
  865. ```{r exp2 rate_score and social comparison}
  866. exp2_rate.feedback_lmm = lmerTest::lmer(rate_score ~ age+gender+FeedCode+
  867. (1+FeedCode|SubID),
  868. data = exp2_interval_efficacy, REML=FALSE)
  869. exp2_rate.feedback_anova = anova(exp2_rate.feedback_lmm)
  870. exp2_rate.feedback_anova
  871. summary(exp2_rate.feedback_lmm)
  873. exp2_rate.feedback_std_ci = model_parameters(
  874. exp2_rate.feedback_lmm,
  875. standardize = "refit",
  876. ci_method = "wald",
  877. ci = 0.95,
  878. effects = "all",
  879. iterations = 1000,
  880. summary = getOption("parameters_mixed_summary", FALSE),
  881. digits = 3)
  883. exp2_rate.feedback_std_ci
  884. ```
  886. ```{r exp2 rate_score and efficacy}
  887. exp2_efficacy.rate_lmm= lmerTest::lmer(rate_score ~age+gender+mbased_efficacy+
  888. (1+mbased_efficacy|SubID),
  889. data = exp2_interval_learn, REML=FALSE)
  890. exp2_efficacy.rate_anova = anova(exp2_efficacy.rate_lmm)
  891. exp2_efficacy.rate_anova
  892. summary(exp2_efficacy.rate_lmm)
  895. exp2_efficacy.rate_std_ci = model_parameters(
  896. exp2_efficacy.rate_lmm,
  897. standardize = "refit",
  898. ci_method = "wald",
  899. ci = 0.95,
  900. effects = "all",
  901. iterations = 1000,
  902. summary = getOption("parameters_mixed_summary", FALSE),
  903. digits = 3)
  905. exp2_efficacy.rate_std_ci
  906. ```
  909. ```{r rate_score and efficacy plot}
  910. exp2_sub_rate_efficacy = exp2_interval_efficacy %>%
  911. group_by(RoundID) %>%
  912. dplyr::summarise(
  913. efficacy_mean = round(mean(mbased_efficacy, na.rm = TRUE), 3),
  914. efficacy_sd = round(sd(mbased_efficacy, na.rm = TRUE), 3),
  915. efficacy_se = efficacy_sd / sqrt(32),
  916. Rate_mean = round(mean(EfficacyProbeRespLin, na.rm = TRUE),3),
  917. Rate_sd = round(sd(EfficacyProbeRespLin, na.rm = TRUE), 3),
  918. Rate_se = Rate_sd / sqrt(34))%>%
  919. drop_na(Rate_mean)
  921. exp2_rate.efficacy_plot = ggplot(data = exp2_sub_rate_efficacy, aes(x = RoundID)) +
  922. # Plot the "rate" line and confidence interval
  923. geom_line(aes(y = Rate_mean, color = "rate"), size = 1) +
  924. geom_ribbon(aes(ymin =Rate_mean - Rate_se,
  925. ymax = Rate_mean + Rate_se, fill = "rate"),
  926. alpha = 0.2) +
  927. # Plot the "CRPS" line and confidence interval
  928. geom_line(aes(y = efficacy_mean, color = "efficacy"), size = 1) +
  929. geom_ribbon(aes(ymin = efficacy_mean - efficacy_se,
  930. ymax = efficacy_mean + efficacy_se, fill = "efficacy"),
  931. alpha = 0.2) +
  932. # Labels and scales
  933. labs(x = "RoundID", y = "Value", color = "Legend", fill = "Legend") +
  934. scale_x_continuous(breaks = seq(0, max(exp2_sub_rate_efficacy$RoundID), by = 24)) +
  935. scale_color_manual(values = c("rate" = "#808080", "efficacy" = "#749857")) + # Specify colors
  936. scale_fill_manual(values = c("rate" = "#A9A9A9", "efficacy" = "#749857")) + # Specify fill colors
  937. ylim(0.25, 0.75) +
  938. theme_minimal() +
  939. theme(legend.position = "top")+gg.side # Place the legend at the top
  941. exp2_rate.efficacy_plot
  942. ```
  945. ```{r exp2 efficacy and social comparison}
  946. exp2_efficacy.feedback_lmm = lmerTest::lmer(mbased_efficacy ~ age+gender+FeedCode+
  947. (1+FeedCode|SubID),
  948. data = exp2_interval_efficacy, REML=FALSE)
  949. exp2_efficacy.feedback_anova = anova(exp2_efficacy.feedback_lmm)
  950. exp2_efficacy.feedback_anova
  951. summary(exp2_efficacy.feedback_lmm)
  953. exp2_efficacy.feedback_std_ci = model_parameters(
  954. exp2_efficacy.feedback_lmm,
  955. standardize = "refit",
  956. ci_method = "wald",
  957. ci = 0.95,
  958. effects = "all",
  959. iterations = 1000,
  960. summary = getOption("parameters_mixed_summary", FALSE),
  961. digits = 3)
  963. exp2_efficacy.feedback_std_ci
  964. ```
  967. ##Exp2 behavior analysis
  968. ### Exp2 CRPS analysis
  969. ```{r exp2 CRPS and social comparison}
  970. exp2_feedback.CRPS = exp2_interval_CRPS %>%
  971. group_by(SubID) %>%
  972. dplyr::filter(FeedPrev == "1"|FeedPrev == "2") %>%
  973. dplyr::mutate(FeedPrev = factor(recode(FeedPrev,"1" = "Downward", "2" = "Upward"),
  974. levels = c("Downward","Upward")))
  976. contrasts(exp2_feedback.CRPS$FeedPrev) = contr.sum
  977. contrasts(exp2_feedback.CRPS$gender) = contr.sum
  979. exp2_feedback.CRPS_lmm = lmerTest::lmer(CRPS~age+gender+FeedPrev+meanCongruency+
  980. (1+FeedPrev+meanCongruency|SubID), exp2_feedback.CRPS, REML=FALSE)
  981. exp2_feedback.CRPS_anova = anova(exp2_feedback.CRPS_lmm)
  982. exp2_feedback.CRPS_anova
  983. summary(exp2_feedback.CRPS_lmm)
  985. exp2_feedback.CRPS_std_ci = model_parameters(
  986. exp2_feedback.CRPS_lmm,
  987. standardize = "refit",
  988. ci_method = "profile",
  989. ci = 0.95,
  990. digits = 3)
  992. exp2_feedback.CRPS_std_ci
  993. ```
  996. ```{r exp2 CRPS and efficacy}
  997. exp2_efficacy.CRPS = exp2_interval_learn %>%
  998. group_by(SubID) %>%
  999. dplyr::group_by(SubID,FeedPrev) %>%
  1000. filter(CRPS > mean(CRPS) - 3 * sd(CRPS),
  1001. CRPS < mean(CRPS) + 3 * sd(CRPS))
  1004. contrasts(exp2_efficacy.CRPS$gender)<-contr.sum
  1006. exp2_efficacy.CRPS_lmm =lmerTest::lmer(CRPS ~ age+gender+mbased_efficacy_prev+meanCongruency+
  1007. (1+mbased_efficacy_prev+meanCongruency|SubID),
  1008. data=exp2_efficacy.CRPS, REML=FALSE)
  1010. summary(exp2_efficacy.CRPS_lmm)
  1011. exp2_efficacy.CRPS_anova = anova(exp2_efficacy.CRPS_lmm)
  1012. exp2_efficacy.CRPS_anova
  1015. exp2_efficacy.CRPS_std_ci = model_parameters(
  1016. exp2_efficacy.CRPS_lmm,
  1017. standardize = "refit",
  1018. ci_method = "profile",
  1019. ci = 0.95,
  1020. digits = 3)
  1022. exp2_efficacy.CRPS_std_ci
  1023. ```
  1027. ### Exp2 RT analysis
  1028. ```{r exp2 RT and social comparison}
  1029. exp2_feedback.RT = exp2_interval_RT %>%
  1030. group_by(SubID) %>%
  1031. filter(FeedPrev == "1"|FeedPrev == "2") %>%
  1032. dplyr::mutate(FeedPrev = factor(recode(FeedPrev,"1" = "Downward", "2" = "Upward"),
  1033. levels = c("Downward","Upward")))
  1035. contrasts(exp2_feedback.RT$FeedPrev) = contr.sum
  1036. contrasts(exp2_feedback.RT$gender) = contr.sum
  1038. exp2_feedback.RT_lmm = lmerTest::lmer(meanRT~age+gender+FeedPrev+meanCongruency+
  1039. (1+FeedPrev+meanCongruency|SubID), exp2_feedback.RT, REML=FALSE)
  1040. exp2_feedback.RT_anova = anova(exp2_feedback.RT_lmm)
  1041. exp2_feedback.RT_anova
  1042. summary(exp2_feedback.RT_lmm)
  1044. exp2_feedback.RT_std_ci = model_parameters(
  1045. exp2_feedback.RT_lmm,
  1046. standardize = "refit",
  1047. df_method = "satterthwaite",
  1048. ci_method = "wald",
  1049. ci = 0.95,
  1050. effects = "all",
  1051. iterations = 1000,
  1052. summary = getOption("parameters_mixed_summary", FALSE),
  1053. digits = 3)
  1055. exp2_feedback.RT_std_ci
  1056. ```
  1058. ```{r exp2 RT and efficacy}
  1059. exp2_efficacy.RT = exp2_interval_learn %>%
  1060. group_by(SubID) %>%
  1061. dplyr::group_by(SubID,FeedPrev) %>%
  1062. filter(meanRT > mean(meanRT) - 3 * sd(meanRT),
  1063. meanRT < mean(meanRT) + 3 * sd(meanRT))
  1066. contrasts(exp2_efficacy.RT$gender) = contr.sum
  1068. exp2_efficacy.RT_lmm = lmerTest::lmer(meanRT ~ age+gender+mbased_efficacy_prev+meanCongruency+
  1069. (1+mbased_efficacy_prev+meanCongruency|SubID),
  1070. data=exp2_efficacy.RT, REML=FALSE)
  1072. exp2_efficacy.RT_anova = anova(exp2_efficacy.RT_lmm)
  1073. exp2_efficacy.RT_anova
  1074. summary(exp2_efficacy.RT_lmm)
  1076. exp2_efficacy.RT_std_ci = model_parameters(
  1077. exp2_efficacy.RT_lmm,
  1078. standardize = "refit",
  1079. df_method = "satterthwaite",
  1080. ci_method = "wald",
  1081. ci = 0.95,
  1082. effects = "all",
  1083. iterations = 1000,
  1084. summary = getOption("parameters_mixed_summary", FALSE),
  1085. digits = 3)
  1087. exp2_efficacy.RT_std_ci
  1088. ```
  1091. ###Exp1 LMM power
  1092. ```{r exp2 efficacy and behavior power}
  1093. ###rate & feedack###
  1094. exp2_rate.feedback_lmm = lmerTest::lmer(rate_score ~ age+gender+FeedCode+
  1095. (1+FeedCode|SubID),
  1096. exp2_interval_efficacy, REML=FALSE)
  1097. exp2_rate.feedback_anova = anova(exp2_rate.feedback_lmm)
  1098. exp2_rate.feedback_anova
  1099. summary(exp2_rate.feedback_lmm)
  1101. model_rate_fb_exp2=powerSim(exp2_rate.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  1104. ###efficacy & rate###
  1105. exp2_efficacy.rate_lmm = lmerTest::lmer(rate_score ~ age+gender+ mbased_efficacy+
  1106. (1+mbased_efficacy|SubID),
  1107. exp2_interval_efficacy, REML=FALSE)
  1108. exp2_efficacy.rate_anova = anova(exp2_efficacy.rate_lmm)
  1109. exp2_efficacy.rate_anova
  1110. summary(exp2_efficacy.rate_lmm)
  1112. model_rate_eff_exp2=powerSim(exp2_efficacy.rate_lmm,fixed("mbased_efficacy", "t"),nsim=1000)
  1114. ###efficacy & feedback###
  1115. exp2_efficacy.feedback_lmm = lmerTest::lmer(mbased_efficacy ~ age+gender+FeedCode+
  1116. (1+FeedCode|SubID),exp2_interval_efficacy, REML=FALSE)
  1117. exp2_efficacy.feedback_anova = anova(exp2_efficacy.feedback_lmm)
  1118. exp2_efficacy.feedback_anova
  1119. summary(exp2_efficacy.feedback_lmm)
  1121. model_eff_fb_exp2=powerSim(exp2_efficacy.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  1123. # ###CRPS & feedback###
  1124. exp2_feedback.CRPS_lmm = lmerTest::lmer(CRPS ~ age+gender+FeedPrev+meanCongruency+
  1125. (1+FeedPrev+meanCongruency|SubID),exp2_feedback.CRPS, REML=FALSE)
  1126. exp2_feedback.CRPS_anova = anova(exp2_feedback.CRPS_lmm)
  1127. exp2_feedback.CRPS_anova
  1128. summary(exp2_feedback.CRPS_lmm)
  1130. model_CRPS_fb_exp2=powerSim(exp2_feedback.CRPS_lmm,fixed("FeedPrev", "f"),nsim=1000)
  1132. ###CRPS & efficacy###
  1133. exp2_efficacy.CRPS_lmm =lmerTest::lmer(CRPS ~ age+gender+mbased_efficacy_prev+meanCongruency+
  1134. (1+mbased_efficacy_prev+meanCongruency|SubID),
  1135. exp2_efficacy.CRPS, REML=FALSE)
  1136. exp2_efficacy.CRPS_anova = anova(exp2_efficacy.CRPS_lmm)
  1137. exp2_efficacy.CRPS_anova
  1138. summary(exp2_efficacy.CRPS_lmm)
  1140. model_CRPS_eff_exp2=powerSim(exp2_efficacy.CRPS_lmm,
  1141. fixed("mbased_efficacy_prev", "t"),nsim=1000)
  1143. # ###RT & feedback###
  1144. exp2_feedback.RT_lmm = lmerTest::lmer(meanRT~age+gender+FeedPrev+meanCongruency+
  1145. (1+FeedPrev+meanCongruency|SubID),
  1146. exp2_feedback.RT, REML=FALSE)
  1147. exp2_feedback.RT_anova = anova(exp2_feedback.RT_lmm)
  1148. exp2_feedback.RT_anova
  1149. summary(exp2_feedback.RT_lmm)
  1151. model_RT_fb_exp2=powerSim(exp2_feedback.RT_lmm,fixed("FeedPrev", "f"),nsim=1000)
  1153. ###RT & efficacy###
  1154. exp2_efficacy.RT_lmm = lmerTest::lmer(meanRT ~ age+gender+mbased_efficacy_prev+meanCongruency+
  1155. (1+mbased_efficacy_prev+meanCongruency|SubID),
  1156. data = exp2_efficacy.RT, REML=FALSE)
  1157. exp2_efficacy.RT_anova = anova(exp2_efficacy.RT_lmm)
  1158. exp2_efficacy.RT_anova
  1159. summary(exp2_efficacy.RT_lmm)
  1161. model_RT_eff_exp2=powerSim(exp2_efficacy.RT_lmm,
  1162. fixed("mbased_efficacy_prev", "t"),nsim=1000)
  1163. ```
  1166. # Exp2 ERP analysis
  1167. ```{r import exp2 ERPS data}
  1168. ###ERP data##
  1169. cue_ERP = read.csv(".../ERP/Cue_alltrial.csv", header = T) %>%
  1170. dplyr::select(SubID,RoundID,CNV)
  1172. feed_ERP = read.csv(".../ERP/Feed_alltrial.csv", header = T) %>%
  1173. dplyr::select(SubID,RoundID,RewP,P3,LPP)
  1175. ###ERSP data##
  1176. cue_tfa_beta=read_excel(".../Data/ERP/TFAcue_db.xlsx") %>%
  1177. dplyr::select(SubID,RoundID,cuebeta)
  1179. ###combine data###
  1180. exp2_feedERPS_combined = merge(exp2_interval_learn, feed_ERP,by =c("SubID","RoundID"), all = TRUE) %>%
  1181. arrange(SubID,RoundID)
  1183. exp2_cue_ERSP = merge(cue_ERP,cue_tfa_beta,by =c("SubID","RoundID"), all = TRUE)
  1184. exp2_cueERSP_combined = merge(exp2_interval_learn, exp2_cue_ERSP,by =c("SubID","RoundID"), all = TRUE)%>%
  1185. arrange(SubID,RoundID)
  1187. ```
  1190. ```{r import exp2 ERPS data}
  1191. exp2_feedERPS_combined$age= as.numeric(exp2_feedERPS_combined$age)
  1192. exp2_feedERPS_combined$gender= as.factor(exp2_feedERPS_combined$gender)
  1193. exp2_feedERPS_combined$RewP= as.numeric(exp2_feedERPS_combined$RewP)
  1194. exp2_feedERPS_combined$P3= as.numeric(exp2_feedERPS_combined$P3)
  1195. exp2_feedERPS_combined$LPP= as.numeric(exp2_feedERPS_combined$LPP)
  1197. exp2_cueERSP_combined$age= as.numeric(exp2_cueERSP_combined$age)
  1198. exp2_cueERSP_combined$mbased_efficacy_prev = as.numeric(exp2_cueERSP_combined$mbased_efficacy_prev)
  1199. exp2_cueERSP_combined$gender= as.factor(exp2_cueERSP_combined$gender)
  1200. exp2_cueERSP_combined$CNV= as.numeric(exp2_cueERSP_combined$CNV)
  1201. exp2_cueERSP_combined$cuebeta= as.numeric(exp2_cueERSP_combined$cuebeta)
  1203. exp2_interval_feedERPS = exp2_feedERPS_combined %>%
  1204. filter(FeedCode == "1"|FeedCode== "2") %>%
  1205. mutate(FeedCode = factor(FeedCode, levels = c("1", "2"),
  1206. labels=c("Downward","Upward")))
  1208. contrasts(exp2_interval_feedERPS$FeedCode) = contr.sum
  1209. contrasts(exp2_interval_feedERPS$gender) = contr.sum
  1211. exp2_interval_cueERPS = exp2_cueERSP_combined %>%
  1212. filter(FeedPrev == "1"|FeedPrev== "2") %>%
  1213. mutate(FeedPrev = factor(FeedPrev, levels = c("1", "2"),
  1214. labels=c("Downward","Upward")))
  1216. contrasts(exp2_interval_cueERPS$FeedPrev) = contr.sum
  1217. contrasts(exp2_interval_cueERPS$gender) = contr.sum
  1219. contrasts(exp2_cueERSP_combined$FeedPrev) = contr.sum
  1220. contrasts(exp2_cueERSP_combined$gender) = contr.sum
  1222. ```
  1223. ### Exp2 Feedback phase analysis
  1224. ```{r RewP and social comparison lmm}
  1225. exp2_RewP.feedback_lmm = lmerTest::lmer(RewP~age+gender+FeedCode+
  1226. (1+FeedCode|SubID),exp2_interval_feedERPS, REML=FALSE)
  1227. exp2_RewP.feedback_anova = anova(exp2_RewP.feedback_lmm)
  1228. exp2_RewP.feedback_anova
  1229. summary(exp2_RewP.feedback_lmm )
  1232. exp2_RewP.feedback_std_ci = model_parameters(
  1233. exp2_RewP.feedback_lmm,
  1234. standardize = "refit",
  1235. df_method = "satterthwaite",
  1236. ci_method = "wald",
  1237. ci = 0.95,
  1238. effects = "all",
  1239. iterations = 1000,
  1240. summary = getOption("parameters_mixed_summary", FALSE),
  1241. digits = 3)
  1243. exp2_RewP.feedback_std_ci
  1244. ```
  1247. ```{r RewP and social comparison des}
  1248. exp2_RewP_feedback = exp2_interval_feedERPS %>%
  1249. dplyr::group_by(SubID) %>%
  1250. drop_na(RewP)%>%
  1251. dplyr::group_by(FeedCode) %>%
  1252. dplyr::summarise(RewP_mean = round(mean(RewP),3),
  1253. RewP_sd = round(sd(RewP),3),
  1254. RewP_se=RewP_sd/sqrt(34))
  1255. exp2_RewP_feedback
  1256. ```
  1258. ```{r RewP and social comparison plot}
  1259. exp2_sub_RewP_feedback = exp2_interval_ERPS %>%
  1260. dplyr::group_by(SubID) %>%
  1261. drop_na(RewP)%>%
  1262. dplyr::group_by(SubID, FeedCode) %>%
  1263. dplyr::summarise(RewP_mean = round(mean(RewP),3),
  1264. RewP_sd = round(sd(RewP),3),
  1265. RewP_se=RewP_sd/sqrt(34))
  1268. exp2_RewP.feedback_plot= ggbarplot(exp2_sub_RewP_feedback, x = "FeedCode",
  1269. y = "RewP_mean",
  1270. alpha = 0.6, ylab= "RewP",
  1271. xlab = "Social Comparison",
  1272. color = "black",fill = "FeedCode",
  1273. add = c("mean_se","jitter"),
  1274. add.params = list(color = "FeedCode"),
  1275. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  1276. legend = 'right', position = position_dodge(width =0.4)) +
  1277. scale_y_continuous(breaks = seq(0, 25, by = 5), expand = c(0, 0)) +
  1278. coord_cartesian(ylim = c(0, 28))+
  1279. theme(axis.text = element_text(size = 16, family = "serif"),
  1280. axis.title = element_text(size = 20, family = "serif"))
  1282. exp2_RewP.feedback_plot = ggpar(exp2_RewP.feedback_plot, legend = 'right') +
  1283. theme(legend.text = element_text(size = 12, family = "serif"),
  1284. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  1286. exp2_RewP.feedback_plot
  1288. ```
  1291. ```{r P3 and social comparison lmm}
  1292. exp2_P3.feedback_lmm = lmerTest::lmer(P3~age+gender+FeedCode+
  1293. (1+FeedCode|SubID),exp2_interval_feedERPS, REML=FALSE)
  1294. exp2_P3.efficacy.fd_anova = anova(exp2_P3.feedback_lmm)
  1295. exp2_P3.efficacy.fd_anova
  1296. summary(exp2_P3.feedback_lmm)
  1298. exp2_P3.feedback_std_ci = model_parameters(
  1299. exp2_P3.feedback_lmm,
  1300. standardize = "refit",
  1301. df_method = "satterthwaite",
  1302. ci_method = "wald",
  1303. ci = 0.95,
  1304. effects = "all",
  1305. iterations = 1000,
  1306. summary = getOption("parameters_mixed_summary", FALSE),
  1307. digits = 3)
  1309. exp2_P3.feedback_std_ci
  1310. ```
  1313. ```{r P3 and social comparison des}
  1314. exp2_P3_feedback = exp2_interval_ERPS %>%
  1315. dplyr::group_by(SubID) %>%
  1316. drop_na(P3)%>%
  1317. dplyr::group_by(FeedCode) %>%
  1318. dplyr::summarise(P3_mean = round(mean(P3),3),
  1319. P3_sd = round(sd(P3),3),
  1320. P3_se=P3_sd/sqrt(34))
  1321. exp2_P3_feedback
  1322. ```
  1324. ```{r P3 and social comparison plot}
  1325. exp2_sub_P3_feedback = exp2_interval_ERPS %>%
  1326. dplyr::group_by(SubID) %>%
  1327. drop_na(P3)%>%
  1328. dplyr::group_by(SubID, FeedCode) %>%
  1329. dplyr::summarise(P3_mean = round(mean(P3),3),
  1330. P3_sd = round(sd(P3),3),
  1331. P3_se=P3_sd/sqrt(34))
  1334. exp2_P3.feedback_plot= ggbarplot(exp2_sub_P3_feedback, x = "FeedCode",
  1335. y = "P3_mean",
  1336. alpha = 0.6, ylab= "P3",
  1337. xlab = "Social Comparison",
  1338. color = "black",fill = "FeedCode",
  1339. add = c("mean_se","jitter"),
  1340. add.params = list(color = "FeedCode"),
  1341. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  1342. legend = 'right', position = position_dodge(width =0.4)) +
  1343. scale_y_continuous(breaks = seq(0, 25, by = 5), expand = c(0, 0)) +
  1344. coord_cartesian(ylim = c(0, 27))+
  1345. theme(axis.text = element_text(size = 16, family = "serif"),
  1346. axis.title = element_text(size = 20, family = "serif"))
  1348. exp2_P3.feedback_plot = ggpar(exp2_P3.feedback_plot, legend = 'right') +
  1349. theme(legend.text = element_text(size = 12, family = "serif"),
  1350. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  1352. exp2_P3.feedback_plot
  1354. ```
  1357. ```{r LPP and social comparison lmm}
  1358. exp2_LPP.feedback_lmm = lmerTest::lmer(LPP~age+gender+FeedCode+
  1359. (1+FeedCode|SubID),
  1360. exp2_interval_feedERPS, REML=FALSE)
  1361. exp2_LPP.feedback_anova = anova(exp2_LPP.feedback_lmm)
  1362. exp2_LPP.feedback_anova
  1363. summary(exp2_LPP.feedback_lmm)
  1365. exp2_LPP.feedback_std_ci = model_parameters(
  1366. exp2_LPP.feedback_lmm,
  1367. standardize = "refit",
  1368. df_method = "satterthwaite",
  1369. ci_method = "wald",
  1370. ci = 0.95,
  1371. effects = "all",
  1372. iterations = 1000,
  1373. summary = getOption("parameters_mixed_summary", FALSE),
  1374. digits = 3)
  1376. exp2_LPP.feedback_std_ci
  1377. ```
  1380. ```{r LPP and social comparison des}
  1381. exp2_LPP_feedback = exp2_interval_feedERPS %>%
  1382. dplyr::group_by(SubID) %>%
  1383. drop_na(LPP)%>%
  1384. dplyr::group_by(FeedCode) %>%
  1385. dplyr::summarise(LPP_mean = round(mean(LPP),3),
  1386. LPP_sd = round(sd(LPP),3),
  1387. LPP_se=LPP_sd/sqrt(34))
  1388. exp2_LPP_feedback
  1389. ```
  1391. ```{r LPP and social comparison plot}
  1392. exp2_sub_LPP_feedback = exp2_interval_ERPS %>%
  1393. dplyr::group_by(SubID) %>%
  1394. drop_na(LPP)%>%
  1395. dplyr::group_by(SubID, FeedCode) %>%
  1396. dplyr::summarise(LPP_mean = round(mean(LPP),3),
  1397. LPP_sd = round(sd(LPP),3),
  1398. LPP_se=LPP_sd/sqrt(34))
  1401. exp2_LPP.feedback_plot= ggbarplot(exp2_sub_LPP_feedback, x = "FeedCode",
  1402. y = "LPP_mean",
  1403. alpha = 0.6, ylab= "LPP",
  1404. xlab = "Social Comparison",
  1405. color = "black",fill = "FeedCode",
  1406. add = c("mean_se","jitter"),
  1407. add.params = list(color = "FeedCode"),
  1408. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  1409. legend = 'right', position = position_dodge(width =0.4)) +
  1410. scale_y_continuous(breaks = seq(-4, 10, by = 2), expand = c(0, 0)) +
  1411. coord_cartesian(ylim = c(-5, 11))+
  1412. theme(axis.text = element_text(size = 16, family = "serif"),
  1413. axis.title = element_text(size = 20, family = "serif"))
  1415. exp2_LPP.feedback_plot = ggpar(exp2_LPP.feedback_plot, legend = 'right') +
  1416. theme(legend.text = element_text(size = 12, family = "serif"),
  1417. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  1419. exp2_LPP.feedback_plot
  1421. ```
  1423. ```{r feedback power}
  1424. exp2_RewP.feedback_lmm = lmerTest::lmer(RewP~age+gender+FeedCode+
  1425. (1+FeedCode|SubID),exp2_interval_feedERPS, REML=FALSE)
  1426. summary(exp2_RewP.feedback_lmm )
  1427. modelRewP=powerSim(exp2_RewP.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  1429. exp2_P3.feedback_lmm = lmerTest::lmer(P3~age+gender+FeedCode+
  1430. (1+FeedCode|SubID),exp2_interval_ERPS, REML=FALSE)
  1432. summary(exp2_P3.feedback_lmm)
  1434. modelP3=powerSim(exp2_P3.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  1436. exp2_LPP.feedback_lmm = lmerTest::lmer(LPP~age+gender+FeedCode+
  1437. (1+FeedCode|SubID),
  1438. exp2_interval_ERPS, REML=FALSE)
  1439. summary(exp2_LPP.feedback_lmm)
  1441. modelLPP=powerSim(exp2_LPP.feedback_lmm,fixed("FeedCode", "f"),nsim=1000)
  1442. ```
  1446. ##Exp2 Cue pahse analysis
  1447. ```{r CNV and social comparison}
  1448. exp2_CNV.feedback_lmm = lmerTest::lmer(CNV~age+gender+FeedPrev+
  1449. (1+FeedPrev|SubID),
  1450. exp2_interval_cueERPS, REML=FALSE)
  1451. exp2_CNV.feedback_anova = anova(exp2_CNV.feedback_lmm)
  1452. exp2_CNV.feedback_anova
  1453. summary(exp2_CNV.feedback_lmm)
  1456. exp2_CNV.feedback_std_ci = model_parameters(
  1457. exp2_CNV.feedback_lmm,
  1458. standardize = "refit",
  1459. df_method = "satterthwaite",
  1460. ci_method = "wald",
  1461. ci = 0.95,
  1462. effects = "all",
  1463. iterations = 1000,
  1464. summary = getOption("parameters_mixed_summary", FALSE),
  1465. digits = 3)
  1467. exp2_CNV.feedback_std_ci
  1468. ```
  1471. ```{r CNV and social comparison des}
  1472. exp2_CNV_feedback = exp2_interval_cueERPS %>%
  1473. dplyr::group_by(SubID) %>%
  1474. drop_na(CNV)%>%
  1475. dplyr::group_by(FeedPrev) %>%
  1476. dplyr::summarise(CNV_mean = round(mean(CNV),3),
  1477. CNV_sd = round(sd(CNV),3),
  1478. CNV_se = CNV_sd/sqrt(34))
  1479. exp2_CNV_feedback
  1480. ```
  1482. ```{r CNV and social comparison plot}
  1483. exp2_sub_CNV_feedback = exp2_interval_cueERPS %>%
  1484. dplyr::group_by(SubID) %>%
  1485. drop_na(CNV)%>%
  1486. dplyr::group_by(SubID, FeedPrev) %>%
  1487. dplyr::summarise(CNV_mean = round(mean(CNV),3),
  1488. CNV_sd = round(sd(CNV),3),
  1489. CNV_se = CNV_sd/sqrt(34))
  1491. exp2_CNV.feedback_plot= ggbarplot(exp2_sub_CNV_feedback, x = "FeedPrev",
  1492. y = "CNV_mean", alpha = 0.6,
  1493. ylab= "CNV", xlab = "Social Comparison",
  1494. color = "black",fill = "FeedPrev",
  1495. add = c("mean_se","jitter"),
  1496. add.params = list(color = "FeedPrev"),
  1497. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  1498. legend = 'right', position = position_dodge(width =0.4))+
  1499. scale_y_continuous(breaks = seq(-8, 8, by = 2), expand = c(0, 0)) +
  1500. coord_cartesian(ylim = c(-7.5, 7.5))+
  1501. theme(axis.text = element_text(size = 16, family = "serif"),
  1502. axis.title = element_text(size = 20, family = "serif"))
  1504. exp2_CNV.feedback_plot = ggpar(exp2_CNV.feedback_plot, legend = 'right') +
  1505. theme(legend.text = element_text(size = 12, family = "serif"),
  1506. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  1507. exp2_CNV.feedback_plot
  1508. ```
  1510. ```{r CNV and efficacy}
  1511. exp2_CNV.efficacy_lmm = lmerTest::lmer(CNV~gender+mbased_efficacy_prev+
  1512. (1+mbased_efficacy_prev|SubID),
  1513. exp2_interval_cueERPS, REML=FALSE)
  1514. exp2_CNV.efficacy_anova = anova(exp2_CNV.efficacy_lmm)
  1515. exp2_CNV.efficacy_anova
  1516. summary(exp2_CNV.efficacy_lmm)
  1518. exp2_CNV.efficacy_std_ci = model_parameters(
  1519. exp2_CNV.efficacy_lmm,
  1520. standardize = "refit",
  1521. df_method = "satterthwaite",
  1522. ci_method = "wald",
  1523. ci = 0.95,
  1524. effects = "all",
  1525. iterations = 1000,
  1526. summary = getOption("parameters_mixed_summary", FALSE),
  1527. digits = 3)
  1529. exp2_CNV.efficacy_std_ci
  1530. ```
  1532. ```{r CNV and efficacy plot}
  1533. eff_df3 <- Effect(c("mbased_efficacy_prev"), exp2_CNV.efficacy_lmm,
  1534. xlevels = list(mbased_efficacy_prev = seq(min(exp2_interval_cueERPS$mbased_efficacy_prev, na.rm = TRUE),
  1535. max(exp2_interval_cueERPS$mbased_efficacy_prev, na.rm = TRUE), 0.1)))
  1536. exp2.CNV.efficacy <- as.data.frame(eff_df3)
  1537. head(exp2.CNV.efficacy)
  1539. len3 <- length(exp2.CNV.efficacy$mbased_efficacy_prev)
  1541. exp2.plot.CNV.efficacy <-
  1542. ggplot() +
  1543. geom_line(data = exp2.CNV.efficacy,
  1544. aes(x = mbased_efficacy_prev, y = fit, color = "S"),
  1545. color = "#F8766D", size = 1, linetype = 1) +
  1546. geom_ribbon(data = exp2.CNV.efficacy,
  1547. aes(x = mbased_efficacy_prev, ymax = fit + se, ymin = fit - se),
  1548. fill = "#F8766D", alpha = 0.1, inherit.aes = FALSE) +
  1549. gg.side +
  1550. scale_y_continuous(breaks = scales::breaks_width(1))
  1552. xmax <- max(exp2.CNV.efficacy$mbased_efficacy_prev, na.rm = TRUE)
  1554. exp2.plot.CNV.efficacy <- exp2.plot.CNV.efficacy +
  1555. coord_cartesian(xlim = c(0, xmax))
  1558. exp2.plot.CNV.efficacy
  1559. ```
  1563. ```{r CNV and effort}
  1564. exp2_CNV.CRPS_lmm = lmerTest::lmer(CRPS~age+gender+mbased_efficacy_prev+CNV+
  1565. (1+mbased_efficacy_prev+CNV|SubID),
  1566. exp2_interval_cueERPS, REML=FALSE)
  1567. exp2_CNV.CRPS_anova = anova(exp2_CNV.CRPS_lmm)
  1568. exp2_CNV.CRPS_anova
  1569. summary(exp2_CNV.CRPS_lmm)
  1571. exp2_CNV.CRPS_std_ci = model_parameters(
  1572. exp2_CNV.CRPS_lmm,
  1573. standardize = "refit",
  1574. df_method = "satterthwaite",
  1575. ci_method = "wald",
  1576. ci = 0.95,
  1577. effects = "all",
  1578. iterations = 1000,
  1579. summary = getOption("parameters_mixed_summary", FALSE),
  1580. digits = 3)
  1582. exp2_CNV.CRPS_std_ci
  1583. ```
  1585. ```{r Cuebeta and social comparison}
  1586. exp2_Cuebeta.fd_lmm = lmerTest::lmer(cuebeta~age+gender+FeedPrev+(1+FeedPrev|SubID),
  1587. exp2_interval_cueERPS, REML=FALSE)
  1588. exp2_Cuebeta.fd_anova = anova(exp2_Cuebeta.fd_lmm)
  1589. exp2_Cuebeta.fd_anova
  1590. summary(exp2_Cuebeta.fd_lmm)
  1592. exp2_Cuebeta.fd_std_ci = model_parameters(
  1593. exp2_Cuebeta.fd_lmm,
  1594. standardize = "refit",
  1595. df_method = "satterthwaite",
  1596. ci_method = "wald",
  1597. ci = 0.95,
  1598. effects = "all",
  1599. iterations = 1000,
  1600. summary = getOption("parameters_mixed_summary", FALSE),
  1601. digits = 3)
  1603. exp2_Cuebeta.fd_std_ci
  1604. ```
  1607. ```{r Cuebeta and social comparison des}
  1608. exp2_Cuebeta_feedback = exp2_interval_cueERPS %>%
  1609. dplyr::group_by(SubID) %>%
  1610. drop_na(cuebeta)%>%
  1611. dplyr::group_by(FeedPrev) %>%
  1612. dplyr::summarise(cuebeta_mean = round(mean(cuebeta),3),
  1613. cuebeta_sd = round(sd(cuebeta),3),
  1614. cuebeta_se=cuebeta_sd/sqrt(34))
  1615. exp2_Cuebeta_feedback
  1616. ```
  1618. ```{r Cuebeta and social comparison plot}
  1619. exp2_sub_cue_feedback = exp2_interval_cueERPS %>%
  1620. dplyr::group_by(SubID) %>%
  1621. drop_na(cuebeta)%>%
  1622. dplyr::group_by(SubID,FeedPrev) %>%
  1623. dplyr::summarise(cuebeta_mean = round(mean(cuebeta),3),
  1624. cuebeta_sd = round(sd(cuebeta),3),
  1625. cuebeta_se=cuebeta_sd/sqrt(34))
  1627. exp2_cuebeta.feedback.plot = ggbarplot(exp2_sub_cue_feedback, x = "FeedPrev",
  1628. y = "cuebeta_mean", alpha = 0.6,
  1629. ylab= "cuebeta", xlab = "Social Comparison",
  1630. color = "black",fill = "FeedPrev",
  1631. add = c("mean_se","jitter"),
  1632. add.params = list(color = "FeedPrev"),
  1633. palette =c("#1b7c3d","#2b6a99"),width = 0.4,
  1634. legend = 'right', position = position_dodge(width =0.4)) +
  1635. scale_y_continuous(
  1636. breaks = seq(-1.5, 1.5, by = 0.5),
  1637. labels = scales::number_format(accuracy = 0.1)
  1638. ) +
  1639. coord_cartesian(ylim = c(-1.5, 1.5))+
  1640. theme(axis.text = element_text(size = 16, family = "serif"),
  1641. axis.title = element_text(size = 20, family = "serif"))
  1643. exp2_cuebeta.feedback.plot = ggpar(exp2_cuebeta.feedback.plot, legend = 'right') +
  1644. theme(legend.text = element_text(size = 12, family = "serif"),
  1645. legend.title = element_text(size = 14, family = "serif", face = "bold"))
  1647. exp2_cuebeta.feedback.plot
  1649. ```
  1653. ```{r Cuebeta and efficacy}
  1654. exp2_Cuebeta.efficacy_lmm = lmerTest::lmer(cuebeta~age+gender+mbased_efficacy_prev+
  1655. (1+mbased_efficacy_prev|SubID),exp2_interval_cueERPS, REML=FALSE)
  1656. exp2_Cuebeta.efficacy_anova = anova(exp2_Cuebeta.efficacy_lmm)
  1657. exp2_Cuebeta.efficacy_anova
  1658. summary(exp2_Cuebeta.efficacy_lmm)
  1660. exp2_Cuebeta.efficacy_std_ci = model_parameters(
  1661. exp2_Cuebeta.efficacy_lmm,
  1662. standardize = "refit",
  1663. df_method = "satterthwaite",
  1664. ci_method = "wald",
  1665. ci = 0.95,
  1666. effects = "all",
  1667. iterations = 1000,
  1668. summary = getOption("parameters_mixed_summary", FALSE),
  1669. digits = 3)
  1671. exp2_Cuebeta.efficacy_std_ci
  1672. ```
  1674. ```{r Cuebeta and efficacy plot}
  1675. eff_df4 <- Effect(c("mbased_efficacy_prev"), exp2_Cuebeta.efficacy_lmm,
  1676. xlevels = list(mbased_efficacy_prev = seq(min(exp2_interval_cueERPS$mbased_efficacy_prev, na.rm = TRUE),
  1677. max(exp2_interval_cueERPS$mbased_efficacy_prev, na.rm = TRUE), 0.1)))
  1678. exp2.Cuebeta.efficacy <- as.data.frame(eff_df4)
  1679. head(exp2.Cuebeta.efficacy)
  1680. len4<-length(exp2_interval_cueERPS$mbased_efficacy_prev)
  1681. exp2.plot.Cuebeta.efficacy <- ggplot()+
  1683. geom_line(data=exp2.Cuebeta.efficacy, aes(x=mbased_efficacy_prev, y=fit,color="S"),color="#F8766D",size=1,linetype=1)+
  1684. geom_ribbon(data=exp2.Cuebeta.efficacy, aes(x=mbased_efficacy_prev, max = fit + se, min = fit- se),
  1685. fill = "#F8766D",alpha=0.1, inherit.aes = FALSE)+gg.side
  1687. xmax <- max(exp2.Cuebeta.efficacy$mbased_efficacy_prev, na.rm = TRUE)
  1689. exp2.plot.Cuebeta.efficacy <- exp2.plot.Cuebeta.efficacy +
  1690. coord_cartesian(xlim = c(0, xmax))
  1693. exp2.plot.Cuebeta.efficacy
  1695. ```
  1700. ```{r cue power}
  1701. exp2_CNV.feedback_lmm = lmerTest::lmer(CNV~age+gender+FeedPrev+
  1702. (1+FeedPrev|SubID),
  1703. exp2_interval_cueERPS, REML=FALSE)
  1704. summary(exp2_CNV.feedback_lmm)
  1705. modelCNV=powerSim(exp2_CNV.feedback_lmm,fixed("FeedPrev", "f"),nsim=1000)
  1707. exp2_CNV.efficacy_lmm = lmerTest::lmer(CNV~age+gender+mbased_efficacy_prev+
  1708. (1+mbased_efficacy_prev|SubID),
  1709. exp2_interval_cueERPS, REML=FALSE)
  1710. summary(exp2_CNV.efficacy_lmm)
  1711. model_eff_CNV=powerSim(exp2_CNV.efficacy_lmm,fixed("mbased_efficacy_prev", "f"),nsim=1000)
  1713. exp2_Cuebeta.fd_lmm = lmerTest::lmer(cuebeta~age+gender+FeedPrev+(1+FeedPrev|SubID),
  1714. exp2_interval_cueERPS, REML=FALSE)
  1715. summary(exp2_Cuebeta.fd_lmm)
  1716. modelCuebeta=powerSim(exp2_Cuebeta.fd_lmm,fixed("FeedPrev", "f"),nsim=1000)
  1718. exp2_Cuebeta.efficacy_lmm = lmerTest::lmer(cuebeta~age+gender+mbased_efficacy_prev+
  1719. (1+mbased_efficacy_prev|SubID),
  1720. exp2_interval_cueERPS, REML=FALSE)
  1721. summary(exp2_Cuebeta.efficacy_lmm)
  1722. model_eff_Cuebeta=powerSim(exp2_Cuebeta.efficacy_lmm,fixed("mbased_efficacy_prev", "f"),nsim=1000)
  1723. ```
  1725. # Exp1 and Exp2 hddm
  1726. ## Exp1 hddm analysis
  1727. ```{r exp1 social comaprison hddm,echo = FALSE}
  1728. ###model1 social comparison feedback ###
  1729. exp1_feedback.hddm = read.csv('.../HDDM/exp1_feedback_cong_traces7000.csv')
  1731. exp1_feedback.traces = exp1_feedback.hddm %>%
  1732. dplyr::select(v_Intercept,v_Feed.T.up.,
  1733. a_Intercept,a_Feed.T.up.)%>%
  1734. dplyr::rename(v_upward = v_Intercept,
  1735. v_diff = v_Feed.T.up.,
  1736. a_upward = a_Intercept,
  1737. a_diff = a_Feed.T.up.)%>%
  1738. mutate(v_downward = v_diff+ v_upward,
  1739. a_downward = a_diff + a_upward)
  1741. # Add the posterior probabilities
  1742. exp1_feedback_hypotheses = list(
  1743. t1 = hypothesis(exp1_feedback.traces, "v_downward > v_upward"),
  1744. t2 = hypothesis(exp1_feedback.traces, "a_downward < a_upward"))
  1746. exp1_feedback.hddm.results = data.frame(
  1747. test = c("v_downward > v_upward",
  1748. "a_downward < a_upward"),
  1749. Post.Prob = sapply(exp1_feedback_hypotheses, function(x) x$hypothesis$Post.Prob),
  1750. estimate = sapply(exp1_feedback_hypotheses, function(x) x$hypothesis$Estimate),
  1751. ci.lower = sapply(exp1_feedback_hypotheses, function(x) x$hypothesis$CI.Lower),
  1752. ci.upper = sapply(exp1_feedback_hypotheses, function(x) x$hypothesis$CI.Upper),
  1753. est.error = sapply(exp1_feedback_hypotheses, function(x) x$hypothesis$Est.Error))
  1755. exp1_feedback.hddm.results
  1756. ```
  1758. ```{r exp1 social comaprison hddm plot}
  1759. exp1_feedback.plot = bind_rows(
  1760. exp1_feedback.traces %>%
  1761. mutate(diff = v_diff,
  1762. Category = 'v')%>%
  1763. dplyr::select(diff,Category),
  1765. exp1_feedback.traces %>%
  1766. mutate(diff = a_diff,
  1767. Category = 'a')%>%
  1768. dplyr::select(diff,Category),
  1769. )
  1773. exp1_feedback.hddm=ggplot(exp1_feedback.plot, aes(x = diff, fill = Category)) +
  1774. geom_density(alpha = 0.4, aes(color = Category)) +
  1775. scale_fill_manual(values = c('#24A669', '#BF6B82')) +
  1776. scale_color_manual(values = c('#24A669', '#BF6B82'))+
  1777. labs(x = 'Social', y = 'Density')+gg.side
  1779. exp1_feedback.hddm
  1781. ```
  1783. ```{r exp1 efficacy hddm,echo = FALSE}
  1784. ###model2 efficacy ###
  1785. exp1_efficacy.hddm = read.csv(".../HDDM/exp1_efficacy_traces.csv")
  1787. exp1_efficacy.traces<-exp1_efficacy.hddm %>%
  1788. dplyr::select(v_Intercept,
  1789. v_mbased_efficacy_prev,
  1790. a_Intercept,
  1791. a_mbased_efficacy_prev)
  1793. # Add the posterior probabilities
  1794. exp1_efficacy_hypotheses = list(
  1795. t1 = hypothesis(exp1_efficacy.traces, "v_Intercept < 0"),
  1796. t2 = hypothesis(exp1_efficacy.traces, "v_mbased_efficacy_prev > 0"),
  1797. t3 = hypothesis(exp1_efficacy.traces, "a_Intercept > 0"),
  1798. t4 = hypothesis(exp1_efficacy.traces, "a_mbased_efficacy_prev < 0"))
  1800. exp1_efficacy.hddm.results = data.frame(
  1801. test = c("v_Intercept < 0",
  1802. "v_mbased_efficacy_prev > 0",
  1803. "a_Intercept > 0",
  1804. "a_mbased_efficacy_prev < 0"),
  1805. Post.Prob = sapply(exp1_efficacy_hypotheses, function(x) x$hypothesis$Post.Prob),
  1806. estimate = sapply(exp1_efficacy_hypotheses, function(x) x$hypothesis$Estimate),
  1807. ci.lower = sapply(exp1_efficacy_hypotheses, function(x) x$hypothesis$CI.Lower),
  1808. ci.upper = sapply(exp1_efficacy_hypotheses, function(x) x$hypothesis$CI.Upper),
  1809. est.error = sapply(exp1_efficacy_hypotheses, function(x) x$hypothesis$Est.Error))
  1811. exp1_efficacy.hddm.results
  1812. ```
  1814. ```{r exp1 efficacy hddm plot}
  1815. exp1_efficacy.plot = bind_rows(
  1816. exp1_efficacy.traces %>%
  1817. dplyr::mutate(Intercept = v_Intercept,
  1818. mbased_efficacy = v_mbased_efficacy_prev,
  1819. Category = 'v')%>%
  1820. dplyr::select(Intercept,mbased_efficacy,Category),
  1822. exp1_efficacy.traces %>%
  1823. dplyr::mutate(Intercept = a_Intercept,
  1824. mbased_efficacy = a_mbased_efficacy_prev,
  1825. Category = 'a')%>%
  1826. dplyr::select(Intercept,mbased_efficacy,Category),
  1827. )
  1829. exp1_efficacy.hddm=ggplot(exp1_efficacy.plot, aes(x = mbased_efficacy, fill = Category)) +
  1830. geom_density(alpha = 0.4, aes(color = Category)) +
  1831. scale_fill_manual(values = c('#24A669', '#BF6B82')) +
  1832. scale_color_manual(values = c('#24A669', '#BF6B82')) +
  1833. labs(x = 'Efficacy', y = 'Density')+gg.side
  1835. exp1_efficacy.hddm
  1837. ```
  1839. ## Exp2 hddm analysis
  1840. ```{r exp2 social comparison hddm}
  1841. ###model1 social comparison feedback ###
  1842. exp2_feedback.hddm = read.csv(".../HDDM/exp2_feedback_cong_traces7000.csv")
  1844. exp2_feedback.traces = exp2_feedback.hddm %>%
  1845. dplyr::select(v_Intercept,v_Feed.T.down.,
  1846. a_Intercept,a_Feed.T.down.)%>%
  1847. dplyr::rename(v_diff = v_Feed.T.down.,
  1848. v_downward = v_Intercept,
  1849. a_diff = a_Feed.T.down.,
  1850. a_downward = a_Intercept)%>%
  1851. mutate(v_upward = v_downward + v_diff,
  1852. a_upward = a_downward + a_diff)
  1855. # Add the posterior probabilities
  1856. exp2_feedback_hypotheses = list(
  1857. t1 = hypothesis(exp2_feedback.traces, "v_downward > v_upward"),
  1858. t2 = hypothesis(exp2_feedback.traces, "a_downward < a_upward"))
  1860. exp2_feedback.hddm.results = data.frame(
  1861. test = c("v_downward > v_upward",
  1862. "a_downward < a_upward"),
  1863. Post.Prob = sapply(exp2_feedback_hypotheses, function(x) x$hypothesis$Post.Prob),
  1864. estimate = sapply(exp2_feedback_hypotheses, function(x) x$hypothesis$Estimate),
  1865. ci.lower = sapply(exp2_feedback_hypotheses, function(x) x$hypothesis$CI.Lower),
  1866. ci.upper = sapply(exp2_feedback_hypotheses, function(x) x$hypothesis$CI.Upper),
  1867. est.error = sapply(exp2_feedback_hypotheses, function(x) x$hypothesis$Est.Error))
  1869. exp2_feedback.hddm.results
  1870. ```
  1872. ```{r exp2 efficacy hddm}
  1873. ###model1 efficacy ###
  1874. exp2_efficacy.hddm = read.csv(".../HDDM/exp2_efficacy_traces.csv")
  1876. exp2_efficacy.traces = exp2_efficacy.hddm %>%
  1877. dplyr::select(v_Intercept,
  1878. v_mbased_efficacy_prev,
  1879. a_Intercept,
  1880. a_mbased_efficacy_prev)
  1882. # Add the posterior probabilities
  1883. exp2_efficacy_hypotheses = list(
  1884. t1 = hypothesis(exp2_efficacy.traces, "v_Intercept < 0"),
  1885. t2 = hypothesis(exp2_efficacy.traces, "v_mbased_efficacy_prev > 0"),
  1886. t3 = hypothesis(exp2_efficacy.traces, "a_Intercept > 0"),
  1887. t4 = hypothesis(exp2_efficacy.traces, "a_mbased_efficacy_prev < 0"))
  1889. exp2_efficacy.hddm.results = data.frame(
  1890. test = c("v_Intercept < 0",
  1891. "v_mbased_efficacy_prev > 0",
  1892. "a_Intercept > 0",
  1893. "a_mbased_efficacy_prev < 0"),
  1894. Post.Prob = sapply(exp2_efficacy_hypotheses, function(x) x$hypothesis$Post.Prob),
  1895. estimate = sapply(exp2_efficacy_hypotheses, function(x) x$hypothesis$Estimate),
  1896. ci.lower = sapply(exp2_efficacy_hypotheses, function(x) x$hypothesis$CI.Lower),
  1897. ci.upper = sapply(exp2_efficacy_hypotheses, function(x) x$hypothesis$CI.Upper),
  1898. est.error = sapply(exp2_efficacy_hypotheses, function(x) x$hypothesis$Est.Error))
  1900. exp2_efficacy.hddm.results
  1901. ```
  1904. ```{r LOOIC plot}
  1905. df1 <- tribble(
  1906. ~Model, ~dLOOIC,
  1907. "Intercept", 552.194,
  1908. "1LR", 287.183,
  1909. "2LR", 0
  1910. ) %>%
  1911. mutate(Model = factor(Model, levels = c("2LR","1LR","Intercept")))
  1913. p1 <- ggplot(df1, aes(x = Model, y = dLOOIC)) +
  1914. geom_col(width = 0.7) +
  1915. geom_text(aes(label = sprintf("%.1f", dLOOIC)),
  1916. vjust = -0.3, size = 3) +
  1917. labs(x = NULL, y = "ΔLOOIC (best model = 0; lower is better)") +
  1918. theme_classic() +
  1919. coord_cartesian(ylim = c(0, max(df1$dLOOIC) * 1.15))+gg.side
  1921. p1
  1922. ```

analysis_HB_scripts.Rmd, no license · at the source

Overview

Authors: Jiarui Dong1,2, Yachao Rong3, Shengjie Ma1,2, Yang Xu1,2, Ping Wei1,2
ORCID iDs: Yang Xu, Ping Wei
  1. School of Psychology, Capital Normal University,Beijing, China
  2. Beijing Key Laboratory of Learning and Cognition, Capital Normal University,Beijing, China
  3. Faculty of Education, Henan Normal University,Xinxiang, China
Institutions: Capital Normal University (China); Henan Normal University (China)
Journal: Communications biology, volume 9, issue 1, article 984
Dates: received 30 September 2025; accepted 28 April 2026; published online 9 May 2026
Type: Research article · Language: English
License: CC BY-NC-ND
Identifiers: DOI 10.1038/s42003-026-10229-5 · PMID 42106486 · PMCID PMC13381887 · OpenAlex W7160691893
Open access: gold, a free copy (OpenAlex)
Status: code verified
Categories: EEG (modality), human (organism)
Methods: Spectral & time-frequency, Statistics, Smoothing, state filtering, decompositions, Preprocessing, Evoked potentials, fMRI & imaging, Physiology & signal measures
Keywords: Human behaviour, Social neuroscience
MeSH: Social Comparison*, Computer Simulation, Contingent Negative Variation, Electroencephalography, Female, Humans, Male, Reward, Self Efficacy (* major topic)
Topic: Neural and Behavioral Psychology Studies (Cognitive Neuroscience, Neuroscience), according to OpenAlex
Funding: National Natural Science Foundation of China (National Science Foundation of China) (32471105, 31971030); Youth Beijing Scholar Program of Beijing Government
Citations: cited by 1 paper (Europe PMC); 74 references in the paper

Abstract

The abstract is not reproduced here: the paper's license (CC BY-NC-ND) does not allow it. Read it in the paper, at the publisher or on Europe PMC.

Repository

Its files are read in the Code ↔ Paper reader above.

OSF vpr9e

License: none: the authors keep all their rights
State: the link answers, verified on 28 September 2026
Evidence: files inventoried
Languages: R (1)
Size: 20 files, 1 script
Software Heritage: not checked
Found in: “Code availability”
Holds: 1 notebook
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: afex (1 file), brms (1 file), easystats (1 file), emmeans (1 file), ggplot2 (1 file), ggpubr (1 file), lmerTest (1 file), psych (1 file), tidyverse (1 file)
Availability: 1 check, the latest on 28 September 2026: the link answers (HTTP 200)
  • 28 September 2026: the link answers (HTTP 200)
1 file
At the source:

Code availability statement

The paper has a code availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

  • it points to the authors' code: OSF vpr9e

Read it in the paper: doi.org/10.1038/s42003-026-10229-5.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 1 repository of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 1 script, each with its path and the digest of its content;
  • no match between paragraphs and code yet;
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

No dataset and no data link were found in the paper.

Data availability statement

The paper has a data availability statement. Its license (CC BY-NC-ND) does not allow reproducing it here; in short, from what the harvester recognized in it:

  • no repository, dataset or request procedure was recognized in it

Read it in the paper: doi.org/10.1038/s42003-026-10229-5.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 28 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 5 authors, 2 keywords, 9 MeSH terms, 2 funders, 71 references.

Cite

This paper

Dong, J., Rong, Y., Ma, S., Xu, Y., & Wei, P. (2026). The neurocomputational mechanisms underlying the impact of social comparison on effort investment. Communications biology, 9(1), 984. https://doi.org/10.1038/s42003-026-10229-5

BibTeX

@article{dong2026neurocomputational,
author = {Dong, Jiarui and Rong, Yachao and Ma, Shengjie and Xu, Yang and Wei, Ping},
title = {{The neurocomputational mechanisms underlying the impact of social comparison on effort investment}},
journal = {Communications biology},
year = {2026},
month = may,
volume = {9},
number = {1},
pages = {984},
publisher = {Nature Publishing Group},
issn = {2399-3642},
doi = {10.1038/s42003-026-10229-5},
url = {https://doi.org/10.1038/s42003-026-10229-5},
pmid = {42106486},
pmcid = {PMC13381887}
}

RIS

TY - JOUR
AU - Dong, Jiarui
AU - Rong, Yachao
AU - Ma, Shengjie
AU - Xu, Yang
AU - Wei, Ping
TI - The neurocomputational mechanisms underlying the impact of social comparison on effort investment
T2 - Communications biology
J2 - Commun Biol
PY - 2026
DA - 2026/05/09
VL - 9
IS - 1
SP - 984
SN - 2399-3642
PB - Nature Publishing Group
DO - 10.1038/s42003-026-10229-5
UR - https://doi.org/10.1038/s42003-026-10229-5
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s42003-026-10229-5",
"type": "article-journal",
"title": "The neurocomputational mechanisms underlying the impact of social comparison on effort investment",
"container-title": "Communications biology",
"author": [
{
"family": "Dong",
"given": "Jiarui"
},
{
"family": "Rong",
"given": "Yachao"
},
{
"family": "Ma",
"given": "Shengjie"
},
{
"family": "Xu",
"given": "Yang"
},
{
"family": "Wei",
"given": "Ping"
}
],
"container-title-short": "Commun Biol",
"volume": "9",
"issue": "1",
"page": "984",
"DOI": "10.1038/s42003-026-10229-5",
"PMID": "42106486",
"PMCID": "PMC13381887",
"ISSN": "2399-3642",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s42003-026-10229-5",
"language": "en",
"issued": {
"date-parts": [
[
2026,
5,
9
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1371/journal.pbio.3003979 [code]
Impaired midfrontal‑motor theta phase synchronization characterizes maladaptive motivational behavior in people with obsessive‑compulsive disorder.
Journal: PLoS biology
In common: afex, psych, emmeans, 3 other tools, EEG, 2 references
[2] doi:10.1073/pnas.2603114123 [code]
The human hippocampus can pattern separate memories by meaning.
Journal: Proceedings of the National Academy of Sciences of the United States of America
In common: afex, psych, easystats, 5 other tools
[3] doi:10.1016/j.neuroimage.2026.122115 [code]
Midfrontal theta power relates to response speeding following frustrative nonreward.
Journal: NeuroImage
In common: psych, easystats, emmeans, 4 other tools, EEG, 1 reference
[4] doi:10.7554/elife.103566 [code]
Effort produces after-effects costly for others but valued for self.
Journal: eLife
In common: psych, emmeans, lmerTest, 3 other tools, EEG, 2 references
[5] doi:10.1073/pnas.2606871123 [code]
Oxytocin modulates the neurocomputational mechanisms engaged in learning rank relationships in social networks.
Journal: Proceedings of the National Academy of Sciences of the United States of America
In common: psych, easystats, emmeans, 4 other tools, 1 reference
[6] doi:10.1038/s41398-026-04341-7 [code]
Aberrant insula activity to negative and reduced learning from positive feedback underlie maladaptive self-beliefs in depression.
Journal: Translational psychiatry
In common: psych, tidyverse, 5 references
[7] doi:10.1371/journal.pone.0353990 [code]
Positive mood enhances accessibility of unrelated concepts in the first language but not in the foreign language.
Journal: PloS one
In common: afex, psych, emmeans, 4 other tools, EEG
[8] doi:10.1038/s41598-026-58046-4 [code]
Dissecting the interplay of model-based control, impulsivity and compulsivity on self-control in daily life.
Journal: Scientific reports
In common: psych, easystats, emmeans, 2 other tools, EEG, 2 references
[9] doi:10.1126/sciadv.aec9291 [code]
Computational mechanisms of perception in autism revealed using games inspired by rodent operant tasks.
Journal: Science advances
In common: afex, psych, easystats, 4 other tools
[10] doi:10.1093/cercor/bhag113 [code]
Long-term reliability and stability of parameterized resting state EEG: evidence from a five-year follow-up.
Journal: Cerebral cortex (New York, N.Y. : 1991)
In common: psych, easystats, emmeans, 4 other tools, EEG

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.