decision_making.tex

% Decision Making Models Template
% Topics: Expected utility, prospect theory, drift-diffusion, Bayesian decisions, reinforcement learning
% Style: Computational neuroscience report with behavioral modeling

\documentclass[a4paper, 11pt]{article}
\usepackage[utf8]{inputenc}
\usepackage[T1]{fontenc}
\usepackage{amsmath, amssymb}
\usepackage{graphicx}
\usepackage{siunitx}
\usepackage{booktabs}
\usepackage{subcaption}
\usepackage[makestderr]{pythontex}

% Theorem environments
\newtheorem{definition}{Definition}[section]
\newtheorem{theorem}{Theorem}[section]
\newtheorem{example}{Example}[section]
\newtheorem{remark}{Remark}[section]

\title{Computational Models of Human Decision Making:\\
From Rational Choice to Bounded Rationality}
\author{Cognitive Neuroscience Laboratory}
\date{\today}

\begin{document}
\maketitle

\begin{abstract}
This report presents a comprehensive computational analysis of human decision-making processes,
spanning classical normative theories to modern descriptive models. We examine expected utility
theory as the rational benchmark, prospect theory's account of systematic deviations from
rationality, drift-diffusion models for response time distributions, Bayesian frameworks for
decision making under uncertainty, and reinforcement learning models of value-based choice.
Computational simulations reveal the parametric signatures that distinguish these frameworks
and their ability to capture key phenomena including loss aversion, probability weighting,
speed-accuracy tradeoffs, and learning dynamics.
\end{abstract}

\section{Introduction}

Decision making is a fundamental cognitive process that has been studied across economics,
psychology, and neuroscience. Classical economic theory posits that humans are rational agents
who maximize expected utility, but decades of empirical research have revealed systematic
deviations from this normative ideal. Modern computational approaches seek to characterize
both the mechanisms underlying choice behavior and the real-time dynamics of the decision process.

\begin{definition}[Decision Problem]
A decision problem consists of a set of alternatives $\mathcal{A}$, a set of possible outcomes
$\mathcal{O}$, and a probability distribution $p(o|a)$ mapping actions to outcomes. The decision
maker selects an action $a^* \in \mathcal{A}$ to optimize some objective function.
\end{definition}

\section{Expected Utility Theory}

\subsection{Theoretical Foundation}

Expected utility theory (EUT), formalized by von Neumann and Morgenstern, assumes that decision
makers assign utilities $u(o)$ to outcomes and choose actions that maximize expected utility.

\begin{definition}[Expected Utility]
For a risky prospect that yields outcome $o_i$ with probability $p_i$, the expected utility is:
\begin{equation}
EU = \sum_{i=1}^{n} p_i \cdot u(o_i)
\end{equation}
where $u: \mathcal{O} \rightarrow \mathbb{R}$ is a utility function.
\end{definition}

\begin{theorem}[Risk Aversion]
A decision maker is risk-averse if $u$ is concave ($u'' < 0$), risk-neutral if $u$ is linear,
and risk-seeking if $u$ is convex ($u'' > 0$).
\end{theorem}

\subsection{Computational Analysis}

We now examine the implications of different utility functions for choice behavior under risk.
The curvature of the utility function determines risk preferences: concave functions exhibit
diminishing marginal utility, leading to risk aversion. We simulate decisions over monetary
gambles with varying utility function parameters to demonstrate these effects.

\begin{pycode}
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import minimize, fsolve
from scipy.stats import norm, gamma
from scipy.integrate import odeint
import matplotlib.patches as mpatches

np.random.seed(42)

# Expected Utility Theory
def power_utility(x, alpha=0.5):
    """Power utility function with risk aversion parameter alpha"""
    if alpha == 1.0:
        return x
    return np.sign(x) * np.abs(x)**alpha

def expected_utility(outcomes, probabilities, alpha=0.5):
    """Calculate expected utility for a prospect"""
    utilities = power_utility(outcomes, alpha)
    return np.sum(probabilities * utilities)

# Simulate choice behavior under different risk attitudes
wealth = 100.0
stake = 50.0
win_prob = 0.5

# Range of risk aversion parameters
alphas = np.linspace(0.3, 1.5, 100)
certainty_equivalents = []

for alpha in alphas:
    # Gamble: win or lose $50 with equal probability
    outcomes = np.array([wealth + stake, wealth - stake])
    probs = np.array([win_prob, 1 - win_prob])
    eu_gamble = expected_utility(outcomes, probs, alpha)

    # Find certainty equivalent: CE such that u(wealth + CE) = EU(gamble)
    def objective(ce):
        return power_utility(wealth + ce, alpha) - eu_gamble

    ce = fsolve(objective, 0)[0]
    certainty_equivalents.append(ce)

certainty_equivalents = np.array(certainty_equivalents)

# Risk premium: expected value - certainty equivalent
expected_value = stake * win_prob - stake * (1 - win_prob)  # = 0 for fair gamble
risk_premium = expected_value - certainty_equivalents

# Create first figure
fig1, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))

# Plot 1: Utility functions with different risk attitudes
x_wealth = np.linspace(0, 200, 500)
ax1.plot(x_wealth, power_utility(x_wealth, 0.5), 'b-', linewidth=2.5, label=r'$\alpha=0.5$ (Risk-averse)')
ax1.plot(x_wealth, power_utility(x_wealth, 1.0), 'g--', linewidth=2.5, label=r'$\alpha=1.0$ (Risk-neutral)')
ax1.plot(x_wealth, power_utility(x_wealth, 1.5), 'r-.', linewidth=2.5, label=r'$\alpha=1.5$ (Risk-seeking)')
ax1.axvline(wealth, color='gray', linestyle=':', alpha=0.5)
ax1.axvline(wealth + stake, color='orange', linestyle=':', alpha=0.5)
ax1.axvline(wealth - stake, color='purple', linestyle=':', alpha=0.5)
ax1.set_xlabel('Wealth (\$)', fontsize=11)
ax1.set_ylabel('Utility $u(x)$', fontsize=11)
ax1.set_title('Expected Utility Theory: Risk Attitudes', fontsize=12, fontweight='bold')
ax1.legend(fontsize=10)
ax1.grid(True, alpha=0.3)

# Plot 2: Risk premium as function of risk aversion
ax2.plot(alphas, risk_premium, 'b-', linewidth=2.5)
ax2.axhline(0, color='black', linestyle='--', alpha=0.5)
ax2.axvline(1.0, color='gray', linestyle=':', alpha=0.5, label='Risk-neutral')
ax2.fill_between(alphas, risk_premium, 0, where=(alphas < 1.0), alpha=0.3, color='red', label='Risk premium (averse)')
ax2.fill_between(alphas, risk_premium, 0, where=(alphas > 1.0), alpha=0.3, color='green', label='Risk premium (seeking)')
ax2.set_xlabel(r'Risk Aversion Parameter $\alpha$', fontsize=11)
ax2.set_ylabel('Risk Premium (\$)', fontsize=11)
ax2.set_title('Risk Premium vs. Risk Attitude', fontsize=12, fontweight='bold')
ax2.legend(fontsize=10)
ax2.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('decision_making_eut.pdf', dpi=150, bbox_inches='tight')
plt.close()

