Introducing Two Compartment exchange Model

.pre-commit-config.yaml

@@ -1,9 +1,9 @@
 repos:
   - repo: https://github.com/astral-sh/ruff-pre-commit
-    rev: v0.4.8  # Use the latest version
+    rev: v0.4.8 # Use the latest version
     hooks:
       - id: ruff
-        args: [--fix]  # Optional: to enable lint fixes
+        args: [--fix] # Optional: to enable lint fixes
       - id: ruff-format
   - repo: https://github.com/pre-commit/pre-commit-hooks
     rev: v4.6.0

docs/examples/tissue/plot_two_cxm.py

@@ -0,0 +1,39 @@
+"""
+=============
+The Two Compartment Exchange Model
+=============
+
+Simulating tissue concentrations from two compartment models with different settings.
+"""
+
+import matplotlib.pyplot as plt
+
+# %%
+# Import necessary packages
+import numpy as np
+import osipi
+
+# %%
+# Generate Parker AIF with default settings.
+
+# Define time points in units of seconds - in this case we use a time
+# resolution of 1 sec and a total duration of 6 minutes.
+t = np.arange(0, 6 * 60, 1, dtype=float)
+
+# Create an AIF with default settings
+ca = osipi.aif_parker(t)
+
+# %%
+# Plot the tissue concentrations for an extracellular volume fraction
+# of 0.2, plasma volume fraction of 0.3, extraction fraction of 0.15
+# and flow rate of 0.2 ml/min
+Ps = 0.05  # Permeability surface area product in ml/min
+Fp = 10  # Flow rate in ml/min
+Ve = 0.2  # Extracellular volume fraction
+Vp = 0.01  # Plasma volume fraction
+ct = osipi.two_compartment_exchange_model(t, ca, Fp=Fp, Ps=Ps, Ve=Ve, Vp=Vp)
+plt.plot(t, ct, "b-", label=f" Fp = {Fp},Ps = {Ps}, Ve = {Ve}, Vp = {Vp}")
+plt.xlabel("Time (sec)")
+plt.ylabel("Tissue concentration (mM)")
+plt.legend()
+plt.show()

docs/references/models/tissue_models/2cxm.md

@@ -0,0 +1,17 @@
+# osipi.Two Compartment Exchange Model
+
+::: osipi.two_compartment_exchange_model
+
+
+## Example using `osipi.two_compartment_exchange_model`
+
+<figure class="mkd-glr-thumbcontainer" tooltip="Simulating tissue concentrations from two compartment exchange model with different settings.">
+    <img alt="The Two Compartment Exchange Model" src="..\..\..\..\generated\gallery\tissue\images\thumb\mkd_glr_plot_two_cxm_thumb.png" />
+    <figcaption class="caption">
+        <span class="caption-text">
+            <a class="reference internal" href="..\..\..\..\generated\gallery\tissue\plot_two_cxm">
+                <span class="std std-ref">The Two Compartment Exchange Model</span>
+            </a>
+        </span>
+    </figcaption>
+</figure>

docs/references/models/tissue_models/index.md

@@ -1,3 +1,4 @@
 
 - [tofts](tofts.md)
 - [extended_tofts](extended_tofts.md)
+- [two_cxm](2cxm.md)

mkdocs.yml

@@ -75,6 +75,7 @@ nav:
               - references/models/tissue_models/index.md
               - osipi.tofts: references/models/tissue_models/tofts.md
               - osipi.extended_tofts: references/models/tissue_models/extended_tofts.md
+              - osipi.two_compartment_exchange: references/models/tissue_models/2cxm.md
 
   - Examples: generated/gallery
 

src/osipi/__init__.py

