OSCR

Targeted stool metabolomics suggests exploratory catecholamine- and tryptophan-linked metabolic features in autism spectrum disorder.

Code ↔ Paper

4 matches 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 4 matches
  1. [1] § Materials and methods › Statistical analysis › Machine learning ↔ statistical_analysis.Rmd, lines 953–1036 · score 0.77 · Brier score, predicted probabilities, positive class, Calibration, LOESS, smoothed
  2. [2] § Results › Discriminative performance of metabolite features in classifying ASD ↔ statistical_analysis.Rmd, lines 491–587 · score 0.62 · cross validated, positive class, random forest, predicting ASD, ROC, AUC
  3. [3] § Results › Model calibration performance ↔ statistical_analysis.Rmd, lines 953–1036 · score 0.62 · Brier score, predicted probabilities, calibration, LOESS, smoothed, curves
  4. [4] § Materials and methods › Statistical analysis › Machine learning ↔ statistical_analysis.Rmd, lines 491–587 · score 0.60 · Random forest classifiers, cross validation, seed, trees, workflow, fold

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,040 lines · 30 KB · MIT · 4 matches

  1. ---
  2. title: "Targeted ASD Stool Metabolomics Analysis"
  3. author: "Kevin Liu"
  4. date: "`r Sys.Date()`"
  5. output:
  6. pdf_document:
  7. toc: true
  8. keep_tex: true
  9. ---
  11. ```{r setup, include=FALSE, message=FALSE}
  12. knitr::opts_chunk$set(echo = TRUE)
  13. options(scipen = 999)
  15. library(tidyverse)
  16. library(egg)
  17. library(rstatix)
  18. library(tidymodels)
  19. library(yardstick)
  20. library(ranger)
  21. library(pROC)
  22. library(patchwork)
  23. library(scales)
  24. library(broom)
  25. library(ggpubr)
  27. theme_set(theme_article() +
  28. theme(aspect.ratio = 1))
  29. ```
  31. ```{r read_data, message=FALSE}
  32. clin_data = read_csv("data/clin_data.csv")
  34. avg_data = read_csv("data/data_averaged.csv") %>%
  35. left_join(clin_data, by = "subject_id") %>%
  36. relocate(subject_id, group, sex, age) %>%
  37. mutate(group = factor(group, levels = c("NAC", "ASD")),
  38. sex = factor(sex, levels = c("Female", "Male")))
  40. norm_data = read_csv("data/data_normalized.csv") %>%
  41. left_join(clin_data, by = "subject_id") %>%
  42. relocate(subject_id, group, sex, age) %>%
  43. mutate(group = factor(group, levels = c("NAC", "ASD")),
  44. sex = factor(sex, levels = c("Female", "Male")))
  46. name_map = read_csv("data/name_map.csv", na = c("", "NaN"))
  48. levels(avg_data$group)
  49. levels(norm_data$group)
  50. ```
  52. ```{r descriptive_stats}
  53. # Figure S1
  54. avg_data %>%
  55. ggplot(aes(x = age, fill = group)) +
  56. geom_histogram(data = subset(avg_data, group == "NAC"),
  57. aes(y = after_stat(count)), binwidth = 1, boundary = 0, closed = "left") +
  58. geom_histogram(data = subset(avg_data, group == "ASD"),
  59. aes(y = -after_stat(count)), binwidth = 1, boundary = 0, closed = "left") +
  60. scale_y_continuous(labels = abs, limits = c(-10, 10)) +
  61. scale_fill_manual(values = c("NAC" = "#21918c", "ASD" = "#fde725")) +
  62. labs(x = "Age (years)", y = "Number of Subjects", fill = "Group") +
  63. theme(legend.position = "bottom") +
  64. geom_hline(yintercept = 0, linewidth = 0.25) +
  65. coord_flip()
  67. age_tbl = avg_data %>%
  68. group_by(group) %>%
  69. summarise(
  70. value = sprintf(
  71. "%.2f ± %.2f, [%.2f, %.2f]",
  72. mean(age, na.rm = TRUE),
  73. sd(age, na.rm = TRUE),
  74. min(age, na.rm = TRUE),
  75. max(age, na.rm = TRUE)
  76. ),
  77. .groups = "drop"
  78. ) %>%
  79. pivot_wider(names_from = group, values_from = value) %>%
  80. mutate(
  81. Characteristic = "Age (mean ± SD, [min, max])")
  83. sex_tbl = avg_data %>%
  84. count(group, sex) %>%
  85. group_by(group) %>%
  86. mutate(value = sprintf("%d (%.1f%%)", n, 100 * n / sum(n))) %>%
  87. ungroup() %>%
  88. select(group, sex, value) %>%
  89. pivot_wider(names_from = group, values_from = value) %>%
  90. mutate(Characteristic = as.character(sex)) %>%
  91. select(-sex)
  93. wilcox.test(age ~ group, data = avg_data)
  94. fisher.test(table(avg_data$group, avg_data$sex))
  96. # Table 1
  97. bind_rows(age_tbl, sex_tbl) %>%
  98. select(Characteristic, NAC, ASD) %>%
  99. write.csv("output_data/1_demographics.csv")
  100. ```
  102. ```{r get_mol_summary}
  103. num_cols = avg_data %>%
  104. select(where(is.numeric), -subject_id, -age) %>%
  105. names()
  107. # Table 2
  108. avg_data %>%
  109. select(group, all_of(num_cols)) %>%
  110. pivot_longer(
  111. cols = -group,
  112. names_to = "variable",
  113. values_to = "value"
  114. ) %>%
  115. group_by(group, variable) %>%
  116. summarise(
  117. n = n(),
  118. n_zero = sum(value == 0, na.rm = TRUE),
  119. pct_zero = n_zero / n * 100,
  120. n_nonzero = sum(value != 0 & !is.na(value)),
  121. # stats excluding zeros
  122. mean = mean(if_else(value == 0, NA_real_, value), na.rm = TRUE),
  123. sd = sd(if_else(value == 0, NA_real_, value), na.rm = TRUE),
  124. min = min(if_else(value == 0, NA_real_, value), na.rm = TRUE),
  125. max = max(if_else(value == 0, NA_real_, value), na.rm = TRUE),
  126. median = median(if_else(value == 0, NA_real_, value), na.rm = TRUE),
  127. .groups = "drop"
  128. ) %>%
  129. ungroup() %>%
  130. mutate(across(where(is.numeric), \(x) round(x, 2))) %>%
  131. mutate(range = str_c("[", min, ", ", max, "]"),
  132. mean_sd = str_c(mean, " (", sd, ")"),
  133. zeros = str_c(n_zero, "/", n, " (", pct_zero, "%)")) %>%
  134. select(group, variable, mean_sd, range, median, zeros) %>%
  135. write_csv("output_data/2_mol_summary.csv")
  136. ```
  139. ```{r viz_data, fig.dpi=300, fig.width=8}
  140. # Figure 1
  141. avg_data %>%
  142. pivot_longer(cols = starts_with("mol_"), names_to = "molecule") %>%
  143. left_join(name_map %>% select(mol_name, compound_acronym), by = c("molecule" = "mol_name")) %>%
  144. mutate(molecule = compound_acronym) %>%
  145. select(-compound_acronym) %>%
  146. ggplot(aes(x = group, y = value, fill = group)) +
  147. geom_violin() +
  148. geom_jitter(width = 0.2, alpha = 0.5) +
  149. facet_wrap(~ molecule, scales = "free_y", ncol = 6) +
  150. labs(x = "", y = "Concentration (nM)",
  151. fill = "Group") +
  152. theme(legend.position = "bottom", aspect.ratio = 1,
  153. strip.text = element_text(size = 6)) +
  154. scale_fill_manual(values = c("#21918c", "#fde725"))
  155. ```
  157. ```{r corr_matrix, fig.width=10}
  158. num_vars = c("age", grep("^mol_", names(avg_data), value = TRUE))
  160. name_lookup = setNames(c(name_map$compound_acronym, "Age"),
  161. c(name_map$mol_name, "age"))
  163. axis_levels = num_vars
  165. corr_lower_by_group = function(df_group, group_label) {
  166. mat = df_group %>%
  167. select(all_of(num_vars)) %>%
  168. as.matrix()
  170. rc = Hmisc::rcorr(mat, type = "spearman")
  171. R = rc$r
  172. P = rc$P
  174. vars = colnames(R)
  175. idx = which(lower.tri(R), arr.ind = TRUE)
  177. tibble(
  178. Group = group_label,
  179. Var1 = vars[idx[,1]],
  180. Var2 = vars[idx[,2]],
  181. r = as.numeric(R[idx]),
  182. p = as.numeric(P[idx])
  183. )
  184. }
  186. cor_df = bind_rows(
  187. corr_lower_by_group(filter(avg_data, group == "ASD"), "ASD"),
  188. corr_lower_by_group(filter(avg_data, group == "NAC"), "NAC")
  189. ) %>%
  190. group_by(Group) %>%
  191. mutate(
  192. p_adj = p.adjust(p, method = "fdr"), # FDR within each group
  193. p_cap = pmax(p_adj, 1e-300), # avoid -Inf for -log10
  194. Var1 = factor(Var1, levels = axis_levels),
  195. Var2 = factor(Var2, levels = axis_levels),
  196. Significant = !is.na(p_adj) & p_adj < 0.05
  197. ) %>%
  198. ungroup() %>%
  199. mutate(
  200. Var1_label = recode(Var1, !!!name_lookup),
  201. Var2_label = recode(Var2, !!!name_lookup)
  202. )
  204. # Figure 5
  205. cor_df %>%
  206. ggplot(aes(x = Var1_label, y = Var2_label)) +
  207. geom_point(aes(fill = r, size = -log10(p_cap)),
  208. shape = 21, color = "grey60", stroke = 0.1, na.rm = TRUE) +
  209. geom_point(data = filter(cor_df, Significant),
  210. aes(size = -log10(p_cap)),
  211. shape = 21, fill = NA, color = "black", stroke = 0.5, na.rm = TRUE) +
  212. scale_fill_gradient2(low = "blue", mid = "white", high = "red",
  213. midpoint = 0, limits = c(-1, 1),
  214. oob = scales::squish, name = "Spearman r") +
  215. scale_size(range = c(1.5, 6), name = NULL, guide = "none") +
  216. coord_fixed() +
  217. facet_wrap(~ Group, ncol = 2) +
  218. theme_light() +
  219. theme(
  220. axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1),
  221. legend.position = "right"
  222. ) +
  223. labs(x = NULL, y = NULL)
  225. # Table S5
  226. cor_df %>%
  227. filter(p_adj < 0.05) %>%
  228. select(Group, Var1_label, Var2_label, r, p, p_adj) %>%
  229. mutate(
  230. r = signif(r, 3),
  231. p = formatC(p, format = "e", digits = 2),
  232. p_adj = formatC(p_adj, format = "e", digits = 2)
  233. ) %>%
  234. write_csv("output_data/s5_groupwise_sigif_corr.csv")
  235. ```
  237. ```{r corr_scatterplots_all_significant_patchwork, fig.dpi=300, fig.width=10, fig.height=12}
  238. # Supplementary Figure 2
  239. # Scatterplots for all FDR-significant within-group correlations
  241. cor_df_labeled = cor_df %>%
  242. select(-any_of(c("Var1_label", "Var2_label"))) %>%
  243. left_join(
  244. name_map %>% select(mol_name, compound_acronym),
  245. by = c("Var1" = "mol_name")
  246. ) %>%
  247. rename(Var1_label = compound_acronym) %>%
  248. left_join(
  249. name_map %>% select(mol_name, compound_acronym),
  250. by = c("Var2" = "mol_name")
  251. ) %>%
  252. rename(Var2_label = compound_acronym) %>%
  253. mutate(
  254. Var1 = as.character(Var1),
  255. Var2 = as.character(Var2),
  256. Var1_label = if_else(Var1 == "age", "Age", Var1_label),
  257. Var2_label = if_else(Var2 == "age", "Age", Var2_label)
  258. )
  260. panel_order = tribble(
  261. ~Group, ~Var1, ~Var2, ~order,
  262. "ASD", "mol_trp", "mol_niacin", 1,
  263. "NAC", "mol_trp", "mol_kyn", 2,
  264. "ASD", "mol_trp", "mol_gaba", 3,
  265. "ASD", "mol_trp", "mol_ne", 4,
  266. "ASD", "mol_trp", "mol_aea", 5,
  267. "NAC", "mol_trp", "mol_aea", 6,
  268. "ASD", "mol_ne", "mol_aea", 7,
  269. "ASD", "mol_niacin", "mol_gaba", 8,
  270. "ASD", "mol_3nt", "mol_ua", 9,
  271. "ASD", "mol_lca", "mol_dca", 10,
  272. "NAC", "mol_lca", "mol_dca", 11
  273. )
  275. sig_corr_stats = cor_df_labeled %>%
  276. filter(Significant == TRUE) %>%
  277. inner_join(panel_order, by = c("Group", "Var1", "Var2")) %>%
  278. arrange(order) %>%
  279. mutate(
  280. q_label = case_when(
  281. p_adj < 0.001 ~ "q < 0.001",
  282. TRUE ~ str_c("q = ", sprintf("%.3f", p_adj))
  283. ),
  284. stat_label = str_c(
  285. "Spearman r = ", sprintf("%.2f", r),
  286. "\nFDR ", q_label
  287. ),
  288. panel_title = str_c(Group, ": ", Var2_label, " vs ", Var1_label)
  289. )
  291. plot_sig_corr_pair = function(group_label, var1, var2, var1_label, var2_label,
  292. stat_label, panel_title) {
  294. plot_df = avg_data %>%
  295. filter(group == group_label) %>%
  296. transmute(
  297. Group = group,
  298. x = .data[[var2]],
  299. y = .data[[var1]]
  300. ) %>%
  301. filter(x > 0, y > 0)
  303. label_df = plot_df %>%
  304. summarise(
  305. x = 10^(min(log10(x), na.rm = TRUE) + 0.05 * diff(range(log10(x), na.rm = TRUE))),
  306. y = 10^(max(log10(y), na.rm = TRUE) - 0.05 * diff(range(log10(y), na.rm = TRUE)))
  307. ) %>%
  308. mutate(label = stat_label)
  310. ggplot(plot_df, aes(x = x, y = y)) +
  311. geom_point(
  312. aes(fill = Group),
  313. shape = 21,
  314. color = "black",
  315. stroke = 0.2,
  316. size = 1.8,
  317. alpha = 0.8
  318. ) +
  319. geom_smooth(
  320. method = "loess",
  321. se = FALSE,
  322. linewidth = 0.5,
  323. color = "black") +
  324. geom_text(
  325. data = label_df,
  326. aes(x = x, y = y, label = label),
  327. inherit.aes = FALSE,
  328. hjust = 0,
  329. vjust = 1,
  330. size = 2.5,
  331. lineheight = 0.95
  332. ) +
  333. scale_x_log10() +
  334. scale_y_log10() +
  335. scale_fill_manual(values = c("NAC" = "#21918c", "ASD" = "#fde725")) +
  336. labs(
  337. title = panel_title,
  338. x = str_c(var2_label, " (nM, log10)"),
  339. y = str_c(var1_label, " (nM, log10)"),
  340. fill = "Group"
  341. ) +
  342. theme(
  343. plot.title = element_text(size = 8),
  344. axis.title = element_text(size = 7),
  345. axis.text = element_text(size = 6),
  346. aspect.ratio = 1
  347. )
  348. }
  350. sig_corr_plots = pmap(
  351. sig_corr_stats %>%
  352. select(
  353. group_label = Group,
  354. var1 = Var1,
  355. var2 = Var2,
  356. var1_label = Var1_label,
  357. var2_label = Var2_label,
  358. stat_label,
  359. panel_title
  360. ),
  361. plot_sig_corr_pair
  362. )
  364. wrap_plots(sig_corr_plots, ncol = 3, guides = "collect") &
  365. theme(
  366. legend.position = "bottom",
  367. legend.title = element_text(size = 9),
  368. legend.text = element_text(size = 9)
  369. )
  370. ```
  372. ```{r univar_tests}
  373. # norm_data %>%
  374. # pivot_longer(cols = starts_with("mol_"), names_to = "metabolite") %>%
  375. # group_by(metabolite) %>%
  376. # wilcox_test(value ~ group) %>%
  377. # adjust_pvalue(method = "fdr") %>%
  378. # arrange(p) %>%
  379. # left_join(name_map %>% select(mol_name, compound_name), by = c("metabolite" = "mol_name"))
  380. ```
  382. ```{r adj_univar_tests}
  383. adj_univar_lm = norm_data %>%
  384. pivot_longer(cols = starts_with("mol_"), names_to = "molecule", values_to = "value") %>%
  385. group_by(molecule) %>%
  386. do(tidy(lm(value ~ group + age + sex, data = .))) %>%
  387. ungroup() %>%
  388. filter(term == "groupASD") %>%
  389. mutate(p_adj = p.adjust(p.value, method = "fdr")) %>%
  390. arrange(p.value) %>%
  391. select(molecule, estimate, std.error, statistic, p.value, p_adj) %>%
  392. left_join(name_map %>% select(mol_name, compound_acronym), by = c("molecule" = "mol_name")) %>%
  393. mutate(molecule = compound_acronym) %>%
  394. select(-compound_acronym)
  395. ```
  397. ```{r log2fc, fig.dpi=300}
  398. log2fc_result = norm_data %>%
  399. pivot_longer(cols = starts_with("mol_"), names_to = "molecule", values_to = "value") %>%
  400. group_by(group, molecule) %>%
  401. summarise(mean_value = mean(value, na.rm = TRUE), .groups = "drop") %>%
  402. pivot_wider(names_from = group, values_from = mean_value) %>%
  403. mutate(
  404. log2_fc = (ASD - NAC) * log2(10) # Convert log10 difference to log2 scale
  405. ) %>%
  406. arrange(desc(log2_fc)) %>%
  407. left_join(name_map %>% select(mol_name, compound_acronym), by = c("molecule" = "mol_name")) %>%
  408. mutate(molecule = compound_acronym) %>%
  409. select(-compound_acronym)
  411. # Figure 2
  412. log2fc_result %>%
  413. mutate(molecule = fct_rev(fct_inorder(molecule))) %>%
  414. ggplot(aes(x = molecule, y = log2_fc, fill = log2_fc > 0)) +
  415. geom_col() +
  416. coord_flip() +
  417. scale_fill_manual(values = c("TRUE" = "#fde725", "FALSE" = "#21918c"),
  418. labels = c("NAC > ASD", "ASD > NAC")) +
  419. labs(
  420. title = "Log2 Fold Change (ASD vs. NAC)",
  421. x = NULL,
  422. y = "Log2 Fold Change",
  423. fill = "Direction"
  424. ) +
  425. ylim(-2.25, 2.25) +
  426. theme(legend.position = "bottom") +
  427. geom_hline(yintercept = 1, linetype = "dashed", color = "gray") +
  428. geom_hline(yintercept = -1, linetype = "dashed", color = "gray")
  430. # Table S2
  431. log2fc_result %>%
  432. left_join(adj_univar_lm, by = "molecule") %>%
  433. mutate(across(where(is.numeric), \(x) round(x, 4))) %>%
  434. arrange(p.value) %>%
  435. left_join(name_map %>% select(compound_acronym, compound_name), by = c("molecule" = "compound_acronym")) %>%
  436. mutate(molecule = str_c(compound_name, " (", molecule, ")")) %>%
  437. select(-compound_name) %>%
  438. write_csv("output_data/s2_log2fc_adj_univar.csv")
  439. ```
  441. ```{r ne_sensitivity, message=FALSE, fig.width=5, fig.height=4}
  442. # NE sensitivity check
  444. # # NE distribution by group
  445. # ggplot(avg_data, aes(x = group, y = mol_ne, color = group)) +
  446. # geom_boxplot(outlier.shape = NA) +
  447. # geom_jitter(width = 0.1, alpha = 0.7, size = 2) +
  448. # labs(x = NULL, y = "Norepinephrine (nM)")
  450. # Top NE values
  451. ne_top = avg_data %>%
  452. arrange(desc(mol_ne)) %>%
  453. select(subject_id, group, sex, age, mol_ne) %>%
  454. slice_head(n = 10)
  456. ne_top
  458. # IQR-based outlier threshold
  459. Q1_ne = quantile(avg_data$mol_ne, 0.25, na.rm = TRUE)
  460. Q3_ne = quantile(avg_data$mol_ne, 0.75, na.rm = TRUE)
  461. IQR_ne = Q3_ne - Q1_ne
  462. upper_bound_ne = Q3_ne + 1.5 * IQR_ne
  464. upper_bound_ne
  466. avg_data_ne = avg_data %>%
  467. mutate(ne_outlier = mol_ne > upper_bound_ne)
  469. # Full-data test
  470. wilcox.test(mol_ne ~ group, data = avg_data_ne)
  472. # Outlier-removed test
  473. avg_data_ne_no_outliers = avg_data_ne %>%
  474. filter(!ne_outlier)
  476. wilcox.test(mol_ne ~ group, data = avg_data_ne_no_outliers)
  478. # Mean and median by group
  479. avg_data_ne %>%
  480. group_by(group) %>%
  481. summarise(
  482. mean_ne = mean(mol_ne, na.rm = TRUE),
  483. median_ne = median(mol_ne, na.rm = TRUE),
  484. min_ne = min(mol_ne, na.rm = TRUE),
  485. max_ne = max(mol_ne, na.rm = TRUE),
  486. .groups = "drop"
  487. )
  488. ```
  491. ```{r univar_rf}
  492. rf_base_data = avg_data %>%
  493. mutate(group = factor(group, levels = c("NAC", "ASD"))) # ASD is the positive class
  495. set.seed(123)
  496. folds = vfold_cv(rf_base_data, v = 5, strata = group)
  498. mol_names = names(rf_base_data)[str_starts(names(rf_base_data), "mol_")]
  500. get_cv_rf_auc = function(mol) {
  501. df = rf_base_data %>% select(group, !!sym(mol))
  503. rec = recipe(group ~ ., data = df) %>%
  504. step_zv(all_predictors()) # Remove zero-variance predictors (just in case)
  506. rf_spec = rand_forest(trees = 500) %>%
  507. set_mode("classification") %>%
  508. set_engine("ranger", importance = "none")
  510. wf = workflow() %>%
  511. add_recipe(rec) %>%
  512. add_model(rf_spec)
  514. res = fit_resamples(
  515. wf,
  516. resamples = folds,
  517. metrics = metric_set(roc_auc),
  518. control = control_resamples(save_pred = TRUE, verbose = FALSE)
  519. )
  521. auc_val = collect_metrics(res) %>%
  522. filter(.metric == "roc_auc") %>%
  523. pull(mean)
  525. tibble(molecule = mol, auc = auc_val)
  526. }
  528. rf_auc_all = map_dfr(mol_names, get_cv_rf_auc) %>%
  529. arrange(desc(auc))
  531. # rf_auc_all %>%
  532. # mutate(molecule = fct_reorder(molecule, auc)) %>%
  533. # ggplot(aes(x = molecule, y = auc)) +
  534. # geom_col(fill = "#21918c") +
  535. # coord_flip() +
  536. # labs(
  537. # title = "Random Forest Classification Performance by Molecule",
  538. # y = "Cross-Validated AUC",
  539. # x = NULL
  540. # ) +
  541. # geom_hline(yintercept = 0.5, linetype = "dashed", color = "gray")
  543. ind_auc_results = map_dfr(mol_names, function(mol) {
  544. df = avg_data %>% select(group, !!sym(mol))
  546. rec = recipe(group ~ ., data = df) %>%
  547. step_zv(all_predictors())
  549. wf = workflow() %>%
  550. add_recipe(rec) %>%
  551. add_model(rand_forest(trees = 500) %>%
  552. set_mode("classification") %>%
  553. set_engine("ranger"))
  555. res = fit_resamples(
  556. wf, resamples = folds,
  557. metrics = metric_set(roc_auc),
  558. control = control_resamples(save_pred = TRUE)
  559. )
  561. preds = collect_predictions(res)
  562. roc_obj = roc(response = preds$group, predictor = preds$.pred_ASD,
  563. levels = c("NAC", "ASD"), direction = "<", ci = TRUE, boot.n = 2000)
  565. tibble(
  566. metabolite = mol,
  567. auc = as.numeric(auc(roc_obj)),
  568. ci_low = as.numeric(ci.auc(roc_obj)[1]),
  569. ci_high = as.numeric(ci.auc(roc_obj)[3])
  570. )
  571. }) %>%
  572. arrange(desc(auc)) %>%
  573. mutate(
  574. auc_ci = sprintf("%.3f [%.3f-%.3f]", auc, ci_low, ci_high)
  575. ) %>%
  576. select(metabolite, auc_ci) %>%
  577. left_join(name_map %>% select(mol_name, compound_name), by = c("metabolite" = "mol_name")) %>%
  578. mutate(metabolite = compound_name) %>%
  579. select(-compound_name)
  581. # Table S3
  582. ind_auc_results %>%
  583. left_join(name_map %>% select(compound_acronym, compound_name), by = c("metabolite" = "compound_name")) %>%
  584. mutate(metabolite = str_c(metabolite, " (", compound_acronym, ")")) %>%
  585. select(-compound_acronym) %>%
  586. write_csv("output_data/s3_univar_rf_auc.csv")
  587. ```
  589. ```{r univar_roc}
  590. top_molecule = rf_auc_all %>%
  591. slice_max(auc, n = 3) %>%
  592. pull(molecule)
  594. rf_top_data = avg_data %>%
  595. mutate(group = factor(group, levels = c("NAC", "ASD"))) %>%
  596. select(group, all_of(top_molecule))
  598. rec = recipe(group ~ ., data = rf_top_data) %>%
  599. step_zv(all_predictors())
  601. rf_spec = rand_forest(trees = 500) %>%
  602. set_mode("classification") %>%
  603. set_engine("ranger")
  605. wf = workflow() %>%
  606. add_recipe(rec) %>%
  607. add_model(rf_spec)
  609. set.seed(123)
  610. folds = vfold_cv(rf_top_data, v = 5, strata = group)
  612. res = fit_resamples(
  613. wf,
  614. resamples = folds,
  615. metrics = metric_set(roc_auc),
  616. control = control_resamples(save_pred = TRUE)
  617. )
  619. preds = collect_predictions(res)
  621. roc_obj = roc(
  622. response = preds$group,
  623. predictor = preds$.pred_ASD,
  624. levels = c("NAC", "ASD"),
  625. direction = "<",
  626. ci = TRUE,
  627. ci.alpha = 0.95,
  628. boot.n = 2000,
  629. stratified = FALSE
  630. )
  632. roc_df = tibble(
  633. specificity = rev(roc_obj$specificities),
  634. sensitivity = rev(roc_obj$sensitivities)
  635. )
  637. auc_val = auc(roc_obj)
  638. ci_vals = ci.auc(roc_obj)
  640. # ggplot(roc_df, aes(x = 1 - specificity, y = sensitivity)) +
  641. # geom_line(color = "#21918c", linewidth = 1.2) +
  642. # geom_abline(linetype = "dashed", color = "gray") +
  643. # annotate(
  644. # "text",
  645. # x = 0.25,
  646. # y = 0.25,
  647. # label = sprintf("AUC = %.3f [%.3f–%.3f]", auc_val, ci_vals[1], ci_vals[3]),
  648. # hjust = 0,
  649. # size = 5
  650. # ) +
  651. # labs(
  652. # title = paste0("ROC Curve with 95% CI for ", paste0(top_molecule, collapse = ", ")),
  653. # x = "1 - Specificity",
  654. # y = "Sensitivity"
  655. # )
  656. ```
  658. ```{r multivar_rf_all, fig.dpi=300, fig.height=10}
  659. set.seed(123)
  661. folds = vfold_cv(avg_data, v = 5, strata = group)
  663. top_mols = rf_auc_all %>%
  664. slice_max(auc, n = 3) %>%
  665. pull(molecule)
  667. model_definitions = tribble(
  668. ~model_name, ~formula, ~panel_group, ~color,
  669. "Age", group ~ age, "Demographics", "#21908CFF",
  670. "Sex", group ~ sex, "Demographics", "#440154FF",
  671. "Age + Sex", group ~ age + sex, "Demographics", "#FDE725FF",
  672. "BH4", group ~ mol_bh4, "Metabolites", "#35B779FF",
  673. "KYN", group ~ mol_kyn, "Metabolites", "#31688EFF",
  674. "GABA", group ~ mol_gaba, "Metabolites", "#440154FF",
  675. "BH4 + KYN + GABA", reformulate(top_mols, response = "group"), "Metabolites", "#26828EFF",
  676. "Age + BH4 + KYN + GABA", reformulate(c("age", top_mols), response = "group"), "Combined", "#5DC863FF",
  677. "Sex + BH4 + KYN + GABA", reformulate(c("sex", top_mols), response = "group"), "Combined", "#3B528BFF",
  678. "Age + Sex + BH4 + KYN + GABA", reformulate(c("age", "sex", top_mols), response = "group"), "Combined", "#481567FF"
  679. )
  681. get_rf_roc = function(formula, model_label) {
  682. rec = recipe(formula, data = avg_data) %>%
  683. step_zv(all_predictors()) %>%
  684. step_dummy(all_nominal_predictors())
  686. wf = workflow() %>%
  687. add_model(rand_forest(trees = 500, mode = "classification") %>%
  688. set_engine("ranger")) %>%
  689. add_recipe(rec)
  691. res = fit_resamples(
  692. wf, resamples = folds,
  693. metrics = metric_set(roc_auc),
  694. control = control_resamples(save_pred = TRUE)
  695. )
  697. preds = collect_predictions(res)
  699. roc_obj = roc(
  700. response = preds$group,
  701. predictor = preds$.pred_ASD,
  702. levels = c("NAC", "ASD"),
  703. direction = "<",
  704. ci = TRUE, boot.n = 2000
  705. )
  707. tibble(
  708. model_name = model_label,
  709. specificity = rev(roc_obj$specificities),
  710. sensitivity = rev(roc_obj$sensitivities),
  711. auc = as.numeric(auc(roc_obj)),
  712. ci_low = as.numeric(ci.auc(roc_obj)[1]),
  713. ci_high = as.numeric(ci.auc(roc_obj)[3])
  714. )
  715. }
  717. roc_results = model_definitions %>%
  718. mutate(roc = map2(formula, model_name, get_rf_roc)) %>%
  719. unnest(roc, names_sep = "_") %>%
  720. rename(
  721. specificity = roc_specificity,
  722. sensitivity = roc_sensitivity,
  723. auc = roc_auc,
  724. ci_low = roc_ci_low,
  725. ci_high = roc_ci_high
  726. )
  728. roc_labeled = roc_results %>%
  729. group_by(model_name) %>%
  730. mutate(
  731. label = sprintf("%s\nAUC = %.3f [%.3f-%.3f]", model_name, first(auc), first(ci_low), first(ci_high))
  732. ) %>%
  733. ungroup()
  735. label_order = roc_labeled %>%
  736. distinct(model_name, label) %>%
  737. arrange(factor(model_name, levels = model_definitions$model_name)) %>%
  738. pull(label)
  740. roc_labeled = roc_labeled %>%
  741. mutate(label = factor(label, levels = label_order))
  743. plot_roc_group = function(data, panel) {
  744. df = filter(data, panel_group == panel)
  746. ggplot(df, aes(x = 1 - specificity, y = sensitivity, color = label)) +
  747. geom_line(linewidth = 1) +
  748. geom_abline(linetype = "dashed", color = "gray") +
  749. scale_color_manual(values = setNames(df$color, df$label)) +
  750. labs(
  751. x = "1-Specificity",
  752. y = "Sensitivity",
  753. color = NULL
  754. ) +
  755. coord_equal()
  756. }
  758. p_demo = plot_roc_group(roc_labeled, "Demographics")
  759. p_metab = plot_roc_group(roc_labeled, "Metabolites")
  760. p_combined = plot_roc_group(roc_labeled, "Combined")
  762. # Figure 3
  763. p_demo / p_metab / p_combined +
  764. plot_annotation(tag_levels = "A", tag_prefix = "(", tag_suffix = ")")
  765. ```
  767. ```{r multivar_rf_metrics}
  768. set.seed(123)
  769. folds = vfold_cv(avg_data, v = 5, strata = group)
  771. compute_model_metrics = function(formula, label) {
  772. rec = recipe(formula, data = avg_data) %>%
  773. step_zv(all_predictors()) %>%
  774. step_dummy(all_nominal_predictors())
  776. rf_spec = rand_forest(trees = 500) %>%
  777. set_mode("classification") %>%
  778. set_engine("ranger")
  780. wf = workflow() %>%
  781. add_model(rf_spec) %>%
  782. add_recipe(rec)
  784. res = fit_resamples(
  785. wf,
  786. resamples = folds,
  787. metrics = metric_set(roc_auc, accuracy, sens, spec, precision, recall, f_meas),
  788. control = control_resamples(save_pred = TRUE)
  789. )
  791. preds = collect_predictions(res)
  793. roc_obj = roc(
  794. response = preds$group,
  795. predictor = preds$.pred_ASD,
  796. levels = c("NAC", "ASD"),
  797. direction = "<",
  798. ci = TRUE,
  799. boot.n = 2000
  800. )
  802. auc_val = as.numeric(auc(roc_obj))
  803. ci_vals = as.numeric(ci.auc(roc_obj))
  805. metric_vals = preds %>%
  806. summarise(
  807. accuracy = accuracy_vec(group, .pred_class),
  808. sens = sens_vec(group, .pred_class, event_level = "second"),
  809. spec = spec_vec(group, .pred_class, event_level = "second"),
  810. precision = precision_vec(group, .pred_class, event_level = "second"),
  811. recall = recall_vec(group, .pred_class, event_level = "second"),
  812. f_meas = f_meas_vec(group, .pred_class, event_level = "second")
  813. )
  815. tibble(
  816. model = label,
  817. auc = auc_val,
  818. auc_ci_low = ci_vals[1],
  819. auc_ci_high = ci_vals[3]
  820. ) %>%
  821. bind_cols(metric_vals)
  822. }
  824. formulas = list(
  825. "Age" = group ~ age,
  826. "Sex" = group ~ sex,
  827. "Age + Sex" = group ~ age + sex,
  828. "BH4" = group ~ mol_bh4,
  829. "KYN" = group ~ mol_kyn,
  830. "GABA" = group ~ mol_gaba,
  831. "BH4 + KYN + GABA" = reformulate(top_mols, response = "group"),
  832. "Sex + BH4 + KYN + GABA" = reformulate(c("sex", top_mols), response = "group"),
  833. "Age + BH4 + KYN + GABA" = reformulate(c("age", top_mols), response = "group"),
  834. "Age + Sex + BH4 + KYN + GABA" = reformulate(c("age", "sex", top_mols), response = "group")
  835. )
  837. metrics_all = imap_dfr(formulas, compute_model_metrics) %>%
  838. mutate(across(where(is.numeric), ~ round(.x, 3)),
  839. auc = paste0(auc, " [", auc_ci_low, "-", auc_ci_high, "]")) %>%
  840. select(-auc_ci_low, -auc_ci_high)
  841. ```
  843. ```{r rf_model_comparison}
  844. # Model comparison to assess whether adding more metabolites
  845. # improves classification performance beyond the top 3 features.
  847. set.seed(123)
  848. folds = vfold_cv(avg_data, v = 5, strata = group)
  850. # Select metabolite sets
  851. top3_mols = rf_auc_all %>%
  852. slice_max(auc, n = 3, with_ties = FALSE) %>%
  853. pull(molecule)
  855. top5_mols = rf_auc_all %>%
  856. slice_max(auc, n = 5, with_ties = FALSE) %>%
  857. pull(molecule)
  859. all_mols = names(avg_data)[str_starts(names(avg_data), "mol_")]
  861. feature_labels = list(
  862. "Top 3 metabolites" = str_c(top3_mols, collapse = ", "),
  863. "Top 3 + Trp" = str_c(c(top3_mols, "mol_trp"), collapse = ", "),
  864. "Top 5 metabolites" = str_c(top5_mols, collapse = ", "),
  865. "All 18 metabolites" = "All metabolites"
  866. )
  868. # Reuse the same modeling framework as in Table S4
  869. compute_model_metrics_compare = function(formula, label) {
  870. rec = recipe(formula, data = avg_data) %>%
  871. step_zv(all_predictors()) %>%
  872. step_dummy(all_nominal_predictors())
  874. rf_spec = rand_forest(trees = 500) %>%
  875. set_mode("classification") %>%
  876. set_engine("ranger")
  878. wf = workflow() %>%
  879. add_model(rf_spec) %>%
  880. add_recipe(rec)
  882. res = fit_resamples(
  883. wf,
  884. resamples = folds,
  885. metrics = metric_set(roc_auc),
  886. control = control_resamples(save_pred = TRUE)
  887. )
  889. preds = collect_predictions(res) %>%
  890. mutate(
  891. .pred_class = factor(
  892. if_else(.pred_ASD >= 0.5, "ASD", "NAC"),
  893. levels = c("NAC", "ASD")
  894. )
  895. )
  897. roc_obj = roc(
  898. response = preds$group,
  899. predictor = preds$.pred_ASD,
  900. levels = c("NAC", "ASD"),
  901. direction = "<",
  902. ci = TRUE,
  903. boot.n = 2000
  904. )
  906. auc_val = as.numeric(auc(roc_obj))
  907. ci_vals = as.numeric(ci.auc(roc_obj))
  909. metric_vals = preds %>%
  910. summarise(
  911. accuracy = accuracy_vec(group, .pred_class),
  912. sens = sens_vec(group, .pred_class, event_level = "second"),
  913. spec = spec_vec(group, .pred_class, event_level = "second"),
  914. precision = precision_vec(group, .pred_class, event_level = "second"),
  915. f_meas = f_meas_vec(group, .pred_class, event_level = "second")
  916. )
  918. tibble(
  919. model = label,
  920. auc = auc_val,
  921. auc_ci_low = ci_vals[1],
  922. auc_ci_high = ci_vals[3]
  923. ) %>%
  924. bind_cols(metric_vals)
  925. }
  927. # Define comparison models
  928. compare_formulas = list(
  929. "Top 3 metabolites" = reformulate(top3_mols, response = "group"),
  930. "Top 3 + Trp" = reformulate(c(top3_mols, "mol_trp"), response = "group"),
  931. "Top 5 metabolites" = reformulate(top5_mols, response = "group"),
  932. "All 18 metabolites" = reformulate(all_mols, response = "group")
  933. )
  935. # Run comparison
  936. rf_model_compare = imap_dfr(compare_formulas, compute_model_metrics_compare) %>%
  937. mutate(
  938. features = feature_labels[model] %>% unlist()
  939. ) %>%
  940. mutate(
  941. across(where(is.numeric), \(x) round(x, 3)),
  942. auc = paste0(auc, " [", auc_ci_low, "-", auc_ci_high, "]")
  943. ) %>%
  944. select(model, features, everything(), -auc_ci_low, -auc_ci_high)
  946. rf_model_compare
  948. # Optional export for supplement
  949. rf_model_compare %>%
  950. write_csv("output_data/s6_rf_model_comparison.csv")
  951. ```
  953. ```{r multivar_rf_calibration, fig.dpi=300, fig.height=10}
  954. get_calibration_metrics = function(formula, model_label) {
  955. rec = recipe(formula, data = avg_data) %>%
  956. step_zv(all_predictors()) %>%
  957. step_dummy(all_nominal_predictors())
  959. wf = workflow() %>%
  960. add_model(rand_forest(trees = 500, mode = "classification") %>% set_engine("ranger")) %>%
  961. add_recipe(rec)
  963. res = fit_resamples(
  964. wf, resamples = folds,
  965. control = control_resamples(save_pred = TRUE)
  966. )
  968. preds = collect_predictions(res) %>%
  969. mutate(group = factor(group, levels = c("NAC", "ASD"))) # ASD = positive class
  971. brier = mean((if_else(preds$group == "ASD", 1, 0) - preds$.pred_ASD)^2)
  973. cal_curve_raw = preds %>%
  974. transmute(
  975. pred = .pred_ASD,
  976. obs = as.integer(group == "ASD"),
  977. model_name = model_label
  978. )
  980. list(
  981. brier = tibble(model_name = model_label, brier_score = brier),
  982. curve_raw = cal_curve_raw
  983. )
  984. }
  986. set.seed(123)
  987. cal_results = model_definitions %>%
  988. mutate(cal = map2(formula, model_name, get_calibration_metrics))
  990. brier_all = map_dfr(cal_results$cal, "brier")
  992. curve_all_raw = map_dfr(cal_results$cal, "curve_raw")
  994. curve_all_annotated = curve_all_raw %>%
  995. left_join(model_definitions %>% select(model_name, panel_group, color), by = "model_name") %>%
  996. left_join(brier_all, by = "model_name") %>%
  997. mutate(
  998. label = sprintf("%s\nBrier = %.3f", model_name, brier_score)
  999. )
  1001. label_order_brier = model_definitions %>%
  1002. left_join(brier_all, by = "model_name") %>%
  1003. mutate(label = sprintf("%s\nBrier = %.3f", model_name, brier_score)) %>%
  1004. pull(label)
  1006. curve_all_annotated = curve_all_annotated %>%
  1007. mutate(label = factor(label, levels = label_order_brier))
  1009. plot_calibration_panel = function(data, panel_name) {
  1010. df = filter(data, panel_group == panel_name)
  1012. ggplot(df, aes(x = pred, y = obs, color = label)) +
  1013. geom_smooth(method = "loess", formula = y ~ x, se = FALSE, span = 1, linewidth = 1) +
  1014. geom_abline(slope = 1, intercept = 0, linetype = "dashed", color = "gray50") +
  1015. scale_color_manual(values = setNames(df$color, df$label)) +
  1016. labs(
  1017. x = "Predicted Probability (ASD)",
  1018. y = "Observed Proportion (ASD)",
  1019. color = NULL
  1020. ) +
  1021. coord_fixed(xlim = c(0, 1), ylim = c(0, 1))
  1022. }
  1024. p_cal_demo = plot_calibration_panel(curve_all_annotated, "Demographics")
  1025. p_cal_metab = plot_calibration_panel(curve_all_annotated, "Metabolites")
  1026. p_cal_combined = plot_calibration_panel(curve_all_annotated, "Combined")
  1028. # Figure 4
  1029. p_cal_demo / p_cal_metab / p_cal_combined +
  1030. plot_annotation(tag_levels = "A", tag_prefix = "(", tag_suffix = ")")
  1032. # Table S4
  1033. metrics_all %>%
  1034. left_join(brier_all %>% mutate(brier_score = signif(brier_score, 3)), by = c("model" = "model_name")) %>%
  1035. write_csv("output_data/s4_multivar_rf_perf_brier_scores.csv")
  1036. ```
  1038. ```{r session-info}
  1039. sessionInfo()
  1040. ```

