21:null 22:["$","title",null,{"children":"Two-compartment pharmacokinetic model with intermittent dialysis | Sami Safadi | Observable"}] 27:Tc11,# Two-compartment pharmacokinetic model with intermittent dialysis This notebook simulates how the concentration of a drug (or toxin) changes in a patient who receives **intermittent hemodialysis**. The body is split into two compartments: - **Central compartment (blood)**: volume V_c, concentration C_c. Drug is removed from here, by body elimination and by the dialyzer. - **Tissue (peripheral) compartment**: volume V_t, concentration C_t. It exchanges drug with the blood but is not cleared directly. The sliders control every parameter, and the plot updates automatically. ## Model The ODEs are written out in the next cell. In words: 1. **Central compartment.** Concentration falls through *body elimination* (CL_body) and *dialysis* (CL_dial, only while a session is running). It also falls as drug moves into tissue at rate k₁₂, and rises as drug returns from tissue at rate k₂₁. The factor V_t/V_c converts the tissue concentration into an equivalent amount per unit of blood volume. 2. **Tissue compartment.** Drug flows in from blood at rate k₁₂ (scaled by V_c/V_t) and flows back out at rate k₂₁. 3. **Dialysis schedule.** The dialyzer clearance acts as a square wave. It is on for the first t_s hours of every interval T, for a total of N sessions, and off otherwise. Because the tissue compartment refills the blood after a session, the blood concentration **rebounds** once dialysis stops. This is a well-known effect for drugs that distribute into tissue. ## Parameters | Symbol | Control | Meaning | |---|---|---| | V_c, V_t | Vc, Vt | Central and tissue volumes (L) | | CL_body | CL_body | Body elimination clearance (L/hr) | | CL_dial | CL_dial | Dialyzer clearance while a session is running (L/hr) | | k₁₂, k₂₁ | k12, k21 | Blood→tissue and tissue→blood transfer rates (1/hr) | | t_s | Session length | Duration of each session (hr) | | T | Interval | Time between session starts (hr) | | N | Number of sessions | Total sessions delivered | | C_c(0), C_t(0) | Initial Cc, Ct | Starting concentrations (mg/L) | Default values: V_c = 10, V_t = 30, CL_body = 2, CL_dial = 15, k₁₂ = 0.5, k₂₁ = 0.3, 3 sessions of 4 hr every 48 hr, simulated for 144 hr, starting at 100 mg/L in both compartments. ## Numerical method The ODEs are integrated with a fixed-step **4th-order Runge–Kutta** scheme (step 0.01 hr). The dialysis on/off state is evaluated at the midpoint of each step, so the sharp session boundaries are not skipped. This is the JavaScript equivalent of the original Python `solve_ivp` code, which used `max_step` for the same reason. ## Notebook structure 1. **Model and solver:** the `simulate` function. 2. **Controls:** parameter sliders and display toggles. 3. **Plot:** concentration over time, with the dialysis sessions shaded.28:T83b,const p = view(Inputs.form( { Vc: Inputs.range([1, 50], {label: "Central volume Vc (L)", value: 10, step: 0.5}), Vt: Inputs.range([1, 150], {label: "Tissue volume Vt (L)", value: 30, step: 0.5}), CL_body: Inputs.range([0, 10], {label: "Body elimination (L/hr)", value: 2, step: 0.1}), k12: Inputs.range([0, 2], {label: "k12 blood→tissue (1/hr)", value: 0.5, step: 0.01}), k21: Inputs.range([0, 2], {label: "k21 tissue→blood (1/hr)", value: 0.3, step: 0.01}), CL_dial: Inputs.range([0, 30], {label: "Dialyzer clearance (L/hr)", value: 15, step: 0.5}), t_session: Inputs.range([0.5, 12], {label: "Session length (hr)", value: 4, step: 0.5}), t_interval: Inputs.range([12, 96], {label: "Interval (hr)", value: 48, step: 1}), num_sessions: Inputs.range([0, 10], {label: "Number of sessions", value: 3, step: 1}), t_max: Inputs.range([24, 480], {label: "Simulation length (hr)", value: 144, step: 6}), C0c: Inputs.range([0, 200], {label: "Initial Cc (mg/L)", value: 100, step: 1}), C0t: Inputs.range([0, 200], {label: "Initial Ct (mg/L)", value: 100, step: 1}), logScale: Inputs.toggle({label: "Log y-axis", value: false}), showTissue: Inputs.toggle({label: "Show tissue curve", value: true}) }, { template: (inputs) => { const card = (title, color, keys) => htl.html`