@@ -7,7 +7,8 @@
 
 from ._tissue import (
     tofts,
-    extended_tofts
+    extended_tofts,
+    two_compartment_exchange_model
 )
 
 from ._signal import (

src/osipi/_tissue.py

@@ -2,6 +2,7 @@
 
 import numpy as np
 from scipy.interpolate import interp1d
+from scipy.signal import convolve
 
 from ._convolution import exp_conv
 
@@ -332,3 +333,84 @@ def extended_tofts(
                 ct = ct_func(t)
 
     return ct
+
+
+def two_compartment_exchange_model(
+    t: np.ndarray,
+    ca: np.ndarray,
+    Fp: float,
+    Ps: float,
+    Ve: float,
+    Vp: float,
+    Ta: float = 30.0,
+) -> np.ndarray:
+    """Two compartment exchange model
+
+    Tracer flows from the AIF to the blood plasma compartment; two-way leakage
+    between the plasma and extracellular compartments(EES) is permitted.
+
+    Args:
+        t: 1D array of times(s).
+        ca: Arterial concentrations in mM for each time point in t.
+        Fp: Blood plasma flow rate into a unit tissue volume in ml/min.
+        Ps: Permeability surface area product in ml/min.
+        Ve: Extracellular volume fraction.
+        Vp: Plasma volume fraction.
+        Ta: Arterial delay time, i.e., difference in onset time between tissue curve and AIF in units of sec.
+
+    Returns:
+        Ct: Tissue concentrations in mM
+
+    References:
+        - Lexicon url: https://osipi.github.io/OSIPI_CAPLEX/perfusionModels/#indicator-kinetic-models
+        - Lexicon code: M.IC1.009
+        - OSIPI name: Two Compartment Exchange Model
+        - Adapted from contributions by: MJT_UoEdinburgh_UK
+
+    """
+    if not np.allclose(np.diff(t), np.diff(t)[0]):
+        warnings.warn(
+            ("Non-uniform time spacing detected. Time array may be" " resampled."),
+            stacklevel=2,
+        )
+
+    # Convert units
+    fp_per_s = Fp / (60.0 * 100.0)
+    ps_per_s = Ps / 60.0
+
+    # Calculate the impulse response function
+    v = Ve + Vp
+
+    # Mean transit time
+    T = v / fp_per_s
+    tc = Vp / fp_per_s
+    te = Ve / ps_per_s
+
+    sig_p = ((T + te) + np.sqrt((T + te) ** 2 - 4 * tc * te)) / (2 * tc * te)
+    sig_n = ((T + tc) - np.sqrt((T + tc) ** 2 - 4 * tc * te)) / (2 * tc * te)
+
+    # Calculate the impulse response function for the plasma compartment and EES
+
+    irf_cp = (
+        Vp
+        * sig_p
+        * sig_n
+        * ((1 - te * sig_n) * np.exp(-t * sig_n) + (te * sig_p - 1.0) * np.exp(-t * sig_p))
+        / (sig_p - sig_n)
+    )
+
+    irf_ce = Ve * sig_p * sig_n * (np.exp(-t * sig_n) - np.exp(-t * sig_p)) / (sig_p - sig_n)
+
+    irf_cp[[0]] /= 2
+    irf_ce[[0]] /= 2
+
+    dt = np.min(np.diff(t))
+
+    # get concentration in plasma and EES
+
+    Cp = dt * convolve(ca, irf_cp, mode="full", method="auto")[: len(t)]
+    Ce = dt * convolve(ca, irf_ce, mode="full", method="auto")[: len(t)]
+
+    # get tissue concentration
+    Ct = Cp + Ce
+    return Ct

tests/test_tissue.py

@@ -124,8 +124,13 @@ def test_tissue_extended_tofts():
     assert np.allclose(ct_conv, ca * 0.3, rtol=1e-4, atol=1e-3)
 
 
+def test_tissue_two_compartment_exchange_model():
+    pass
+
+
 if __name__ == "__main__":
     test_tissue_tofts()
     test_tissue_extended_tofts()
+    test_tissue_two_compartment_exchange_model()
 
     print("All tissue concentration model tests passed!!")