Bayesian Workflow for Generative Modeling in Computational Psychiatry

Alexander J. Hess; Sandra Iglesias; Laura Köchli; Stephanie Marino; Matthias Müller-Schrader; Lionel Rigoux; Christoph Mathys; Olivia K. Harrison; Jakob Heinzle; Stefan Frässle; Klaas Enno Stephan

doi:10.5334/cpsy.116

Figures & Tables

The speed-incentivised reward learning (SPIRL) task. A shows the trial structure of the SPIRL task protocol. A yellow and a green fractal were presented on every trial together with a time bar indicating the remaining time of the 1.7s long response window. The participants had to predict on each trial, which fractal would be rewarded monetarily. After the response window, the trial outcome was revealed (reward/no reward) concurrently with the start of the new response window of the next trial. In B, the probability of reward for one of the two fractals over the entire 160 trials is displayed (black line). The individual trial outcomes for this fractal are indicated by black dots (1 = reward, 0 = no reward). The reward probabilities of the two fractals were complementary (summing to 1) across the entire task. The colour shadings represent different phases during the task (white: stable; blue: high volatility, grey: highly unpredictable).

Bayesian workflow for generative modelling in Computational Psychiatry. The general steps of Bayesian workflow are indicated by grey boxes with labels. These include the specification of a model space, prior specification, model inversion and validation of computation, model comparison as well as model evaluation. The flowchart represents how the concrete analyses steps of our application map onto the general framework of Bayesian workflow. Above the dashed line are analysis steps that involve the pilot data set as well as synthetic data generated using parameter values sampled from the priors. Below the dashed line are analysis steps including the main data set and synthetic data generated using parameter values sampled from the posteriors. Filled red boxes are steps of the pipeline that involve synthetic data. The green boxes highlight analysis steps that specifically refer to our research aim 1 (comparison of RT model families) and research aim 2 (assessment of individual RT model parameters).

Graphical model representation and example belief trajectories. In the left part of the figure, a schematic representation of the generative model of the 3-level eHGF for binary inputs (perceptual model) is presented on top. Below, the response data modalities are visualised (response model). Shaded circles represent known quantities (inputs shaded black, response data shaded red). Unshaded circles represent estimated time-independent parameters (black circles) and time-varying states with trial indices as superscript. Dashed lines indicate how the response model depends on the inferential process in the perceptual model. Solid lines indicate generative processes. Dark turquoise lines indicate the probabilistic network on trial k. Light turquoise lines indicate the network at other points in time. On the right side of the figure, average belief trajectories of the perceptual model are shown in blue in the upper four panels. Note that for illustration of trajectories and responses, we used M1 (for more details on M1, please refer to Table 1). The four belief trajectories represent the average trajectory over participants (*N_main* = 59) of the states that enter the log RT GLM of M1. In the lower two panels, average response data over participants are shown. For the binary response modality, the red line represents average binary prediction ( ${\bar{y}}_{bin}$ ), the blue line represents the average belief about the probability of the outcome ( ${\hat{μ}}_{1}$ ) according to M1. For the continuous response modality, the red line represents average log-transformed RTs and the average predicted log RTs by M1 are shown in blue.

Table 1

Model Space. All seven models in our model space are composed of a perceptual (Prc) model and an observation (Obs) or response model. The perceptual model as well as the binary part of the response model is held constant across all seven models. The equations of the log RT GLMs (continuous part of the response model) of M1–M4 (family of informed RT models) contain belief trajectories of the perceptual model as regressors. The update equations for the perceptual model (eHGF) are listed in the analysis plan. M5–M7 (uninformed RT family) predict log RTs independent of the perceptual model.