# Calculate key statistics for reporting
alpha_moderate = 0.5
eu_gamble_moderate = expected_utility(
    np.array([wealth + stake, wealth - stake]),
    np.array([win_prob, 1 - win_prob]),
    alpha_moderate
)
ce_moderate = fsolve(lambda ce: power_utility(wealth + ce, alpha_moderate) - eu_gamble_moderate, 0)[0]
rp_moderate = expected_value - ce_moderate
\end{pycode}

\begin{figure}[htbp]
\centering
\includegraphics[width=\textwidth]{decision_making_eut.pdf}
\caption{Expected utility theory predictions for risky choice. Left panel shows utility functions
with varying risk attitudes: concave (risk-averse, $\alpha=0.5$), linear (risk-neutral, $\alpha=1.0$),
and convex (risk-seeking, $\alpha=1.5$). Vertical lines indicate current wealth (\$100) and potential
outcomes (\$150, \$50) from a 50-50 gamble. Right panel displays the risk premium as a function of
the risk aversion parameter $\alpha$, demonstrating that risk-averse agents ($\alpha < 1$) demand
compensation to accept fair gambles, while risk-seeking agents ($\alpha > 1$) accept unfavorable
gambles. For $\alpha=0.5$, the risk premium is approximately \$\py{f"{abs(rp_moderate):.2f}"},
indicating the agent requires the gamble's expected value to exceed this amount.}
\label{fig:eut}
\end{figure}

The analysis demonstrates that a moderately risk-averse agent ($\alpha = \py{f"{alpha_moderate}"}$)
has a certainty equivalent of \$\py{f"{ce_moderate:.2f}"} for the gamble, yielding a risk premium
of \$\py{f"{abs(rp_moderate):.2f}"}. This quantifies the intuition that risk-averse individuals
prefer certain outcomes over risky prospects with the same expected value.

\section{Prospect Theory}

\subsection{Value Function and Loss Aversion}

Kahneman and Tversky's prospect theory revolutionized the study of decision making by documenting
systematic violations of expected utility theory. The theory proposes that people evaluate outcomes
relative to a reference point and exhibit loss aversion: losses loom larger than gains.

\begin{definition}[Prospect Theory Value Function]
The value function is defined as:
\begin{equation}
v(x) = \begin{cases}
x^\alpha & \text{if } x \geq 0 \\
-\lambda |x|^\beta & \text{if } x < 0
\end{cases}
\end{equation}
where $\alpha, \beta \in (0,1)$ capture diminishing sensitivity, and $\lambda > 1$ represents
loss aversion (typical value: $\lambda \approx 2.25$).
\end{definition}

\subsection{Probability Weighting}

Prospect theory also incorporates nonlinear probability weighting, whereby small probabilities
are overweighted and large probabilities are underweighted.

\begin{definition}[Probability Weighting Function]
The one-parameter Prelec weighting function is:
\begin{equation}
w(p) = \exp\{-(-\ln p)^\gamma\}
\end{equation}
where $\gamma < 1$ produces the characteristic inverse-S shape.
\end{definition}

\subsection{Computational Simulation}

We now implement the full prospect theory model to analyze how loss aversion and probability
weighting jointly affect decision making. The model predicts specific phenomena such as the
fourfold pattern of risk attitudes and preference reversals that cannot be explained by
expected utility theory.

\begin{pycode}
# Prospect Theory Implementation
def prospect_value_function(x, alpha=0.88, beta=0.88, lam=2.25):
    """Kahneman-Tversky value function with loss aversion"""
    if isinstance(x, np.ndarray):
        v = np.zeros_like(x)
        v[x >= 0] = x[x >= 0]**alpha
        v[x < 0] = -lam * np.abs(x[x < 0])**beta
        return v
    else:
        if x >= 0:
            return x**alpha
        else:
            return -lam * np.abs(x)**beta

def prelec_weighting(p, gamma=0.61):
    """Prelec probability weighting function"""
    return np.exp(-(-np.log(p))**gamma)

def prospect_theory_value(outcomes, probabilities, alpha=0.88, beta=0.88, lam=2.25, gamma=0.61):
    """Calculate prospect theory value for a gamble"""
    values = prospect_value_function(outcomes, alpha, beta, lam)
    weights = prelec_weighting(probabilities, gamma)
    return np.sum(weights * values)

# Parameters from Tversky & Kahneman (1992)
alpha_tk = 0.88
beta_tk = 0.88
lambda_tk = 2.25
gamma_tk = 0.61

# Create value function plot
x_range = np.linspace(-100, 100, 500)
v_values = prospect_value_function(x_range, alpha_tk, beta_tk, lambda_tk)

fig2, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(14, 11))

# Plot 1: Value function
ax1.plot(x_range, v_values, 'b-', linewidth=2.5)
ax1.axhline(0, color='black', linestyle='-', linewidth=0.5)
ax1.axvline(0, color='black', linestyle='-', linewidth=0.5)
# Add reference lines showing loss aversion
x_ref = 50
v_gain = prospect_value_function(x_ref, alpha_tk, beta_tk, lambda_tk)
v_loss = prospect_value_function(-x_ref, alpha_tk, beta_tk, lambda_tk)
ax1.plot([x_ref, x_ref], [0, v_gain], 'g--', alpha=0.5)
ax1.plot([-x_ref, -x_ref], [v_loss, 0], 'r--', alpha=0.5)
ax1.plot([0, x_ref], [v_gain, v_gain], 'g--', alpha=0.5)
ax1.plot([0, -x_ref], [v_loss, v_loss], 'r--', alpha=0.5)
ax1.text(x_ref+5, v_gain/2, f'Gain: {v_gain:.1f}', fontsize=9, color='green')
ax1.text(-x_ref-25, v_loss/2, f'Loss: {v_loss:.1f}', fontsize=9, color='red')
ax1.set_xlabel('Outcome (Relative to Reference Point)', fontsize=11)
ax1.set_ylabel('Subjective Value $v(x)$', fontsize=11)
ax1.set_title(f'Prospect Theory Value Function ($\\lambda={lambda_tk}$)', fontsize=12, fontweight='bold')
ax1.grid(True, alpha=0.3)

# Plot 2: Probability weighting
p_range = np.linspace(0.01, 0.99, 100)
w_values = prelec_weighting(p_range, gamma_tk)
ax2.plot(p_range, w_values, 'b-', linewidth=2.5, label=f'$\\gamma={gamma_tk}$')
ax2.plot([0, 1], [0, 1], 'k--', linewidth=1, label='Identity (no weighting)')
ax2.fill_between(p_range, p_range, w_values, where=(w_values > p_range),
                  alpha=0.3, color='red', label='Overweighting')
ax2.fill_between(p_range, p_range, w_values, where=(w_values < p_range),
                  alpha=0.3, color='blue', label='Underweighting')