statistical_analysis.Rmd at commit f07a3a4, under MIT · at the source

Overview

Authors: Kevin Liu1, Huixi Li2, Shaohan Zhang3, Muya Xi1, Junru Zhu2,4, Jonathan Chen3, William Xu3, Alexander Xie1, Alexandros Makriyannis2, Jason J. Guo2,4, Xue-Jun Kong1,3,5
  1. Fetal-Neonatal Neuroimaging and Developmental Science Center, Division of Newborn Medicine, Boston Children’s Hospital, Boston, MA, United States
  2. Department of Pharmaceutical Sciences, Center for Drug Discovery, Northeastern University, Boston, MA, United States
  3. Athinoula A. Martinos Center for Biomedical Imaging, Massachusetts General Hospital, Charlestown, MA, United States
  4. Department of Chemistry and Chemical Biology, Barnett Institute for Chemical and Biological Analysis, Northeastern University, Boston, MA, United States
  5. Department of Medicine and Psychiatry, Beth Israel Deaconess Medical Center, Boston, MA, United States
Institutions: Boston Children's Hospital (United States); Northeastern University (United States); Massachusetts General Hospital (United States); Athinoula A. Martinos Center for Biomedical Imaging (United States); Beth Israel Deaconess Medical Center (United States)
Journal: Frontiers in neuroscience, volume 20, article 1858005
Dates: received 16 April 2026; accepted 22 June 2026; published online 3 July 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.3389/fnins.2026.1858005 · PMID 42490967 · PMCID PMC13375968 · OpenAlex W7167221654
Open access: gold, a free copy (OpenAlex)
Status: code verified
Categories: genetics / omics (modality), autism (population), clinical / translational (subfield)
Methods: Connectivity, Statistics, Machine learning
Keywords: autism spectrum disorder, biomarkers, catecholamines, gut-brain axis, machine learning, microbiome, stool metabolomics, tryptophan metabolism
Topic: Autism Spectrum Disorder Research (Cognitive Neuroscience, Neuroscience), according to OpenAlex
Funding: Boston Children's Hospital (92436); massachusetts general hospital (233263)
Citations: not cited yet (Europe PMC); 62 references in the paper

