From 2188bf59f76189469096d19da3db2b0d14eab42a Mon Sep 17 00:00:00 2001 From: Mohamed Nasser Date: Sun, 15 Sep 2024 19:11:22 +0300 Subject: [PATCH] Added use case and corrected units --- docs/examples/tissue/plot_two_cxm.py | 14 ++++++---- docs/references/models/tissue_models/2cxm.md | 3 ++- src/osipi/__init__.py | 2 +- src/osipi/_tissue.py | 27 +++++++++++++------- tests/test_tissue.py | 26 +++++++++++-------- 5 files changed, 46 insertions(+), 26 deletions(-) diff --git a/docs/examples/tissue/plot_two_cxm.py b/docs/examples/tissue/plot_two_cxm.py index e6b7d2a..d5572d2 100644 --- a/docs/examples/tissue/plot_two_cxm.py +++ b/docs/examples/tissue/plot_two_cxm.py @@ -18,7 +18,7 @@ # 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) +t = np.arange(0, 6 * 60, 1, dtype=float) # Create an AIF with default settings ca = osipi.aif_parker(t) @@ -27,12 +27,16 @@ # 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 -E = 0.15 # Extraction fraction -Fp = 0.2 # Flow rate in ml/min +E = 0.3 # Extraction fraction +Fp = 0.1 # Flow rate in ml/min Ve = 0.2 # Extracellular volume fraction -Vp = 0.3 # Plasma volume fraction -ct = osipi.two_cxm(t, ca, E=E, Fp=Fp, Ve=Ve, Vp=Vp) +Vp = 0.1 # Plasma volume fraction +ct = osipi.two_compartment_exchange_model(t, ca, E=E, Fp=Fp, Ve=Ve, Vp=Vp) plt.plot(t, ct, "b-", label=f"E = {E}, Fp = {Fp}, Ve = {Ve}, Vp = {Vp}") +t = np.arange(0, 6 * 60, 1, dtype=float) +ct = osipi.two_compartment_exchange_model(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.2) +plt.plot(t, ct, "r-", label="E = 0.15, Fp = 0.2, Ve = 0.2, Vp = 2") plt.xlabel("Time (sec)") plt.ylabel("Tissue concentration (mM)") +plt.legend() plt.show() diff --git a/docs/references/models/tissue_models/2cxm.md b/docs/references/models/tissue_models/2cxm.md index 38bc95f..0b3d636 100644 --- a/docs/references/models/tissue_models/2cxm.md +++ b/docs/references/models/tissue_models/2cxm.md @@ -1,6 +1,7 @@ # osipi.Two Compartment Exchange Model -::: osipi.two_cxm +::: osipi.two_compartment_exchange_model + ## Example using `osipi.two_cxm` diff --git a/src/osipi/__init__.py b/src/osipi/__init__.py index 06fdc38..ffe955b 100644 --- a/src/osipi/__init__.py +++ b/src/osipi/__init__.py @@ -8,7 +8,7 @@ from ._tissue import ( tofts, extended_tofts, - two_cxm + two_compartment_exchange_model ) from ._signal import ( diff --git a/src/osipi/_tissue.py b/src/osipi/_tissue.py index 220fac4..8441e3e 100755 --- a/src/osipi/_tissue.py +++ b/src/osipi/_tissue.py @@ -334,7 +334,7 @@ def extended_tofts( return ct -def two_cxm( +def two_compartment_exchange_model( t: np.ndarray, ca: np.ndarray, E: float, @@ -342,9 +342,8 @@ def two_cxm( Ve: float, Vp: float, Ta: float = 30.0, - discretization_method: str = "conv", ) -> np.ndarray: - """The 2CX model allows bi-directional exchange of indicator between vascular and + """The 2 compartment exchange model allows bi-directional exchange of indicator between vascular and extra vascular extracellular compartments. Indicator is assumed to be well mixed within each compartment. @@ -352,8 +351,8 @@ def two_cxm( t (np.ndarray): array of time points in units of sec.[OSIPI code Q.GE1.004] ca (np.ndarray): Arterial concentrations in mM for each time point in t.[OSIPI code Q.IC1.001] - E (float): Extraction fraction in units of mL/min/100g. - Fp (float): Plasma flow in units of mL/min/100g. [OSIPI code Q.PH1.003] + E (float): Extraction fraction in units of mL/min/100mL. + Fp (float): Plasma flow in units of mL/min/100mL. [OSIPI code Q.PH1.003] Ve (float): Relative volume fraction of the extracellular extravascular compartment (e). [OSIPI code Q.PH1.001.[e]] Vp (float): Relative volyme fraction of the plasma @@ -395,15 +394,25 @@ def two_cxm( Calculate tissue concentrations and plot: - >>> E = 0.15 # in units of mL/min/100g - >>> Fp = 0.2 # in units of mL/min/100g - >>> Ve = 0.2 # takes values from 0 to 1 + >>> E = 0.15 # in units of mL/min/100mL + >>> Fp = 5 # in units of mL/min/100mL + >>> Ve = 0.1 # takes values from 0 to 1 >>> Vp = 0.3 # takes values from 0 to 1 - >>> ct = osipi.two_cxm(t, ca, E, Fp, Ve, Vp) + >>> ct = osipi.two_compartment_exchange_model(t, ca, E, Fp, Ve, Vp) >>> plt.plot(t, ca, "r", t, ct, "b") """ + if Vp == 0: + Ktrans = E * Fp + + ct = tofts(t, ca, Ktrans, Ve) + return ct + + # Convert units + t /= 60.0 + Ta /= 60.0 + 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 diff --git a/tests/test_tissue.py b/tests/test_tissue.py index 09f8a98..e7f9537 100755 --- a/tests/test_tissue.py +++ b/tests/test_tissue.py @@ -124,35 +124,41 @@ def test_tissue_extended_tofts(): assert np.allclose(ct_conv, ca * 0.3, rtol=1e-4, atol=1e-3) -def test_tissue_2compartment_model(): +def test_tissue_two_compartment_exchange_model(): # 1. Basic operation of the function - test that the peak tissue + # concentration is less than the peak AIF t = np.arange(0, 6 * 60, 1, dtype=float) ca = osipi.aif_parker(t) - ct = osipi.two_cxm(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.3) + ct = osipi.two_compartment_exchange_model(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.3) assert np.round(np.max(ct)) < np.round(np.max(ca)) # 2. Basic operation of the function - test with non-uniform spacing of + # time array t = np.geomspace(1, 6 * 60 + 1, num=360) - 1 ca = osipi.aif_parker(t) - ct = osipi.two_cxm(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.3) + ct = osipi.two_compartment_exchange_model(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.3) assert np.round(np.max(ct)) < np.round(np.max(ca)) # 3. The offset option - test that the tissue concentration is shifted - t = np.arange(0, 6 * 60, 1) + # from the AIF by the specified offset time + t = np.arange(0, 6 * 60, 1, dtype=float) ca = osipi.aif_parker(t) - ct = osipi.two_cxm(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.3, Ta=60.0) + ct = osipi.two_compartment_exchange_model(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.1, Ta=60.0) assert (np.min(np.where(ct > 0.0)) - np.min(np.where(ca > 0.0)) - 1) * 1 == 60.0 - # 4. Test that the area under the ct curve is approximately the sum of - # the extracellular volume and the plasma volume - t = np.arange(0, 6 * 60, 1) + # 4. Test specific use cases + + # Test when Vp is 0 the output matches tofts model + t = np.arange(0, 6 * 60, 1, dtype=float) ca = osipi.aif_parker(t) - ct = osipi.two_cxm(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0.3) - assert math.isclose(np.trapz(ct, t) / np.trapz(ca, t), 0.2 + 0.3, abs_tol=1e-1) + ct = osipi.two_compartment_exchange_model(t, ca, E=0.15, Fp=0.2, Ve=0.2, Vp=0) + ct_tofts = osipi.tofts(t, ca, Ktrans=0.03, ve=0.2) + assert np.allclose(ct, ct_tofts, rtol=1e-4, atol=1e-3) if __name__ == "__main__": test_tissue_tofts() test_tissue_extended_tofts() + test_tissue_two_compartment_exchange_model() print("All tissue concentration model tests passed!!")