ax2.set_xlabel('Objective Probability $p$', fontsize=11)
ax2.set_ylabel('Decision Weight $w(p)$', fontsize=11)
ax2.set_title('Prelec Probability Weighting Function', fontsize=12, fontweight='bold')
ax2.legend(fontsize=9)
ax2.grid(True, alpha=0.3)

# Plot 3: Fourfold pattern of risk attitudes
# Gains domain
gain_amounts = np.linspace(10, 100, 20)
p_low = 0.01
p_high = 0.90

ce_gain_low = []
ce_gain_high = []

for gain in gain_amounts:
    # Low probability gain
    pt_val_low = prospect_theory_value(np.array([gain, 0]), np.array([p_low, 1-p_low]),
                                        alpha_tk, beta_tk, lambda_tk, gamma_tk)
    # Find CE
    ce_low = fsolve(lambda x: prospect_value_function(x, alpha_tk, beta_tk, lambda_tk) - pt_val_low, gain*p_low)[0]
    ce_gain_low.append(ce_low)

    # High probability gain
    pt_val_high = prospect_theory_value(np.array([gain, 0]), np.array([p_high, 1-p_high]),
                                         alpha_tk, beta_tk, lambda_tk, gamma_tk)
    ce_high = fsolve(lambda x: prospect_value_function(x, alpha_tk, beta_tk, lambda_tk) - pt_val_high, gain*p_high)[0]
    ce_gain_high.append(ce_high)

expected_gain_low = gain_amounts * p_low
expected_gain_high = gain_amounts * p_high

ax3.plot(expected_gain_low, ce_gain_low, 'b-', linewidth=2.5, label=f'Low prob ({p_low})')
ax3.plot(expected_gain_high, ce_gain_high, 'g-', linewidth=2.5, label=f'High prob ({p_high})')
ax3.plot([0, 100], [0, 100], 'k--', linewidth=1, label='Risk-neutral')
ax3.set_xlabel('Expected Value', fontsize=11)
ax3.set_ylabel('Certainty Equivalent', fontsize=11)
ax3.set_title('Fourfold Pattern: Gains Domain', fontsize=12, fontweight='bold')
ax3.legend(fontsize=10)
ax3.grid(True, alpha=0.3)

# Plot 4: Loss aversion magnitude
lambda_values = np.linspace(1.0, 3.5, 50)
mixed_gamble_acceptance = []

for lam in lambda_values:
    # Mixed gamble: 50% chance to win $100, 50% chance to lose $100
    outcomes = np.array([100, -100])
    probs = np.array([0.5, 0.5])
    pt_val = prospect_theory_value(outcomes, probs, alpha_tk, beta_tk, lam, gamma_tk)
    mixed_gamble_acceptance.append(pt_val)

mixed_gamble_acceptance = np.array(mixed_gamble_acceptance)

ax4.plot(lambda_values, mixed_gamble_acceptance, 'b-', linewidth=2.5)
ax4.axhline(0, color='red', linestyle='--', linewidth=2, label='Indifference threshold')
ax4.axvline(lambda_tk, color='gray', linestyle=':', alpha=0.7, label=f'Typical $\\lambda={lambda_tk}$')
ax4.fill_between(lambda_values, mixed_gamble_acceptance, 0,
                  where=(mixed_gamble_acceptance > 0), alpha=0.3, color='green', label='Accept')
ax4.fill_between(lambda_values, mixed_gamble_acceptance, 0,
                  where=(mixed_gamble_acceptance < 0), alpha=0.3, color='red', label='Reject')
ax4.set_xlabel('Loss Aversion Parameter $\\lambda$', fontsize=11)
ax4.set_ylabel('Prospect Value', fontsize=11)
ax4.set_title('Mixed Gamble Acceptance vs. Loss Aversion', fontsize=12, fontweight='bold')
ax4.legend(fontsize=10)
ax4.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('decision_making_prospect.pdf', dpi=150, bbox_inches='tight')
plt.close()

# Calculate key statistics
loss_aversion_ratio = abs(v_loss / v_gain)
overweighting_01 = prelec_weighting(0.01, gamma_tk) / 0.01
underweighting_90 = prelec_weighting(0.90, gamma_tk) / 0.90
\end{pycode}

\begin{figure}[htbp]
\centering
\includegraphics[width=\textwidth]{decision_making_prospect.pdf}
\caption{Prospect theory components and predictions. (a) Value function showing loss aversion:
losses loom larger than equivalent gains by a factor of $\lambda=\py{f"{lambda_tk}"}$, resulting
in a loss aversion ratio of \py{f"{loss_aversion_ratio:.2f}"}. The function is steeper for losses
than gains, demonstrating diminishing sensitivity in both domains. (b) Prelec probability weighting
function ($\gamma=\py{f"{gamma_tk}"}$) exhibits inverse-S shape: small probabilities (e.g., 0.01)
are overweighted by a factor of \py{f"{overweighting_01:.2f}"}, while large probabilities (e.g., 0.90)
are underweighted to \py{f"{underweighting_90:.3f}"} of their objective value. (c) Fourfold pattern
in gains domain: risk-seeking for low-probability gains (lottery tickets) and risk-averse for
high-probability gains (insurance). (d) Mixed gamble acceptance threshold: agents with typical loss
aversion ($\lambda \approx \py{f"{lambda_tk}"}$) reject 50-50 mixed gambles even when expected value
is zero, requiring gains to be approximately twice as large as losses to accept.}
\label{fig:prospect}
\end{figure}

The prospect theory analysis reveals that with parameters $\alpha = \py{f"{alpha_tk}"}$,
$\beta = \py{f"{beta_tk}"}$, $\lambda = \py{f"{lambda_tk}"}$, and $\gamma = \py{f"{gamma_tk}"}$,
losses are experienced as \py{f"{loss_aversion_ratio:.2f}"} times more impactful than equivalent
gains. This loss aversion explains the rejection of symmetric mixed gambles and the endowment effect.

\section{Drift-Diffusion Model}

\subsection{Sequential Sampling Framework}

The drift-diffusion model (DDM) provides a mechanistic account of two-alternative forced choice
decisions, including both choice and response time distributions. Evidence accumulates over time
as a noisy process until reaching a decision boundary.

\begin{definition}[Drift-Diffusion Process]
The accumulated evidence $x(t)$ evolves according to:
\begin{equation}
dx = \mu \, dt + \sigma \, dW
\end{equation}
where $\mu$ is the drift rate (evidence quality), $\sigma$ is diffusion noise, and $dW$ is a
Wiener increment. Decision occurs when $x(t)$ reaches boundary $\pm a$.
\end{definition}

\begin{theorem}[Mean Decision Time]
For symmetric boundaries at $\pm a$ with starting point $z=0$, the mean decision time is:
\begin{equation}
\mathbb{E}[T] = \frac{a}{\mu} \tanh\left(\frac{a\mu}{\sigma^2}\right) + T_{er}
\end{equation}
where $T_{er}$ is non-decision time (encoding and response execution).
\end{theorem}