Abstract

Background: Gut-brain axis dysregulation and microbiome-linked metabolic alterations have been implicated in autism spectrum disorder (ASD), but the contribution of gut-derived neuroactive metabolites remains incompletely characterized.

Methods: We conducted a cross-sectional case-control study of 59 participants (32 ASD, 27 controls) and quantified 18 stool metabolites related to catecholamine synthesis, inhibitory neurotransmission, and tryptophan-linked NAD+-precursor metabolism using targeted liquid chromatography-tandem mass spectrometry. Group differences were assessed using fold-change analysis and linear models adjusted for age and sex. Random forest models evaluated classification performance, and within-group Spearman correlations were used to examine metabolic relationships.

Results: Norepinephrine showed the largest increase in ASD, whereas dopamine and tetrahydrobiopterin exhibited nominal group differences that did not remain significant after correction for multiple testing. A three-metabolite panel comprising tetrahydrobiopterin, γ-aminobutyric acid, and kynurenine showed exploratory discrimination between groups (area under the receiver operating characteristic curve = 0.750, 95% confidence interval 0.622–0.878), but this performance requires external validation. Correlation analysis revealed conserved bile acid coupling in both groups. In controls, tryptophan was positively associated with kynurenine, whereas this relationship was not observed in ASD. Instead, ASD samples showed broader associations between tryptophan and metabolites linked to neurotransmission and NAD+-precursor metabolism.