MODEL	PRC MODEL	OBS MODEL (CONTINUOUS)	OBS MODEL (BINARY)
M1	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0} + β_{1} S^{(k - 1)} + β_{2} {\hat{σ}}_{1}^{(k)} + β_{3} {\hat{σ}}_{2}^{(k)} + β_{4} e^{{\hat{μ}}_{3}^{(k)}}, Σ)$ with $S^{(k - 1)} \overset{def}{=} {\begin{matrix} - \log_{2} ({\hat{μ}}_{1}^{(k - 1)}) if u^{(k - 1)} = 1 \\ - \log_{2} (1 - {\hat{μ}}_{1}^{(k - 1)}) if u^{(k - 1)} = 0 \end{matrix}$	Unit-square sigmoid
M2	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0} + β_{1} σ_{2}^{(k)} + β_{2} {\hat{μ}}_{3}^{(k)}, Σ)$	Unit-square sigmoid
M3	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0} + β_{1} {\hat{σ}}_{1}^{(k)} + β_{2} σ_{2}^{(k)} + β_{3} e^{(κ_{2} μ_{3}^{(k - 1)} + ω_{2})}, Σ)$	Unit-square sigmoid
M4	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0} + β_{1} σ_{2}^{(k)} + β_{2} e^{(κ_{2} μ_{3}^{(k - 1)} + ω_{2})}, Σ)$	Unit-square sigmoid
M5	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0}, Σ)$	Unit-square sigmoid
M6	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0} + β_{1} \frac{k}{K}, Σ)$	Unit-square sigmoid
M7	eHGF	$\log (y_{rt}^{(k)}) ~ N (β_{0}^{(k)}, Σ)$ with $β_{0}^{(k)} \overset{def}{=} {\begin{matrix} a, if y_{bin}^{(k - 1)} = u^{(k - 1)} \\ b, if y_{bin}^{(k - 1)} \neq u^{(k - 1)} \end{matrix}$	Unit-square sigmoid

Binary responses and continuous log-transformed response times. In A, the red line represents a trial-wise summary of percentage of incorrect responses (inverted for true probabilities 0.2 indicated by the black line) over all participants (*N_main* = 59) and the average absolute prediction error about the trial outcome ( $| {\bar{δ}}_{1} |$ ) of M1 in red. Black dots represent the reward of fractal A on each trial (1 = reward, 0 = no reward) and the black line shows the underlying probability structure of the task. B shows the histogram of log-transformed RTs over all participants in ms. The histogram of residuals of log RT model fits obtained by M1 are visualised in C.

Prior configurations of M1. A shows the *empirical prior* densities for each free parameter of M1 (solid line) as estimated using MAP estimates (black dots) obtained from a separate pilot data set (*N_pilot* = 20) using the *initial priors* (dashed lines). A detailed description of M1 can be found in Table 1 and in the main text. Prior predictive distributions under the *empirical priors* of M1 are displayed for both response data modalities. In B on the left, belief trajectories about the outcome ( ${\hat{μ}}_{1}$ ) at the lowest level of the eHGF are displayed in blue with the thick blue line representing the resulting belief trajectory using the *empirical prior* mean parameter values and the yellow line representing the average simulated binary response (*N_sim* = 100). The green vertical lines indicate trials with minimum (dark green), median (green), and maximum (light green) variance across the simulated trajectories as illustrated on the right. In B on the right, histograms of simulated ${\hat{μ}}_{1}$ values are presented in green for trials with maximum, median and minimum ${\hat{μ}}_{1}$ variance (from light to dark) across simulated trajectories. See the traces on the left for an indication where these trials are located within the trajectory. In C on the left, simulated log RT data are shown in blue, the dashed black lines are the boundaries given by the length of the response window in each trial. The thick blue line represents the average simulated log RT trajectory (*N_sim* = 100). In C on the right, histograms of simulated log RT values are presented in green for trials with maximum, median and minimum log RT variance (from light to dark) across simulated trajectories.

Validation of Computation. A shows results from family-level recovery analysis comparing Ef (left) and XP (right) values with true family frequencies. B depicts results from individual model level recovery analysis. 7 × 7 confusion matrices for LME winner frequencies (left) and PXP scores (right) are shown with data generating models on the y-axis and recovered models on the x-axis. C visualizes simulated vs. estimated values of all free parameters of M1 (parameter recovery). Correlation between simulated and estimated parameter values is indicated using Pearson correlation coefficients r. In D, recoverability of the ω₃ parameter of the 3-level eHGF for binary inputs is visualised for all seven models in the model space including Pearson correlation coefficients r. Models from the informed RT family (M1–M4) show consistently better ω₃ recovery compared to models of the uninformed RT family (M5–M7).

Hypothesis testing. A shows the results of the family level RFX BMS (Efs on the left, XPs on the right) with the informed RT family clearly outperforming the uninformed RT model family. B displays Efs (left) and PXPs (right) resulting from individual model level RFX BMS. M1 can be identified as the clear winning model. C shows raincloud plots of the MAP estimates of the M1 GLM regressors (generated using the RainCloudPlots library). Fine black vertical lines indicate the *initial prior* mean values (i.e. 0) and a black star indicates significantly different prior and posterior means of β₂ which scales the influence of ${\hat{σ}}_{1}^{(k)}$ on log RTs (one-sample t-test, p < 0.001).

Posterior predictive checks for M1. Data from three participants of the main data set are shown. Participants were chosen according to goodness of model fit of M1, i.e. participant 21 with a high log likelihood value, participant 29 showing average goodness of fit and participant 5 showing the worst fit. A displays adjusted correctness of binary responses for these participants in red. Blue circles are the simulated adjusted correctness values resulting from sampled parameter values of the subject-specific posteriors of M1. The blue probability densities are the estimated posterior predictive densities based on the samples drawn from the posteriors (*N_ppc* = 100) using kernel density estimation as implemented in the RainCloudPlots library. In B, we show empirical log RT trajectories of the three participants in red. Fine blue lines are simulated log RT trajectories resulting from sampled parameter values of the subject-specific posteriors and the thick blue line represents the predicted log RT when using the MAP estimates of M1 for each participant to generate synthetic RT data. The green vertical lines indicate trials with minimum (dark green), median (green), and maximum (light green) variance across the simulated trajectories as illustrated on the right. The histograms on the right visualise the distribution of synthetic log RT data generated by simulating from the subject-specific posteriors for trials with maximum, median, and minimum variance (from light to dark green) across all simulated trajectories. See the traces on the left for an indication where these trials are located within the trajectory.

BDA	Bayesian data analysis
BFGS	Broyden-Fletcher-Goldfarb-Shanno
BMS	Bayesian Model Selection
CI	Confidence Interval
CP	Computational Psychiatry
Ef	Expected posterior frequency
eHGF	3-level enhanced Hierarchical Gaussian Filter for binary inputs
FFX	fixed-effects
GLM	General Linear Model
GRW	Gaussian random walks
HGF	Hierarchical Gaussian Filter
LME	Log model evidence
MAP	Maximum a posteriori
MCMC	Markov Chain Monte Carlo
PXP	Protected exceedance probability
RFX	Random effects
RL	Reinforcement Learning
RT	Response times
SBC	Simulation-based Calibration
SPIRL	Speed-incentivised associative reward learning
TAPAS	Translational Algorithms for Psychiatry Advancing Science
TN	Translational Neuromodeling
VB	Variational Bayes
XP	Exceedance probability

Bayesian Workflow for Generative Modeling in Computational Psychiatry

Figures & Tables

Figure 1

Figure 2

Figure 3

Table 1

Figure 4

Figure 5

Figure 6

Figure 7

Figure 8

Paradigm

My account