\subsection{Simulation of Decision Dynamics}

We simulate the drift-diffusion process using Euler-Maruyama discretization to generate response
time distributions and examine the speed-accuracy tradeoff. The model predicts that increasing
boundary separation increases accuracy but slows responses, while drift rate affects both.

\begin{pycode}
# Drift-Diffusion Model Simulation
def simulate_ddm_trial(drift_rate, boundary, noise=1.0, dt=0.001, non_decision_time=0.3):
    """Simulate a single DDM trial"""
    evidence = 0
    t = 0
    max_time = 10.0  # Prevent infinite loops

    while abs(evidence) < boundary and t < max_time:
        evidence += drift_rate * dt + noise * np.sqrt(dt) * np.random.randn()
        t += dt

    choice = 1 if evidence >= boundary else -1
    reaction_time = t + non_decision_time
    return choice, reaction_time

def simulate_ddm_experiment(drift_rate, boundary, n_trials=1000, noise=1.0, non_decision_time=0.3):
    """Simulate multiple DDM trials"""
    choices = []
    rts = []

    for _ in range(n_trials):
        choice, rt = simulate_ddm_trial(drift_rate, boundary, noise, non_decision_time=non_decision_time)
        choices.append(choice)
        rts.append(rt)

    return np.array(choices), np.array(rts)

# Simulate DDM with different parameters
drift_rates = [0.1, 0.3, 0.5]
boundaries = [0.8, 1.2, 1.6]
n_trials = 2000

fig3 = plt.figure(figsize=(14, 10))

# Drift rate manipulation
ax1 = plt.subplot(2, 3, 1)
for drift in drift_rates:
    choices, rts = simulate_ddm_experiment(drift, boundary=1.2, n_trials=n_trials)
    correct = choices[choices == 1]
    correct_rts = rts[choices == 1]
    ax1.hist(correct_rts, bins=50, alpha=0.6, density=True, label=f'$\\mu={drift}$')

ax1.set_xlabel('Response Time (s)', fontsize=11)
ax1.set_ylabel('Density', fontsize=11)
ax1.set_title('RT Distributions: Drift Rate Effect', fontsize=12, fontweight='bold')
ax1.legend(fontsize=10)
ax1.set_xlim(0, 5)

# Boundary manipulation
ax2 = plt.subplot(2, 3, 2)
for bound in boundaries:
    choices, rts = simulate_ddm_experiment(drift_rate=0.3, boundary=bound, n_trials=n_trials)
    correct = choices[choices == 1]
    correct_rts = rts[choices == 1]
    ax2.hist(correct_rts, bins=50, alpha=0.6, density=True, label=f'$a={bound}$')

ax2.set_xlabel('Response Time (s)', fontsize=11)
ax2.set_ylabel('Density', fontsize=11)
ax2.set_title('RT Distributions: Boundary Effect', fontsize=12, fontweight='bold')
ax2.legend(fontsize=10)
ax2.set_xlim(0, 5)

# Speed-accuracy tradeoff
ax3 = plt.subplot(2, 3, 3)
boundaries_sac = np.linspace(0.5, 2.0, 10)
accuracy_sac = []
mean_rt_sac = []

for bound in boundaries_sac:
    choices, rts = simulate_ddm_experiment(drift_rate=0.3, boundary=bound, n_trials=1000)
    accuracy_sac.append(np.mean(choices == 1))
    mean_rt_sac.append(np.mean(rts))

ax3.plot(mean_rt_sac, accuracy_sac, 'bo-', linewidth=2, markersize=6)
ax3.set_xlabel('Mean Response Time (s)', fontsize=11)
ax3.set_ylabel('Accuracy', fontsize=11)
ax3.set_title('Speed-Accuracy Tradeoff', fontsize=12, fontweight='bold')
ax3.grid(True, alpha=0.3)

# Sample trajectories
ax4 = plt.subplot(2, 3, 4)
drift_demo = 0.4
bound_demo = 1.5
n_trajectories = 20

for _ in range(n_trajectories):
    evidence_path = [0]
    t_path = [0]
    evidence = 0
    t = 0
    dt = 0.01

    while abs(evidence) < bound_demo and t < 5:
        evidence += drift_demo * dt + 1.0 * np.sqrt(dt) * np.random.randn()
        t += dt
        evidence_path.append(evidence)
        t_path.append(t)

    color = 'green' if evidence >= bound_demo else 'red'
    alpha = 0.3
    ax4.plot(t_path, evidence_path, color=color, alpha=alpha, linewidth=1)

ax4.axhline(bound_demo, color='blue', linestyle='--', linewidth=2, label='Upper boundary')
ax4.axhline(-bound_demo, color='orange', linestyle='--', linewidth=2, label='Lower boundary')
ax4.axhline(0, color='black', linestyle='-', linewidth=0.5)
ax4.set_xlabel('Time (s)', fontsize=11)
ax4.set_ylabel('Evidence $x(t)$', fontsize=11)
ax4.set_title(f'Sample Trajectories ($\\mu={drift_demo}$, $a={bound_demo}$)', fontsize=12, fontweight='bold')
ax4.legend(fontsize=10)
ax4.set_xlim(0, 5)
ax4.grid(True, alpha=0.3)

# Quantile probability plot
ax5 = plt.subplot(2, 3, 5)
choices_qp, rts_qp = simulate_ddm_experiment(drift_rate=0.4, boundary=1.2, n_trials=5000)
correct_rts = rts_qp[choices_qp == 1]
error_rts = rts_qp[choices_qp == -1]

quantiles = [0.1, 0.3, 0.5, 0.7, 0.9]
correct_quantiles = np.quantile(correct_rts, quantiles)
error_quantiles = np.quantile(error_rts, quantiles) if len(error_rts) > 10 else np.zeros(len(quantiles))

ax5.plot(quantiles, correct_quantiles, 'go-', linewidth=2, markersize=8, label='Correct')
if len(error_rts) > 10:
    ax5.plot(quantiles, error_quantiles, 'ro-', linewidth=2, markersize=8, label='Error')
ax5.set_xlabel('Quantile', fontsize=11)
ax5.set_ylabel('Response Time (s)', fontsize=11)
ax5.set_title('Quantile-Probability Plot', fontsize=12, fontweight='bold')
ax5.legend(fontsize=10)
ax5.grid(True, alpha=0.3)

# Drift rate vs accuracy/RT
ax6 = plt.subplot(2, 3, 6)
drift_range = np.linspace(0.05, 0.8, 15)
acc_by_drift = []
rt_by_drift = []

for drift in drift_range:
    choices_d, rts_d = simulate_ddm_experiment(drift, boundary=1.2, n_trials=500)
    acc_by_drift.append(np.mean(choices_d == 1))
    rt_by_drift.append(np.mean(rts_d))

ax6_twin = ax6.twinx()
line1 = ax6.plot(drift_range, acc_by_drift, 'b-o', linewidth=2, markersize=5, label='Accuracy')
line2 = ax6_twin.plot(drift_range, rt_by_drift, 'r-s', linewidth=2, markersize=5, label='Mean RT')