Conclusion: Stool metabolite profiling revealed altered organization of tryptophan- and catecholamine-linked metabolic associations in ASD and identified a small metabolite panel with exploratory discriminative potential. These findings provide a foundation for future studies examining gut-derived neuroactive metabolites in ASD and their relationship to gut-brain axis biology.

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 4 matches between paragraphs and lines of code.

kevinliu-bmb/targeted-asd-stool-metabolomics

License: MIT
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: f07a3a45dec086833da3a29b2dcc84d758b69361, 4 July 2026
Languages: R (1)
Size: 4 files, 1 script
Software Heritage: not archived
Found in: the text, “Statistical analysis”
Holds: README, license file, 1 notebook
Not found: CITATION.cff, environment file, tests, continuous integration, documentation
Tools: broom (1 file), ggpubr (1 file), patchwork (1 file), pROC (1 file), rstatix (1 file), tidyverse (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
3 files
At the source: github.com/kevinliu-bmb/targeted-asd-stool-metabolomics

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;
  • 4 matches 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 statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

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, pages, dates, 11 authors, 8 keywords, 2 funders, 61 references.

Cite

This paper

Liu, K., Li, H., Zhang, S., Xi, M., Zhu, J., Chen, J., Xu, W., Xie, A., Makriyannis, A., Guo, J. J., & Kong, X.-J. (2026). Targeted stool metabolomics suggests exploratory catecholamine- and tryptophan-linked metabolic features in autism spectrum disorder. Frontiers in neuroscience, 20, 1858005. https://doi.org/10.3389/fnins.2026.1858005

BibTeX

@article{liu2026targeted,
author = {Liu, Kevin and Li, Huixi and Zhang, Shaohan and Xi, Muya and Zhu, Junru and Chen, Jonathan and Xu, William and Xie, Alexander and Makriyannis, Alexandros and Guo, Jason J. and Kong, Xue-Jun},
title = {{Targeted stool metabolomics suggests exploratory catecholamine- and tryptophan-linked metabolic features in autism spectrum disorder}},
journal = {Frontiers in neuroscience},
year = {2026},
month = jul,
volume = {20},
pages = {1858005},
publisher = {Frontiers Media SA},
issn = {1662-4548},
doi = {10.3389/fnins.2026.1858005},
url = {https://doi.org/10.3389/fnins.2026.1858005},
pmid = {42490967},
pmcid = {PMC13375968}
}

RIS

TY - JOUR
AU - Liu, Kevin
AU - Li, Huixi
AU - Zhang, Shaohan
AU - Xi, Muya
AU - Zhu, Junru
AU - Chen, Jonathan
AU - Xu, William
AU - Xie, Alexander
AU - Makriyannis, Alexandros
AU - Guo, Jason J.
AU - Kong, Xue-Jun
TI - Targeted stool metabolomics suggests exploratory catecholamine- and tryptophan-linked metabolic features in autism spectrum disorder
T2 - Frontiers in neuroscience
J2 - Front Neurosci
PY - 2026
DA - 2026/07/03
VL - 20
SP - 1858005
SN - 1662-4548
PB - Frontiers Media SA
DO - 10.3389/fnins.2026.1858005
UR - https://doi.org/10.3389/fnins.2026.1858005
LA - en
ER -

CSL-JSON

{
"id": "10.3389/fnins.2026.1858005",
"type": "article-journal",
"title": "Targeted stool metabolomics suggests exploratory catecholamine- and tryptophan-linked metabolic features in autism spectrum disorder",
"container-title": "Frontiers in neuroscience",
"author": [
{
"family": "Liu",
"given": "Kevin"
},
{
"family": "Li",
"given": "Huixi"
},
{
"family": "Zhang",
"given": "Shaohan"
},
{
"family": "Xi",
"given": "Muya"
},
{
"family": "Zhu",
"given": "Junru"
},
{
"family": "Chen",
"given": "Jonathan"
},
{
"family": "Xu",
"given": "William"
},
{
"family": "Xie",
"given": "Alexander"
},
{
"family": "Makriyannis",
"given": "Alexandros"
},
{
"family": "Guo",
"given": "Jason J."
},
{
"family": "Kong",
"given": "Xue-Jun"
}
],
"container-title-short": "Front Neurosci",
"volume": "20",
"page": "1858005",
"DOI": "10.3389/fnins.2026.1858005",
"PMID": "42490967",
"PMCID": "PMC13375968",
"ISSN": "1662-4548",
"publisher": "Frontiers Media SA",
"URL": "https://doi.org/10.3389/fnins.2026.1858005",
"language": "en",
"issued": {
"date-parts": [
[
2026,
7,
3
]
]
}
}

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.1126/sciadv.aec9291 [code]
Computational mechanisms of perception in autism revealed using games inspired by rodent operant tasks.
Journal: Science advances
In common: pROC, rstatix, broom, 3 other tools, autism
[2] doi:10.1093/bioinformatics/btag592 [code]
Network-based stratification of allele-specific expression reveals patient subgroups in Huntington's disease.
Journal: Bioinformatics (Oxford, England)
In common: pROC, rstatix, broom, 3 other tools, genetics / omics
[3] doi:10.1186/s12888-026-08178-8 [code]
Machine learning model to identify gut microbiome-derived metabolites as potential biomarkers of autism spectrum disorder: a pilot study.
Journal: BMC psychiatry
In common: autism, clinical / translational, genetics / omics, 4 references
[4] doi:10.1016/j.isci.2026.115657 [code]
Integration of machine learning to develop a disulfidptosis model for predicting glioma prognosis, immunotherapy response, and drug.
Journal: iScience
In common: pROC, broom, ggpubr, 2 other tools, clinical / translational
[5] doi:10.64898/2026.05.08.26348885 [code]
Insights from nine nights of self-applied, low-density sleep EEG during sleep restriction therapy: a proof-of-concept evaluation
Journal: medRxiv (preprint)
In common: rstatix, broom, ggpubr, 2 other tools, clinical / translational
[6] doi:10.1093/braincomms/fcag236 [code]
Dynamic, state-dependent characteristics of cognitive fluctuations in Lewy body dementia: a magnetoencephalography study.
Journal: Brain communications
In common: pROC, broom, ggpubr, 2 other tools
[7] doi:10.1073/pnas.2609132123 [code]
A human lysosomal storage disorder toolkit for decoding proteome landscapes in cortical-like and dopaminergic-like induced neurons.
Journal: Proceedings of the National Academy of Sciences of the United States of America
In common: rstatix, broom, ggpubr, 2 other tools, genetics / omics
[8] doi:10.1002/trc2.70257 [code]
Blood DNA methylation signature of cognitive reserve moderates the association between CSF tau pathology and memory in prodromal Alzheimer's disease.
Journal: Alzheimer's & dementia (New York, N. Y.)
In common: rstatix, broom, ggpubr, 2 other tools, genetics / omics
[9] doi:10.1093/neuonc/noag128 [code]
Spatially-resolved single-cell imaging of melanoma brain metastases identifies localized immune patterns predictive of immune checkpoint blockade response.
Journal: Neuro-oncology
In common: pROC, rstatix, patchwork, 1 other tool, clinical / translational, genetics / omics
[10] doi:10.1093/braincomms/fcag272 [code]
Strenuous physical activity is associated with a younger age of amyotrophic lateral sclerosis onset in two independent cohorts.
Journal: Brain communications
In common: pROC, rstatix, broom, 1 other tool, clinical / translational

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.