OSCR

Obesity is associated with greater variability of reward signals in the nucleus accumbens.

Code ↔ Paper

1 match between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 1 match
  1. [1] § Methods › Data analysis › Behavioral data ↔ TUE002_EAT_trial_based_plots_share.R, lines 25–81 · score 0.50 · reward magnitude, subBED, age, uncertain, food, BMI

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 · 1,101 lines · 41 KB · no license · 1 match

  1. library(ggplot2)
  2. library(foreign)
  3. library(MASS)
  4. library(cowplot)
  5. library(viridis)
  6. library(readxl)
  7. library(tidyverse)
  8. library(ggdist)
  9. library(dplyr)
  10. library(ggridges)
  11. library(brms)
  14. theme_set(theme_cowplot(font_size=12))
  16. dT <- read.csv("/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/data/TUE002_trial_cue_responses.csv")
  17. dT_fb <- read.csv("/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/data/TUE002_trial_feedback_responses_more_mvmt.csv")
  19. dQ <- read_excel("../data/TUE002_Sample_allS1_inclFEV.xlsx")
  20. d <- read_excel("/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/data/TUE002_S1_incl_diagnosis.xlsx")
  23. d$rBMI <- rank(d$BMI_S1)
  25. #Prep TFEQ for plotting
  26. dQ <-
  27. dQ %>%
  28. mutate(
  29. bTFEQ_Dis = ntile(FEV_disinh, n=4))
  31. dQ <-
  32. dQ %>%
  33. mutate(bTFEQ_Dis = if_else(bTFEQ_Dis>2, bTFEQ_Dis - 1, bTFEQ_Dis))
  35. dQ$fTFEQ_Dis <- factor(dQ$bTFEQ_Dis, labels = c("Low", "Medium", "High"))
  37. dT <- merge(x = dT, y = d, by = "ID", all.x=TRUE)
  39. dT$fBinge <- factor(dT$Group, labels = c('No BE', 'subBED', 'BED'))
  40. dT$fBMI <- cut(dT$BMI_S1, breaks = c(-100,25,30,'Inf'), labels = c("Normal", "Overweight", "Obese"))
  41. dT$fMoney <- factor(dT$Money, labels = c("Food","Money"))
  42. dT$RewMag <- factor(dT$RewardMagnitude, labels = c("Low","High"))
  43. dT$fUnc <- factor(dT$Uncertainty, labels = c("Certain","Uncertain"))
  45. dT$cBMI = dT$BMI_S1 - mean(dT$BMI_S1, na.rm = TRUE)
  46. dT$zBMI = scale(dT$BMI_S1)
  47. dT$cAge = dT$Age - mean(dT$Age, na.rm = TRUE)
  49. dT <- merge(x = dT, y = dQ[ , c("ID","FEV_disinh","bTFEQ_Dis","fTFEQ_Dis")], by = "ID", all.x=TRUE)
  52. dT$cRewM <- dT$RewardMagnitude
  54. dT$cTFEQ_Dis <- dT$FEV_disinh - mean(dT$FEV_disinh)
  56. dT <- left_join(dT,d_SD %>% select(ID, Uncertainty,std_Cue_aMTL,std_Cue_pMTL,std_Cue_toMTL))
  59. dT_fb <- merge(x = dT_fb, y = d, by = "ID", all.x=TRUE)
  61. dT_fb$fBinge <- factor(dT_fb$Group, labels = c('No BE', 'subBED', 'BED'))
  62. dT_fb$fBMI <- cut(dT_fb$BMI_S1, breaks = c(-100,25,30,'Inf'), labels = c("Normal", "Overweight", "Obese"))
  63. dT_fb$fMoney <- factor(dT_fb$Money, labels = c("Food","Money"))
  64. dT_fb$RewMag <- factor(dT_fb$Reward.Magnitude, labels = c("Low","High"))
  65. dT_fb$fUnc <- factor(dT_fb$Uncertainty, labels = c("Certain","Uncertain"))
  67. dT_fb$cBMI = dT_fb$BMI_S1 - mean(dT_fb$BMI_S1, na.rm = TRUE)
  68. dT_fb$zBMI = scale(dT_fb$BMI_S1)
  69. dT_fb$cAge = dT$Age - mean(dT_fb$Age, na.rm = TRUE)
  71. dT_fb <- merge(x = dT_fb, y = dQ[ , c("ID","FEV_disinh","bTFEQ_Dis","fTFEQ_Dis")], by = "ID", all.x=TRUE)
  74. dT_fb$cRewM <- (dT_fb$Reward.Magnitude-5.5)/4.5
  76. dT_fb$cTFEQ_Dis <- dT_fb$FEV_disinh - mean(dT_fb$FEV_disinh)
  78. dT_fb <- left_join(dT_fb,d_SD %>% select(ID, Uncertainty,std_Cue_aMTL,std_Cue_pMTL,std_Cue_toMTL))
  80. dT_S1_Unc <- read_csv("/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/data/TUE002_grEAT_S1_resid_Unc.csv")
  81. dT_S1_Cer <- read_csv("/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/data/TUE002_grEAT_S1_resid_Cer.csv")
  85. dT_S1_Cer$fBinge <- factor(dT_S1_Cer$fBinge)
  86. dT_S1_Cer$fBinge <- fct_relevel(dT_S1_Cer$fBinge, "No BE")
  87. dT_S1_Cer$fMoney <- factor(dT_S1_Cer$cRewT, labels = c("Food","Money"))
  88. dT_S1_Cer$zBMI = scale(dT_S1_Cer$BMI_S1)
  89. dT_S1_Cer <-
  90. dT_S1_Cer %>%
  91. mutate(bTFEQ_Dis = ntile(FEV_disinh, n=4)) %>%
  92. mutate(bTFEQ_Dis = if_else(bTFEQ_Dis>2, bTFEQ_Dis - 1, bTFEQ_Dis))
  94. dT_S1_Cer$fTFEQ_Dis <- factor(dT_S1_Cer$bTFEQ_Dis, labels = c("Low", "Medium", "High"))
  96. dT_S1_Unc$fBinge <- factor(dT_S1_Unc$fBinge)
  97. dT_S1_Unc$fBinge <- fct_relevel(dT_S1_Unc$fBinge, "No BE","subBED","BED")
  98. dT_S1_Unc$fMoney <- factor(dT_S1_Unc$cRewT, labels = c("Food","Money"))
  99. dT_S1_Unc$zBMI = scale(dT_S1_Unc$BMI_S1)
  101. dT_S1_Unc <-
  102. dT_S1_Unc %>%
  103. mutate(bTFEQ_Dis = ntile(FEV_disinh, n=4)) %>%
  104. mutate(bTFEQ_Dis = if_else(bTFEQ_Dis>2, bTFEQ_Dis - 1, bTFEQ_Dis))
  106. dT_S1_Unc$fTFEQ_Dis <- factor(dT_S1_Unc$bTFEQ_Dis, labels = c("Low", "Medium", "High"))
  108. dT_S1_Unc %>%
  109. dplyr::select(ID,fBinge,BMI_S1,FEV_disinh,Age) %>%
  110. unique()-> demo
  113. dT_S1_Cer$cTFEQ_Dis <- dT_S1_Cer$FEV_disinh - mean(dT_S1_Cer$FEV_disinh)
  114. dT_S1_Unc$cTFEQ_Dis <- dT_S1_Unc$FEV_disinh - mean(dT_S1_Unc$FEV_disinh)
  116. dT_S1_Cer$certain <- 1
  117. dT_S1_Unc$certain <- 0
  119. dT_S1 <- rbind(dT_S1_Cer,dT_S1_Unc)
  124. #Plots
  125. #Figure 3
  126. pF1a <-
  127. ggplot(aes(x = fBinge,y = Cue_NAcc, fill = fBinge), data=dT) +
  128. stat_halfeye(size = 5, alpha = 0.8) +
  129. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#86c8f7','#4a87af','#004c6d')) +
  130. coord_cartesian(ylim = c(-0.52, 0.52)) +
  131. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  132. axis.text.x = element_text(size = 11.0), legend.position = 'none',
  133. strip.text.x = element_text(margin = margin(0.15,0,0.15,0, "cm")),
  134. plot.title = element_text(size = 13, hjust = 0.5),
  135. plot.subtitle = element_text(size = 12, hjust = 0.5))+
  136. facet_grid(fMoney ~ .) +
  137. geom_hline(yintercept = 0, color = 'grey60', linewidth = 0.5) +
  138. ggtitle("NAcc cue responses", "Split by reward type") +
  139. ylab(label = 'Trial residuals') +
  140. xlab(label = 'Group')
  142. pF1b <-
  143. ggplot(aes(x = fBMI,y = Cue_NAcc, fill = fBMI), data=dT) +
  144. stat_halfeye(size = 5, alpha = 0.8) +
  145. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#E3D0D8','#aea3b0','#827081')) +
  146. scale_x_discrete(labels = c("Normal", "Over-\nweight","Obese")) +
  147. coord_cartesian(ylim = c(-0.52, 0.52)) +
  148. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  149. axis.text.x = element_text(size = 11.0), legend.position = 'none',
  150. strip.text.x = element_text(margin = margin(0.10,0,0.10,0, "cm")))+
  151. facet_grid(fMoney ~ .) +
  152. geom_hline(yintercept = 0, color = 'grey60', linewidth = 0.5) +
  153. ylab(label = 'Trial residuals') +
  154. xlab(label = 'BMI cat')
  156. pF1c <-
  157. ggplot(aes(x = fTFEQ_Dis,y = Cue_NAcc, fill = fTFEQ_Dis), data=dT) +
  158. stat_halfeye(size = 5, alpha = 0.8) +
  159. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#2AFC98','#16C172','#214F4B')) +
  160. coord_cartesian(ylim = c(-0.52, 0.52)) +
  161. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  162. axis.text.x = element_text(size = 11.0), legend.position = 'none',
  163. strip.text.x = element_text(margin = margin(0.10,0,0.10,0, "cm")))+
  164. facet_grid(fMoney ~ .) +
  165. geom_hline(yintercept = 0, color = 'grey60', linewidth = 0.5) +
  166. ylab(label = 'Trial residuals') +
  167. xlab(label = 'TFEQ Disinhibition')
  169. pF1d <-
  170. ggplot(aes(y = as.factor(rBMI),x = Cue_NAcc, fill = fBinge), data=dT) +
  171. #stat_halfeye(size = 5, alpha = 0.8) +
  172. geom_density_ridges(rel_min_height = 0.005, scale = 7, alpha = 0.75, color = 'grey80') +
  173. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#86c8f7','#4a87af','#004c6d')) +
  174. coord_cartesian(xlim = c(-0.52, 0.52)) +
  175. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  176. axis.text.x = element_text(size = 11.0), axis.text.y=element_blank(),
  177. axis.ticks.y=element_blank(), strip.text.x = element_text(margin = margin(0.10,0,0.10,0, "cm")),
  178. plot.title = element_text(size = 13, hjust = 0.5, margin=margin(0,0,40,0)),
  179. legend.position = 'none', legend.key.width = unit(1.5, 'cm'),) +
  180. facet_grid(. ~ fBinge) +
  181. geom_vline(xintercept = 0, color = 'grey60', linewidth = 0.5) +
  182. ggtitle('Variability of NAcc cue responses', '') +
  183. xlab(label = 'Trial residuals') +
  184. ylab(label = 'ID [ranked by BMI]')
  192. p1x <- plot_grid(pF1a, pF1b, pF1c, labels = c("b", "c", "d"), label_size = 12,
  193. ncol=1, rel_heights = c(1.2, 1, 1), align = "v", axis = 1)
  195. p1 <- plot_grid(pF1d, p1x, labels = c("a", ""), label_size = 12,
  196. ncol=2, rel_widths = c(1.25, 1), align = "v", axis = 1)
  200. ggsave("../../Plots/TUE002_EAT_NAccVAR.png",
  201. plot = p1, height = 8, width = 6.5, units = "in", dpi = 600, bg = "white")
  206. #Figure 4 dlpfc
  208. ci_data_Binge <- dT %>%
  209. group_by(fBinge, fMoney) %>%
  210. median_qi(Cue_DLPFC, .width = 0.95) %>%
  211. filter(fBinge=="No BE")
  213. ci_data_BMI <- dT %>%
  214. group_by(fBMI, fMoney) %>%
  215. median_qi(Cue_DLPFC, .width = 0.95) %>%
  216. filter(fBMI=="Normal")
  218. ci_data_TFEQ <- dT %>%
  219. group_by(fTFEQ_Dis, fMoney) %>%
  220. median_qi(Cue_DLPFC, .width = 0.95) %>%
  221. filter(fTFEQ_Dis=="Low")
  224. pF2a <-
  225. ggplot(aes(x = fBinge,y = Cue_DLPFC, fill = fBinge), data=dT) +
  226. stat_halfeye(size = 5, alpha = 0.8) +
  227. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#86c8f7','#4a87af','#004c6d')) +
  228. geom_hline(data = ci_data_Binge, aes(yintercept = .lower), color = "black", linetype = "dashed", alpha = .5) +
  229. geom_hline(data = ci_data_Binge, aes(yintercept = .upper), color = "black", linetype = "dashed", alpha = .5) +
  230. coord_cartesian(ylim = c(-0.7, 0.7)) +
  231. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  232. axis.text.x = element_text(size = 12.0), strip.text.x = element_text(margin = margin(0.15,0,0.15,0, "cm")))+
  233. facet_grid(. ~ fMoney) +
  234. geom_hline(yintercept = 0, color = 'grey60', linewidth = 0.5) +
  235. ylab(label = 'Variability\ncue responses') +
  236. xlab(label = 'Group')
  238. pF2b <-
  239. ggplot(aes(x = fBMI,y = Cue_DLPFC, fill = fBMI), data=dT) +
  240. stat_halfeye(size = 5, alpha = 0.8) +
  241. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#E3D0D8','#aea3b0','#827081')) +
  242. geom_hline(data = ci_data_BMI, aes(yintercept = .lower), color = "black", linetype = "dashed", alpha = .5) +
  243. geom_hline(data = ci_data_BMI, aes(yintercept = .upper), color = "black", linetype = "dashed", alpha = .5) +
  244. coord_cartesian(ylim = c(-0.7, 0.7)) +
  245. scale_x_discrete(labels = c("Normal", "Over-\nweight","Obese")) +
  246. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  247. axis.text.x = element_text(size = 12.0), strip.text.x = element_text(margin = margin(0.15,0,0.15,0, "cm")))+
  248. facet_grid(. ~ fMoney) +
  249. geom_hline(yintercept = 0, color = 'grey60', linewidth = 0.5) +
  250. ylab(label = 'Variability\ncue responses') +
  251. xlab(label = 'BMI category')
  253. pF2c <-
  254. ggplot(aes(x = fTFEQ_Dis,y = Cue_DLPFC, fill = fTFEQ_Dis), data=dT) +
  255. stat_halfeye(size = 5, alpha = 0.8) +
  256. scale_fill_manual(guide = guide_legend(title="Group"),values = c('#2AFC98','#16C172','#214F4B')) +
  257. geom_hline(data = ci_data_TFEQ, aes(yintercept = .lower), color = "black", linetype = "dashed", alpha = .5) +
  258. geom_hline(data = ci_data_TFEQ, aes(yintercept = .upper), color = "black", linetype = "dashed", alpha = .5) +
  259. coord_cartesian(ylim = c(-0.7, 0.7)) +
  260. theme(text = element_text(face = 'bold',size = 12.0),axis.text = element_text(face = 'plain',size = 12.0),
  261. axis.text.x = element_text(size = 12.0), strip.text.x = element_text(margin = margin(0.15,0,0.15,0, "cm")))+
  262. facet_grid(. ~ fMoney) +
  263. geom_hline(yintercept = 0, color = 'grey60', linewidth = 0.5) +
  264. ylab(label = 'Variability\ncue responses') +
  265. xlab(label = 'TFEQ Disinhibition')
  267. p2x <- plot_grid(pF2a, pF2b, pF2c, labels = c("a", "b", "c"), label_size = 12,
  268. ncol=1, rel_heights = c(1.2, 1, 1), align = "v", axis = 1)
  271. ggsave("../../Plots/TUE002_EAT_DLPFCVAR.png",
  272. plot = p2x, height = 8, width = 6.5, units = "in", dpi = 600, bg = "white")
  277. #brms
  278. b1.10a <-
  279. brm(data = dT,
  280. family = gaussian,
  281. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  282. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  283. prior = c(prior(normal(0, 0.2), class = Intercept),
  285. prior(normal(0, 1), class = Intercept, dpar = sigma),
  286. prior(normal(0, 0.2), class = b, dpar = sigma),
  287. prior(exponential(1), class = sd, dpar = sigma),
  289. prior(lkj(2), class = cor)),
  290. iter = 5000, warmup = 2000, chains = 4, cores = 4,
  291. seed = 14,
  292. sample_prior = TRUE,
  293. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc")
  295. b1.10b <-
  296. brm(data = dT,
  297. family = gaussian,
  298. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  299. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  300. prior = c(prior(normal(0, 0.2), class = Intercept),
  302. prior(normal(0, 1), class = Intercept, dpar = sigma),
  303. prior(normal(0, 0.2), class = b, dpar = sigma),
  304. prior(exponential(1), class = sd, dpar = sigma),
  306. prior(lkj(2), class = cor)),
  307. iter = 5000, warmup = 2000, chains = 4, cores = 4,
  308. seed = 14,
  309. sample_prior = TRUE,
  310. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_aMTL_DVARS")
  315. b1.10bb <-
  316. brm(data = dT,
  317. family = gaussian,
  318. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  319. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + std_Cue_aMTL + log_FD + (1 + fMoney + cRewM |i| ID)),
  320. prior = c(prior(normal(0, 0.2), class = Intercept),
  322. prior(normal(0, 1), class = Intercept, dpar = sigma),
  323. prior(normal(0, 0.2), class = b, dpar = sigma),
  324. prior(exponential(1), class = sd, dpar = sigma),
  326. prior(lkj(2), class = cor)),
  327. iter = 5000, warmup = 2000, chains = 4, cores = 4,
  328. seed = 14,
  329. sample_prior = TRUE,
  330. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_moremvmt_aMTL_logFD")
  333. b1.10da<-
  334. brm(data = dT_fb,
  335. family = gaussian,
  336. bf(Feedback_NAcc ~ 1 + (1 |i| ID),
  337. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  338. prior = c(prior(normal(0, 0.2), class = Intercept),
  340. prior(normal(0, 1), class = Intercept, dpar = sigma),
  341. prior(normal(0, 0.2), class = b, dpar = sigma),
  342. prior(exponential(1), class = sd, dpar = sigma),
  344. prior(lkj(2), class = cor)),
  345. iter = 5000, warmup = 2000, chains = 4, cores = 4,
  346. seed = 14,
  347. sample_prior = TRUE,
  348. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Feedback_NAcc_moremvmt_aMTL")
  351. b1.10d<-
  352. brm(data = dT_fb,
  353. family = gaussian,
  354. bf(Feedback_NAcc ~ 1 + (1 |i| ID),
  355. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + std_Cue_aMTL + log_FD + (1 + fMoney + cRewM |i| ID)),
  356. prior = c(prior(normal(0, 0.2), class = Intercept),
  358. prior(normal(0, 1), class = Intercept, dpar = sigma),
  359. prior(normal(0, 0.2), class = b, dpar = sigma),
  360. prior(exponential(1), class = sd, dpar = sigma),
  362. prior(lkj(2), class = cor)),
  363. iter = 5000, warmup = 2000, chains = 4, cores = 4,
  364. seed = 14,
  365. sample_prior = TRUE,
  366. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Feedback_NAcc_moremvmt_aMTL_logFD")
  372. dT_Unc <- filter(dT, Uncertainty ==1)
  373. dT_Cer <- filter(dT, Uncertainty ==0)
  375. b1.2a <-
  376. brm(data = dT_Unc,
  377. family = gaussian,
  378. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  379. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  380. prior = c(prior(normal(0, 0.2), class = Intercept),
  382. prior(normal(0, 1), class = Intercept, dpar = sigma),
  383. prior(normal(0, 0.2), class = b, dpar = sigma),
  384. prior(exponential(1), class = sd, dpar = sigma),
  386. prior(lkj(2), class = cor)),
  387. iter = 4000, warmup = 800, chains = 4, cores = 4,
  388. seed = 14,
  389. sample_prior = TRUE,
  390. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_Unc")
  394. b1.2a <-
  395. brm(data = dT_Unc,
  396. family = gaussian,
  397. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  398. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  399. prior = c(prior(normal(0, 0.2), class = Intercept),
  401. prior(normal(0, 1), class = Intercept, dpar = sigma),
  402. prior(normal(0, 0.2), class = b, dpar = sigma),
  403. prior(exponential(1), class = sd, dpar = sigma),
  405. prior(lkj(2), class = cor)),
  406. iter = 4000, warmup = 800, chains = 4, cores = 4,
  407. seed = 14,
  408. sample_prior = TRUE,
  409. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_Unc_aMTL")
  414. b1.3 <-
  415. brm(data = dT_Cer,
  416. family = gaussian,
  417. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  418. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + scale(SD_resRelEff) + scale(SD_resR_Want) + (1 + fMoney + cRewM |i| ID)),
  419. prior = c(prior(normal(0, 0.2), class = Intercept),
  421. prior(normal(0, 1), class = Intercept, dpar = sigma),
  422. prior(normal(0, 0.2), class = b, dpar = sigma),
  423. prior(exponential(1), class = sd, dpar = sigma),
  425. prior(lkj(2), class = cor)),
  426. iter = 4000, warmup = 800, chains = 4, cores = 4,
  427. seed = 14,
  428. sample_prior = TRUE,
  429. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_Cer_SDz")
  431. b1.7 <-
  432. brm(data = dT_fb,
  433. family = gaussian,
  434. bf(Feedback_NAcc ~ 1 + (1 |i| ID),
  435. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  436. prior = c(prior(normal(0, 0.2), class = Intercept),
  438. prior(normal(0, 1), class = Intercept, dpar = sigma),
  439. prior(normal(0, 0.2), class = b, dpar = sigma),
  440. prior(exponential(1), class = sd, dpar = sigma),
  442. prior(lkj(2), class = cor)),
  443. iter = 4000, warmup = 800, chains = 4, cores = 4,
  444. seed = 14,
  445. sample_prior = TRUE,
  446. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Feedback_NAcc_fDis")
  449. b1.7a <-
  450. brm(data = dT_fb,
  451. family = gaussian,
  452. bf(Feedback_NAcc ~ 1 + (1 |i| ID),
  453. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM + std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  454. prior = c(prior(normal(0, 0.2), class = Intercept),
  456. prior(normal(0, 1), class = Intercept, dpar = sigma),
  457. prior(normal(0, 0.2), class = b, dpar = sigma),
  458. prior(exponential(1), class = sd, dpar = sigma),
  460. prior(lkj(2), class = cor)),
  461. iter = 4000, warmup = 800, chains = 4, cores = 4,
  462. seed = 14,
  463. sample_prior = TRUE,
  464. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Feedback_NAcc_fDis_aMTL")
  466. b1.7b <-
  467. brm(data = dT_fb,
  468. family = gaussian,
  469. bf(Feedback_NAcc ~ 1 + (1 |i| ID),
  470. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM + std_Cue_aMTL + logFD + (1 + fMoney + cRewM |i| ID)),
  471. prior = c(prior(normal(0, 0.2), class = Intercept),
  473. prior(normal(0, 1), class = Intercept, dpar = sigma),
  474. prior(normal(0, 0.2), class = b, dpar = sigma),
  475. prior(exponential(1), class = sd, dpar = sigma),
  477. prior(lkj(2), class = cor)),
  478. iter = 4000, warmup = 800, chains = 4, cores = 4,
  479. seed = 14,
  480. sample_prior = TRUE,
  481. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Feedback_NAcc_fDis_aMTL_logFD")
  484. b1.8 <-
  485. brm(data = dT_Unc,
  486. family = gaussian,
  487. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  488. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM +std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  489. prior = c(prior(normal(0, 0.2), class = Intercept),
  491. prior(normal(0, 1), class = Intercept, dpar = sigma),
  492. prior(normal(0, 0.2), class = b, dpar = sigma),
  493. prior(exponential(1), class = sd, dpar = sigma),
  495. prior(lkj(2), class = cor)),
  496. iter = 4000, warmup = 800, chains = 4, cores = 4,
  497. seed = 14,
  498. sample_prior = TRUE,
  499. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_fDis_Unc_aMTL")
  501. print(summary(b1.8), digits = 3)
  503. b1.9 <-
  504. brm(data = dT_Cer,
  505. family = gaussian,
  506. bf(Cue_NAcc ~ 1 + (1 |i| ID),
  507. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM +std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  508. prior = c(prior(normal(0, 0.2), class = Intercept),
  510. prior(normal(0, 1), class = Intercept, dpar = sigma),
  511. prior(normal(0, 0.2), class = b, dpar = sigma),
  512. prior(exponential(1), class = sd, dpar = sigma),
  514. prior(lkj(2), class = cor)),
  515. iter = 4000, warmup = 800, chains = 4, cores = 4,
  516. seed = 14,
  517. sample_prior = TRUE,
  518. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_NAcc_fDis_Cer_aMTL")
  520. print(summary(b1.9), digits = 3)
  522. b2.0 <-
  523. brm(data = dT,
  524. family = gaussian,
  525. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  526. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + log_FD + (1 + fMoney + cRewM |i| ID)),
  527. prior = c(prior(normal(0, 0.2), class = Intercept),
  529. prior(normal(0, 1), class = Intercept, dpar = sigma),
  530. prior(normal(0, 0.2), class = b, dpar = sigma),
  531. prior(exponential(1), class = sd, dpar = sigma),
  533. prior(lkj(2), class = cor)),
  534. iter = 4000, warmup = 800, chains = 4, cores = 4,
  535. seed = 14,
  536. sample_prior = TRUE,
  537. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_aMTL")
  539. b2.0a <-
  540. brm(data = dT,
  541. family = gaussian,
  542. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  543. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + lstd_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  544. prior = c(prior(normal(0, 0.2), class = Intercept),
  546. prior(normal(0, 1), class = Intercept, dpar = sigma),
  547. prior(normal(0, 0.2), class = b, dpar = sigma),
  548. prior(exponential(1), class = sd, dpar = sigma),
  550. prior(lkj(2), class = cor)),
  551. iter = 4000, warmup = 800, chains = 4, cores = 4,
  552. seed = 14,
  553. sample_prior = TRUE,
  554. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_aMTL")
  556. b2.0b <-
  557. brm(data = dT,
  558. family = gaussian,
  559. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  560. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + log_FD +std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  561. prior = c(prior(normal(0, 0.2), class = Intercept),
  563. prior(normal(0, 1), class = Intercept, dpar = sigma),
  564. prior(normal(0, 0.2), class = b, dpar = sigma),
  565. prior(exponential(1), class = sd, dpar = sigma),
  567. prior(lkj(2), class = cor)),
  568. iter = 4000, warmup = 800, chains = 4, cores = 4,
  569. seed = 14,
  570. sample_prior = TRUE,
  571. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_aMTL_logFD")
  573. print(summary(b2.0), digits = 3)
  576. b2.1 <-
  577. brm(data = dT,
  578. family = gaussian,
  579. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  580. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + scale(SD_resRelEff) + scale(SD_resR_Want) + (1 + fMoney + cRewM |i| ID)),
  581. prior = c(prior(normal(0, 0.2), class = Intercept),
  583. prior(normal(0, 1), class = Intercept, dpar = sigma),
  584. prior(normal(0, 0.2), class = b, dpar = sigma),
  585. prior(exponential(1), class = sd, dpar = sigma),
  587. prior(lkj(2), class = cor)),
  588. iter = 4000, warmup = 800, chains = 4, cores = 4,
  589. seed = 14,
  590. sample_prior = TRUE,
  591. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_SDz")
  593. print(summary(b2.1), digits = 3)
  595. b2.2 <-
  596. brm(data = dT_Unc,
  597. family = gaussian,
  598. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  599. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  600. prior = c(prior(normal(0, 0.2), class = Intercept),
  602. prior(normal(0, 1), class = Intercept, dpar = sigma),
  603. prior(normal(0, 0.2), class = b, dpar = sigma),
  604. prior(exponential(1), class = sd, dpar = sigma),
  606. prior(lkj(2), class = cor)),
  607. iter = 4000, warmup = 800, chains = 4, cores = 4,
  608. seed = 14,
  609. sample_prior = TRUE,
  610. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_Unc")
  612. print(summary(b2.2), digits = 3)
  614. b2.3 <-
  615. brm(data = dT_Cer,
  616. family = gaussian,
  617. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  618. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  619. prior = c(prior(normal(0, 0.2), class = Intercept),
  621. prior(normal(0, 1), class = Intercept, dpar = sigma),
  622. prior(normal(0, 0.2), class = b, dpar = sigma),
  623. prior(exponential(1), class = sd, dpar = sigma),
  625. prior(lkj(2), class = cor)),
  626. iter = 4000, warmup = 800, chains = 4, cores = 4,
  627. seed = 14,
  628. sample_prior = TRUE,
  629. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_Cer")
  631. print(summary(b2.3), digits = 3)
  633. b2.7 <-
  634. brm(data = dT,
  635. family = gaussian,
  636. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  637. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  638. prior = c(prior(normal(0, 0.2), class = Intercept),
  640. prior(normal(0, 1), class = Intercept, dpar = sigma),
  641. prior(normal(0, 0.2), class = b, dpar = sigma),
  642. prior(exponential(1), class = sd, dpar = sigma),
  644. prior(lkj(2), class = cor)),
  645. iter = 4000, warmup = 800, chains = 4, cores = 4,
  646. seed = 14,
  647. sample_prior = TRUE,
  648. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_fDisL")
  650. b2.7a <-
  651. brm(data = dT,
  652. family = gaussian,
  653. bf(Cue_DLPFC ~ 1 + (1 |i| ID),
  654. sigma ~ 1 + (fTFEQ_Dis + zBMI) * fMoney * cRewM +std_Cue_aMTL + (1 + fMoney + cRewM |i| ID)),
  655. prior = c(prior(normal(0, 0.2), class = Intercept),
  657. prior(normal(0, 1), class = Intercept, dpar = sigma),
  658. prior(normal(0, 0.2), class = b, dpar = sigma),
  659. prior(exponential(1), class = sd, dpar = sigma),
  661. prior(lkj(2), class = cor)),
  662. iter = 4000, warmup = 800, chains = 4, cores = 4,
  663. seed = 14,
  664. sample_prior = TRUE,
  665. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_DLPFC_fDis_aMTL")
  667. print(summary(b2.7), digits = 3)
  669. b3.1 <-
  670. brm(data = dT,
  671. family = gaussian,
  672. bf(Cue_toMTL ~ 1 + (1 |i| ID),
  673. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  674. prior = c(prior(normal(0, 0.2), class = Intercept),
  676. prior(normal(0, 1), class = Intercept, dpar = sigma),
  677. prior(normal(0, 0.2), class = b, dpar = sigma),
  678. prior(exponential(1), class = sd, dpar = sigma),
  680. prior(lkj(2), class = cor)),
  681. iter = 4000, warmup = 800, chains = 4, cores = 4,
  682. seed = 14,
  683. sample_prior = TRUE,
  684. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_toMTL")
  686. print(summary(b3.1), digits = 3)
  688. b4.1 <-
  689. brm(data = dT,
  690. family = gaussian,
  691. bf(Cue_aMTL ~ 1 + (1 |i| ID),
  692. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  693. prior = c(prior(normal(0, 0.2), class = Intercept),
  695. prior(normal(0, 1), class = Intercept, dpar = sigma),
  696. prior(normal(0, 0.2), class = b, dpar = sigma),
  697. prior(exponential(1), class = sd, dpar = sigma),
  699. prior(lkj(2), class = cor)),
  700. iter = 4000, warmup = 800, chains = 4, cores = 4,
  701. seed = 14,
  702. sample_prior = TRUE,
  703. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_aMTL")
  705. print(summary(b4.1), digits = 3)
  707. b5.1 <-
  708. brm(data = dT,
  709. family = gaussian,
  710. bf(Cue_pMTL ~ 1 + (1 |i| ID),
  711. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  712. prior = c(prior(normal(0, 0.2), class = Intercept),
  714. prior(normal(0, 1), class = Intercept, dpar = sigma),
  715. prior(normal(0, 0.2), class = b, dpar = sigma),
  716. prior(exponential(1), class = sd, dpar = sigma),
  718. prior(lkj(2), class = cor)),
  719. iter = 4000, warmup = 800, chains = 4, cores = 4,
  720. seed = 14,
  721. sample_prior = TRUE,
  722. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_pMTL")
  724. print(summary(b5.1), digits = 3)
  727. b3.2 <-
  728. brm(data = dT_Cer,
  729. family = gaussian,
  730. bf(Cue_toMTL ~ 1 + (1 |i| ID),
  731. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  732. prior = c(prior(normal(0, 0.2), class = Intercept),
  734. prior(normal(0, 1), class = Intercept, dpar = sigma),
  735. prior(normal(0, 0.2), class = b, dpar = sigma),
  736. prior(exponential(1), class = sd, dpar = sigma),
  738. prior(lkj(2), class = cor)),
  739. iter = 4000, warmup = 800, chains = 4, cores = 4,
  740. seed = 14,
  741. sample_prior = TRUE,
  742. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_toMTL_Cer")
  744. print(summary(b3.2), digits = 3)
  746. b4.2 <-
  747. brm(data = dT_Cer,
  748. family = gaussian,
  749. bf(Cue_aMTL ~ 1 + (1 |i| ID),
  750. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  751. prior = c(prior(normal(0, 0.2), class = Intercept),
  753. prior(normal(0, 1), class = Intercept, dpar = sigma),
  754. prior(normal(0, 0.2), class = b, dpar = sigma),
  755. prior(exponential(1), class = sd, dpar = sigma),
  757. prior(lkj(2), class = cor)),
  758. iter = 4000, warmup = 800, chains = 4, cores = 4,
  759. seed = 14,
  760. sample_prior = TRUE,
  761. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_aMTL_Cer")
  763. print(summary(b4.2), digits = 3)
  765. b5.2 <-
  766. brm(data = dT_Cer,
  767. family = gaussian,
  768. bf(Cue_pMTL ~ 1 + (1 |i| ID),
  769. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  770. prior = c(prior(normal(0, 0.2), class = Intercept),
  772. prior(normal(0, 1), class = Intercept, dpar = sigma),
  773. prior(normal(0, 0.2), class = b, dpar = sigma),
  774. prior(exponential(1), class = sd, dpar = sigma),
  776. prior(lkj(2), class = cor)),
  777. iter = 4000, warmup = 800, chains = 4, cores = 4,
  778. seed = 14,
  779. sample_prior = TRUE,
  780. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_pMTL_Cer")
  782. print(summary(b5.2), digits = 3)
  785. b3.3 <-
  786. brm(data = dT_Unc,
  787. family = gaussian,
  788. bf(Cue_toMTL ~ 1 + (1 |i| ID),
  789. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  790. prior = c(prior(normal(0, 0.2), class = Intercept),
  792. prior(normal(0, 1), class = Intercept, dpar = sigma),
  793. prior(normal(0, 0.2), class = b, dpar = sigma),
  794. prior(exponential(1), class = sd, dpar = sigma),
  796. prior(lkj(2), class = cor)),
  797. iter = 4000, warmup = 800, chains = 4, cores = 4,
  798. seed = 14,
  799. sample_prior = TRUE,
  800. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_toMTL_Unc")
  802. print(summary(b3.3), digits = 3)
  804. b4.3 <-
  805. brm(data = dT_Unc,
  806. family = gaussian,
  807. bf(Cue_aMTL ~ 1 + (1 |i| ID),
  808. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  809. prior = c(prior(normal(0, 0.2), class = Intercept),
  811. prior(normal(0, 1), class = Intercept, dpar = sigma),
  812. prior(normal(0, 0.2), class = b, dpar = sigma),
  813. prior(exponential(1), class = sd, dpar = sigma),
  815. prior(lkj(2), class = cor)),
  816. iter = 4000, warmup = 800, chains = 4, cores = 4,
  817. seed = 14,
  818. sample_prior = TRUE,
  819. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_aMTL_Unc")
  821. print(summary(b4.3), digits = 3)
  823. b5.3 <-
  824. brm(data = dT_Unc,
  825. family = gaussian,
  826. bf(Cue_pMTL ~ 1 + (1 |i| ID),
  827. sigma ~ 1 + (fBinge + zBMI) * fMoney * cRewM + (1 + fMoney + cRewM |i| ID)),
  828. prior = c(prior(normal(0, 0.2), class = Intercept),
  830. prior(normal(0, 1), class = Intercept, dpar = sigma),
  831. prior(normal(0, 0.2), class = b, dpar = sigma),
  832. prior(exponential(1), class = sd, dpar = sigma),
  834. prior(lkj(2), class = cor)),
  835. iter = 4000, warmup = 800, chains = 4, cores = 4,
  836. seed = 14,
  837. sample_prior = TRUE,
  838. file = "/mnt/big_data/SynologyDrive/Paper/drafts/TUE002_EATVAR/analysis/brms_AK/b1.Cue_pMTL_Unc")
  840. print(summary(b5.3), digits = 3)
  843. #Behavioral data analysis variance ratio tests
  846. dEAT_S1U <- dT_S1_Unc
  849. n_boot <- 1000
  850. unique_ids <- unique(dEAT_S1U$ID)
  851. conditions <- unique(dEAT_S1U$Money)
  853. group_id_list <- dEAT_S1U %>%
  854. distinct(ID, fBinge) %>%
  855. group_split(fBinge)
  856. #Wanting
  857. # Step 1: Compute observed variance ratios (Group C vs A)
  858. observed_variances <- dEAT_S1U %>%
  859. group_by(fBinge, Money) %>%
  860. summarise(Variance = var(resR_Want), .groups = "drop")
  862. observed_ratios <- observed_variances %>%
  863. filter(fBinge %in% c("No BE", "subBED")) %>%
  864. pivot_wider(names_from = fBinge, values_from = Variance) %>%
  865. mutate(VarRatio_CvsA = subBED/`No BE`)
  867. observed_combined_variances <- dEAT_S1U %>%
  868. filter(fBinge %in% c("No BE", "subBED")) %>%
  869. group_by(fBinge) %>%
  870. summarise(Variance = var(resR_Want), .groups = "drop")
  872. observed_combined_ratio <- observed_combined_variances %>%
  873. pivot_wider(names_from = fBinge, values_from = Variance) %>%
  874. mutate(VarRatio_CvsA = subBED / `No BE`) %>%
  875. pull(VarRatio_CvsA)
  877. # Step 2: Bootstrapping
  878. #bootstrap_var_results <- vector("list", n_boot)
  879. variance_ratios <- data.frame()
  881. set.seed(123)
  883. for (i in 1:n_boot) {
  884. sampled_ids <- group_id_list %>%
  885. lapply(function(group_df) {
  886. sample(group_df$ID, size = nrow(group_df), replace = TRUE)
  887. }) %>%
  888. unlist()
  890. boot_data <- dEAT_S1U %>%
  891. semi_join(data.frame(ID = sampled_ids), by = "ID")
  893. var_summary <- boot_data %>%
  894. group_by(fBinge, Money) %>%
  895. summarise(Variance = var(resR_Want), .groups = "drop")
  897. # Calculate ratios for each Money
  898. for (cond in conditions) {
  899. var_C <- var_summary %>% filter(fBinge == "subBED", Money == cond) %>% pull(Variance)
  900. var_A <- var_summary %>% filter(fBinge == "No BE", Money == cond) %>% pull(Variance)
  902. var_ratio <- if (length(var_C) == 1 && length(var_A) == 1 && var_A != 0) var_C / var_A else NA
  904. variance_ratios <- rbind(variance_ratios, data.frame(
  905. Bootstrap = i,
  906. Money = cond,
  907. VarRatio_CvsA = var_ratio
  908. ))
  911. }
  913. boot_combined_var <- boot_data %>%
  914. filter(fBinge %in% c("No BE", "subBED")) %>%
  915. group_by(fBinge) %>%
  916. summarise(Variance = var(resR_Want), .groups = "drop")
  918. if (all(c("No BE", "subBED") %in% boot_combined_var$fBinge)) {
  919. var_C <- boot_combined_var %>% filter(fBinge == "subBED") %>% pull(Variance)
  920. var_A <- boot_combined_var %>% filter(fBinge == "No BE") %>% pull(Variance)
  921. var_ratio_combined <- ifelse(var_A != 0, var_C / var_A, NA)
  922. } else {
  923. var_ratio_combined <- NA
  924. }
  926. # Store combined ratio
  927. variance_ratios <- rbind(variance_ratios, data.frame(
  928. Bootstrap = i,
  929. Money = "Combined",
  930. VarRatio_CvsA = var_ratio_combined
  931. ))
  932. }
  934. # Step 3: Significance testing (2-sided p-values)
  936. # Function to compute 2-sided p-value
  937. two_sided_p <- function(boot_dist, obs_val) {
  938. mean(boot_dist <= 1, na.rm = TRUE)
  939. }
  942. # Individual condition p-values
  943. # Final p-value computation
  944. p_values <- variance_ratios %>%
  945. group_by(Money) %>%
  946. mutate(
  947. Observed = case_when(
  948. Money == "Combined" ~ observed_combined_ratio,
  949. TRUE ~ observed_ratios$VarRatio_CvsA[match(Money, observed_ratios$Money)]
  950. ),
  951. P_value = two_sided_p(VarRatio_CvsA, Observed)
  952. ) %>%
  953. select(Money,P_value) %>%
  954. unique()
  956. print(p_values)
  958. #effort
  961. n_boot <- 1000
  962. unique_ids <- unique(dEAT_S1U$ID)
  963. conditions <- unique(dEAT_S1U$Money)
  965. group_id_list <- dEAT_S1U %>%
  966. distinct(ID, fBinge) %>%
  967. group_split(fBinge)
  969. # Step 1: Compute observed variance ratios (Group C vs A)
  970. observed_variances <- dEAT_S1U %>%
  971. group_by(fBinge, Money) %>%
  972. summarise(Variance = var(resRelEff), .groups = "drop")
  974. observed_ratios <- observed_variances %>%
  975. filter(fBinge %in% c("No BE", "subBED")) %>%
  976. pivot_wider(names_from = fBinge, values_from = Variance) %>%
  977. mutate(VarRatio_CvsA = subBED/`No BE`)
  979. observed_combined_variances <- dEAT_S1U %>%
  980. filter(fBinge %in% c("No BE", "subBED")) %>%
  981. group_by(fBinge) %>%
  982. summarise(Variance = var(resRelEff), .groups = "drop")
  984. observed_combined_ratio <- observed_combined_variances %>%
  985. pivot_wider(names_from = fBinge, values_from = Variance) %>%
  986. mutate(VarRatio_CvsA = subBED / `No BE`) %>%
  987. pull(VarRatio_CvsA)
  989. # Step 2: Bootstrapping
  990. #bootstrap_var_results <- vector("list", n_boot)
  991. variance_ratios <- data.frame()
  993. set.seed(108)
  995. for (i in 1:n_boot) {
  996. sampled_ids <- group_id_list %>%
  997. lapply(function(group_df) {
  998. sample(group_df$ID, size = nrow(group_df), replace = TRUE)
  999. }) %>%
  1000. unlist()
  1002. boot_data <- dEAT_S1U %>%
  1003. semi_join(data.frame(ID = sampled_ids), by = "ID")
  1005. var_summary <- boot_data %>%
  1006. group_by(fBinge, Money) %>%
  1007. summarise(Variance = var(resRelEff), .groups = "drop")
  1009. # Calculate ratios for each Money
  1010. for (cond in conditions) {
  1011. var_C <- var_summary %>% filter(fBinge == "subBED", Money == cond) %>% pull(Variance)
  1012. var_A <- var_summary %>% filter(fBinge == "No BE", Money == cond) %>% pull(Variance)
  1014. var_ratio <- if (length(var_C) == 1 && length(var_A) == 1 && var_A != 0) var_C / var_A else NA
  1016. variance_ratios <- rbind(variance_ratios, data.frame(
  1017. Bootstrap = i,
  1018. Money = cond,
  1019. VarRatio_CvsA = var_ratio
  1020. ))
  1023. }
  1025. boot_combined_var <- boot_data %>%
  1026. filter(fBinge %in% c("No BE", "subBED")) %>%
  1027. group_by(fBinge) %>%
  1028. summarise(Variance = var(resRelEff), .groups = "drop")
  1030. if (all(c("No BE", "subBED") %in% boot_combined_var$fBinge)) {
  1031. var_C <- boot_combined_var %>% filter(fBinge == "subBED") %>% pull(Variance)
  1032. var_A <- boot_combined_var %>% filter(fBinge == "No BE") %>% pull(Variance)
  1033. var_ratio_combined <- ifelse(var_A != 0, var_C / var_A, NA)
  1034. } else {
  1035. var_ratio_combined <- NA
  1036. }
  1038. # Store combined ratio
  1039. variance_ratios <- rbind(variance_ratios, data.frame(
  1040. Bootstrap = i,
  1041. Money = "Combined",
  1042. VarRatio_CvsA = var_ratio_combined
  1043. ))
  1044. }
  1046. # Step 3: Significance testing (2-sided p-values)
  1048. # Function to compute 2-sided p-value
  1049. two_sided_p <- function(boot_dist, obs_val) {
  1050. mean(boot_dist <= 1, na.rm = TRUE)
  1051. }
  1053. two_sided_p <- function(boot_dist, obs_val) {
  1054. mean(abs(boot_dist - obs) >= abs(obs_val - 1), na.rm = TRUE)
  1055. }
  1057. # Individual condition p-values
  1058. # Final p-value computation
  1059. p_values <- variance_ratios %>%
  1060. group_by(Money) %>%
  1061. mutate(
  1062. Observed = case_when(
  1063. Money == "Combined" ~ observed_combined_ratio,
  1064. TRUE ~ observed_ratios$VarRatio_CvsA[match(Money, observed_ratios$Money)]
  1065. ),
  1066. P_value = two_sided_p(VarRatio_CvsA, Observed)
  1067. ) %>%
  1068. select(Money,P_value) %>%
  1069. unique()
  1071. print(p_values)
  1073. #bootstrap correlation with BMI
  1075. correlation_fn <- function(data, indices) {
  1076. d <- data[indices, ]
  1077. return(cor(d$cBMI, d$SD_resRelEff, method = "pearson"))
  1078. }
  1080. # Run bootstrap with 1000 resamples
  1081. results_resRelEff <- boot(data = dAgg_S1U, statistic = correlation_fn, R = 1000)
  1083. # View results
  1084. print(results_resRelEff)
  1086. # 95% Confidence interval
  1087. boot.ci(results_resRelEff, type = "perc")
  1089. correlation_fn <- function(data, indices) {
  1090. d <- data[indices, ]
  1091. return(cor(d$cBMI, d$SD_resR_Want, method = "pearson"))
  1092. }
  1094. # Run bootstrap with 1000 resamples
  1095. results_resRWant <- boot(data = dAgg_S1U, statistic = correlation_fn, R = 1000)
  1097. # View results
  1098. print(results_resRWant)
  1100. # 95% Confidence interval
  1101. boot.ci(results_resRWant, type = "perc")

TUE002_EAT_trial_based_plots_share.R, no license · at the source

Overview

Authors: Mechteld M. van den Hoek Ostende1,2, Anne Kühnel3,4, Monja P. Neuser1, Thomas Dresler1,4,5, Jennifer Svaldi2,4, Nils B. Kroemer1,3,4,6
ORCID iDs: Anne Kühnel, Thomas Dresler, Nils B. Kroemer
  1. Department of Psychiatry and Psychotherapy, Tübingen Center for Mental Health, University of Tübingen,Tübingen, Germany
  2. Department of Psychology, Tübingen Center for Mental Health, University of Tübingen,Tübingen, Germany
  3. Section of Medical Psychology, Department of Psychiatry and Psychotherapy, University Hospital Bonn, University of Bonn,Bonn, Germany
  4. German Center for Mental Health (DZPG), partner site Tübingen,Tübingen, Germany
  5. LEAD Graduate School & Research Network, University of Tübingen,Tübingen, Germany
  6. German Center for Diabetes Research (DZD),Neuherberg, Germany
Institutions: University of Tübingen (Germany); University Hospital Bonn (Germany); Deutsches Zentrum für Psychische Gesundheit (Germany); German Center for Diabetes Research (Germany)
Journal: Translational psychiatry, volume 16, issue 1, article 370
Dates: received 4 November 2025; accepted 8 June 2026; published online 20 July 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1038/s41398-026-04172-6 · PMID 42477311 · PMCID PMC13385738 · OpenAlex W4414077146
Open access: gold, a free copy (OpenAlex)
Preprint: doi:10.1101/2025.09.04.25335092
Status: code verified
Categories: human (organism), other condition (population), clinical / translational (subfield)
Methods: Connectivity, Statistics, fMRI & imaging
Keywords: Human behaviour, Neuroscience, Psychiatric disorders, Diagnostic markers
MeSH: Binge-Eating Disorder*, Nucleus Accumbens*, Obesity*, Reward*, Adult, Cues, Female, Humans, Magnetic Resonance Imaging, Male (* major topic)
Topic: Adipose Tissue and Metabolism (Physiology, Medicine), according to OpenAlex
Funding: Deutsche Forschungsgemeinschaft (DFG) (KR 4555/7-1, KR 4555/9-1, KR4555/10_1)
Citations: cited by 1 paper (Europe PMC); 84 references in the paper

Abstract

Binge eating disorder (BED) is characterized by repeated episodes of binge eating accompanied by a loss of control. Although the neurobiological underpinnings of binge eating (BE) episodes are not fully understood, there are indications that variability in nucleus accumbens (NAcc) responses could lead to increased variability in food intake. Here, we assessed whether BED is associated with higher intra-individual variability in behavioral and neuroimaging indices of reward responses. To this end, patients with BED (n = 35, MBMI = 33.2 kg/m2 ± 6.8), participants with subsyndromal BED (n = 21, MBMI = 29.0 kg/m2 ± 7.2), and individuals without symptoms of binge eating (n = 23, MBMI = 32.3 kg/m2 ± 6.5) completed an effort allocation task with concurrent functional magnetic resonance imaging. In line with our hypothesis, we found that patients with BED had higher variability in subjective wanting ratings of food (F34,21 = 1.48, pboot = 0.024), but not effort exertion (F34,21 = 1.13, pboot = 0.30). Crucially, trial-by-trial variability in NAcc responses during the presentation of cues was associated with a higher BMI (b = 0.11, 95%CI [0.03, 0.19], BF10 = 11.1) and disinhibited eating (b = 0.19, 95%CI [0.01, 0.36], BF10 = 4.0) across groups, whereas NAcc variability was only marginally elevated in patients with BED (b = 0.12, 95%CI [−0.04, 0.29], BF10 = 1.2, P > 0|data = 88%). Our results support the idea that BMI and disinhibited eating are associated with more variable NAcc responses, which may contribute to the symptoms of BED. However, this association is only weakly indicative of clinical severity of BED.

Reproduced under the paper's license (CC BY), from the paper cited above.

Repository

Its files are read in the Code ↔ Paper reader above, with 1 match between paragraphs and lines of code.

OSF tewpn

License: none: the authors keep all their rights
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Languages: R (1)
Size: 5 files, 1 script
Software Heritage: not checked
Found in: “Code availability”
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: brms (1 file), cowplot (1 file), ggplot2 (1 file), tidyverse (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers (HTTP 200)
  • 27 September 2026: the link answers (HTTP 200)
1 file
At the source: osf.io/tewpn/?view_only=e35bb7a5f9574cc4b8faa54ed2979d01

Code availability

Trial-wise data (estimated cue responses and behavior) and analysis code are available at https://osf.io/tewpn/?view_only=e35bb7a5f9574cc4b8faa54ed2979d01.

Reproduced under the paper's license (CC BY), from the paper cited above.

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;
  • 1 match between paragraphs of the paper and lines of the code (method lexical-v1);
  • 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

Trial-wise data (estimated cue responses and behavior) is available at https://osf.io/tewpn/?view_only=e35bb7a5f9574cc4b8faa54ed2979d01.

Reproduced under the paper's license (CC BY), from the paper cited above.

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, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 6 authors, 4 keywords, 10 MeSH terms, 1 funder, 80 references.

Cite

This paper

van den Hoek Ostende, M. M., Kühnel, A., Neuser, M. P., Dresler, T., Svaldi, J., & Kroemer, N. B. (2026). Obesity is associated with greater variability of reward signals in the nucleus accumbens. Translational psychiatry, 16(1), 370. https://doi.org/10.1038/s41398-026-04172-6

BibTeX

@article{vandenhoekostende2026obesity,
author = {van den Hoek Ostende, Mechteld M. and Kühnel, Anne and Neuser, Monja P. and Dresler, Thomas and Svaldi, Jennifer and Kroemer, Nils B.},
title = {{Obesity is associated with greater variability of reward signals in the nucleus accumbens}},
journal = {Translational psychiatry},
year = {2026},
month = jul,
volume = {16},
number = {1},
pages = {370},
publisher = {Nature Publishing Group},
issn = {2158-3188},
doi = {10.1038/s41398-026-04172-6},
url = {https://doi.org/10.1038/s41398-026-04172-6},
pmid = {42477311},
pmcid = {PMC13385738}
}

RIS

TY - JOUR
AU - van den Hoek Ostende, Mechteld M.
AU - Kühnel, Anne
AU - Neuser, Monja P.
AU - Dresler, Thomas
AU - Svaldi, Jennifer
AU - Kroemer, Nils B.
TI - Obesity is associated with greater variability of reward signals in the nucleus accumbens
T2 - Translational psychiatry
J2 - Transl Psychiatry
PY - 2026
DA - 2026/07/20
VL - 16
IS - 1
SP - 370
SN - 2158-3188
PB - Nature Publishing Group
DO - 10.1038/s41398-026-04172-6
UR - https://doi.org/10.1038/s41398-026-04172-6
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s41398-026-04172-6",
"type": "article-journal",
"title": "Obesity is associated with greater variability of reward signals in the nucleus accumbens",
"container-title": "Translational psychiatry",
"author": [
{
"family": "van den Hoek Ostende",
"given": "Mechteld M."
},
{
"family": "Kühnel",
"given": "Anne"
},
{
"family": "Neuser",
"given": "Monja P."
},
{
"family": "Dresler",
"given": "Thomas"
},
{
"family": "Svaldi",
"given": "Jennifer"
},
{
"family": "Kroemer",
"given": "Nils B."
}
],
"container-title-short": "Transl Psychiatry",
"volume": "16",
"issue": "1",
"page": "370",
"DOI": "10.1038/s41398-026-04172-6",
"PMID": "42477311",
"PMCID": "PMC13385738",
"ISSN": "2158-3188",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s41398-026-04172-6",
"language": "en",
"issued": {
"date-parts": [
[
2026,
7,
20
]
]
}
}

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.1111/infa.70114 [code]
Statistics in Motion: Does the Infant Motor System Predict Actions Based on Their Transitional Probability?
Journal: Infancy : the official journal of the International Society on Infant Studies
In common: brms, cowplot, tidyverse, 3 references
[2] doi:10.1093/cercor/bhag032 [code]
Pulvinar-posterior superior temporal sulcus connectivity contributes to non-conscious emotion processing in affective blindsight.
Journal: Cerebral cortex (New York, N.Y. : 1991)
In common: brms, cowplot, ggplot2, 1 other tool, other condition, 1 reference
[3] doi:10.1162/imag.a.1258 [code]
Non-specific increase in alpha power during a neurofeedback session targeting its downregulation.
Journal: Imaging neuroscience (Cambridge, Mass.)
In common: brms, cowplot, ggplot2, 1 other tool, 1 reference
[4] doi:10.1038/s44271-026-00431-w [code]
Alpha power increases spontaneously during a neurofeedback session.
Journal: Communications psychology
In common: brms, cowplot, ggplot2, 1 other tool, 1 reference
[5] doi:10.1038/s41598-026-53424-4 [code]
Cognitive control networks causally support implicit emotion regulation: evidence from dlPFC stimulation and directed functional connectivity.
Journal: Scientific reports
In common: brms, cowplot, ggplot2, 1 other tool, 1 reference
[6] doi:10.3758/s13415-026-01417-1 [code]
Computational signatures of exertion and rest underlie moment-to-moment dynamics of subjective perceptions of effort and fatigue.
Journal: Cognitive, affective & behavioral neuroscience
In common: 4 references
[7] doi:10.1162/nol.a.266 [code]
The Role of the Anterior Temporal Lobe in Reading: An HD-tDCS Study.
Journal: Neurobiology of language (Cambridge, Mass.)
In common: brms, ggplot2, tidyverse, 2 references
[8] doi:10.1162/opmi.a.372 [code]
Broadening the Agent Preference Hypothesis Through Experiencers: Eye-Tracking and EEG Evidence of Proto-Agents and Proto-Patients.
Journal: Open mind : discoveries in cognitive science
In common: brms, ggplot2, tidyverse, clinical / translational, 1 reference
[9] doi:10.1111/psyp.70285 [code]
Spontaneous Modulation of Alpha Power During a Neurofeedback Session Without Instructions.
Journal: Psychophysiology
In common: brms, cowplot, tidyverse, 1 reference
[10] doi:10.1016/j.cell.2026.05.047 [code]
An emergent disease-associated motor neuron state precedes cell death in ALS.
Journal: Cell
In common: cowplot, ggplot2, tidyverse, other condition, 2 references

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.