ax6.set_xlabel('Drift Rate $\\mu$', fontsize=11)
ax6.set_ylabel('Accuracy', fontsize=11, color='blue')
ax6_twin.set_ylabel('Mean RT (s)', fontsize=11, color='red')
ax6.set_title('Drift Rate Effects', fontsize=12, fontweight='bold')
ax6.tick_params(axis='y', labelcolor='blue')
ax6_twin.tick_params(axis='y', labelcolor='red')
ax6.grid(True, alpha=0.3)

lines = line1 + line2
labels = [l.get_label() for l in lines]
ax6.legend(lines, labels, fontsize=10, loc='center right')

plt.tight_layout()
plt.savefig('decision_making_ddm.pdf', dpi=150, bbox_inches='tight')
plt.close()

# Calculate key DDM statistics
drift_key = 0.4
boundary_key = 1.2
choices_key, rts_key = simulate_ddm_experiment(drift_key, boundary_key, n_trials=5000)
accuracy_key = np.mean(choices_key == 1)
mean_rt_key = np.mean(rts_key)
std_rt_key = np.std(rts_key)
median_rt_key = np.median(rts_key)
\end{pycode}

\begin{figure}[htbp]
\centering
\includegraphics[width=\textwidth]{decision_making_ddm.pdf}
\caption{Drift-diffusion model simulations of two-alternative forced choice. (a) Response time
distributions for different drift rates ($\mu$): higher evidence quality produces faster, more
accurate decisions with reduced RT variability. (b) Boundary separation effects: increasing the
decision threshold ($a$) from 0.8 to 1.6 produces slower but more conservative responses with
broader RT distributions. (c) Speed-accuracy tradeoff: adjusting boundary separation produces
the characteristic inverted-U relationship between mean RT and accuracy. (d) Twenty sample
evidence accumulation trajectories (green = correct, red = error) showing the stochastic nature
of the decision process. (e) Quantile-probability plot showing that correct responses are
systematically faster than errors at all quantiles. (f) Drift rate effects on both accuracy
and mean RT: with $a=\py{f"{boundary_key}"}$, a drift rate of $\mu=\py{f"{drift_key}"}$ yields
accuracy = \py{f"{accuracy_key:.3f}"} and mean RT = \py{f"{mean_rt_key:.3f}"} s (SD = \py{f"{std_rt_key:.3f}"} s).}
\label{fig:ddm}
\end{figure}

The drift-diffusion analysis with $\mu = \py{f"{drift_key}"}$ and $a = \py{f"{boundary_key}"}$
produces mean RT = \py{f"{mean_rt_key:.2f}"} s (median = \py{f"{median_rt_key:.2f}"} s) with
accuracy = \py{f"{accuracy_key:.1%}"}. The model successfully captures empirical RT distributions,
including their characteristic positive skew and the speed-accuracy tradeoff.

\section{Bayesian Decision Making}

\subsection{Decision Making Under Uncertainty}

Bayesian decision theory provides a normative framework for optimal decision making when beliefs
must be updated from noisy evidence. The decision maker maintains a posterior distribution over
states and chooses actions to minimize expected loss.

\begin{definition}[Bayes' Rule]
Given evidence $e$ and hypotheses $h_i$, the posterior probability is:
\begin{equation}
P(h_i | e) = \frac{P(e | h_i) P(h_i)}{\sum_j P(e | h_j) P(h_j)}
\end{equation}
\end{definition}

\begin{theorem}[Optimal Bayesian Decision]
The action $a^*$ that minimizes expected loss $L(a, h)$ is:
\begin{equation}
a^* = \arg\min_a \sum_h L(a, h) P(h | e)
\end{equation}
\end{theorem}

\subsection{Sequential Evidence Accumulation}

We simulate a medical diagnosis scenario where a physician must decide whether a patient has
a disease based on sequential test results. The Bayesian framework prescribes how beliefs should
be updated and when sufficient evidence has accumulated to make a decision.

\begin{pycode}
# Bayesian Decision Making
def bayesian_update(prior, likelihood_pos, likelihood_neg, evidence):
    """Update beliefs using Bayes' rule"""
    if evidence == 1:  # Positive test
        posterior_unnorm = prior * likelihood_pos
    else:  # Negative test
        posterior_unnorm = prior * likelihood_neg

    return posterior_unnorm / np.sum(posterior_unnorm)

# Medical diagnosis scenario
n_tests = 20
sensitivity = 0.85  # P(test+ | disease+)
specificity = 0.90  # P(test- | disease-)
prior_prob_disease = 0.10

# True state: patient has disease
true_state = 1

# Generate test results
np.random.seed(123)
test_results = []
for _ in range(n_tests):
    if true_state == 1:
        test_results.append(1 if np.random.rand() < sensitivity else 0)
    else:
        test_results.append(0 if np.random.rand() < specificity else 1)

# Sequential belief updating
beliefs = [prior_prob_disease]
for test in test_results:
    # Prior over disease states
    prior = np.array([1 - beliefs[-1], beliefs[-1]])

    # Likelihood: P(test result | disease state)
    if test == 1:
        likelihood = np.array([1 - specificity, sensitivity])
    else:
        likelihood = np.array([specificity, 1 - sensitivity])

    # Bayesian update
    posterior = bayesian_update(prior, likelihood, likelihood, test)
    beliefs.append(posterior[1])

beliefs = np.array(beliefs)

# Create Bayesian decision figure
fig4, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(14, 10))

# Plot 1: Sequential belief updating
ax1.plot(range(len(beliefs)), beliefs, 'b-o', linewidth=2, markersize=4)
ax1.axhline(0.5, color='red', linestyle='--', linewidth=2, label='Decision threshold')
ax1.axhline(prior_prob_disease, color='gray', linestyle=':', alpha=0.7, label='Prior')
ax1.fill_between(range(len(beliefs)), 0.5, beliefs, where=(beliefs > 0.5),
                  alpha=0.3, color='green', label='Diagnose disease')
ax1.fill_between(range(len(beliefs)), 0.5, beliefs, where=(beliefs <= 0.5),
                  alpha=0.3, color='blue', label='Diagnose healthy')
ax1.set_xlabel('Number of Tests', fontsize=11)
ax1.set_ylabel('P(Disease | Evidence)', fontsize=11)
ax1.set_title('Sequential Bayesian Belief Updating', fontsize=12, fontweight='bold')
ax1.legend(fontsize=9)
ax1.grid(True, alpha=0.3)
ax1.set_ylim(0, 1)

# Plot 2: Likelihood ratio and evidence strength
log_likelihood_ratios = []
for test in test_results:
    if test == 1:
        lr = sensitivity / (1 - specificity)
    else:
        lr = (1 - sensitivity) / specificity
    log_likelihood_ratios.append(np.log(lr))

cumulative_log_lr = np.cumsum(log_likelihood_ratios)

