# Tripartite Synapse - Multi-Scale Computational Model # ===================================================== # Presynaptic + Postsynaptic perspectives, fully integrated. # # Change log: # ORIG - present from the original document # NEW - added in the missing-behavior integration pass # DET - deterministic Ca2+-driven vesicle release # NKA - explicit Na/K-ATPase V_pre decay and ATP cost # POST-ATP - postsynaptic Ca2+ dynamics and ATP loop # FIX - corrections applied in this pass: # * NT_released_this_window accumulator (was missing entirely) # * k_rec_fast / k_rec_slow converted to /s, recruitment moved to Loop 2 # * dt_slow_s added # * mGluR now reads NT_released_this_window (not NT_cleft snapshot) # * IP3 now reads NT_released_this_window (not cleared_NT residual) # * wave_active flag + conversion_efficiency boost on astrocyte wave # * CDI rise gated to spike window only # # Clock structure: # Loop 1 - dt = 1 ms (Ca2+, vesicle release, traces, postsynaptic fast) # Loop 2 - dt = 1000 ms (astrocyte clearance, eCB, mGluR, recruitment) # Loop 3 - dt = 60000 ms (glutamine shuttle, metabolic health) # # ======================================================================= # THREE CLOSED LOOPS # ======================================================================= # # PRESYNAPTIC: # NT loop : release (ms) -> cleft -> astrocyte clearance (s) -> # glutamine shuttle (min) -> RP refill -> RRP -> release # Ca2+ loop : VGCC influx (ms) -> Tr_Ca -> recruitment speed (s) -> # eCB retrograde from post (s) -> VGCC suppression # ATP loop : NKA + pump costs (ms) -> ATP_demand (min) -> ATP_level -> # pump_scale -> Ca2+ clearance rate -> CDI recovery # # POSTSYNAPTIC: # NT detection loop : NT_cleft -> AMPA -> V_post -> desensitization -> # reduces next response # Ca2+ coincidence : NMDA (NT + V_post) -> Ca_post -> eCB -> pre brake # ATP loop : NKA + PMCA costs (ms) -> ATP_demand_post (min) -> # ATP_level_post -> pump_scale_post -> Ca_post clearance # # SHARED: # eCB_level : post synthesises -> pre reads (retrograde brake) # NT_cleft : pre releases -> post detects -> astrocyte clears # Glucose_level : astrocyte supplies both sides from same budget # # ======================================================================= # METABOLIC SILENCING CASCADE (presynaptic) # ======================================================================= # [CASCADE 1] HIGH FIRING -> VESICLE DEPLETION (~seconds) # release rate >> recruitment rate -> N_RRP -> 0 # [CASCADE 2] HIGH FIRING -> ATP DEPLETION (~minutes) # NKA + PMCA + docking demand > glucose-driven supply # [CASCADE 3] LOW ATP -> PUMP FAILURE # pump_scale = Hill(ATP_level) -> cleared_PMCA/SERCA fall # [CASCADE 4] PUMP FAILURE -> RESIDUAL Ca2+ STAYS HIGH # Ca_micro persists between spikes # [CASCADE 5] RESIDUAL Ca2+ -> CDI LOCKS VGCCs SHUT # CDI rise (spike only) + recovery blocked by Ca2+ -> CDI -> 1 # [CASCADE 6] SYNAPSE SILENCES (excitotoxicity protection) # effective_conductance = N_VGCC*(1-eCB)*(1-CDI)*(1-mGluR*alpha) # -> 0; NCX auto-reset when drive stops # # POSTSYNAPTIC ATP CASCADE (no CDI equivalent -> dangerous): # [POST-ATP 1] HIGH V_post + NMDA -> ATP_demand_post rises # [POST-ATP 2] ATP_level_post falls -> pump_scale_post falls # [POST-ATP 3] Ca_post clearance slows -> Ca_post stays elevated # [POST-ATP 4] Ca_post > eCB_threshold without real coincidence # -> false retrograde signal suppresses presynapse # [POST-ATP 5] Critically low ATP_post -> runaway Ca_post -> excitotoxicity # ======================================================================= import numpy as np # ----------------------------------------------------------------------- # CLOCK # ----------------------------------------------------------------------- dt = 1.0 # ms dt_slow = 1000.0 # ms dt_meta = 60_000.0 # ms High_Freq_Multiplier = int(dt_slow / dt) # 1000 Metabolic_Multiplier = int(dt_meta / dt) # 60000 dt_s = dt / 1000.0 # 0.001 s/step - for /s rate constants in Loop 1 dt_slow_s = dt_slow / 1000.0 # 1.0 s/step - for /s rate constants in Loop 2 # ----------------------------------------------------------------------- # PRESYNAPTIC PARAMETERS # ----------------------------------------------------------------------- # -- Voltage / membrane -- tau_V_pre = 2.0 # ms - AP waveform decay (Na/K-ATPase recharge) V_pre_peak = 1.0 # a.u. - normalised AP peak V_rest = 0.0 # a.u. - resting potential V_pre_voltage = -10.0 # mV - driving force for compute_flux NKA_cost_per_AP = 0.002 # ATP units per AP (dominant drain at high rates) # -- Ca2+ influx & buffering -- N_VGCC = 100 # number of VGCCs (ceiling of effective_conductance) k_flux = 0.05 # Ca2+ influx per open channel per unit driving force B_total = 1.0 # total buffer capacity (normalised) tau_buffer_rebind = 200.0 # ms - buffer recharge time constant # -- Ca2+ clearance (/ms constants) -- k_PMCA = 0.03 # ATP-dependent primary pump k_NCX = 0.10 # ATP-independent floor k_SERCA = 0.01 # ATP-dependent ER pump ATP_half = 0.3 # Hill half-saturation for presynaptic pumps ATP_cost_PMCA = 0.0005 # ATP per unit Ca2+ extruded by PMCA ATP_cost_SERCA = 0.0002 # ATP per unit Ca2+ pumped into ER ATP_cost_docking = 0.001 # ATP per vesicle docked (RP->RRP) # -- Deterministic release (Hill + NT suppression) -- k_rel = 0.5 # max releasable fraction of RRP per spike KD_rel = 1.0 # half-saturation [Ca2+] n_rel = 4 # Hill cooperativity (synaptotagmin-1) NT_suppression_weight = 0.3 # max NT_cleft brake on release fraction NT_suppression_sat = 50.0 # NT_cleft level that saturates suppression # -- CDI -- k_CDI_rise = 0.8 # /s - CDI build rate (applied * dt_s, spike only) Ca_micro_saturation = 2.0 # normalisation ceiling for CDI recovery k_CDI_rec = 0.015 # /s - CDI de-inactivation rate (applied * dt_s) # -- Vesicle pools -- Max_RRP = 20 Max_RP = 200 # -- Calcium trace -- tau_Tr_Ca = 1000.0 # ms T_high = 0.6 # Tr_Ca threshold -> fast recruitment T_low = 0.2 # Tr_Ca threshold -> slow recruitment # -- RP->RRP recruitment (/s, runs in Loop 2) -- k_rec_fast = 5.0 # /s - fast recruitment (at Tr_Ca > T_high) k_rec_slow = 0.5 # /s - slow recruitment (at Tr_Ca < T_low) # -- NT accumulator for Loop 2 signals -- NT_window_sat = 40.0 # vesicles/s that saturates mGluR and IP3 # at 20 Hz releasing ~2/spike = 40/s # -- eCB retrograde brake -- tau_eCB_rise = 2000.0 tau_eCB_decay = 10_000.0 eCB_threshold = 0.7 # Ca_post level that triggers eCB synthesis # -- mGluR presynaptic autoreceptor -- Km_mGluR = 0.5 tau_mGluR = 2000.0 # ms alpha_mGluR = 0.4 # max fractional VGCC suppression # -- Astrocyte / IP3 -- tau_IP3 = 3000.0 # ms IP3_threshold = 0.8 wave_boost = 0.2 # conversion_efficiency boost when wave fires tau_wave_decay = 2 # metabolic cycles before boost decays back # -- Glutamine shuttle -- conversion_efficiency_base = 0.8 # -- NT cleft -- tau_NT_decay = 5.0 # ms # ----------------------------------------------------------------------- # POSTSYNAPTIC PARAMETERS # ----------------------------------------------------------------------- # -- NMDA coincidence detection -- k_NMDA = 0.08 # Ca_post influx per unit NT * (1 - Mg_block) per ms V_NMDA_half = 0.3 # V_post at which Mg block is 50% lifted # -- Ca_post clearance -- k_Ca_post_clear = 0.05 # /ms - ATP-dependent PMCA in spine k_Ca_post_NCX = 0.02 # /ms - ATP-independent NCX floor ATP_half_post = 0.3 # Hill half-saturation for postsynaptic pumps # -- Postsynaptic ATP costs -- NKA_cost_per_bAP_post = 0.002 # ATP per unit V_post per s (continuous) ATP_cost_Ca_post_pump = 0.0005 # ATP per unit Ca_post cleared ATP_demand_scale_post = 50.0 # normalisation (same as presynaptic) # -- Receptor desensitization -- tau_membrane = 20.0 # ms tau_desensitization = 500.0 # ms # ----------------------------------------------------------------------- # HELPER FUNCTIONS # ----------------------------------------------------------------------- def compute_flux(conductance, voltage): return k_flux * conductance * abs(voltage) def deterministic_release(N_RRP, Ca_micro, NT_cleft): # Hill equation: Ca2+ sensor cooperativity (synaptotagmin-1, n=4) Ca_n = Ca_micro ** n_rel release_frac = k_rel * Ca_n / (Ca_n + KD_rel ** n_rel) # NT suppression: physical crowding + fast local autoreceptors NT_norm = min(1.0, NT_cleft / NT_suppression_sat) release_frac = release_frac * (1.0 - NT_suppression_weight * NT_norm) release_frac = np.clip(release_frac, 0.0, 1.0) return max(0.0, release_frac * N_RRP) def map_trace_to_speed(Tr_Ca): # Returns /s recruitment rate based on Tr_Ca level if Tr_Ca > T_high: return k_rec_fast elif Tr_Ca < T_low: return k_rec_slow else: t = (Tr_Ca - T_low) / (T_high - T_low) return k_rec_slow + t * (k_rec_fast - k_rec_slow) def compute_pump_atp_factor(atp, atp_half): # Hill function: ATP gates pump speed (shared by pre and post) return (atp ** 2) / (atp ** 2 + atp_half ** 2) def compute_EPSP(receptor_conductance): return receptor_conductance * 0.1 def compute_astrocyte_metabolic_health(Glucose_level, ATP_demand_accumulated, demand_scale=50.0): # Converts glucose supply and accumulated demand into ATP_level (0->1) # and conversion_efficiency (0->1). Both sides use this function with # their own demand accumulators but the same Glucose_level — shared # metabolic vulnerability. health = np.clip(Glucose_level - ATP_demand_accumulated / demand_scale, 0.0, 1.0) return health, health # (conversion_efficiency, ATP_level) def trigger_slow_astrocyte_calcium_wave(): # Placeholder - gliotransmitter release over ~10 s pass # ----------------------------------------------------------------------- # STATE VARIABLES # ----------------------------------------------------------------------- # -- Presynaptic membrane -- V_pre_state = 0.0 # -- Presynaptic Ca2+ -- Ca_micro = 0.0 Ca_ER = 0.5 Ca_buffer_bound = 0.0 B_free = B_total # -- CDI -- CDI_factor = 0.0 # -- Vesicle pools -- N_RRP = 15.0 N_RP = 150.0 # -- Calcium trace -- Tr_Ca = 0.0 # -- NT cleft -- NT_cleft = 0.0 # -- NT accumulator for slow signals -- # FIX: this was missing. Accumulates every ms in Loop 1, # consumed by mGluR and IP3 in Loop 2, reset each second. NT_released_this_window = 0.0 # -- Postsynaptic membrane + receptors -- V_post = 0.0 receptor_conductance = 0.0 Desensitization_level = 0.0 V_post_history = [] # -- Postsynaptic Ca2+ (spine compartment) -- Ca_post = 0.0 # Driven by NMDA coincidence (NT + V_post). Cleared by PMCA (ATP-gated) # and NCX (always). Drives eCB synthesis. No CDI equivalent -> # elevated Ca_post under ATP failure has no self-limiting mechanism. # -- Retrograde / autoreceptor -- eCB_level = 0.0 mGluR_activation = 0.0 # -- Astrocyte -- IP3 = 0.0 wave_active = 0 # countdown: cycles remaining of wave boost Glutamine_pool = 50.0 # -- Presynaptic ATP -- ATP_level = 1.0 ATP_demand = 0.0 conversion_efficiency = conversion_efficiency_base Glucose_level = 1.0 # set < 1.0 to engage metabolic silencing # -- Postsynaptic ATP -- ATP_level_post = 1.0 # separate pool; same glucose budget as presynaptic ATP_demand_post = 0.0 # accumulates from NKA (V_post) and PMCA (Ca_post) # ----------------------------------------------------------------------- # MAIN SIMULATION LOOP # ----------------------------------------------------------------------- def run_simulation(spike_train, total_steps, bAP_train=None): """ spike_train : list of int - presynaptic AP timestep indices total_steps : int bAP_train : list of int - postsynaptic bAP timestep indices (optional) if None, no bAPs are delivered """ global V_pre_state global Ca_micro, Ca_ER, Ca_buffer_bound, B_free global CDI_factor global N_RRP, N_RP, Tr_Ca, NT_cleft, NT_released_this_window global V_post, receptor_conductance, Desensitization_level, V_post_history global Ca_post global eCB_level, mGluR_activation global IP3, wave_active, Glutamine_pool global ATP_level, ATP_demand, conversion_efficiency, Glucose_level global ATP_level_post, ATP_demand_post log = {k: [] for k in [ "V_pre_state", "Ca_micro", "Ca_ER", "CDI_factor", "B_free", "N_RRP", "N_RP", "Tr_Ca", "NT_cleft", "V_post", "Ca_post", "eCB_level", "mGluR_activation", "released_NT", "ATP_level", "ATP_demand", "ATP_level_post", "ATP_demand_post", ]} spike_set = set(spike_train) bAP_set = set(bAP_train) if bAP_train else set() for step in range(total_steps): # ============================================================== # LOOP 1 — HIGH-FREQUENCY (dt = 1 ms) # ============================================================== V_pre = 1 if step in spike_set else 0 bAP = 1 if step in bAP_set else 0 released_NT = 0.0 # -- 1A. PRESYNAPTIC MEMBRANE / Na-K-ATPase ------------------- # AP fires: membrane jumps to peak, then decays with tau_V_pre. # Ca2+ influx uses V_pre_state (continuous) not binary V_pre, # giving a temporal influx profile that tapers as membrane repolarises. if V_pre == 1: V_pre_state = V_pre_peak ATP_demand += NKA_cost_per_AP # dominant presynaptic ATP cost V_pre_state += (V_rest - V_pre_state) * dt / tau_V_pre # -- 1B. PRESYNAPTIC Ca2+ INFLUX ------------------------------ # Three multiplicative brakes on effective_conductance: # eCB_level : retrograde brake from postsynapse (Loop 2) # CDI_factor : Ca2+-dependent inactivation (below) # mGluR_activation : autoreceptor brake (Loop 2) effective_conductance = ( N_VGCC * (1.0 - eCB_level) * (1.0 - CDI_factor) * (1.0 - mGluR_activation * alpha_mGluR) ) raw_influx = compute_flux(effective_conductance, V_pre_state) # Buffer proteins capture a fraction immediately (fast sponge). # B_free -> 0 during sustained bursting -> capture_fraction -> 0 # -> full raw_influx enters Ca_micro (CASCADE 4 acceleration). capture_fraction = B_free / B_total captured = raw_influx * capture_fraction B_free = max(0.0, B_free - captured) Ca_buffer_bound += captured Ca_micro += (raw_influx - captured) # -- 1C. VESICLE RELEASE -------------------------------------- # Deterministic: Hill Ca2+ sensor * NT suppression * N_RRP. # Runs every ms that Ca_micro > 0 (release profile follows Ca2+ # transient, not locked to spike flag). if N_RRP > 0 and Ca_micro > 0: released_NT = deterministic_release(N_RRP, Ca_micro, NT_cleft) released_NT = min(released_NT, N_RRP) N_RRP -= released_NT NT_cleft += released_NT # FIX: accumulate for Loop 2 mGluR and IP3 signals. # This is the only correct way to feed slow signals from fast # events — snapshot of NT_cleft at Loop 2 time would be ~0 # because passive diffusion has already cleared it. NT_released_this_window += released_NT # Passive NT diffusion out of cleft each ms. NT_cleft *= (1.0 - dt / tau_NT_decay) NT_cleft = max(0.0, NT_cleft) # -- 1D. PRESYNAPTIC Ca2+ CLEARANCE --------------------------- # pump_scale: Hill(ATP_level) — bridges Loop 3 ATP to Loop 1 clearance. # NCX is ATP-independent (floor); PMCA and SERCA are ATP-gated. pump_scale = compute_pump_atp_factor(ATP_level, ATP_half) cleared_PMCA = k_PMCA * Ca_micro * pump_scale cleared_NCX = k_NCX * Ca_micro cleared_SERCA = k_SERCA * Ca_micro * pump_scale Ca_micro -= (cleared_PMCA + cleared_NCX + cleared_SERCA) Ca_micro = max(0.0, Ca_micro) Ca_ER += cleared_SERCA ATP_demand += ATP_cost_PMCA * cleared_PMCA ATP_demand += ATP_cost_SERCA * cleared_SERCA # Buffer recharge: bound Ca2+ slowly re-releases back to cytosol. # During pump failure this sustains Ca_micro elevation (CASCADE 4). rebind_flux = Ca_buffer_bound * dt / tau_buffer_rebind Ca_micro += rebind_flux Ca_buffer_bound = max(0.0, Ca_buffer_bound - rebind_flux) B_free = B_total - Ca_buffer_bound # -- 1E. CDI — RISE (spike only) AND RECOVERY (every ms) ------ # RISE: Ca2+ entering through open channels inactivates them locally. # Gated to spike window — requires channels to be open. # (Running every ms was wrong: CDI needs Ca2+ flowing through # the channel, not ambient cytosolic Ca2+.) if V_pre == 1: CDI_factor += k_CDI_rise * Ca_micro * dt_s # RECOVERY: continuous, suppressed when Ca_micro is high. # Self-locking: pump failure -> Ca_micro high -> recovery ~0 # -> CDI_factor -> 1 -> effective_conductance -> 0 (CASCADE 5-6). CDI_recovery_rate = k_CDI_rec * (1.0 - Ca_micro / Ca_micro_saturation) CDI_factor = np.clip(CDI_factor - CDI_recovery_rate * dt_s, 0.0, 1.0) # -- 1F. CALCIUM TRACE ---------------------------------------- # Leaky integrator — integrates full Ca2+ waveform every ms # including inter-spike clearance. Drives Loop 2 recruitment speed. Tr_Ca = Tr_Ca + (Ca_micro - Tr_Ca / tau_Tr_Ca) * dt # -- 1G. POSTSYNAPTIC: NT DETECTION & AMPA -------------------- # Desensitization reduces effective NT — sustained NT exposure # progressively silences receptors (postsynaptic equivalent of CDI). effective_NT = released_NT * (1.0 - Desensitization_level) receptor_conductance += effective_NT * 0.05 receptor_conductance *= (1.0 - dt / tau_membrane) V_post += compute_EPSP(receptor_conductance) - (V_post / tau_membrane) * dt V_post = max(0.0, V_post) Desensitization_level += NT_cleft * 0.001 * dt Desensitization_level -= (Desensitization_level / tau_desensitization) * dt Desensitization_level = np.clip(Desensitization_level, 0.0, 1.0) V_post_history.append(V_post) if len(V_post_history) > 5000: V_post_history.pop(0) # -- 1H. POSTSYNAPTIC: NMDA COINCIDENCE DETECTION ------------- # Ca_post enters only when BOTH conditions hold simultaneously: # (1) NT_cleft > 0 — ligand gate (glutamate present) # (2) V_post elevated — voltage gate (Mg2+ block lifted) # Mg block removal is a sigmoid of V_post. # bAP (backpropagating AP) boosts V_post further, enabling # full NMDA opening and larger Ca_post surge. V_post_effective = V_post + (bAP * 0.5) # bAP adds depolarisation Mg_block_removal = V_post_effective / (V_post_effective + V_NMDA_half) NMDA_Ca_influx = k_NMDA * NT_cleft * Mg_block_removal Ca_post += NMDA_Ca_influx # Postsynaptic NKA: membrane recharge cost proportional to V_post. # [POST-ATP 1] Dominant postsynaptic ATP drain at high activity. ATP_demand_post += NKA_cost_per_bAP_post * V_post * dt_s # -- 1I. POSTSYNAPTIC: Ca_post CLEARANCE ---------------------- # pump_scale_post: Hill(ATP_level_post) — same structure as presynaptic. # NCX is ATP-independent floor (enables auto-reset after ATP recovery). # [POST-ATP 3] When pump_scale_post falls, Ca_post stays elevated -> # eCB threshold crossed without genuine coincidence -> false retrograde. pump_scale_post = compute_pump_atp_factor(ATP_level_post, ATP_half_post) cleared_Ca_post_pump = k_Ca_post_clear * Ca_post * pump_scale_post cleared_Ca_post_NCX = k_Ca_post_NCX * Ca_post Ca_post -= (cleared_Ca_post_pump + cleared_Ca_post_NCX) Ca_post = max(0.0, Ca_post) # [POST-ATP 2] ATP cost of postsynaptic PMCA. ATP_demand_post += ATP_cost_Ca_post_pump * cleared_Ca_post_pump # -- RECORD --------------------------------------------------- log["V_pre_state"].append(V_pre_state) log["Ca_micro"].append(Ca_micro) log["Ca_ER"].append(Ca_ER) log["CDI_factor"].append(CDI_factor) log["B_free"].append(B_free) log["N_RRP"].append(N_RRP) log["N_RP"].append(N_RP) log["Tr_Ca"].append(Tr_Ca) log["NT_cleft"].append(NT_cleft) log["V_post"].append(V_post) log["Ca_post"].append(Ca_post) log["eCB_level"].append(eCB_level) log["mGluR_activation"].append(mGluR_activation) log["released_NT"].append(released_NT) log["ATP_level"].append(ATP_level) log["ATP_demand"].append(ATP_demand) log["ATP_level_post"].append(ATP_level_post) log["ATP_demand_post"].append(ATP_demand_post) # ============================================================== # LOOP 2 — SLOW / ASTROCYTE (dt_slow = 1 s) # ============================================================== if (step % High_Freq_Multiplier) == 0: # Astrocyte EAAT clearance — active NT removal from cleft. cleared_NT = NT_cleft * 0.3 NT_cleft = max(0.0, NT_cleft - cleared_NT) # FIX: IP3 integrates NT_released_this_window (total release # since last Loop 2), not the post-diffusion NT_cleft residual # which is ~0 by the time Loop 2 runs. IP3 += NT_released_this_window - (IP3 / tau_IP3) * dt_slow IP3 = max(0.0, IP3) if IP3 > IP3_threshold: trigger_slow_astrocyte_calcium_wave() # FIX: wave boosts conversion_efficiency in the next mins cycle. # The astrocyte responds to heavy load by upregulating its # recycling machinery — shipping more glutamine back to the # presynapse. Boost decays over tau_wave_decay metabolic cycles. wave_active = tau_wave_decay # FIX: mGluR reads NT_released_this_window (accumulated release # load), not NT_cleft snapshot. NT_cleft is ~0 at Loop 2 time # due to diffusion; the accumulator correctly represents the # burst load the autoreceptor has sensed during this window. NT_window_norm = min(1.0, NT_released_this_window / NT_window_sat) mGluR_target = NT_window_norm mGluR_activation += (mGluR_target - mGluR_activation) * (dt_slow / tau_mGluR) mGluR_activation = np.clip(mGluR_activation, 0.0, 1.0) # FIX: reset accumulator for next window. NT_released_this_window = 0.0 # eCB retrograde synthesis: now driven by Ca_post (spine Ca2+), # not V_post_history. The actual eCB synthesis in the spine is # triggered by Ca2+-dependent enzymes (DAGL, PLC), not voltage. # Under normal conditions Ca_post only rises with coincidence. # Under POST-ATP failure Ca_post stays elevated without genuine # coincidence -> false retrograde signal (POST-ATP 4). recent_Ca_post = (np.mean(log["Ca_post"][-2000:]) if len(log["Ca_post"]) >= 2000 else (np.mean(log["Ca_post"]) if log["Ca_post"] else 0.0)) eCB_signal = max(0.0, recent_Ca_post - eCB_threshold) if eCB_signal > 0: eCB_level += eCB_signal * (dt_slow / tau_eCB_rise) else: eCB_level -= eCB_level * (dt_slow / tau_eCB_decay) eCB_level = np.clip(eCB_level, 0.0, 1.0) # FIX: RP->RRP recruitment moved here from Loop 1. # Biological timescale: vesicle docking and priming take seconds, # not milliseconds. k_rec_fast/slow are /s; * dt_slow_s = 1.0 s # gives dimensionless per-step fraction — no hidden unit scaling. current_recruitment_rate = map_trace_to_speed(Tr_Ca) # /s refill_amount = (current_recruitment_rate * dt_slow_s * N_RP * (Max_RRP - N_RRP) / Max_RRP) refill_amount = max(0.0, refill_amount) refill_amount = min(refill_amount, N_RP) N_RRP = min(N_RRP + refill_amount, Max_RRP) N_RP = max(0.0, N_RP - refill_amount) ATP_demand += ATP_cost_docking * refill_amount # ============================================================== # LOOP 3 — METABOLIC (dt_meta = 1 min) # ============================================================== if (step % Metabolic_Multiplier) == 0: # Presynaptic ATP: glucose supply minus accumulated demand. conversion_efficiency, ATP_level = compute_astrocyte_metabolic_health( Glucose_level, ATP_demand ) ATP_demand = 0.0 # FIX: wave boost applied to conversion_efficiency. # Astrocyte calcium wave (triggered by high IP3) upregulates # glutamine synthetase -> faster NT recycling -> more RP refill. # Boost decays over tau_wave_decay cycles. if wave_active > 0: conversion_efficiency = min(1.0, conversion_efficiency + wave_boost) wave_active -= 1 # Glutamine shuttle: astrocyte converts cleared NT to glutamine, # presynapse repackages it into vesicles -> N_RP replenished. refill_RP = Glutamine_pool * conversion_efficiency N_RP = min(Max_RP, N_RP + refill_RP) Glutamine_pool = max(0.0, Glutamine_pool - refill_RP) # Postsynaptic ATP: same glucose budget, own demand accumulator. # Both sides draw from Glucose_level -> shared metabolic vulnerability. # Presynaptic silence reduces NT -> less NMDA -> less Ca_post -> # less ATP_demand_post: presynaptic protection indirectly # protects the postsynapse. _, ATP_level_post = compute_astrocyte_metabolic_health( Glucose_level, ATP_demand_post, ATP_demand_scale_post ) ATP_demand_post = 0.0 return log # ----------------------------------------------------------------------- # EXAMPLE USAGE # ----------------------------------------------------------------------- if __name__ == "__main__": import matplotlib.pyplot as plt total_steps = 10_000 # 10 seconds # Presynaptic 20 Hz burst for 2 s. spike_train = list(range(0, 2000, 50)) # Postsynaptic bAPs coincident with every 5th presynaptic spike # (simulates partial coincidence for NMDA activation). bAP_train = list(range(0, 2000, 250)) results = run_simulation(spike_train, total_steps, bAP_train) t = np.arange(total_steps) * dt fig, axes = plt.subplots(8, 1, figsize=(12, 18), sharex=True) fig.suptitle("Tripartite Synapse — Presynaptic + Postsynaptic", fontsize=13) axes[0].plot(t, results["V_pre_state"], color="slateblue", lw=0.8) axes[0].set_ylabel("V_pre") axes[0].set_title("Presynaptic membrane (AP waveform)", fontsize=9, loc="left") axes[1].plot(t, results["Ca_micro"], color="darkorange", lw=0.8) axes[1].set_ylabel("[Ca2+] pre") axes[1].set_title("CASCADE 4 — presynaptic Ca2+", fontsize=9, loc="left") axes[2].plot(t, results["CDI_factor"], color="firebrick", lw=0.8, label="CDI") axes[2].plot(t, results["B_free"], color="steelblue", lw=0.8, label="Buffer free") axes[2].set_ylabel("CDI / Buffer") axes[2].set_title("CASCADE 5 — CDI lock-out", fontsize=9, loc="left") axes[2].legend(fontsize=8) axes[3].plot(t, results["N_RRP"], color="teal", lw=0.8, label="RRP") axes[3].plot(t, results["N_RP"], color="purple", lw=0.8, label="RP") axes[3].set_ylabel("Vesicles") axes[3].set_title("CASCADE 1 — vesicle depletion", fontsize=9, loc="left") axes[3].legend(fontsize=8) axes[4].plot(t, results["NT_cleft"], color="darkgreen", lw=0.8, label="NT cleft") axes[4].plot(t, results["mGluR_activation"], color="saddlebrown", lw=0.8, label="mGluR") axes[4].plot(t, results["eCB_level"], color="crimson", lw=0.8, label="eCB") axes[4].set_ylabel("Cleft / Feedback") axes[4].set_title("CASCADE 6 — three brakes on conductance", fontsize=9, loc="left") axes[4].legend(fontsize=8) axes[5].plot(t, results["V_post"], color="navy", lw=0.8, label="V_post") axes[5].plot(t, results["Ca_post"], color="coral", lw=0.8, label="Ca_post (spine)") axes[5].set_ylabel("Postsynaptic") axes[5].set_title("Postsynaptic potential + NMDA spine Ca2+", fontsize=9, loc="left") axes[5].legend(fontsize=8) axes[6].plot(t, results["ATP_level"], color="goldenrod", lw=0.8, label="ATP pre") axes[6].plot(t, results["ATP_level_post"], color="darkorange", lw=0.8, label="ATP post") axes[6].set_ylabel("ATP level") axes[6].set_title("CASCADE 2 / POST-ATP — presynaptic and postsynaptic ATP", fontsize=9, loc="left") axes[6].legend(fontsize=8) axes[7].plot(t, results["ATP_demand"], color="tomato", lw=0.8, label="demand pre") axes[7].plot(t, results["ATP_demand_post"], color="orangered", lw=0.8, label="demand post") axes[7].set_ylabel("ATP demand") axes[7].set_title("Accumulated ATP demand (resets each min cycle)", fontsize=9, loc="left") axes[7].set_xlabel("Time (ms)") axes[7].legend(fontsize=8) plt.tight_layout() plt.savefig("./synapse_simulation.png", dpi=150) plt.close() print("Done.")