${title}${keys.map(k => inputs[k])}

`; return htl.html`

${[ card("Patient (distribution & elimination)", "#4269d0", ["Vc", "Vt", "CL_body", "k12", "k21"]), card("Dialysis schedule", "#efb118", ["CL_dial", "t_session", "t_interval", "num_sessions", "t_max"]), card("Initial state", "#3ca951", ["C0c", "C0t"]), card("Display", "#9498a0", ["logScale", "showTissue"]) ]}

`; } } ));29:T41f,const sim = simulate(p); const longData = sim.flatMap((d) => [ {t: d.t, conc: d.Cc, compartment: "Central (blood)"}, ...(p.showTissue ? [{t: d.t, conc: d.Ct, compartment: "Tissue"}] : []) ]); const sessionStarts = d3.range(p.num_sessions).map((i) => i * p.t_interval); const interactivePlot = display(Plot.plot({ title: "Drug concentration with intermittent dialysis", subtitle: `${p.num_sessions} × ${p.t_session} h sessions every ${p.t_interval} h (shaded); Vc = ${p.Vc} L, Vt = ${p.Vt} L`, x: {label: "Time (hr)"}, y: { label: "Concentration (mg/L)", type: p.logScale ? "log" : "linear", grid: true, domain: p.logScale ? [Math.max(0.1, d3.min(longData, (d) => d.conc) * 0.8), d3.max(longData, (d) => d.conc) * 1.1] : undefined }, color: {legend: true, domain: ["Central (blood)", "Tissue"]}, marks: [ Plot.rectX(sessionStarts, {x1: (d) => d, x2: (d) => d + p.t_session, fill: "gray", fillOpacity: 0.2}), Plot.lineY(longData, {x: "t", y: "conc", stroke: "compartment", tip: true}), Plot.ruleY([0]) ] }));2a:T410,function simulate({Vc, Vt, CL_body, CL_dial, k12, k21, t_session, t_interval, num_sessions, t_max, C0c, C0t}, dt = 0.01) { const f = (Cc, Ct, on) => [ -(CL_body / Vc) * Cc - (on ? CL_dial / Vc : 0) * Cc - k12 * Cc + k21 * (Vt / Vc) * Ct, k12 * (Vc / Vt) * Cc - k21 * Ct ]; let Cc = C0c, Ct = C0t; const out = []; const n = Math.round(t_max / dt); for (let i = 0; i <= n; i++) { const t = i * dt; if (i % 10 === 0) out.push({t, Cc, Ct, dialyzing: false}); // midpoint of the step decides the dialysis switch state const tm = t + dt / 2; const on = (tm % t_interval) <= t_session && Math.floor(tm / t_interval) < num_sessions; if (i % 10 === 0) out[out.length - 1].dialyzing = on; const k1 = f(Cc, Ct, on); const k2 = f(Cc + dt/2*k1[0], Ct + dt/2*k1[1], on); const k3 = f(Cc + dt/2*k2[0], Ct + dt/2*k2[1], on); const k4 = f(Cc + dt*k3[0], Ct + dt*k3[1], on); Cc += dt/6*(k1[0] + 2*k2[0] + 2*k3[0] + k4[0]); Ct += dt/6*(k1[1] + 2*k2[1] + 2*k3[1] + k4[1]); } return out; }2b:[] 25:{"id":"96d1ae348e31a347","slug":"two-compartment-pharmacokinetic-model-with-intermittent-dialysis","title":"Two-compartment pharmacokinetic model with intermittent dialysis","thumbnail":"2b56aa2beaeb4ca1e7fbd128ac6e7f15e93aa426817090ddf81470cfb6e51dc0","access_level":"public","create_time":"$D2026-10-01T18:28:05.966Z","update_time":"$D2026-10-01T19:25:29.113Z","view_time":null,"first_public_time":"$D2026-10-01T19:18:32.544Z","has_unpublished_edits":null,"fork":false,"stars":0,"starred":false,"owner":{"type":"team","id":"1d804b626bd08da7","login":"samisaf","name":"Sami Safadi","avatar_url":"https://avatars-next.observableusercontent.com/a2adb3491825954cb8c600a2e043624c461158aa29ae90abb5d3fab14240f3de","tier":"starter_2024"},"version":43,"publish_version":43,"role":null,"authored":false,"authors":[{"type":"individual","id":"a139fe285b2a1eaf","login":"samisaf","name":"Sami Safadi","avatar_url":"https://avatars-next.observableusercontent.com/a2adb3491825954cb8c600a2e043624c461158aa29ae90abb5d3fab14240f3de","tier":null}],"members":[{"type":"team","id":"1d804b626bd08da7","login":"samisaf","name":"Sami Safadi","avatar_url":"https://avatars-next.observableusercontent.com/a2adb3491825954cb8c600a2e043624c461158aa29ae90abb5d3fab14240f3de","tier":"starter_2024","role":"editor"}],"comments":[],"custom_thumbnail":false,"trashed":false,"lapsed":null,"subscribed":null,"forks":0,"collection":null,"type":"notebook","body":{"titling":"auto","title":"Two-compartment pharmacokinetic model with intermittent dialysis","cells":[{"id":1,"value":"$27","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":5,"value":"\\begin{aligned}\n\\frac{dC_c}{dt} &= -\\frac{CL_{body}}{V_c}\\,C_c \\;-\\; \\frac{CL_{dial}(t)}{V_c}\\,C_c \\;-\\; k_{12}\\,C_c \\;+\\; k_{21}\\,\\frac{V_t}{V_c}\\,C_t \\\\[1.2em]\n\\frac{dC_t}{dt} &= k_{12}\\,\\frac{V_c}{V_t}\\,C_c \\;-\\; k_{21}\\,C_t \\\\[1.2em]\nCL_{dial}(t) &= \\begin{cases} CL_{dial} & \\text{if } (t \\bmod T) \\le t_s \\ \\text{ and } \\ \\lfloor t/T \\rfloor < N \\\\ 0 & \\text{otherwise} \\end{cases} \\\\[1.2em]\nC_c(0) &= C_{c,0}, \\qquad C_t(0) = C_{t,0}\n\\end{aligned}","pinned":false,"hidden":false,"mode":"tex","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":10,"value":"## Interactive Plot","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":3,"value":"$28","pinned":true,"hidden":false,"mode":"js","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":12,"value":"$29","pinned":true,"hidden":false,"mode":"js","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":8,"value":"## Example 2: a drug with a high volume of distribution\n\nDrugs such as digoxin distribute extensively into tissue (V_t ≫ V_c). The settings below use V_c = 15 L and V_t = 500 L. k₂₁ is chosen so that k₁₂·V_c = k₂₁·V_t, which means both compartments start at equilibrium. These values are illustrative, not fitted to a specific drug.\n\nDialysis can only clear the small blood pool. Blood concentration drops during a session, then **rebounds** as the large tissue reservoir refills the blood. Little drug is removed overall.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":9,"value":"const highVd = {\n Vc: 15, Vt: 500, CL_body: 3, CL_dial: 15,\n k12: 0.6, k21: 0.018,\n t_session: 4, t_interval: 48, num_sessions: 3, t_max: 144,\n C0c: 100, C0t: 100\n};\nconst highVdData = simulate(highVd).flatMap((d) => [\n {t: d.t, conc: d.Cc, compartment: \"Central (blood)\"},\n {t: d.t, conc: d.Ct, compartment: \"Tissue\"}\n]);\nconst highVdPlot = display(Plot.plot({\n title: \"High volume of distribution: strong post-dialysis rebound\",\n subtitle: `${highVd.num_sessions} × ${highVd.t_session} h sessions every ${highVd.t_interval} h (shaded); Vc = ${highVd.Vc} L, Vt = ${highVd.Vt} L`,\n x: {label: \"Time (hr)\"},\n y: {label: \"Concentration (mg/L)\", grid: true},\n color: {legend: true, domain: [\"Central (blood)\", \"Tissue\"]},\n marks: [\n Plot.rectX(d3.range(highVd.num_sessions).map((i) => i * highVd.t_interval), {x1: (d) => d, x2: (d) => d + highVd.t_session, fill: \"gray\", fillOpacity: 0.2}),\n Plot.lineY(highVdData, {x: \"t\", y: \"conc\", stroke: \"compartment\", tip: true}),\n Plot.ruleY([0])\n ]\n}));","pinned":true,"hidden":false,"mode":"js","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"},{"id":11,"value":"$2a","pinned":true,"hidden":false,"mode":"ts","author":"$undefined","summary":"$undefined","providerMetadata":"$undefined","attachments":"$undefined","output":"$undefined","database":"$undefined"}],"files":"$W2b","stdlib":"2"}} 26:"$25:body"