ax2.plot(range(1, len(cumulative_log_lr)+1), cumulative_log_lr, 'b-o', linewidth=2, markersize=4)
ax2.axhline(0, color='black', linestyle='-', linewidth=0.5)
ax2.axhline(np.log(3), color='green', linestyle='--', alpha=0.7, label='Moderate evidence (3:1)')
ax2.axhline(-np.log(3), color='red', linestyle='--', alpha=0.7, label='Moderate evidence (1:3)')
ax2.set_xlabel('Test Number', fontsize=11)
ax2.set_ylabel('Cumulative Log Likelihood Ratio', fontsize=11)
ax2.set_title('Evidence Accumulation (Log Odds)', fontsize=12, fontweight='bold')
ax2.legend(fontsize=9)
ax2.grid(True, alpha=0.3)

# Plot 3: Optimal decision boundary with costs
costs_fn = 1.0  # False negative (missed disease)
costs_fp_range = np.linspace(0.1, 2.0, 50)
optimal_thresholds = []

for cost_fp in costs_fp_range:
    # Optimal threshold: P(disease) where expected costs are equal
    # cost_fp * (1 - p) = cost_fn * p
    threshold = cost_fp / (cost_fp + costs_fn)
    optimal_thresholds.append(threshold)

ax3.plot(costs_fp_range, optimal_thresholds, 'b-', linewidth=2.5)
ax3.axvline(costs_fn, color='gray', linestyle=':', alpha=0.7, label='Equal costs')
ax3.axhline(0.5, color='gray', linestyle=':', alpha=0.7)
ax3.set_xlabel('Cost of False Positive', fontsize=11)
ax3.set_ylabel('Optimal Decision Threshold', fontsize=11)
ax3.set_title(f'Optimal Threshold vs. Cost Ratio (FN cost = {costs_fn})', fontsize=12, fontweight='bold')
ax3.legend(fontsize=10)
ax3.grid(True, alpha=0.3)

# Plot 4: ROC curve and decision boundaries
specificities = np.linspace(0.01, 0.99, 100)
sensitivities_roc = []

# For different decision thresholds
for spec in specificities:
    # Implied sensitivity at this operating point (simplified model)
    sens = 1 - (1 - spec) * 0.5
    sensitivities_roc.append(max(0, min(1, sens)))

ax4.plot(1 - specificities, sensitivities_roc, 'b-', linewidth=2.5, label='ROC curve')
ax4.plot([0, 1], [0, 1], 'k--', linewidth=1, label='Chance')
ax4.fill_between([0, 1], [0, 1], [1, 1], alpha=0.2, color='green', label='Better than chance')
ax4.scatter([1-specificity], [sensitivity], s=200, c='red', marker='*',
            edgecolors='black', linewidths=2, label=f'Test performance ({sensitivity},{specificity})', zorder=5)
ax4.set_xlabel('False Positive Rate (1 - Specificity)', fontsize=11)
ax4.set_ylabel('True Positive Rate (Sensitivity)', fontsize=11)
ax4.set_title('ROC Curve: Diagnostic Performance', fontsize=12, fontweight='bold')
ax4.legend(fontsize=9)
ax4.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('decision_making_bayesian.pdf', dpi=150, bbox_inches='tight')
plt.close()

# Calculate key statistics
final_belief = beliefs[-1]
n_positive_tests = np.sum(test_results)
final_log_odds = cumulative_log_lr[-1]
\end{pycode}

\begin{figure}[htbp]
\centering
\includegraphics[width=\textwidth]{decision_making_bayesian.pdf}
\caption{Bayesian decision making for medical diagnosis. (a) Sequential belief updating: starting
from a prior of P(disease) = \py{f"{prior_prob_disease}"}, the posterior probability evolves with
each test result (sensitivity = \py{f"{sensitivity}"}, specificity = \py{f"{specificity}"}). After
\py{n_tests} tests with \py{n_positive_tests} positive results, final belief reaches \py{f"{final_belief:.3f}"}.
The decision threshold of 0.5 determines when to diagnose disease. (b) Cumulative log likelihood ratio
tracks total evidence strength: values above log(3) represent moderate evidence for disease. Final log
odds = \py{f"{final_log_odds:.2f}"}, indicating strong evidence. (c) Optimal decision threshold depends
on cost ratio: when false negatives are more costly than false positives, threshold should decrease to
favor diagnosing disease. With equal costs, threshold = 0.5. (d) ROC curve shows tradeoff between true
positive rate (sensitivity) and false positive rate (1-specificity). Test performance (red star) indicates
better-than-chance discrimination, with area under curve representing overall diagnostic accuracy.}
\label{fig:bayesian}
\end{figure}

The Bayesian analysis shows that after \py{n_tests} sequential tests, the posterior probability
of disease is \py{f"{final_belief:.3f}"}, representing a log odds of \py{f"{final_log_odds:.2f}"}.
This demonstrates how weak evidence (single test) accumulates to produce strong diagnostic confidence,
formalizing the intuitive notion that multiple corroborating tests increase certainty.

\section{Reinforcement Learning}

\subsection{Value-Based Decision Making}

Reinforcement learning provides a framework for learning to make decisions through trial and error.
Agents learn value functions that predict expected future reward and use these to guide action selection.

\begin{definition}[Q-Learning]
The action-value function $Q(s, a)$ represents expected return from taking action $a$ in state $s$.
The Q-learning update rule is:
\begin{equation}
Q(s_t, a_t) \leftarrow Q(s_t, a_t) + \alpha \left[ r_t + \gamma \max_{a'} Q(s_{t+1}, a') - Q(s_t, a_t) \right]
\end{equation}
where $\alpha$ is learning rate and $\gamma$ is discount factor.
\end{definition}

\begin{definition}[Softmax Action Selection]
Actions are selected probabilistically according to:
\begin{equation}
P(a|s) = \frac{\exp(\beta Q(s,a))}{\sum_{a'} \exp(\beta Q(s,a'))}
\end{equation}
where $\beta$ is inverse temperature controlling exploration-exploitation.
\end{definition}

\subsection{Multi-Armed Bandit Simulation}

We simulate a classic reinforcement learning task where an agent must learn which of several slot
machines (arms) provides the highest reward. The agent faces the exploration-exploitation dilemma:
whether to exploit known good options or explore alternatives that might be better.

\begin{pycode}
# Reinforcement Learning: Multi-armed bandit
class MultiArmedBandit:
    def __init__(self, n_arms, true_values):
        self.n_arms = n_arms
        self.true_values = true_values

    def pull_arm(self, arm):
        """Pull arm and receive reward"""
        return self.true_values[arm] + np.random.randn() * 0.5

class QLearningAgent:
    def __init__(self, n_arms, learning_rate=0.1, temperature=1.0):
        self.n_arms = n_arms
        self.Q = np.zeros(n_arms)  # Q-values
        self.alpha = learning_rate
        self.beta = temperature
        self.action_counts = np.zeros(n_arms)

    def select_action(self):
        """Softmax action selection"""
        exp_Q = np.exp(self.beta * self.Q)
        probs = exp_Q / np.sum(exp_Q)
        return np.random.choice(self.n_arms, p=probs)

    def update(self, action, reward):
        """Q-learning update"""
        self.Q[action] += self.alpha * (reward - self.Q[action])
        self.action_counts[action] += 1

# Setup bandit task
n_arms = 4
true_values = np.array([1.0, 2.5, 1.5, 0.5])  # True mean rewards
n_trials = 500

# Simulate different learning rates
learning_rates = [0.01, 0.1, 0.5]
temperature = 2.0

fig5, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(14, 10))

for lr in learning_rates:
    bandit = MultiArmedBandit(n_arms, true_values)
    agent = QLearningAgent(n_arms, learning_rate=lr, temperature=temperature)

    Q_history = [agent.Q.copy()]
    rewards_history = []

    for trial in range(n_trials):
        action = agent.select_action()
        reward = bandit.pull_arm(action)
        agent.update(action, reward)

        Q_history.append(agent.Q.copy())
        rewards_history.append(reward)

    Q_history = np.array(Q_history)

    # Plot Q-value learning
    for arm in range(n_arms):
        if lr == 0.1:  # Only plot for one learning rate to avoid clutter
            ax1.plot(Q_history[:, arm], label=f'Arm {arm+1} (true: {true_values[arm]:.1f})', alpha=0.8)

ax1.axhline(0, color='black', linestyle='-', linewidth=0.5)
for arm in range(n_arms):
    ax1.axhline(true_values[arm], color='gray', linestyle=':', alpha=0.5)
ax1.set_xlabel('Trial', fontsize=11)
ax1.set_ylabel('Q-value', fontsize=11)
ax1.set_title(f'Q-Learning Dynamics ($\\alpha={0.1}$, $\\beta={temperature}$)', fontsize=12, fontweight='bold')
ax1.legend(fontsize=9)
ax1.grid(True, alpha=0.3)

# Compare learning rates
ax2_data = []
for lr in learning_rates:
    bandit = MultiArmedBandit(n_arms, true_values)
    agent = QLearningAgent(n_arms, learning_rate=lr, temperature=temperature)

    cumulative_reward = 0
    cumulative_rewards = []

    for trial in range(n_trials):
        action = agent.select_action()
        reward = bandit.pull_arm(action)
        agent.update(action, reward)
        cumulative_reward += reward
        cumulative_rewards.append(cumulative_reward)

    ax2.plot(cumulative_rewards, linewidth=2, label=f'$\\alpha={lr}$')

# Optimal cumulative reward
optimal_reward_per_trial = np.max(true_values)
optimal_cumulative = np.arange(1, n_trials+1) * optimal_reward_per_trial
ax2.plot(optimal_cumulative, 'k--', linewidth=2, label='Optimal')

ax2.set_xlabel('Trial', fontsize=11)
ax2.set_ylabel('Cumulative Reward', fontsize=11)
ax2.set_title('Learning Rate Comparison', fontsize=12, fontweight='bold')
ax2.legend(fontsize=10)
ax2.grid(True, alpha=0.3)

# Temperature effects on exploration
temperatures = [0.5, 2.0, 5.0]
ax3_colors = ['blue', 'green', 'red']

for temp, color in zip(temperatures, ax3_colors):
    bandit = MultiArmedBandit(n_arms, true_values)
    agent = QLearningAgent(n_arms, learning_rate=0.1, temperature=temp)

    action_history = []

    for trial in range(n_trials):
        action = agent.select_action()
        reward = bandit.pull_arm(action)
        agent.update(action, reward)
        action_history.append(action)

    # Plot action selection frequency over time
    window = 50
    action_probs = []
    for i in range(window, n_trials):
        recent_actions = action_history[i-window:i]
        probs = [recent_actions.count(a) / window for a in range(n_arms)]
        action_probs.append(probs)

    action_probs = np.array(action_probs)
    best_arm = np.argmax(true_values)

    ax3.plot(range(window, n_trials), action_probs[:, best_arm],
             color=color, linewidth=2, label=f'$\\beta={temp}$')

ax3.set_xlabel('Trial', fontsize=11)
ax3.set_ylabel('P(Select Best Arm)', fontsize=11)
ax3.set_title('Temperature Effects: Exploration vs. Exploitation', fontsize=12, fontweight='bold')
ax3.legend(fontsize=10)
ax3.grid(True, alpha=0.3)
ax3.set_ylim(0, 1)

# Final Q-values and action selection probabilities
ax4_agent = QLearningAgent(n_arms, learning_rate=0.1, temperature=2.0)
bandit_final = MultiArmedBandit(n_arms, true_values)

for trial in range(n_trials):
    action = ax4_agent.select_action()
    reward = bandit_final.pull_arm(action)
    ax4_agent.update(action, reward)

x_pos = np.arange(n_arms)
width = 0.35

bars1 = ax4.bar(x_pos - width/2, true_values, width, label='True values', color='steelblue', edgecolor='black')
bars2 = ax4.bar(x_pos + width/2, ax4_agent.Q, width, label='Learned Q-values', color='coral', edgecolor='black')

# Add action selection probabilities as text
exp_Q = np.exp(ax4_agent.beta * ax4_agent.Q)
action_probs_final = exp_Q / np.sum(exp_Q)
for i, (bar, prob) in enumerate(zip(bars2, action_probs_final)):
    height = bar.get_height()
    ax4.text(bar.get_x() + bar.get_width()/2., height + 0.1,
             f'P={prob:.2f}', ha='center', va='bottom', fontsize=9)

ax4.set_xlabel('Arm', fontsize=11)
ax4.set_ylabel('Value', fontsize=11)
ax4.set_title(f'Learned Values After {n_trials} Trials', fontsize=12, fontweight='bold')
ax4.set_xticks(x_pos)
ax4.set_xticklabels([f'Arm {i+1}' for i in range(n_arms)])
ax4.legend(fontsize=10)
ax4.grid(True, alpha=0.3, axis='y')

plt.tight_layout()
plt.savefig('decision_making_rl.pdf', dpi=150, bbox_inches='tight')
plt.close()

# Calculate key RL statistics
final_Q = ax4_agent.Q
Q_errors = np.abs(final_Q - true_values)
mean_Q_error = np.mean(Q_errors)
best_arm = np.argmax(true_values)
learned_best = np.argmax(final_Q)
correct_identification = (learned_best == best_arm)
\end{pycode}

\begin{figure}[htbp]
\centering
\includegraphics[width=\textwidth]{decision_making_rl.pdf}
\caption{Reinforcement learning in a four-armed bandit task with true arm values of
[\py{', '.join([f"{v:.1f}" for v in true_values])}]. (a) Q-learning dynamics with $\alpha=0.1$
and $\beta=\py{f"{temperature}"}$: Q-values converge toward true values through reward feedback,
with learning speed depending on action selection frequency. Dotted lines show true values.
(b) Learning rate effects on cumulative reward: faster learning ($\alpha=0.5$) initially
outperforms but may be less stable; slower learning ($\alpha=0.01$) produces smoother but
delayed convergence. Black dashed line shows optimal performance (always selecting best arm).
(c) Temperature parameter controls exploration-exploitation balance: low temperature ($\beta=0.5$)
produces high exploration and slow convergence to optimal arm; high temperature ($\beta=5.0$)
rapidly exploits current best estimate but may get stuck on suboptimal choices. (d) Final learned
Q-values after \py{n_trials} trials (orange) compared to true values (blue), with action selection
probabilities shown above bars. Agent correctly identifies best arm (Arm 2) with Q-value
\py{f"{final_Q[1]:.2f}"} vs. true \py{f"{true_values[1]:.2f}"}, mean absolute error = \py{f"{mean_Q_error:.2f}"}.}
\label{fig:rl}
\end{figure}

The reinforcement learning simulation demonstrates value-based decision making through Q-learning.
After \py{n_trials} trials with $\alpha = 0.1$ and $\beta = \py{f"{temperature}"}$, the agent's
learned Q-values have mean absolute error of \py{f"{mean_Q_error:.3f}"} relative to true values,
and it correctly identifies the best arm (Arm \py{best_arm + 1}) with selection probability
\py{f"{action_probs_final[best_arm]:.2f}"}. This illustrates how experience-based learning
converges to near-optimal policies.

\section{Model Comparison and Integration}

The five computational frameworks presented here address different aspects of decision making:

\begin{itemize}
\item \textbf{Expected Utility Theory} provides the normative benchmark for rational choice under risk
\item \textbf{Prospect Theory} captures systematic deviations including loss aversion and probability distortion
\item \textbf{Drift-Diffusion Models} explain response time distributions and speed-accuracy tradeoffs
\item \textbf{Bayesian Decision Making} prescribes optimal belief updating under uncertainty
\item \textbf{Reinforcement Learning} accounts for value learning through experience
\end{itemize}

\begin{remark}[Integration of Frameworks]
Recent research integrates these approaches: drift-diffusion can be derived from Bayesian sequential
analysis; reinforcement learning models can incorporate prospect theory's asymmetric value function;
neural implementations suggest the brain approximates Bayesian inference through sampling mechanisms.
\end{remark}

\section{Conclusions}

This computational analysis reveals the mechanistic underpinnings of human decision making:

\begin{enumerate}
\item Expected utility with $\alpha = \py{f"{alpha_moderate}"}$ predicts risk premium of
\$\py{f"{abs(rp_moderate):.2f}"} for a 50-50 gamble of $\pm$\$50
\item Prospect theory with $\lambda = \py{f"{lambda_tk}"}$ captures loss aversion, with losses
weighted \py{f"{loss_aversion_ratio:.2f}"} times gains
\item Drift-diffusion with $\mu = \py{f"{drift_key}"}$ and $a = \py{f"{boundary_key}"}$ produces
mean RT = \py{f"{mean_rt_key:.2f}"} s and accuracy = \py{f"{accuracy_key:.1%}"}
\item Bayesian updating with sensitivity = \py{f"{sensitivity}"} and specificity = \py{f"{specificity}"}
achieves posterior confidence of \py{f"{final_belief:.3f}"} after \py{n_tests} tests
\item Q-learning with $\alpha = 0.1$ converges to within \py{f"{mean_Q_error:.3f}"} of true values
after \py{n_trials} trials
\end{enumerate}

These computational models provide quantitative frameworks for understanding both optimal decision
strategies and systematic deviations from rationality observed in human behavior.

\section*{References}

\begin{itemize}
\item von Neumann, J., \& Morgenstern, O. (1944). \textit{Theory of Games and Economic Behavior}. Princeton University Press.
\item Kahneman, D., \& Tversky, A. (1979). Prospect theory: An analysis of decision under risk. \textit{Econometrica}, 47(2), 263-291.
\item Tversky, A., \& Kahneman, D. (1992). Advances in prospect theory: Cumulative representation of uncertainty. \textit{Journal of Risk and Uncertainty}, 5(4), 297-323.
\item Ratcliff, R., \& McKoon, G. (2008). The diffusion decision model: Theory and data for two-choice decision tasks. \textit{Neural Computation}, 20(4), 873-922.
\item Ratcliff, R., Smith, P. L., Brown, S. D., \& McKoon, G. (2016). Diffusion decision model: Current issues and history. \textit{Trends in Cognitive Sciences}, 20(4), 260-281.
\item Bogacz, R., Brown, E., Moehlis, J., Holmes, P., \& Cohen, J. D. (2006). The physics of optimal decision making: A formal analysis of models of performance in two-alternative forced-choice tasks. \textit{Psychological Review}, 113(4), 700-765.
\item Edwards, W., Lindman, H., \& Savage, L. J. (1963). Bayesian statistical inference for psychological research. \textit{Psychological Review}, 70(3), 193-242.
\item Sutton, R. S., \& Barto, A. G. (2018). \textit{Reinforcement Learning: An Introduction} (2nd ed.). MIT Press.
\item Daw, N. D., O'Doherty, J. P., Dayan, P., Seymour, B., \& Dolan, R. J. (2006). Cortical substrates for exploratory decisions in humans. \textit{Nature}, 441(7095), 876-879.
\item Dayan, P., \& Daw, N. D. (2008). Decision theory, reinforcement learning, and the brain. \textit{Cognitive, Affective, \& Behavioral Neuroscience}, 8(4), 429-453.
\item Gold, J. I., \& Shadlen, M. N. (2007). The neural basis of decision making. \textit{Annual Review of Neuroscience}, 30, 535-574.
\item Busemeyer, J. R., \& Diederich, A. (2010). \textit{Cognitive Modeling}. Sage Publications.
\item Prelec, D. (1998). The probability weighting function. \textit{Econometrica}, 66(3), 497-527.
\item Rescorla, R. A., \& Wagner, A. R. (1972). A theory of Pavlovian conditioning: Variations in the effectiveness of reinforcement and nonreinforcement. In A. H. Black \& W. F. Prokasy (Eds.), \textit{Classical Conditioning II: Current Research and Theory} (pp. 64-99). Appleton-Century-Crofts.
\item Glimcher, P. W. (2011). \textit{Foundations of Neuroeconomic Analysis}. Oxford University Press.
\item Lee, M. D., \& Wagenmakers, E. J. (2013). \textit{Bayesian Cognitive Modeling: A Practical Course}. Cambridge University Press.
\item Wald, A. (1947). \textit{Sequential Analysis}. John Wiley \& Sons.
\item Tversky, A., \& Fox, C. R. (1995). Weighing risk and uncertainty. \textit{Psychological Review}, 102(2), 269-283.
\item Forstmann, B. U., Ratcliff, R., \& Wagenmakers, E. J. (2016). Sequential sampling models in cognitive neuroscience: Advantages, applications, and extensions. \textit{Annual Review of Psychology}, 67, 641-666.
\item Niv, Y. (2009). Reinforcement learning in the brain. \textit{Journal of Mathematical Psychology}, 53(3), 139-154.
\end{itemize}

\end{document}