From 25b29d916cafa4cdc1af305a526238c885d27902 Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Wed, 14 Aug 2024 17:29:27 +0200 Subject: [PATCH 01/10] Oliver's approach Dynamic scaling adjusted. Also boundaries adjusted Remove normalisation as signal isnormalized. Remove algorithm-specific tolerance boundaries such that default boundaries are used Take medians of parameters as means give some machine precission error... --- phantoms/MR_XCAT_qMRI/sim_ivim_sig.py | 6 +- tests/IVIMmodels/unit_tests/algorithms.json | 106 +++-------------- tests/IVIMmodels/unit_tests/generic.json | 108 +++++++++--------- tests/IVIMmodels/unit_tests/generic_four.json | 108 +++++++++--------- tests/IVIMmodels/unit_tests/generic_one.json | 108 +++++++++--------- .../IVIMmodels/unit_tests/generic_three.json | 108 +++++++++--------- tests/IVIMmodels/unit_tests/generic_two.json | 108 +++++++++--------- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 39 ++++--- 8 files changed, 312 insertions(+), 379 deletions(-) diff --git a/phantoms/MR_XCAT_qMRI/sim_ivim_sig.py b/phantoms/MR_XCAT_qMRI/sim_ivim_sig.py index df334c35..3db09bc8 100644 --- a/phantoms/MR_XCAT_qMRI/sim_ivim_sig.py +++ b/phantoms/MR_XCAT_qMRI/sim_ivim_sig.py @@ -459,9 +459,9 @@ def parse_bvalues_file(file_path): signals = np.squeeze(voxels[int(voxels.shape[0] * voxel_selector_fraction)]).tolist() generic_data[name] = { 'noise': noise, - 'D': np.mean(Dim[selector], axis=0), - 'f': np.mean(fim[selector], axis=0), - 'Dp': np.mean(Dpim[selector], axis=0), + 'D': np.median(Dim[selector], axis=0), + 'f': np.median(fim[selector], axis=0), + 'Dp': np.median(Dpim[selector], axis=0), 'data': signals } generic_data['config'] = { diff --git a/tests/IVIMmodels/unit_tests/algorithms.json b/tests/IVIMmodels/unit_tests/algorithms.json index a39df9fb..16d35af3 100644 --- a/tests/IVIMmodels/unit_tests/algorithms.json +++ b/tests/IVIMmodels/unit_tests/algorithms.json @@ -16,86 +16,32 @@ ], "IAR_LU_biexp": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "IAR_LU_modified_mix": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "IAR_LU_modified_topopro": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "IAR_LU_segmented_2step": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "IAR_LU_segmented_3step": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "IAR_LU_subtracted": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "ETP_SRI_LinearFitting": { @@ -117,51 +63,31 @@ "Dp": 0 }, "dynamic_rtol": { - "offset": 0, - "noise": 0, - "ratio": 0, - "noiseCrossRatio": 10, - "f": 1, + "offset": 0.001, + "noise": 200, + "scale_f": false, + "f": 2, "D": 1, - "Dp": 2 + "Dp": 3 }, "dynamic_atol": { "offset": 5e-2, - "noise": 0, - "ratio": 0, - "noiseCrossRatio": 10, + "noise": 10, + "scale_f": false, "f": 1, - "D": 1, + "D": 0.001, "Dp": 2 } } }, "OGC_AmsterdamUMC_biexp_segmented": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } }, "OJ_GU_seg": { "tolerances": { - "rtol": { - "f": 5, - "D": 5, - "Dp": 25 - }, - "atol": { - "f": 1e-2, - "D": 1e-2, - "Dp": 1e-1 - } + } } } diff --git a/tests/IVIMmodels/unit_tests/generic.json b/tests/IVIMmodels/unit_tests/generic.json index 223c4173..c1e2477b 100644 --- a/tests/IVIMmodels/unit_tests/generic.json +++ b/tests/IVIMmodels/unit_tests/generic.json @@ -1,8 +1,8 @@ { "Myocardium LV": { "noise": 0.0005, - "D": 0.0023999999999999994, - "f": 0.14999999999999997, + "D": 0.0024, + "f": 0.15, "Dp": 0.08, "data": [ 1.0001231751746, @@ -28,8 +28,8 @@ "myocardium RV": { "noise": 0.0005, "D": 0.0024, - "f": 0.14999999999999997, - "Dp": 0.07999999999999999, + "f": 0.15, + "Dp": 0.08, "data": [ 1.000180079067396, 0.986990380915803, @@ -53,7 +53,7 @@ }, "myocardium ra": { "noise": 0.0005, - "D": 0.0014999999999999998, + "D": 0.0015, "f": 0.07, "Dp": 0.07, "data": [ @@ -79,7 +79,7 @@ }, "Blood RV": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 1.0, "Dp": 0.1, "data": [ @@ -105,9 +105,9 @@ }, "Blood RA": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9996412591901256, 0.9045050114102029, @@ -132,8 +132,8 @@ "muscle": { "noise": 0.0005, "D": 0.00137, - "f": 0.10000000000000002, - "Dp": 0.026299999999999994, + "f": 0.1, + "Dp": 0.0263, "data": [ 1.0000325379134485, 0.9958416765116361, @@ -158,7 +158,7 @@ "Liver": { "noise": 0.0005, "D": 0.0015, - "f": 0.11000000000000007, + "f": 0.11, "Dp": 0.1, "data": [ 0.9998555372369898, @@ -183,9 +183,9 @@ }, "gall bladder": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 0.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9989060357343036, 0.9970876206615882, @@ -209,9 +209,9 @@ }, "esophagus": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 0.9993171743502156, 0.9891133604383114, @@ -235,9 +235,9 @@ }, "esophagus cont": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 1.0001332280689843, 0.9895985875899938, @@ -261,9 +261,9 @@ }, "st wall": { "noise": 0.0005, - "D": 0.0015000000000000002, + "D": 0.0015, "f": 0.3, - "Dp": 0.012000000000000002, + "Dp": 0.012, "data": [ 0.9998513623271699, 0.9948379132892078, @@ -287,7 +287,7 @@ }, "Stomach Contents": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 0.0, "Dp": 0.0, "data": [ @@ -313,9 +313,9 @@ }, "pancreas": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.14999999999999997, - "Dp": 0.009999999999999998, + "D": 0.0013, + "f": 0.15, + "Dp": 0.01, "data": [ 1.0003875681734176, 0.9973129244965385, @@ -340,8 +340,8 @@ "Right kydney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999999, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 0.9999350419508581, 0.9960432237822785, @@ -365,7 +365,7 @@ }, "right kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, "Dp": 0.019, "data": [ @@ -392,8 +392,8 @@ "Left kidney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999998, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 1.0011945854174855, 0.9955396602100857, @@ -417,9 +417,9 @@ }, "left kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, - "Dp": 0.018999999999999996, + "Dp": 0.019, "data": [ 0.999784883468332, 0.9956275637829168, @@ -443,9 +443,9 @@ }, "spleen": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.19999999999999998, - "Dp": 0.030000000000000006, + "D": 0.0013, + "f": 0.2, + "Dp": 0.03, "data": [ 0.9995828860807867, 0.9934063169138365, @@ -469,8 +469,8 @@ }, "spinal cord": { "noise": 0.0005, - "D": 0.0004099999999999999, - "f": 0.17799999999999994, + "D": 0.00041, + "f": 0.178, "Dp": 0.0289, "data": [ 0.9996757942460722, @@ -495,9 +495,9 @@ }, "Bone Marrow": { "noise": 0.0005, - "D": 0.00042999999999999994, + "D": 0.00043, "f": 0.145, - "Dp": 0.05000000000000001, + "Dp": 0.05, "data": [ 1.0007913301979523, 0.9918076927329963, @@ -521,9 +521,9 @@ }, "Artery": { "noise": 0.0005, - "D": 0.002999999999999999, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9997800501422142, 0.9045819629490028, @@ -549,7 +549,7 @@ "noise": 0.0005, "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999999, + "Dp": 0.1, "data": [ 1.0001149665130338, 0.9052158548961062, @@ -574,8 +574,8 @@ "asc lower intestine": { "noise": 0.0005, "D": 0.00131, - "f": 0.6899999999999998, - "Dp": 0.028999999999999995, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9985203088030562, 0.9799920249288326, @@ -599,9 +599,9 @@ }, "trans lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9997303817205065, 0.9794580744935948, @@ -625,9 +625,9 @@ }, "desc lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999998, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9998169952297454, 0.9794581567785913, @@ -651,9 +651,9 @@ }, "small intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.02900000000000001, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.000672478806307, 0.9801710601503958, @@ -678,7 +678,7 @@ "pericardium": { "noise": 0.0005, "D": 0.003, - "f": 0.07000000000000002, + "f": 0.07, "Dp": 0.01, "data": [ 1.0001269336613792, diff --git a/tests/IVIMmodels/unit_tests/generic_four.json b/tests/IVIMmodels/unit_tests/generic_four.json index 1907c958..794c9494 100644 --- a/tests/IVIMmodels/unit_tests/generic_four.json +++ b/tests/IVIMmodels/unit_tests/generic_four.json @@ -1,8 +1,8 @@ { "Myocardium LV": { "noise": 0.0005, - "D": 0.0023999999999999994, - "f": 0.14999999999999997, + "D": 0.0024, + "f": 0.15, "Dp": 0.08, "data": [ 1.0001231751746, @@ -28,8 +28,8 @@ "myocardium RV": { "noise": 0.0005, "D": 0.0024, - "f": 0.14999999999999997, - "Dp": 0.07999999999999999, + "f": 0.15, + "Dp": 0.08, "data": [ 1.000180079067396, 0.986990380915803, @@ -53,7 +53,7 @@ }, "myocardium ra": { "noise": 0.0005, - "D": 0.0014999999999999998, + "D": 0.0015, "f": 0.07, "Dp": 0.07, "data": [ @@ -79,7 +79,7 @@ }, "Blood RV": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 1.0, "Dp": 0.1, "data": [ @@ -105,9 +105,9 @@ }, "Blood RA": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9996412591901256, 0.9045050114102029, @@ -132,8 +132,8 @@ "muscle": { "noise": 0.0005, "D": 0.00137, - "f": 0.10000000000000002, - "Dp": 0.026299999999999994, + "f": 0.1, + "Dp": 0.0263, "data": [ 1.0000325379134485, 0.9958416765116361, @@ -158,7 +158,7 @@ "Liver": { "noise": 0.0005, "D": 0.0015, - "f": 0.11000000000000007, + "f": 0.11, "Dp": 0.1, "data": [ 0.9998555372369898, @@ -183,9 +183,9 @@ }, "gall bladder": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 0.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9989060357343036, 0.9970876206615882, @@ -209,9 +209,9 @@ }, "esophagus": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 0.9993171743502156, 0.9891133604383114, @@ -235,9 +235,9 @@ }, "esophagus cont": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 1.0001332280689843, 0.9895985875899938, @@ -261,9 +261,9 @@ }, "st wall": { "noise": 0.0005, - "D": 0.0015000000000000002, + "D": 0.0015, "f": 0.3, - "Dp": 0.012000000000000002, + "Dp": 0.012, "data": [ 0.9998513623271699, 0.9948379132892078, @@ -287,7 +287,7 @@ }, "Stomach Contents": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 0.0, "Dp": 0.0, "data": [ @@ -313,9 +313,9 @@ }, "pancreas": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.14999999999999997, - "Dp": 0.009999999999999998, + "D": 0.0013, + "f": 0.15, + "Dp": 0.01, "data": [ 1.0003875681734176, 0.9973129244965385, @@ -340,8 +340,8 @@ "Right kydney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999999, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 0.9999350419508581, 0.9960432237822785, @@ -365,7 +365,7 @@ }, "right kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, "Dp": 0.019, "data": [ @@ -392,8 +392,8 @@ "Left kidney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999998, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 1.0011945854174855, 0.9955396602100857, @@ -417,9 +417,9 @@ }, "left kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, - "Dp": 0.018999999999999996, + "Dp": 0.019, "data": [ 0.999784883468332, 0.9956275637829168, @@ -443,9 +443,9 @@ }, "spleen": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.19999999999999998, - "Dp": 0.030000000000000006, + "D": 0.0013, + "f": 0.2, + "Dp": 0.03, "data": [ 0.9995828860807867, 0.9934063169138365, @@ -469,8 +469,8 @@ }, "spinal cord": { "noise": 0.0005, - "D": 0.0004099999999999999, - "f": 0.17799999999999994, + "D": 0.00041, + "f": 0.178, "Dp": 0.0289, "data": [ 0.9996757942460722, @@ -495,9 +495,9 @@ }, "Bone Marrow": { "noise": 0.0005, - "D": 0.00042999999999999994, + "D": 0.00043, "f": 0.145, - "Dp": 0.05000000000000001, + "Dp": 0.05, "data": [ 1.0007913301979523, 0.9918076927329963, @@ -521,9 +521,9 @@ }, "Artery": { "noise": 0.0005, - "D": 0.002999999999999999, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9997800501422142, 0.9045819629490028, @@ -549,7 +549,7 @@ "noise": 0.0005, "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999999, + "Dp": 0.1, "data": [ 1.0001149665130338, 0.9052158548961062, @@ -574,8 +574,8 @@ "asc lower intestine": { "noise": 0.0005, "D": 0.00131, - "f": 0.6899999999999998, - "Dp": 0.028999999999999995, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9985203088030562, 0.9799920249288326, @@ -599,9 +599,9 @@ }, "trans lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9997303817205065, 0.9794580744935948, @@ -625,9 +625,9 @@ }, "desc lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999998, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9998169952297454, 0.9794581567785913, @@ -651,9 +651,9 @@ }, "small intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.02900000000000001, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.000672478806307, 0.9801710601503958, @@ -678,7 +678,7 @@ "pericardium": { "noise": 0.0005, "D": 0.003, - "f": 0.07000000000000002, + "f": 0.07, "Dp": 0.01, "data": [ 1.0001269336613792, diff --git a/tests/IVIMmodels/unit_tests/generic_one.json b/tests/IVIMmodels/unit_tests/generic_one.json index b6b2e3e4..3f19a5bf 100644 --- a/tests/IVIMmodels/unit_tests/generic_one.json +++ b/tests/IVIMmodels/unit_tests/generic_one.json @@ -1,8 +1,8 @@ { "Myocardium LV": { "noise": 0.0005, - "D": 0.0023999999999999994, - "f": 0.14999999999999997, + "D": 0.0024, + "f": 0.15, "Dp": 0.08, "data": [ 1.0004147615612788, @@ -30,8 +30,8 @@ "myocardium RV": { "noise": 0.0005, "D": 0.0024, - "f": 0.14999999999999997, - "Dp": 0.07999999999999999, + "f": 0.15, + "Dp": 0.08, "data": [ 1.000181304748344, 0.9861747700418392, @@ -57,7 +57,7 @@ }, "myocardium ra": { "noise": 0.0005, - "D": 0.0014999999999999998, + "D": 0.0015, "f": 0.07, "Dp": 0.07, "data": [ @@ -85,7 +85,7 @@ }, "Blood RV": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 1.0, "Dp": 0.1, "data": [ @@ -113,9 +113,9 @@ }, "Blood RA": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 1.0004892903199283, 0.9046983073964768, @@ -142,8 +142,8 @@ "muscle": { "noise": 0.0005, "D": 0.00137, - "f": 0.10000000000000002, - "Dp": 0.026299999999999994, + "f": 0.1, + "Dp": 0.0263, "data": [ 0.9997142768189549, 0.996013043803113, @@ -170,7 +170,7 @@ "Liver": { "noise": 0.0005, "D": 0.0015, - "f": 0.11000000000000007, + "f": 0.11, "Dp": 0.1, "data": [ 1.0000569769441432, @@ -197,9 +197,9 @@ }, "gall bladder": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 0.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9995238156474249, 0.997362029165271, @@ -225,9 +225,9 @@ }, "esophagus": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 0.9997002304518158, 0.989154675581577, @@ -253,9 +253,9 @@ }, "esophagus cont": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 1.00006491694293, 0.9904119975513289, @@ -281,9 +281,9 @@ }, "st wall": { "noise": 0.0005, - "D": 0.0015000000000000002, + "D": 0.0015, "f": 0.3, - "Dp": 0.012000000000000002, + "Dp": 0.012, "data": [ 0.9996146504640959, 0.9955650448011039, @@ -309,7 +309,7 @@ }, "Stomach Contents": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 0.0, "Dp": 0.0, "data": [ @@ -337,9 +337,9 @@ }, "pancreas": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.14999999999999997, - "Dp": 0.009999999999999998, + "D": 0.0013, + "f": 0.15, + "Dp": 0.01, "data": [ 0.9999289994471732, 0.996955045473085, @@ -366,8 +366,8 @@ "Right kydney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999999, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 0.9994471646163887, 0.9964587942652097, @@ -393,7 +393,7 @@ }, "right kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, "Dp": 0.019, "data": [ @@ -422,8 +422,8 @@ "Left kidney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999998, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 1.0005861007545922, 0.9958085543879336, @@ -449,9 +449,9 @@ }, "left kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, - "Dp": 0.018999999999999996, + "Dp": 0.019, "data": [ 1.0006968468762572, 0.9952601890647235, @@ -477,9 +477,9 @@ }, "spleen": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.19999999999999998, - "Dp": 0.030000000000000006, + "D": 0.0013, + "f": 0.2, + "Dp": 0.03, "data": [ 1.0004068919271272, 0.9945410903182299, @@ -505,8 +505,8 @@ }, "spinal cord": { "noise": 0.0005, - "D": 0.0004099999999999999, - "f": 0.17799999999999994, + "D": 0.00041, + "f": 0.178, "Dp": 0.0289, "data": [ 1.0003503482616034, @@ -533,9 +533,9 @@ }, "Bone Marrow": { "noise": 0.0005, - "D": 0.00042999999999999994, + "D": 0.00043, "f": 0.145, - "Dp": 0.05000000000000001, + "Dp": 0.05, "data": [ 0.999970318173524, 0.9926415227826851, @@ -561,9 +561,9 @@ }, "Artery": { "noise": 0.0005, - "D": 0.002999999999999999, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 1.0013666527439296, 0.9048201228244678, @@ -591,7 +591,7 @@ "noise": 0.0005, "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999999, + "Dp": 0.1, "data": [ 1.0002465126990085, 0.9053141846538422, @@ -618,8 +618,8 @@ "asc lower intestine": { "noise": 0.0005, "D": 0.00131, - "f": 0.6899999999999998, - "Dp": 0.028999999999999995, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0003100507264253, 0.9808596670792784, @@ -645,9 +645,9 @@ }, "trans lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0010109101657108, 0.9796622963563709, @@ -673,9 +673,9 @@ }, "desc lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999998, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0003251408422316, 0.9793306142345213, @@ -701,9 +701,9 @@ }, "small intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.02900000000000001, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0001139232803467, 0.9795892407472606, @@ -730,7 +730,7 @@ "pericardium": { "noise": 0.0005, "D": 0.003, - "f": 0.07000000000000002, + "f": 0.07, "Dp": 0.01, "data": [ 1.0003177369531795, diff --git a/tests/IVIMmodels/unit_tests/generic_three.json b/tests/IVIMmodels/unit_tests/generic_three.json index 259025e8..6344aced 100644 --- a/tests/IVIMmodels/unit_tests/generic_three.json +++ b/tests/IVIMmodels/unit_tests/generic_three.json @@ -1,8 +1,8 @@ { "Myocardium LV": { "noise": 0.0005, - "D": 0.0023999999999999994, - "f": 0.14999999999999997, + "D": 0.0024, + "f": 0.15, "Dp": 0.08, "data": [ 0.9991809304406896, @@ -27,8 +27,8 @@ "myocardium RV": { "noise": 0.0005, "D": 0.0024, - "f": 0.14999999999999997, - "Dp": 0.07999999999999999, + "f": 0.15, + "Dp": 0.08, "data": [ 1.000186442403044, 0.9866124740349732, @@ -51,7 +51,7 @@ }, "myocardium ra": { "noise": 0.0005, - "D": 0.0014999999999999998, + "D": 0.0015, "f": 0.07, "Dp": 0.07, "data": [ @@ -76,7 +76,7 @@ }, "Blood RV": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 1.0, "Dp": 0.1, "data": [ @@ -101,9 +101,9 @@ }, "Blood RA": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 1.0005368625085327, 0.9049245909086914, @@ -127,8 +127,8 @@ "muscle": { "noise": 0.0005, "D": 0.00137, - "f": 0.10000000000000002, - "Dp": 0.026299999999999994, + "f": 0.1, + "Dp": 0.0263, "data": [ 0.9990590772436959, 0.9960448630641721, @@ -152,7 +152,7 @@ "Liver": { "noise": 0.0005, "D": 0.0015, - "f": 0.11000000000000007, + "f": 0.11, "Dp": 0.1, "data": [ 1.0009183383140396, @@ -176,9 +176,9 @@ }, "gall bladder": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 0.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.999011898784023, 0.9967377786835717, @@ -201,9 +201,9 @@ }, "esophagus": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 1.000171208144279, 0.9886326948608724, @@ -226,9 +226,9 @@ }, "esophagus cont": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 0.9998205519492651, 0.9889688370839962, @@ -251,9 +251,9 @@ }, "st wall": { "noise": 0.0005, - "D": 0.0015000000000000002, + "D": 0.0015, "f": 0.3, - "Dp": 0.012000000000000002, + "Dp": 0.012, "data": [ 1.0001009276028885, 0.9959064821861017, @@ -276,7 +276,7 @@ }, "Stomach Contents": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 0.0, "Dp": 0.0, "data": [ @@ -301,9 +301,9 @@ }, "pancreas": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.14999999999999997, - "Dp": 0.009999999999999998, + "D": 0.0013, + "f": 0.15, + "Dp": 0.01, "data": [ 1.000608445907503, 0.9979376297547163, @@ -327,8 +327,8 @@ "Right kydney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999999, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 0.9999709391013328, 0.9969334386920174, @@ -351,7 +351,7 @@ }, "right kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, "Dp": 0.019, "data": [ @@ -377,8 +377,8 @@ "Left kidney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999998, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 0.9999764343046507, 0.9966946653480023, @@ -401,9 +401,9 @@ }, "left kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, - "Dp": 0.018999999999999996, + "Dp": 0.019, "data": [ 1.0001915373938952, 0.995734805186699, @@ -426,9 +426,9 @@ }, "spleen": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.19999999999999998, - "Dp": 0.030000000000000006, + "D": 0.0013, + "f": 0.2, + "Dp": 0.03, "data": [ 0.9996597541995026, 0.9921927935155999, @@ -451,8 +451,8 @@ }, "spinal cord": { "noise": 0.0005, - "D": 0.0004099999999999999, - "f": 0.17799999999999994, + "D": 0.00041, + "f": 0.178, "Dp": 0.0289, "data": [ 0.999653292219937, @@ -476,9 +476,9 @@ }, "Bone Marrow": { "noise": 0.0005, - "D": 0.00042999999999999994, + "D": 0.00043, "f": 0.145, - "Dp": 0.05000000000000001, + "Dp": 0.05, "data": [ 1.0000172385794917, 0.9929136520496381, @@ -501,9 +501,9 @@ }, "Artery": { "noise": 0.0005, - "D": 0.002999999999999999, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 1.0007900407561108, 0.9053446322852805, @@ -528,7 +528,7 @@ "noise": 0.0005, "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999999, + "Dp": 0.1, "data": [ 1.000522845586714, 0.90536598797916, @@ -552,8 +552,8 @@ "asc lower intestine": { "noise": 0.0005, "D": 0.00131, - "f": 0.6899999999999998, - "Dp": 0.028999999999999995, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9995171229676644, 0.9807828623959373, @@ -576,9 +576,9 @@ }, "trans lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 0.9993466584617889, 0.9795617139251637, @@ -601,9 +601,9 @@ }, "desc lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999998, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0003127229429096, 0.9795841525059689, @@ -626,9 +626,9 @@ }, "small intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.02900000000000001, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 0.999560283259107, 0.9799310281771558, @@ -652,7 +652,7 @@ "pericardium": { "noise": 0.0005, "D": 0.003, - "f": 0.07000000000000002, + "f": 0.07, "Dp": 0.01, "data": [ 0.9999552010727751, diff --git a/tests/IVIMmodels/unit_tests/generic_two.json b/tests/IVIMmodels/unit_tests/generic_two.json index bf3a4338..20c88d58 100644 --- a/tests/IVIMmodels/unit_tests/generic_two.json +++ b/tests/IVIMmodels/unit_tests/generic_two.json @@ -1,8 +1,8 @@ { "Myocardium LV": { "noise": 0.0005, - "D": 0.0023999999999999994, - "f": 0.14999999999999997, + "D": 0.0024, + "f": 0.15, "Dp": 0.08, "data": [ 1.0004147615612788, @@ -30,8 +30,8 @@ "myocardium RV": { "noise": 0.0005, "D": 0.0024, - "f": 0.14999999999999997, - "Dp": 0.07999999999999999, + "f": 0.15, + "Dp": 0.08, "data": [ 1.000181304748344, 0.9861747700418392, @@ -57,7 +57,7 @@ }, "myocardium ra": { "noise": 0.0005, - "D": 0.0014999999999999998, + "D": 0.0015, "f": 0.07, "Dp": 0.07, "data": [ @@ -85,7 +85,7 @@ }, "Blood RV": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 1.0, "Dp": 0.1, "data": [ @@ -113,9 +113,9 @@ }, "Blood RA": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 1.0004892903199283, 0.9046983073964768, @@ -142,8 +142,8 @@ "muscle": { "noise": 0.0005, "D": 0.00137, - "f": 0.10000000000000002, - "Dp": 0.026299999999999994, + "f": 0.1, + "Dp": 0.0263, "data": [ 0.9997142768189549, 0.996013043803113, @@ -170,7 +170,7 @@ "Liver": { "noise": 0.0005, "D": 0.0015, - "f": 0.11000000000000007, + "f": 0.11, "Dp": 0.1, "data": [ 1.0000569769441432, @@ -197,9 +197,9 @@ }, "gall bladder": { "noise": 0.0005, - "D": 0.0029999999999999996, + "D": 0.003, "f": 0.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 0.9995238156474249, 0.997362029165271, @@ -225,9 +225,9 @@ }, "esophagus": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 0.9997002304518158, 0.989154675581577, @@ -253,9 +253,9 @@ }, "esophagus cont": { "noise": 0.0005, - "D": 0.0016700000000000003, - "f": 0.31999999999999995, - "Dp": 0.030000000000000002, + "D": 0.00167, + "f": 0.32, + "Dp": 0.03, "data": [ 1.00006491694293, 0.9904119975513289, @@ -281,9 +281,9 @@ }, "st wall": { "noise": 0.0005, - "D": 0.0015000000000000002, + "D": 0.0015, "f": 0.3, - "Dp": 0.012000000000000002, + "Dp": 0.012, "data": [ 0.9996146504640959, 0.9955650448011039, @@ -309,7 +309,7 @@ }, "Stomach Contents": { "noise": 0.0005, - "D": 0.0030000000000000005, + "D": 0.003, "f": 0.0, "Dp": 0.0, "data": [ @@ -337,9 +337,9 @@ }, "pancreas": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.14999999999999997, - "Dp": 0.009999999999999998, + "D": 0.0013, + "f": 0.15, + "Dp": 0.01, "data": [ 0.9999289994471732, 0.996955045473085, @@ -366,8 +366,8 @@ "Right kydney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999999, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 0.9994471646163887, 0.9964587942652097, @@ -393,7 +393,7 @@ }, "right kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, "Dp": 0.019, "data": [ @@ -422,8 +422,8 @@ "Left kidney cortex": { "noise": 0.0005, "D": 0.00212, - "f": 0.09699999999999998, - "Dp": 0.019999999999999997, + "f": 0.097, + "Dp": 0.02, "data": [ 1.0005861007545922, 0.9958085543879336, @@ -449,9 +449,9 @@ }, "left kidney medulla": { "noise": 0.0005, - "D": 0.0020900000000000003, + "D": 0.00209, "f": 0.158, - "Dp": 0.018999999999999996, + "Dp": 0.019, "data": [ 1.0006968468762572, 0.9952601890647235, @@ -477,9 +477,9 @@ }, "spleen": { "noise": 0.0005, - "D": 0.0012999999999999997, - "f": 0.19999999999999998, - "Dp": 0.030000000000000006, + "D": 0.0013, + "f": 0.2, + "Dp": 0.03, "data": [ 1.0004068919271272, 0.9945410903182299, @@ -505,8 +505,8 @@ }, "spinal cord": { "noise": 0.0005, - "D": 0.0004099999999999999, - "f": 0.17799999999999994, + "D": 0.00041, + "f": 0.178, "Dp": 0.0289, "data": [ 1.0003503482616034, @@ -533,9 +533,9 @@ }, "Bone Marrow": { "noise": 0.0005, - "D": 0.00042999999999999994, + "D": 0.00043, "f": 0.145, - "Dp": 0.05000000000000001, + "Dp": 0.05, "data": [ 0.999970318173524, 0.9926415227826851, @@ -561,9 +561,9 @@ }, "Artery": { "noise": 0.0005, - "D": 0.002999999999999999, + "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999998, + "Dp": 0.1, "data": [ 1.0013666527439296, 0.9048201228244678, @@ -591,7 +591,7 @@ "noise": 0.0005, "D": 0.003, "f": 1.0, - "Dp": 0.09999999999999999, + "Dp": 0.1, "data": [ 1.0002465126990085, 0.9053141846538422, @@ -618,8 +618,8 @@ "asc lower intestine": { "noise": 0.0005, "D": 0.00131, - "f": 0.6899999999999998, - "Dp": 0.028999999999999995, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0003100507264253, 0.9808596670792784, @@ -645,9 +645,9 @@ }, "trans lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0010109101657108, 0.9796622963563709, @@ -673,9 +673,9 @@ }, "desc lower intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999998, - "Dp": 0.029000000000000012, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0003251408422316, 0.9793306142345213, @@ -701,9 +701,9 @@ }, "small intestine": { "noise": 0.0005, - "D": 0.0013099999999999997, - "f": 0.6899999999999997, - "Dp": 0.02900000000000001, + "D": 0.00131, + "f": 0.69, + "Dp": 0.029, "data": [ 1.0001139232803467, 0.9795892407472606, @@ -730,7 +730,7 @@ "pericardium": { "noise": 0.0005, "D": 0.003, - "f": 0.07000000000000002, + "f": 0.07, "Dp": 0.01, "data": [ 1.0003177369531795, diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index f4ecc2b8..5dd08e7b 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -114,24 +114,31 @@ def signal_helper(signal): signal = np.asarray(signal) - signal = np.abs(signal) - signal /= signal[0] - ratio = 1 / signal[0] - return signal, ratio + signal = np.abs(signal) #Oliver: do we need this for unit testing? It will bias the fits by forcing positive signals. + signal /= signal[0] # this scales noise differenly. I have now made sure signal is produced with amplitude 1 if right flags are selected. Then SNR is fixed. Then, this sentence may be redundant + #ratio = 1 / signal[0] #this is 1 per definition (See sentence above) + return signal -def tolerances_helper(tolerances, ratio, noise): +def tolerances_helper(tolerances, data): + epsilon = 0.0001 if "dynamic_rtol" in tolerances: dyn_rtol = tolerances["dynamic_rtol"] - scale = dyn_rtol["offset"] + dyn_rtol["ratio"]*ratio + dyn_rtol["noise"]*noise + dyn_rtol["noiseCrossRatio"]*ratio*noise - tolerances["rtol"] = {"f": scale*dyn_rtol["f"], "D": scale*dyn_rtol["D"], "Dp": scale*dyn_rtol["Dp"]} + scale = dyn_rtol["offset"] + dyn_rtol["noise"]*data["noise"] + if dyn_rtol["scale_f"]: + tolerances["rtol"] = {"f": scale*dyn_rtol["f"], "D": scale*dyn_rtol["D"]/(1-data['f']+epsilon), "Dp": scale*dyn_rtol["Dp"]/(data['f']+epsilon)} + else: + tolerances["rtol"] = {"f": scale * dyn_rtol["f"], "D": scale * dyn_rtol["D"]/0.9,"Dp": scale * dyn_rtol["Dp"]/0.1} else: - tolerances["rtol"] = tolerances.get("rtol", {"f": 5, "D": 5, "Dp": 25}) + tolerances["rtol"] = tolerances.get("rtol", {"f": 0.1, "D": 0.1, "Dp": 0.3}) if "dynamic_atol" in tolerances: dyn_atol = tolerances["dynamic_atol"] - scale = dyn_atol["offset"] + dyn_atol["ratio"]*ratio + dyn_atol["noise"]*noise + dyn_atol["noiseCrossRatio"]*ratio*noise - tolerances["atol"] = {"f": scale*dyn_atol["f"], "D": scale*dyn_atol["D"], "Dp": scale*dyn_atol["Dp"]} + scale = dyn_atol["offset"] + dyn_atol["noise"]*data["noise"] + if dyn_atol["scale_f"]: + tolerances["atol"] = {"f": scale*dyn_atol["f"], "D": scale*dyn_atol["D"]/(1-data['f']+epsilon), "Dp": scale*dyn_atol["Dp"]/(data['f']+epsilon)} + else: + tolerances["atol"] = {"f": scale*dyn_atol["f"], "D": scale*dyn_atol["D"]/0.9, "Dp": scale*dyn_atol["Dp"]/0.1} else: - tolerances["atol"] = tolerances.get("atol", {"f": 1e-2, "D": 1e-2, "Dp": 1e-1}) + tolerances["atol"] = tolerances.get("atol", {"f": 2e-1, "D": 5e-4, "Dp": 4e-2}) return tolerances def data_ivim_fit_saved(): @@ -165,12 +172,12 @@ def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances, mark = pytest.mark.xfail(reason="xfail", strict=xfail["strict"]) request.node.add_marker(mark) fit = OsipiBase(algorithm=algorithm, **kwargs) - signal, ratio = signal_helper(data["data"]) - tolerances = tolerances_helper(tolerances, ratio, data["noise"]) + signal = signal_helper(data["data"]) + tolerances = tolerances_helper(tolerances, data) [f_fit, Dp_fit, D_fit] = fit.osipi_fit(signal, bvals) - npt.assert_allclose(data['f'], f_fit, rtol=tolerances["rtol"]["f"], atol=tolerances["atol"]["f"]) + npt.assert_allclose(f_fit,data['f'], rtol=tolerances["rtol"]["f"], atol=tolerances["atol"]["f"]) if data['f']<0.80: # we need some signal for D to be detected - npt.assert_allclose(data['D'], D_fit, rtol=tolerances["rtol"]["D"], atol=tolerances["atol"]["D"]) + npt.assert_allclose(D_fit,data['D'], rtol=tolerances["rtol"]["D"], atol=tolerances["atol"]["D"]) if data['f']>0.03: #we need some f for D* to be interpretable - npt.assert_allclose(data['Dp'], Dp_fit, rtol=tolerances["rtol"]["Dp"], atol=tolerances["atol"]["Dp"]) + npt.assert_allclose(Dp_fit,data['Dp'], rtol=tolerances["rtol"]["Dp"], atol=tolerances["atol"]["Dp"]) From eea774c4cc45d95475993e462d0615bb3e84f597 Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Tue, 20 Aug 2024 16:02:48 +0200 Subject: [PATCH 02/10] removed expected fail as it passes --- tests/IVIMmodels/unit_tests/algorithms.json | 3 --- 1 file changed, 3 deletions(-) diff --git a/tests/IVIMmodels/unit_tests/algorithms.json b/tests/IVIMmodels/unit_tests/algorithms.json index 16d35af3..bbfe548b 100644 --- a/tests/IVIMmodels/unit_tests/algorithms.json +++ b/tests/IVIMmodels/unit_tests/algorithms.json @@ -45,9 +45,6 @@ } }, "ETP_SRI_LinearFitting": { - "xfail_names": { - "Vein": true - }, "options": { "thresholds": [500] }, From 8d1a4188c70043b6d16572a8e900bbead9e09f1d Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Fri, 3 Jan 2025 17:35:45 +0100 Subject: [PATCH 03/10] changes made as requested --- src/standardized/ETP_SRI_LinearFitting.py | 1 + src/standardized/OJ_GU_seg.py | 2 +- src/standardized/PvH_KB_NKI_IVIMfit.py | 1 + tests/IVIMmodels/unit_tests/algorithms.json | 2 -- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 15 ++++----------- 5 files changed, 7 insertions(+), 14 deletions(-) diff --git a/src/standardized/ETP_SRI_LinearFitting.py b/src/standardized/ETP_SRI_LinearFitting.py index da0eb5eb..7764aa5e 100644 --- a/src/standardized/ETP_SRI_LinearFitting.py +++ b/src/standardized/ETP_SRI_LinearFitting.py @@ -57,6 +57,7 @@ def ivim_fit(self, signals, bvalues=None, linear_fit_option=False, **kwargs): Returns: _type_: _description_ """ + signals[signals<0.0000001]=0.0000001 if bvalues is None: bvalues = self.bvalues diff --git a/src/standardized/OJ_GU_seg.py b/src/standardized/OJ_GU_seg.py index 03a3cd08..7ec4b1d0 100644 --- a/src/standardized/OJ_GU_seg.py +++ b/src/standardized/OJ_GU_seg.py @@ -61,7 +61,7 @@ def ivim_fit(self, signals, bvalues=None): bthr = 200 else: bthr = self.thresholds[0] - + signals[signals<0.00001]=0.00001 fit_results = seg(signals, bvalues, bthr) f = fit_results['f'] diff --git a/src/standardized/PvH_KB_NKI_IVIMfit.py b/src/standardized/PvH_KB_NKI_IVIMfit.py index 3d21ad5c..dd8246b3 100644 --- a/src/standardized/PvH_KB_NKI_IVIMfit.py +++ b/src/standardized/PvH_KB_NKI_IVIMfit.py @@ -52,6 +52,7 @@ def ivim_fit(self, signals, bvalues=None): #bvalues = np.array(bvalues) bvalues = bvalues.tolist() #NKI code expects a list instead of nparray # reshape signal as the NKI code expects a 4D array + signals[signals<0.00001]=0.00001 signals = np.reshape(signals, (1, 1, 1, len(signals))) # assuming that in this test the signals are always single voxel fit_results = self.NKI_algorithm(signals,bvalues) diff --git a/tests/IVIMmodels/unit_tests/algorithms.json b/tests/IVIMmodels/unit_tests/algorithms.json index bbfe548b..c2e9869e 100644 --- a/tests/IVIMmodels/unit_tests/algorithms.json +++ b/tests/IVIMmodels/unit_tests/algorithms.json @@ -62,7 +62,6 @@ "dynamic_rtol": { "offset": 0.001, "noise": 200, - "scale_f": false, "f": 2, "D": 1, "Dp": 3 @@ -70,7 +69,6 @@ "dynamic_atol": { "offset": 5e-2, "noise": 10, - "scale_f": false, "f": 1, "D": 0.001, "Dp": 2 diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index 5dd08e7b..c4f4fbc9 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -114,29 +114,22 @@ def signal_helper(signal): signal = np.asarray(signal) - signal = np.abs(signal) #Oliver: do we need this for unit testing? It will bias the fits by forcing positive signals. + #signal[signal < 0] = 0.00001 signal /= signal[0] # this scales noise differenly. I have now made sure signal is produced with amplitude 1 if right flags are selected. Then SNR is fixed. Then, this sentence may be redundant - #ratio = 1 / signal[0] #this is 1 per definition (See sentence above) return signal def tolerances_helper(tolerances, data): - epsilon = 0.0001 + epsilon = 0.001 if "dynamic_rtol" in tolerances: dyn_rtol = tolerances["dynamic_rtol"] scale = dyn_rtol["offset"] + dyn_rtol["noise"]*data["noise"] - if dyn_rtol["scale_f"]: - tolerances["rtol"] = {"f": scale*dyn_rtol["f"], "D": scale*dyn_rtol["D"]/(1-data['f']+epsilon), "Dp": scale*dyn_rtol["Dp"]/(data['f']+epsilon)} - else: - tolerances["rtol"] = {"f": scale * dyn_rtol["f"], "D": scale * dyn_rtol["D"]/0.9,"Dp": scale * dyn_rtol["Dp"]/0.1} + tolerances["rtol"] = {"f": scale*dyn_rtol["f"], "D": scale*dyn_rtol["D"]/(1-data['f']+epsilon), "Dp": scale*dyn_rtol["Dp"]/(data['f']+epsilon)} else: tolerances["rtol"] = tolerances.get("rtol", {"f": 0.1, "D": 0.1, "Dp": 0.3}) if "dynamic_atol" in tolerances: dyn_atol = tolerances["dynamic_atol"] scale = dyn_atol["offset"] + dyn_atol["noise"]*data["noise"] - if dyn_atol["scale_f"]: - tolerances["atol"] = {"f": scale*dyn_atol["f"], "D": scale*dyn_atol["D"]/(1-data['f']+epsilon), "Dp": scale*dyn_atol["Dp"]/(data['f']+epsilon)} - else: - tolerances["atol"] = {"f": scale*dyn_atol["f"], "D": scale*dyn_atol["D"]/0.9, "Dp": scale*dyn_atol["Dp"]/0.1} + tolerances["atol"] = {"f": scale*dyn_atol["f"], "D": scale*dyn_atol["D"]/(1-data['f']+epsilon), "Dp": scale*dyn_atol["Dp"]/(data['f']+epsilon)} else: tolerances["atol"] = tolerances.get("atol", {"f": 2e-1, "D": 5e-4, "Dp": 4e-2}) return tolerances From 1bcb389781e3257f561ceeecd3fafcf852a369be Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Sat, 4 Jan 2025 15:09:04 +0100 Subject: [PATCH 04/10] Boundary testing and timing Added testing for whether bounds are respected and whether timing of algorithms make sense --- src/standardized/OGC_AmsterdamUMC_biexp.py | 4 +- tests/IVIMmodels/unit_tests/algorithms.json | 47 +++++++++----------- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 47 +++++++++++++++++++- 3 files changed, 70 insertions(+), 28 deletions(-) diff --git a/src/standardized/OGC_AmsterdamUMC_biexp.py b/src/standardized/OGC_AmsterdamUMC_biexp.py index 130ffc42..3aa244ce 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp.py @@ -27,7 +27,7 @@ class OGC_AmsterdamUMC_biexp(OsipiBase): required_initial_guess_optional = True accepted_dimensions = 1 # Not sure how to define this for the number of accepted dimensions. Perhaps like the thresholds, at least and at most? - def __init__(self, bvalues=None, thresholds=None, bounds=([0, 0, 0.005, 0.7],[0.005, 0.7, 0.2, 1.3]), initial_guess=None, fitS0=False): + def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=None, fitS0=False): """ Everything this algorithm requires should be implemented here. Number of segmentation thresholds, bounds, etc. @@ -45,7 +45,7 @@ def initialize(self, bounds, initial_guess, fitS0): else: self.bounds=bounds if initial_guess is None: - self.initial_guess = [0.001, 0.001, 0.01, 1] + self.initial_guess = [0.001, 0.1, 0.01, 1] else: self.initial_guess = initial_guess self.fitS0=fitS0 diff --git a/tests/IVIMmodels/unit_tests/algorithms.json b/tests/IVIMmodels/unit_tests/algorithms.json index c2e9869e..ff06f143 100644 --- a/tests/IVIMmodels/unit_tests/algorithms.json +++ b/tests/IVIMmodels/unit_tests/algorithms.json @@ -15,34 +15,22 @@ "OJ_GU_seg" ], "IAR_LU_biexp": { - "tolerances": { - - } + "skip_bounds": true }, "IAR_LU_modified_mix": { - "tolerances": { - - } + "skip_bounds": true }, "IAR_LU_modified_topopro": { - "tolerances": { - - } + "skip_bounds": true }, "IAR_LU_segmented_2step": { - "tolerances": { - - } + "skip_bounds": true }, "IAR_LU_segmented_3step": { - "tolerances": { - - } + "skip_bounds": true }, "IAR_LU_subtracted": { - "tolerances": { - - } + "skip_bounds": true }, "ETP_SRI_LinearFitting": { "options": { @@ -73,16 +61,25 @@ "D": 0.001, "Dp": 2 } - } + }, + "skip_bounds": true }, "OGC_AmsterdamUMC_biexp_segmented": { - "tolerances": { - - } + "skip_bounds": true }, "OJ_GU_seg": { - "tolerances": { - - } + "test_bounds": false + }, + "OGC_AmsterdamUMC_Bayesian_biexp": { + "test_bounds": false + }, + "MUMC_biexp": { + "test_bounds": false + }, + "PvH_KB_NKI_IVIMfit": { + "test_bounds": false + }, + "OGC_AmsterdamUMC_biexp": { + "test_bounds": true } } diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index c4f4fbc9..9329daab 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -4,6 +4,7 @@ import json import pathlib import os +import time from src.wrappers.OsipiBase import OsipiBase from utilities.data_simulation.GenerateData import GenerateData @@ -158,19 +159,63 @@ def data_ivim_fit_saved(): tolerances = algorithm_dict.get("tolerances", {}) yield name, bvals, data, algorithm, xfail, kwargs, tolerances - @pytest.mark.parametrize("name, bvals, data, algorithm, xfail, kwargs, tolerances", data_ivim_fit_saved()) def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances, request): if xfail["xfail"]: mark = pytest.mark.xfail(reason="xfail", strict=xfail["strict"]) request.node.add_marker(mark) + start_time = time.time() # Record the start time fit = OsipiBase(algorithm=algorithm, **kwargs) signal = signal_helper(data["data"]) tolerances = tolerances_helper(tolerances, data) [f_fit, Dp_fit, D_fit] = fit.osipi_fit(signal, bvals) + elapsed_time = time.time() - start_time # Calculate elapsed time npt.assert_allclose(f_fit,data['f'], rtol=tolerances["rtol"]["f"], atol=tolerances["atol"]["f"]) if data['f']<0.80: # we need some signal for D to be detected npt.assert_allclose(D_fit,data['D'], rtol=tolerances["rtol"]["D"], atol=tolerances["atol"]["D"]) if data['f']>0.03: #we need some f for D* to be interpretable npt.assert_allclose(Dp_fit,data['Dp'], rtol=tolerances["rtol"]["Dp"], atol=tolerances["atol"]["Dp"]) + assert elapsed_time < 2, f"Algorithm {name} took {elapsed_time} seconds, which is longer than 2 second to fit per voxel" #less than 0.5 seconds per voxel + +def bound_input(): + # Find the algorithms from algorithms.json + file = pathlib.Path(__file__) + algorithm_path = file.with_name('algorithms.json') + with algorithm_path.open() as f: + algorithm_information = json.load(f) + + # Load generic test data generated from the included phantom: phantoms/MR_XCAT_qMRI + generic = file.with_name('generic.json') + with generic.open() as f: + all_data = json.load(f) + + algorithms = algorithm_information["algorithms"] + bvals = all_data.pop('config') + bvals = bvals['bvalues'] + for name, data in all_data.items(): + for algorithm in algorithms: + algorithm_dict = algorithm_information.get(algorithm, {}) + xfail = {"xfail": name in algorithm_dict.get("xfail_names", {}), + "strict": algorithm_dict.get("xfail_names", {}).get(name, True)} + kwargs = algorithm_dict.get("options", {}) + tolerances = algorithm_dict.get("tolerances", {}) + test_bounds = algorithm_dict.get("test_bounds", {}) + if test_bounds: + yield name, bvals, data, algorithm, xfail, kwargs, tolerances + + +@pytest.mark.parametrize("name, bvals, data, algorithm, xfail, kwargs, tolerances", bound_input()) +def test_bounds(name, bvals, data, algorithm, xfail, kwargs, tolerances, request): + bounds = ([0.0008, 0.2, 0.01, 1.1], [0.0012, 0.3, 0.02, 1.3]) + if xfail["xfail"]: + mark = pytest.mark.xfail(reason="xfail", strict=xfail["strict"]) + request.node.add_marker(mark) + # deliberately have silly bounds to see whether they are used + fit = OsipiBase(algorithm=algorithm, bounds=bounds, initial_guess = [0.001, 0.25, 0.015, 1.2], **kwargs) + signal = signal_helper(data["data"]) + tolerances = tolerances_helper(tolerances, data) + [f_fit, Dp_fit, D_fit] = fit.osipi_fit(signal, bvals) + assert bounds[0][0] <= D_fit <= bounds[1][0], f"Result {D_fit} out of bounds for data: {name}" + assert bounds[0][1] <= f_fit <= bounds[1][1], f"Result {f_fit} out of bounds for data: {name}" + assert bounds[0][2] <= Dp_fit <= bounds[1][2], f"Result {Dp_fit} out of bounds for data: {name}" From 1d26badb6ab495a533537e4c0f964d3c0b05914b Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Sat, 4 Jan 2025 15:40:12 +0100 Subject: [PATCH 05/10] Test whether initial guesses and fit bounds are reasonable --- tests/IVIMmodels/unit_tests/algorithms.json | 16 +++++----- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 33 ++++++++++++++++++++ 2 files changed, 41 insertions(+), 8 deletions(-) diff --git a/tests/IVIMmodels/unit_tests/algorithms.json b/tests/IVIMmodels/unit_tests/algorithms.json index ff06f143..3046d2ff 100644 --- a/tests/IVIMmodels/unit_tests/algorithms.json +++ b/tests/IVIMmodels/unit_tests/algorithms.json @@ -15,22 +15,22 @@ "OJ_GU_seg" ], "IAR_LU_biexp": { - "skip_bounds": true + "test_bounds": false }, "IAR_LU_modified_mix": { - "skip_bounds": true + "test_bounds": false }, "IAR_LU_modified_topopro": { - "skip_bounds": true + "test_bounds": false }, "IAR_LU_segmented_2step": { - "skip_bounds": true + "test_bounds": false }, "IAR_LU_segmented_3step": { - "skip_bounds": true + "test_bounds": false }, "IAR_LU_subtracted": { - "skip_bounds": true + "test_bounds": false }, "ETP_SRI_LinearFitting": { "options": { @@ -62,10 +62,10 @@ "Dp": 2 } }, - "skip_bounds": true + "test_bounds": false }, "OGC_AmsterdamUMC_biexp_segmented": { - "skip_bounds": true + "test_bounds": false }, "OJ_GU_seg": { "test_bounds": false diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index 9329daab..dc9fc6e6 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -177,6 +177,39 @@ def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances, npt.assert_allclose(Dp_fit,data['Dp'], rtol=tolerances["rtol"]["Dp"], atol=tolerances["atol"]["Dp"]) assert elapsed_time < 2, f"Algorithm {name} took {elapsed_time} seconds, which is longer than 2 second to fit per voxel" #less than 0.5 seconds per voxel + +def algorithms(): + # Find the algorithms from algorithms.json + file = pathlib.Path(__file__) + algorithm_path = file.with_name('algorithms.json') + with algorithm_path.open() as f: + algorithm_information = json.load(f) + algorithms = algorithm_information["algorithms"] + for algorithm in algorithms: + yield algorithm + +@pytest.mark.parametrize("algorithm", algorithms()) +def test_default_bounds_and_initial_guesses(algorithm): + fit = OsipiBase(algorithm=algorithm) + #assert fit.bounds is not None, f"For {algorithm}, there is no default fit boundary" + #assert fit.initial_guess is not None, f"For {algorithm}, there is no default fit initial guess" + if fit.bounds is not None: + assert 0 <= fit.bounds[0][0] <= 0.003, f"For {algorithm}, the default lower bound of D {fit.bounds[0][0]} is unrealistic" + assert 0 <= fit.bounds[1][0] <= 0.01, f"For {algorithm}, the default upper bound of D {fit.bounds[1][0]} is unrealistic" + assert 0 <= fit.bounds[0][1] <= 1, f"For {algorithm}, the default lower bound of f {fit.bounds[0][1]} is unrealistic" + assert 0 <= fit.bounds[1][1] <= 1, f"For {algorithm}, the default upper bound of f {fit.bounds[1][1]} is unrealistic" + assert 0.003 <= fit.bounds[0][2] <= 0.05, f"For {algorithm}, the default lower bound of Ds {fit.bounds[0][2]} is unrealistic" + assert 0.003 <= fit.bounds[1][2] <= 0.5, f"For {algorithm}, the default upper bound of Ds {fit.bounds[1][2]} is unrealistic" + assert 0 <= fit.bounds[0][3] <= 1, f"For {algorithm}, the default lower bound of S {fit.bounds[0][3]} is unrealistic; note data is normaized" + assert 1 <= fit.bounds[1][3] <= 1000, f"For {algorithm}, the default upper bound of S {fit.bounds[1][3]} is unrealistic; note data is normaized" + assert fit.bounds[1][0] <= fit.bounds[0][2], f"For {algorithm}, the default upper bound of D {fit.bounds[1][0]} is higher than lower bound of D* {fit.bounds[0][2]}" + if fit.initial_guess is not None: + assert 0.0008 <= fit.initial_guess[0] <= 0.002, f"For {algorithm}, the default initial guess for D {fit.initial_guess[0]} is unrealistic" + assert 0 <= fit.initial_guess[1] <= 0.5, f"For {algorithm}, the default initial guess for f {fit.initial_guess[1]} is unrealistic" + assert 0.003 <= fit.initial_guess[2] <= 0.1, f"For {algorithm}, the default initial guess for Ds {fit.initial_guess[2]} is unrealistic" + assert 0.9 <= fit.initial_guess[3] <= 1.1, f"For {algorithm}, the default initial guess for S {fit.initial_guess[3]} is unrealistic; note signal is normalized" + + def bound_input(): # Find the algorithms from algorithms.json file = pathlib.Path(__file__) From 889a2c115f281958d7dc03edd172657e03479258 Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Sat, 4 Jan 2025 21:45:17 +0100 Subject: [PATCH 06/10] fixed errors after merge with main Also changed fit_result["D*"] back to fit_result["Dp"], as now both names were used throughout the scripts --- WrapImage/nifti_wrapper.py | 2 +- ...github_and_IVIM_Analysis_using_Python.ipynb | 18 +++++++++--------- src/standardized/ETP_SRI_LinearFitting.py | 4 ++-- src/standardized/IAR_LU_biexp.py | 4 ++-- src/standardized/IAR_LU_modified_mix.py | 2 +- src/standardized/IAR_LU_modified_topopro.py | 2 +- src/standardized/IAR_LU_segmented_2step.py | 2 +- src/standardized/IAR_LU_segmented_3step.py | 2 +- src/standardized/IAR_LU_subtracted.py | 2 +- .../OGC_AmsterdamUMC_Bayesian_biexp.py | 2 +- src/standardized/OGC_AmsterdamUMC_biexp.py | 2 +- .../OGC_AmsterdamUMC_biexp_segmented.py | 2 +- src/standardized/OJ_GU_seg.py | 2 +- src/standardized/PV_MUMC_biexp.py | 2 +- src/standardized/PvH_KB_NKI_IVIMfit.py | 2 +- src/wrappers/OsipiBase.py | 8 ++++---- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 10 +++++----- .../unit_tests/test_ivim_synthetic.py | 4 ++-- 18 files changed, 36 insertions(+), 36 deletions(-) diff --git a/WrapImage/nifti_wrapper.py b/WrapImage/nifti_wrapper.py index 150e6106..2117e57d 100644 --- a/WrapImage/nifti_wrapper.py +++ b/WrapImage/nifti_wrapper.py @@ -111,7 +111,7 @@ def loop_over_first_n_minus_1_dimensions(arr): for idx, view in tqdm(loop_over_first_n_minus_1_dimensions(data), desc=f"{args.algorithm} is fitting", dynamic_ncols=True, total=total_iteration): fit_result = fit.osipi_fit(view, bvals) f_image.append(fit_result["f"]) - Dp_image.append(fit_result["D*"]) + Dp_image.append(fit_result["Dp"]) D_image.append(fit_result["D"]) # Convert lists to NumPy arrays diff --git a/doc/Introduction_to_TF24_IVIM-MRI_CodeCollection_github_and_IVIM_Analysis_using_Python.ipynb b/doc/Introduction_to_TF24_IVIM-MRI_CodeCollection_github_and_IVIM_Analysis_using_Python.ipynb index 63bd908c..a9bf1e09 100644 --- a/doc/Introduction_to_TF24_IVIM-MRI_CodeCollection_github_and_IVIM_Analysis_using_Python.ipynb +++ b/doc/Introduction_to_TF24_IVIM-MRI_CodeCollection_github_and_IVIM_Analysis_using_Python.ipynb @@ -532,7 +532,7 @@ { "data": { "text/plain": [ - "{'f': array(0.04609779), 'D*': array(0.01136011), 'D': array(0.00071134)}" + "{'f': array(0.04609779), 'Dp': array(0.01136011), 'D': array(0.00071134)}" ] }, "execution_count": 15, @@ -569,7 +569,7 @@ { "data": { "text/plain": [ - "{'f': array(0.04611801), 'D*': array(0.0113541), 'D': array(0.0007113)}" + "{'f': array(0.04611801), 'Dp': array(0.0113541), 'D': array(0.0007113)}" ] }, "execution_count": 16, @@ -636,10 +636,10 @@ "#plot the results of algorithm 1\n", "plt.subplot(121)\n", "plt.plot(np.unique(bval),signal_1dir,'x')\n", - "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['D*'])+(1-fit['f'])*np.exp(-np.unique(bval)*fit['D']))\n", - "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['D*']))\n", + "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['Dp'])+(1-fit['f'])*np.exp(-np.unique(bval)*fit['D']))\n", + "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['Dp']))\n", "plt.plot(np.unique(bval),(1-fit['f'])*np.exp(-np.unique(bval)*fit['D']))\n", - "plt.legend(['measured data','model fit','D*','D'])\n", + "plt.legend(['measured data','model fit','Dp','D'])\n", "plt.ylabel('S/S0')\n", "plt.xlabel('b-value [s/mm^2]')\n", "plt.title('algorithm 1')\n", @@ -650,10 +650,10 @@ "#plot the results of algorithm 2\n", "plt.subplot(122)\n", "plt.plot(np.unique(bval),signal_1dir,'x')\n", - "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['D*'])+(1-fit['f'])*np.exp(-np.unique(bval)*fit['D']))\n", - "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['D*']))\n", + "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['Dp'])+(1-fit['f'])*np.exp(-np.unique(bval)*fit['D']))\n", + "plt.plot(np.unique(bval),fit['f']*np.exp(-np.unique(bval)*fit['Dp']))\n", "plt.plot(np.unique(bval),(1-fit['f'])*np.exp(-np.unique(bval)*fit['D']))\n", - "plt.legend(['measured data','model fit','D*','D'])\n", + "plt.legend(['measured data','model fit','Dp','D'])\n", "plt.ylabel('S/S0')\n", "plt.xlabel('b-value [s/mm^2]')\n", "plt.title('algorithm 2')\n" @@ -818,7 +818,7 @@ "data": { "text/plain": [ "{'f': array([0., 0., 0., ..., 0., 0., 0.]),\n", - " 'D*': array([0., 0., 0., ..., 0., 0., 0.]),\n", + " 'Dp': array([0., 0., 0., ..., 0., 0., 0.]),\n", " 'D': array([0., 0., 0., ..., 0., 0., 0.])}" ] }, diff --git a/src/standardized/ETP_SRI_LinearFitting.py b/src/standardized/ETP_SRI_LinearFitting.py index 6d1deee0..a1cc866b 100644 --- a/src/standardized/ETP_SRI_LinearFitting.py +++ b/src/standardized/ETP_SRI_LinearFitting.py @@ -71,14 +71,14 @@ def ivim_fit(self, signals, bvalues=None, linear_fit_option=False, **kwargs): f, Dstar = ETP_object.linear_fit(bvalues, signals) results["f"] = f - results["D*"] = Dstar + results["Dp"] = Dstar return results else: f, D, Dstar = ETP_object.ivim_fit(bvalues, signals) results["f"] = f - results["D*"] = Dstar + results["Dp"] = Dstar results["D"] = D return results diff --git a/src/standardized/IAR_LU_biexp.py b/src/standardized/IAR_LU_biexp.py index 4f4f7231..b34d474c 100644 --- a/src/standardized/IAR_LU_biexp.py +++ b/src/standardized/IAR_LU_biexp.py @@ -78,7 +78,7 @@ def ivim_fit(self, signals, bvalues, **kwargs): results = {} results["f"] = fit_results.model_params[1] - results["D*"] = fit_results.model_params[2] + results["Dp"] = fit_results.model_params[2] results["D"] = fit_results.model_params[3] return results @@ -110,7 +110,7 @@ def ivim_fit_full_volume(self, signals, bvalues, **kwargs): results = {} results["f"] = fit_results.model_params[..., 1] - results["D*"] = fit_results.model_params[..., 2] + results["Dp"] = fit_results.model_params[..., 2] results["D"] = fit_results.model_params[..., 3] return results \ No newline at end of file diff --git a/src/standardized/IAR_LU_modified_mix.py b/src/standardized/IAR_LU_modified_mix.py index 9571a569..fc21560d 100644 --- a/src/standardized/IAR_LU_modified_mix.py +++ b/src/standardized/IAR_LU_modified_mix.py @@ -81,7 +81,7 @@ def ivim_fit(self, signals, bvalues, **kwargs): #D = fit_results.model_params[3] results = {} results["f"] = fit_results.model_params[1] - results["D*"] = fit_results.model_params[2] + results["Dp"] = fit_results.model_params[2] results["D"] = fit_results.model_params[3] return results \ No newline at end of file diff --git a/src/standardized/IAR_LU_modified_topopro.py b/src/standardized/IAR_LU_modified_topopro.py index 9e39fc20..8d4ac0c0 100644 --- a/src/standardized/IAR_LU_modified_topopro.py +++ b/src/standardized/IAR_LU_modified_topopro.py @@ -83,7 +83,7 @@ def ivim_fit(self, signals, bvalues, **kwargs): #return f, Dstar, D results = {} results["f"] = fit_results.model_params[1] - results["D*"] = fit_results.model_params[2] + results["Dp"] = fit_results.model_params[2] results["D"] = fit_results.model_params[3] return results \ No newline at end of file diff --git a/src/standardized/IAR_LU_segmented_2step.py b/src/standardized/IAR_LU_segmented_2step.py index 3fb5cdd8..c4e0b660 100644 --- a/src/standardized/IAR_LU_segmented_2step.py +++ b/src/standardized/IAR_LU_segmented_2step.py @@ -84,7 +84,7 @@ def ivim_fit(self, signals, bvalues, thresholds=None, **kwargs): #return f, Dstar, D results = {} results["f"] = fit_results.model_params[1] - results["D*"] = fit_results.model_params[2] + results["Dp"] = fit_results.model_params[2] results["D"] = fit_results.model_params[3] return results \ No newline at end of file diff --git a/src/standardized/IAR_LU_segmented_3step.py b/src/standardized/IAR_LU_segmented_3step.py index 170b8ee0..2711b2d2 100644 --- a/src/standardized/IAR_LU_segmented_3step.py +++ b/src/standardized/IAR_LU_segmented_3step.py @@ -83,7 +83,7 @@ def ivim_fit(self, signals, bvalues, **kwargs): #return f, Dstar, D results = {} results["f"] = fit_results.model_params[1] - results["D*"] = fit_results.model_params[2] + results["Dp"] = fit_results.model_params[2] results["D"] = fit_results.model_params[3] return results \ No newline at end of file diff --git a/src/standardized/IAR_LU_subtracted.py b/src/standardized/IAR_LU_subtracted.py index d5508f22..2fa8f154 100644 --- a/src/standardized/IAR_LU_subtracted.py +++ b/src/standardized/IAR_LU_subtracted.py @@ -83,7 +83,7 @@ def ivim_fit(self, signals, bvalues, **kwargs): #return f, Dstar, D results = {} results["f"] = fit_results.model_params[1] - results["D*"] = fit_results.model_params[2] + results["Dp"] = fit_results.model_params[2] results["D"] = fit_results.model_params[3] return results \ No newline at end of file diff --git a/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py b/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py index 1b8c8f94..620a651b 100644 --- a/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py +++ b/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py @@ -81,6 +81,6 @@ def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): results = {} results["D"] = fit_results[0] results["f"] = fit_results[1] - results["D*"] = fit_results[2] + results["Dp"] = fit_results[2] return results \ No newline at end of file diff --git a/src/standardized/OGC_AmsterdamUMC_biexp.py b/src/standardized/OGC_AmsterdamUMC_biexp.py index 1516ec5d..55595904 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp.py @@ -71,6 +71,6 @@ def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): results = {} results["D"] = fit_results[0] results["f"] = fit_results[1] - results["D*"] = fit_results[2] + results["Dp"] = fit_results[2] return results \ No newline at end of file diff --git a/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py b/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py index c452f7b5..9153f773 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py @@ -72,6 +72,6 @@ def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): results = {} results["D"] = fit_results[0] results["f"] = fit_results[1] - results["D*"] = fit_results[2] + results["Dp"] = fit_results[2] return results \ No newline at end of file diff --git a/src/standardized/OJ_GU_seg.py b/src/standardized/OJ_GU_seg.py index a62f62cd..eb9e2e2a 100644 --- a/src/standardized/OJ_GU_seg.py +++ b/src/standardized/OJ_GU_seg.py @@ -66,7 +66,7 @@ def ivim_fit(self, signals, bvalues=None): results = {} results["f"] = fit_results['f'] - results["D*"] = fit_results['Dstar'] + results["Dp"] = fit_results['Dstar'] results["D"] = fit_results['D'] return results \ No newline at end of file diff --git a/src/standardized/PV_MUMC_biexp.py b/src/standardized/PV_MUMC_biexp.py index f08ddd1c..10793870 100644 --- a/src/standardized/PV_MUMC_biexp.py +++ b/src/standardized/PV_MUMC_biexp.py @@ -51,7 +51,7 @@ def ivim_fit(self, signals, bvalues=None): results = {} results["f"] = fit_results[1] - results["D*"] = fit_results[2] + results["Dp"] = fit_results[2] results["D"] = fit_results[0] return results diff --git a/src/standardized/PvH_KB_NKI_IVIMfit.py b/src/standardized/PvH_KB_NKI_IVIMfit.py index ea98c093..804e73ea 100644 --- a/src/standardized/PvH_KB_NKI_IVIMfit.py +++ b/src/standardized/PvH_KB_NKI_IVIMfit.py @@ -59,6 +59,6 @@ def ivim_fit(self, signals, bvalues=None): results = {} results["D"] = fit_results[0][0,0,0]/1000 results["f"] = fit_results[1][0,0,0] - results["D*"] = fit_results[2][0,0,0]/1000 + results["Dp"] = fit_results[2][0,0,0]/1000 return results diff --git a/src/wrappers/OsipiBase.py b/src/wrappers/OsipiBase.py index 715f1c4c..3a41cfe8 100644 --- a/src/wrappers/OsipiBase.py +++ b/src/wrappers/OsipiBase.py @@ -89,8 +89,8 @@ def osipi_fit(self, data, bvalues=None, **kwargs): # result_keys is a list of strings of parameter names, e.g. "S0", "f1", "f2", etc. result_keys = self.result_keys else: - # Default is ["f", "D*", "D"] - self.result_keys = ["f", "D*", "D"] + # Default is ["f", "Dp", "D"] + self.result_keys = ["f", "Dp", "D"] results = {} for key in self.result_keys: @@ -132,8 +132,8 @@ def osipi_fit_full_volume(self, data, bvalues=None, **kwargs): # result_keys is a list of strings of parameter names, e.g. "S0", "f1", "f2", etc. result_keys = self.result_keys else: - # Default is ["f", "D*", "D"] - self.result_keys = ["f", "D*", "D"] + # Default is ["f", "Dp", "D"] + self.result_keys = ["f", "Dp", "D"] # Create the results dictionary results = {} diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index ae3381a4..72d9b0bf 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -61,7 +61,7 @@ def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances, if xfail["xfail"]: mark = pytest.mark.xfail(reason="xfail", strict=xfail["strict"]) request.node.add_marker(mark) - signal, ratio = signal_helper(data["data"]) + signal = signal_helper(data["data"]) tolerances = tolerances_helper(tolerances, data) start_time = time.time() # Record the start time fit = OsipiBase(algorithm=algorithm, **kwargs) @@ -73,7 +73,7 @@ def to_list_if_needed(value): "name": name, "algorithm": algorithm, "f_fit": to_list_if_needed(fit_result['f']), - "Dp_fit": to_list_if_needed(fit_result['D*']), + "Dp_fit": to_list_if_needed(fit_result['Dp']), "D_fit": to_list_if_needed(fit_result['D']), "f": to_list_if_needed(data['f']), "Dp": to_list_if_needed(data['Dp']), @@ -158,12 +158,12 @@ def test_bounds(name, bvals, data, algorithm, xfail, kwargs, tolerances, request # deliberately have silly bounds to see whether they are used fit = OsipiBase(algorithm=algorithm, bounds=bounds, initial_guess = [0.001, 0.25, 0.015, 1.2], **kwargs) signal = signal_helper(data["data"]) - tolerances = tolerances_helper(tolerances, data) - fit_results = fit.osipi_fit(signal, bvals) + fit_result = fit.osipi_fit(signal, bvals) assert bounds[0][0] <= fit_result['D'] <= bounds[1][0], f"Result {fit_result['D']} out of bounds for data: {name}" assert bounds[0][1] <= fit_result['f'] <= bounds[1][1], f"Result {fit_result['f']} out of bounds for data: {name}" assert bounds[0][2] <= fit_result['Dp'] <= bounds[1][2], f"Result {fit_result['Dp']} out of bounds for data: {name}" - assert bounds[0][3] <= fit_result['S0'] <= bounds[1][3], f"Result {fit_result['S0']} out of bounds for data: {name}" + # S0 is not returned as argument... + #assert bounds[0][3] <= fit_result['S0'] <= bounds[1][3], f"Result {fit_result['S0']} out of bounds for data: {name}" \ No newline at end of file diff --git a/tests/IVIMmodels/unit_tests/test_ivim_synthetic.py b/tests/IVIMmodels/unit_tests/test_ivim_synthetic.py index 58f20eb5..20741327 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_synthetic.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_synthetic.py @@ -40,9 +40,9 @@ def test_generated(ivim_algorithm, ivim_data, SNR, rtol, atol, fit_count, rician fit_result = fit.osipi_fit(signal, bvals) #, prior_in=prior time_delta += datetime.datetime.now() - start_time if save_file is not None: - save_file.writerow([ivim_algorithm, name, SNR, idx, f, Dp, D, fit_result["f"], fit_result["D*"], fit_result["D"], *signal]) + save_file.writerow([ivim_algorithm, name, SNR, idx, f, Dp, D, fit_result["f"], fit_result["Dp"], fit_result["D"], *signal]) # save_results(save_file, ivim_algorithm, name, SNR, idx, [f, Dp, D], [f_fit, Dp_fit, D_fit]) - npt.assert_allclose([f, Dp, D], [fit_result["f"], fit_result["D*"], fit_result["D"]], rtol, atol) + npt.assert_allclose([f, Dp, D], [fit_result["f"], fit_result["Dp"], fit_result["D"]], rtol, atol) if save_duration_file is not None: save_duration_file.writerow([ivim_algorithm, name, SNR, time_delta/datetime.timedelta(microseconds=1), fit_count]) # save_duration(save_duration_file, ivim_algorithm, name, SNR, time_delta, fit_count) From 66cb0db2c642dadc3c3ae0946cab473286763f58 Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Mon, 6 Jan 2025 16:57:51 +0100 Subject: [PATCH 07/10] bounds test automatically check whether bounds are relevant --- src/original/OGC_AmsterdamUMC/LSQ_fitting.py | 45 ++++++++++--------- src/standardized/ETP_SRI_LinearFitting.py | 5 ++- src/standardized/IAR_LU_biexp.py | 5 ++- src/standardized/IAR_LU_modified_mix.py | 5 ++- src/standardized/IAR_LU_modified_topopro.py | 5 ++- src/standardized/IAR_LU_segmented_2step.py | 5 ++- src/standardized/IAR_LU_segmented_3step.py | 5 ++- src/standardized/IAR_LU_subtracted.py | 5 ++- .../OGC_AmsterdamUMC_Bayesian_biexp.py | 31 ++++++++----- src/standardized/OGC_AmsterdamUMC_biexp.py | 5 ++- .../OGC_AmsterdamUMC_biexp_segmented.py | 11 ++--- src/standardized/OJ_GU_seg.py | 5 ++- src/standardized/PV_MUMC_biexp.py | 4 ++ src/standardized/PvH_KB_NKI_IVIMfit.py | 4 ++ src/wrappers/OsipiBase.py | 3 +- tests/IVIMmodels/unit_tests/algorithms.json | 35 +-------------- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 42 +++++++++-------- 17 files changed, 118 insertions(+), 102 deletions(-) diff --git a/src/original/OGC_AmsterdamUMC/LSQ_fitting.py b/src/original/OGC_AmsterdamUMC/LSQ_fitting.py index 5720d790..b4cec611 100644 --- a/src/original/OGC_AmsterdamUMC/LSQ_fitting.py +++ b/src/original/OGC_AmsterdamUMC/LSQ_fitting.py @@ -124,29 +124,30 @@ def fit_segmented(bvalues, dw_data, bounds=([0, 0, 0.005],[0.005, 0.7, 0.2]), cu :return Dp: Fitted Dp :return S0: Fitted S0 """ - p0 = [p0[0] * 1000, p0[1] * 10, p0[2] * 10, p0[3]] try: # determine high b-values and data for D + dw_data=dw_data/np.mean(dw_data[bvalues==0]) high_b = bvalues[bvalues >= cutoff] high_dw_data = dw_data[bvalues >= cutoff] # correct the bounds. Note that S0 bounds determine the max and min of f - bounds1 = ([bounds[0][0] * 1000., 0.7 - bounds[1][1]], [bounds[1][0] * 1000., 1.3 - bounds[0][ - 1]]) # By bounding S0 like this, we effectively insert the boundaries of f - # fit for S0' and D - params, _ = curve_fit(lambda b, Dt, int: int * np.exp(-b * Dt / 1000), high_b, high_dw_data, - p0=(p0[0], p0[3]-p0[1]/10), + bounds1 = ([bounds[0][0], 0], [bounds[1][0], 10000000000]) + params, _ = curve_fit(lambda b, Dt, int: int * np.exp(-b * Dt ), high_b, high_dw_data, + p0=(p0[0], p0[3]-p0[1]), bounds=bounds1) - Dt, Fp = params[0] / 1000, 1 - params[1] + Dt, Fp = 0+params[0], 1 - params[1] + if Fp < bounds[0][1] : Fp = bounds[0][1] + if Fp > bounds[1][1] : Fp = bounds[1][1] + # remove the diffusion part to only keep the pseudo-diffusion dw_data_remaining = dw_data - (1 - Fp) * np.exp(-bvalues * Dt) - bounds2 = (bounds[0][2]*10, bounds[1][2]*10) + bounds2 = (bounds[0][2], bounds[1][2]) # fit for D* params, _ = curve_fit(lambda b, Dp: Fp * np.exp(-b * Dp), bvalues, dw_data_remaining, p0=(p0[2]), bounds=bounds2) - Dp = params[0] - return Dt, Fp, Dp + Dp = 0+params[0] + return Dt, np.float64(Fp), Dp except: # if fit fails, return zeros - # print('segnetned fit failed') + # print('segmented fit failed') return 0., 0., 0. @@ -235,17 +236,17 @@ def fit_least_squares(bvalues, dw_data, S0_output=False, fitS0=True, try: if not fitS0: # bounds are rescaled such that each parameter changes at roughly the same rate to help fitting. - bounds = ([bounds[0][0] * 1000, bounds[0][1] * 10, bounds[0][2] * 10], + bounds2 = ([bounds[0][0] * 1000, bounds[0][1] * 10, bounds[0][2] * 10], [bounds[1][0] * 1000, bounds[1][1] * 10, bounds[1][2] * 10]) p0=[p0[0]*1000,p0[1]*10,p0[2]*10] - params, _ = curve_fit(ivimN_noS0, bvalues, dw_data, p0=p0, bounds=bounds) + params, _ = curve_fit(ivimN_noS0, bvalues, dw_data, p0=p0, bounds=bounds2) S0 = 1 else: # bounds are rescaled such that each parameter changes at roughly the same rate to help fitting. - bounds = ([bounds[0][0] * 1000, bounds[0][1] * 10, bounds[0][2] * 10, bounds[0][3]], + bounds2 = ([bounds[0][0] * 1000, bounds[0][1] * 10, bounds[0][2] * 10, bounds[0][3]], [bounds[1][0] * 1000, bounds[1][1] * 10, bounds[1][2] * 10, bounds[1][3]]) p0=[p0[0]*1000,p0[1]*10,p0[2]*10,p0[3]] - params, _ = curve_fit(ivimN, bvalues, dw_data, p0=p0, bounds=bounds) + params, _ = curve_fit(ivimN, bvalues, dw_data, p0=p0, bounds=bounds2) S0 = params[3] # correct for the rescaling of parameters Dt, Fp, Dp = params[0] / 1000, params[1] / 10, params[2] / 10 @@ -261,7 +262,7 @@ def fit_least_squares(bvalues, dw_data, S0_output=False, fitS0=True, Dt, Fp, Dp = fit_segmented(bvalues, dw_data, bounds=bounds) return Dt, Fp, Dp, 1 else: - return fit_segmented(bvalues, dw_data) + return fit_segmented(bvalues, dw_data, bounds=bounds) def fit_least_squares_array_tri_exp(bvalues, dw_data, S0_output=True, fitS0=True, njobs=4, @@ -561,19 +562,19 @@ def neg_log_prior(p): Dt, Fp, Dp = p[0], p[1], p[2] # make D* self.bounds[1][i] : fit_results[0] = self.bounds[1][i]-epsilon + fit_results = self.OGC_algorithm(bvalues, signals, self.neg_log_prior, x0=fit_results, fitS0=self.fitS0, bounds=self.bounds) results = {} results["D"] = fit_results[0] diff --git a/src/standardized/OGC_AmsterdamUMC_biexp.py b/src/standardized/OGC_AmsterdamUMC_biexp.py index 55595904..84c7c8cd 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp.py @@ -38,20 +38,23 @@ def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=Non #super(OGC_AmsterdamUMC_biexp, self).__init__(bvalues, bounds, initial_guess, fitS0) super(OGC_AmsterdamUMC_biexp, self).__init__(bvalues=bvalues, bounds=bounds, initial_guess=initial_guess) self.OGC_algorithm = fit_least_squares - #self.initialize(bounds, initial_guess, fitS0) self.fitS0=fitS0 + self.initialize(bounds, initial_guess, fitS0) def initialize(self, bounds, initial_guess, fitS0): if bounds is None: self.bounds=([0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]) else: self.bounds=bounds + self.use_bounds = True if initial_guess is None: self.initial_guess = [0.001, 0.1, 0.01, 1] else: self.initial_guess = initial_guess + self.use_initial_guess = True self.fitS0=fitS0 + def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): """Perform the IVIM fit diff --git a/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py b/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py index 9153f773..07bc6c37 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py @@ -27,7 +27,7 @@ class OGC_AmsterdamUMC_biexp_segmented(OsipiBase): required_initial_guess_optional = True accepted_dimensions = 1 # Not sure how to define this for the number of accepted dimensions. Perhaps like the thresholds, at least and at most? - def __init__(self, bvalues=None, thresholds=150, bounds=([0, 0, 0.005],[0.005, 0.7, 0.2]), initial_guess=[0.001, 0.01, 0.01,1]): + def __init__(self, bvalues=None, thresholds=150, bounds=None, initial_guess=None): """ Everything this algorithm requires should be implemented here. Number of segmentation thresholds, bounds, etc. @@ -42,19 +42,21 @@ def __init__(self, bvalues=None, thresholds=150, bounds=([0, 0, 0.005],[0.005, 0 def initialize(self, bounds, initial_guess, thresholds=300): if bounds is None: - self.bounds=([0, 0, 0.005, 0.7],[0.005, 0.7, 0.2, 1.3]) + self.bounds=([0, 0, 0.005, 0.7],[0.005, 1, 0.1, 1.3]) else: self.bounds=bounds + self.use_bounds=True if initial_guess is None: self.initial_guess = [0.001, 0.1, 0.03, 1] else: self.initial_guess = initial_guess + self.use_initial_guess=True if thresholds is None: self.thresholds = 150 else: self.thresholds = thresholds - def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): + def ivim_fit(self, signals, bvalues, **kwargs): """Perform the IVIM fit Args: @@ -64,8 +66,7 @@ def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): Returns: _type_: _description_ """ - if initial_guess is not None and len(initial_guess) == 4: - self.initial_guess = initial_guess + bvalues=np.array(bvalues) fit_results = self.OGC_algorithm(bvalues, signals, bounds=self.bounds, cutoff=self.thresholds, p0=self.initial_guess) diff --git a/src/standardized/OJ_GU_seg.py b/src/standardized/OJ_GU_seg.py index eb9e2e2a..5ad2f170 100644 --- a/src/standardized/OJ_GU_seg.py +++ b/src/standardized/OJ_GU_seg.py @@ -35,7 +35,10 @@ def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=Non the requirements. """ super(OJ_GU_seg, self).__init__(bvalues, thresholds, bounds, initial_guess) - + if bounds is not None: + print('warning, bounds from wrapper are not (yet) used in this algorithm') + self.use_bounds = False + self.use_initial_guess = False # Check the inputs # Initialize the algorithm diff --git a/src/standardized/PV_MUMC_biexp.py b/src/standardized/PV_MUMC_biexp.py index 10793870..d99e3595 100644 --- a/src/standardized/PV_MUMC_biexp.py +++ b/src/standardized/PV_MUMC_biexp.py @@ -33,6 +33,10 @@ def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=Non """ super(PV_MUMC_biexp, self).__init__(bvalues, None, bounds, None) self.PV_algorithm = fit_least_squares + if bounds is not None: + print('warning, bounds from wrapper are not (yet) used in this algorithm') + self.use_bounds = False + self.use_initial_guess = False def ivim_fit(self, signals, bvalues=None): diff --git a/src/standardized/PvH_KB_NKI_IVIMfit.py b/src/standardized/PvH_KB_NKI_IVIMfit.py index 804e73ea..4c1952c6 100644 --- a/src/standardized/PvH_KB_NKI_IVIMfit.py +++ b/src/standardized/PvH_KB_NKI_IVIMfit.py @@ -37,6 +37,10 @@ def __init__(self, bvalues=None, thresholds=None,bounds=None,initial_guess=None) """ super(PvH_KB_NKI_IVIMfit, self).__init__(bvalues, thresholds,bounds,initial_guess) self.NKI_algorithm = generate_IVIMmaps_standalone + if bounds is not None: + print('warning, bounds from wrapper are not (yet) used in this algorithm') + self.use_bounds = False + self.use_initial_guess = False def ivim_fit(self, signals, bvalues=None): diff --git a/src/wrappers/OsipiBase.py b/src/wrappers/OsipiBase.py index 3a41cfe8..6c5d380e 100644 --- a/src/wrappers/OsipiBase.py +++ b/src/wrappers/OsipiBase.py @@ -14,7 +14,8 @@ def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=Non self.thresholds = np.asarray(thresholds) if thresholds is not None else None self.bounds = np.asarray(bounds) if bounds is not None else None self.initial_guess = np.asarray(initial_guess) if initial_guess is not None else None - + self.use_bounds = True + self.use_initial_guess = True # If the user inputs an algorithm to OsipiBase, it is intereprete as initiating # an algorithm object with that name. if algorithm: diff --git a/tests/IVIMmodels/unit_tests/algorithms.json b/tests/IVIMmodels/unit_tests/algorithms.json index 3046d2ff..51a56500 100644 --- a/tests/IVIMmodels/unit_tests/algorithms.json +++ b/tests/IVIMmodels/unit_tests/algorithms.json @@ -15,22 +15,7 @@ "OJ_GU_seg" ], "IAR_LU_biexp": { - "test_bounds": false - }, - "IAR_LU_modified_mix": { - "test_bounds": false - }, - "IAR_LU_modified_topopro": { - "test_bounds": false - }, - "IAR_LU_segmented_2step": { - "test_bounds": false - }, - "IAR_LU_segmented_3step": { - "test_bounds": false - }, - "IAR_LU_subtracted": { - "test_bounds": false + "fail_first_time": true }, "ETP_SRI_LinearFitting": { "options": { @@ -63,23 +48,5 @@ } }, "test_bounds": false - }, - "OGC_AmsterdamUMC_biexp_segmented": { - "test_bounds": false - }, - "OJ_GU_seg": { - "test_bounds": false - }, - "OGC_AmsterdamUMC_Bayesian_biexp": { - "test_bounds": false - }, - "MUMC_biexp": { - "test_bounds": false - }, - "PvH_KB_NKI_IVIMfit": { - "test_bounds": false - }, - "OGC_AmsterdamUMC_biexp": { - "test_bounds": true } } diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index 72d9b0bf..3a405141 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -47,6 +47,7 @@ def data_ivim_fit_saved(): algorithms = algorithm_information["algorithms"] bvals = all_data.pop('config') bvals = bvals['bvalues'] + first=True for name, data in all_data.items(): for algorithm in algorithms: algorithm_dict = algorithm_information.get(algorithm, {}) @@ -54,10 +55,15 @@ def data_ivim_fit_saved(): "strict": algorithm_dict.get("xfail_names", {}).get(name, True)} kwargs = algorithm_dict.get("options", {}) tolerances = algorithm_dict.get("tolerances", {}) - yield name, bvals, data, algorithm, xfail, kwargs, tolerances - -@pytest.mark.parametrize("name, bvals, data, algorithm, xfail, kwargs, tolerances", data_ivim_fit_saved()) -def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances, request, record_property): + skiptime=False + if first == True: + if algorithm_dict.get("fail_first_time", {}) == True: + skiptime = True + first = False + yield name, bvals, data, algorithm, xfail, kwargs, tolerances, skiptime + +@pytest.mark.parametrize("name, bvals, data, algorithm, xfail, kwargs, tolerances, skiptime", data_ivim_fit_saved()) +def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances,skiptime, request, record_property): if xfail["xfail"]: mark = pytest.mark.xfail(reason="xfail", strict=xfail["strict"]) request.node.add_marker(mark) @@ -87,7 +93,9 @@ def to_list_if_needed(value): npt.assert_allclose(fit_result['D'],data['D'], rtol=tolerances["rtol"]["D"], atol=tolerances["atol"]["D"]) if data['f']>0.03: #we need some f for D* to be interpretable npt.assert_allclose(fit_result['Dp'],data['Dp'], rtol=tolerances["rtol"]["Dp"], atol=tolerances["atol"]["Dp"]) - assert elapsed_time < 2, f"Algorithm {name} took {elapsed_time} seconds, which is longer than 2 second to fit per voxel" #less than 0.5 seconds per voxel + #assert fit_result['D'] < fit_result['Dp'], f"D {fit_result['D']} is larger than D* {fit_result['Dp']} for {name}" + if not skiptime: + assert elapsed_time < 0.5, f"Algorithm {name} took {elapsed_time} seconds, which is longer than 2 second to fit per voxel" #less than 0.5 seconds per voxel def algorithms(): @@ -144,26 +152,22 @@ def bound_input(): "strict": algorithm_dict.get("xfail_names", {}).get(name, True)} kwargs = algorithm_dict.get("options", {}) tolerances = algorithm_dict.get("tolerances", {}) - test_bounds = algorithm_dict.get("test_bounds", {}) - if test_bounds: - yield name, bvals, data, algorithm, xfail, kwargs, tolerances + yield name, bvals, data, algorithm, xfail, kwargs, tolerances @pytest.mark.parametrize("name, bvals, data, algorithm, xfail, kwargs, tolerances", bound_input()) def test_bounds(name, bvals, data, algorithm, xfail, kwargs, tolerances, request): bounds = ([0.0008, 0.2, 0.01, 1.1], [0.0012, 0.3, 0.02, 1.3]) - if xfail["xfail"]: - mark = pytest.mark.xfail(reason="xfail", strict=xfail["strict"]) - request.node.add_marker(mark) # deliberately have silly bounds to see whether they are used fit = OsipiBase(algorithm=algorithm, bounds=bounds, initial_guess = [0.001, 0.25, 0.015, 1.2], **kwargs) - signal = signal_helper(data["data"]) - fit_result = fit.osipi_fit(signal, bvals) - - assert bounds[0][0] <= fit_result['D'] <= bounds[1][0], f"Result {fit_result['D']} out of bounds for data: {name}" - assert bounds[0][1] <= fit_result['f'] <= bounds[1][1], f"Result {fit_result['f']} out of bounds for data: {name}" - assert bounds[0][2] <= fit_result['Dp'] <= bounds[1][2], f"Result {fit_result['Dp']} out of bounds for data: {name}" - # S0 is not returned as argument... - #assert bounds[0][3] <= fit_result['S0'] <= bounds[1][3], f"Result {fit_result['S0']} out of bounds for data: {name}" + if fit.use_bounds: + signal = signal_helper(data["data"]) + fit_result = fit.osipi_fit(signal, bvals) + + assert bounds[0][0] <= fit_result['D'] <= bounds[1][0], f"Result {fit_result['D']} out of bounds for data: {name}" + assert bounds[0][1] <= fit_result['f'] <= bounds[1][1], f"Result {fit_result['f']} out of bounds for data: {name}" + assert bounds[0][2] <= fit_result['Dp'] <= bounds[1][2], f"Result {fit_result['Dp']} out of bounds for data: {name}" + # S0 is not returned as argument... + #assert bounds[0][3] <= fit_result['S0'] <= bounds[1][3], f"Result {fit_result['S0']} out of bounds for data: {name}" \ No newline at end of file From ba9fa9186eeb5137cfca2658d557fefff1308a68 Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Fri, 7 Mar 2025 13:50:10 +0100 Subject: [PATCH 08/10] update --- src/original/OGC_AmsterdamUMC/LSQ_fitting.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/original/OGC_AmsterdamUMC/LSQ_fitting.py b/src/original/OGC_AmsterdamUMC/LSQ_fitting.py index b4cec611..053c1eb4 100644 --- a/src/original/OGC_AmsterdamUMC/LSQ_fitting.py +++ b/src/original/OGC_AmsterdamUMC/LSQ_fitting.py @@ -134,17 +134,17 @@ def fit_segmented(bvalues, dw_data, bounds=([0, 0, 0.005],[0.005, 0.7, 0.2]), cu params, _ = curve_fit(lambda b, Dt, int: int * np.exp(-b * Dt ), high_b, high_dw_data, p0=(p0[0], p0[3]-p0[1]), bounds=bounds1) - Dt, Fp = 0+params[0], 1 - params[1] - if Fp < bounds[0][1] : Fp = bounds[0][1] - if Fp > bounds[1][1] : Fp = bounds[1][1] + Dt, Fp = params[0], 1 - params[1] + if Fp < bounds[0][1] : Fp = np.float64(bounds[0][1]) + if Fp > bounds[1][1] : Fp = np.float64(bounds[1][1]) # remove the diffusion part to only keep the pseudo-diffusion dw_data_remaining = dw_data - (1 - Fp) * np.exp(-bvalues * Dt) bounds2 = (bounds[0][2], bounds[1][2]) # fit for D* params, _ = curve_fit(lambda b, Dp: Fp * np.exp(-b * Dp), bvalues, dw_data_remaining, p0=(p0[2]), bounds=bounds2) - Dp = 0+params[0] - return Dt, np.float64(Fp), Dp + Dp = params[0] + return Dt, Fp, Dp except: # if fit fails, return zeros # print('segmented fit failed') From d785c4d64761e7f1cdda5716d64bfc774064a325 Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Wed, 12 Mar 2025 10:59:41 +0100 Subject: [PATCH 09/10] Automated testing compatible with Wrapper some work on initial guess testing, but did not work robust. --- src/original/OGC_AmsterdamUMC/LSQ_fitting.py | 14 ++++--- .../OGC_AmsterdamUMC_Bayesian_biexp.py | 21 +++++++--- src/standardized/OGC_AmsterdamUMC_biexp.py | 24 ++++++----- .../OGC_AmsterdamUMC_biexp_segmented.py | 25 ++++++++---- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 40 +++++++++++++++++-- 5 files changed, 93 insertions(+), 31 deletions(-) diff --git a/src/original/OGC_AmsterdamUMC/LSQ_fitting.py b/src/original/OGC_AmsterdamUMC/LSQ_fitting.py index 053c1eb4..f0e0e60b 100644 --- a/src/original/OGC_AmsterdamUMC/LSQ_fitting.py +++ b/src/original/OGC_AmsterdamUMC/LSQ_fitting.py @@ -130,6 +130,8 @@ def fit_segmented(bvalues, dw_data, bounds=([0, 0, 0.005],[0.005, 0.7, 0.2]), cu high_b = bvalues[bvalues >= cutoff] high_dw_data = dw_data[bvalues >= cutoff] # correct the bounds. Note that S0 bounds determine the max and min of f + assert bounds[0][1] < p0[1] + assert bounds[1][1] > p0[1] bounds1 = ([bounds[0][0], 0], [bounds[1][0], 10000000000]) params, _ = curve_fit(lambda b, Dt, int: int * np.exp(-b * Dt ), high_b, high_dw_data, p0=(p0[0], p0[3]-p0[1]), @@ -238,15 +240,15 @@ def fit_least_squares(bvalues, dw_data, S0_output=False, fitS0=True, # bounds are rescaled such that each parameter changes at roughly the same rate to help fitting. bounds2 = ([bounds[0][0] * 1000, bounds[0][1] * 10, bounds[0][2] * 10], [bounds[1][0] * 1000, bounds[1][1] * 10, bounds[1][2] * 10]) - p0=[p0[0]*1000,p0[1]*10,p0[2]*10] - params, _ = curve_fit(ivimN_noS0, bvalues, dw_data, p0=p0, bounds=bounds2) + p1=[p0[0]*1000,p0[1]*10,p0[2]*10] + params, _ = curve_fit(ivimN_noS0, bvalues, dw_data, p0=p1, bounds=bounds2) S0 = 1 else: # bounds are rescaled such that each parameter changes at roughly the same rate to help fitting. bounds2 = ([bounds[0][0] * 1000, bounds[0][1] * 10, bounds[0][2] * 10, bounds[0][3]], [bounds[1][0] * 1000, bounds[1][1] * 10, bounds[1][2] * 10, bounds[1][3]]) - p0=[p0[0]*1000,p0[1]*10,p0[2]*10,p0[3]] - params, _ = curve_fit(ivimN, bvalues, dw_data, p0=p0, bounds=bounds2) + p1=[p0[0]*1000,p0[1]*10,p0[2]*10,p0[3]] + params, _ = curve_fit(ivimN, bvalues, dw_data, p0=p1, bounds=bounds2) S0 = params[3] # correct for the rescaling of parameters Dt, Fp, Dp = params[0] / 1000, params[1] / 10, params[2] / 10 @@ -259,10 +261,10 @@ def fit_least_squares(bvalues, dw_data, S0_output=False, fitS0=True, # if fit fails, then do a segmented fit instead # print('lsq fit failed, trying segmented') if S0_output: - Dt, Fp, Dp = fit_segmented(bvalues, dw_data, bounds=bounds) + Dt, Fp, Dp = fit_segmented(bvalues, dw_data, bounds=bounds,p0=p0) return Dt, Fp, Dp, 1 else: - return fit_segmented(bvalues, dw_data, bounds=bounds) + return fit_segmented(bvalues, dw_data, bounds=bounds,p0=p0) def fit_least_squares_array_tri_exp(bvalues, dw_data, S0_output=True, fitS0=True, njobs=4, diff --git a/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py b/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py index 4c5c671a..a248166d 100644 --- a/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py +++ b/src/standardized/OGC_AmsterdamUMC_Bayesian_biexp.py @@ -34,7 +34,7 @@ class OGC_AmsterdamUMC_Bayesian_biexp(OsipiBase): supported_initial_guess = True supported_thresholds = True - def __init__(self, bvalues=None, thresholds=None, bounds=([0, 0, 0.005, 0.7],[0.005, 0.7, 0.2, 1.3]), initial_guess=[0.001, 0.1, 0.01, 1], fitS0=True, prior_in=None): + def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=None, fitS0=True, prior_in=None): """ Everything this algorithm requires should be implemented here. @@ -51,11 +51,20 @@ def __init__(self, bvalues=None, thresholds=None, bounds=([0, 0, 0.005, 0.7],[0. self.initialize(bounds, initial_guess, fitS0, prior_in) self.fit_segmented=fit_segmented - def initialize(self, bounds=([0, 0, 0.005, 0.7],[0.005, 0.7, 0.2, 1.3]), initial_guess=[0.001, 0.1, 0.01, 1], fitS0=True, prior_in=None): - self.bounds=bounds - self.use_bounds = True - self.initial_guess = initial_guess + def initialize(self, bounds=None, initial_guess=None, fitS0=True, prior_in=None): + if bounds is None: + print('warning, no bounds were defined, so default bounds are used of [0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]') + self.bounds=([0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]) + else: + self.bounds=bounds + if initial_guess is None: + print('warning, no initial guesses were defined, so default bounds are used of [0.001, 0.001, 0.01, 1]') + self.initial_guess = [0.001, 0.001, 0.01, 1] + else: + self.initial_guess = initial_guess self.use_initial_guess = True + self.use_bounds = True + self.thresholds = 150 if prior_in is None: print('using a flat prior between bounds') self.neg_log_prior=flat_neg_log_prior([self.bounds[0][0],self.bounds[1][0]],[self.bounds[0][1],self.bounds[1][1]],[self.bounds[0][2],self.bounds[1][2]],[self.bounds[0][3],self.bounds[1][3]]) @@ -81,7 +90,7 @@ def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): epsilon = 0.000001 fit_results = fit_segmented(bvalues, signals, bounds=self.bounds, cutoff=self.thresholds, p0=self.initial_guess) - fit_results=fit_results+(1,) + fit_results=np.array(fit_results+(1,)) for i in range(4): if fit_results[i] < self.bounds[0][i] : fit_results[0] = self.bounds[0][i]+epsilon if fit_results[i] > self.bounds[1][i] : fit_results[0] = self.bounds[1][i]-epsilon diff --git a/src/standardized/OGC_AmsterdamUMC_biexp.py b/src/standardized/OGC_AmsterdamUMC_biexp.py index 901a3596..49a16b27 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp.py @@ -33,7 +33,7 @@ class OGC_AmsterdamUMC_biexp(OsipiBase): supported_initial_guess = True supported_thresholds = False - def __init__(self, bvalues=None, thresholds=None, bounds=([0, 0, 0.005, 0.7],[0.005, 0.7, 0.2, 1.3]), initial_guess=[0.001, 0.1, 0.01, 1], fitS0=True): + def __init__(self, bvalues=None, thresholds=None, bounds=None, initial_guess=None, fitS0=True): """ Everything this algorithm requires should be implemented here. Number of segmentation thresholds, bounds, etc. @@ -48,14 +48,22 @@ def __init__(self, bvalues=None, thresholds=None, bounds=([0, 0, 0.005, 0.7],[0. self.initialize(bounds, initial_guess, fitS0) def initialize(self, bounds, initial_guess, fitS0): - self.bounds=bounds - self.use_bounds = True - self.initial_guess = initial_guess - self.use_initial_guess = True + if bounds is None: + print('warning, no bounds were defined, so default bounds are used of [0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]') + self.bounds=([0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]) + else: + self.bounds=bounds + if initial_guess is None: + print('warning, no initial guesses were defined, so default bounds are used of [0.001, 0.001, 0.01, 1]') + self.initial_guess = [0.001, 0.001, 0.01, 1] + else: + self.initial_guess = initial_guess + self.use_initial_guess = True self.fitS0=fitS0 + self.use_initial_guess = True + self.use_bounds = True - - def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): + def ivim_fit(self, signals, bvalues, **kwargs): """Perform the IVIM fit Args: @@ -66,8 +74,6 @@ def ivim_fit(self, signals, bvalues, initial_guess=None, **kwargs): _type_: _description_ """ - if initial_guess is not None and len(initial_guess) == 4: - self.initial_guess = initial_guess bvalues=np.array(bvalues) fit_results = self.OGC_algorithm(bvalues, signals, p0=self.initial_guess, bounds=self.bounds, fitS0=self.fitS0) diff --git a/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py b/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py index 96db65f7..68f6e5bc 100644 --- a/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py +++ b/src/standardized/OGC_AmsterdamUMC_biexp_segmented.py @@ -33,7 +33,7 @@ class OGC_AmsterdamUMC_biexp_segmented(OsipiBase): supported_initial_guess = True supported_thresholds = True - def __init__(self, bvalues=None, thresholds=150, bounds=([0, 0, 0.005],[0.005, 0.7, 0.2]), initial_guess=[0.001, 0.01, 0.01,1]): + def __init__(self, bvalues=None, thresholds=150, bounds=None, initial_guess=None): """ Everything this algorithm requires should be implemented here. Number of segmentation thresholds, bounds, etc. @@ -47,12 +47,23 @@ def __init__(self, bvalues=None, thresholds=150, bounds=([0, 0, 0.005],[0.005, 0 def initialize(self, bounds, initial_guess, thresholds): - self.bounds=bounds - self.use_bounds=True - self.initial_guess = initial_guess - self.use_initial_guess=True - self.thresholds = thresholds - + if bounds is None: + print('warning, no bounds were defined, so default bounds are used of [0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]') + self.bounds=([0, 0, 0.005, 0.7],[0.005, 1.0, 0.2, 1.3]) + else: + self.bounds=bounds + if initial_guess is None: + print('warning, no initial guesses were defined, so default bounds are used of [0.001, 0.001, 0.01, 1]') + self.initial_guess = [0.001, 0.001, 0.01, 1] + else: + self.initial_guess = initial_guess + self.use_initial_guess = True + self.use_bounds = True + if thresholds is None: + self.thresholds = 150 + print('warning, no thresholds were defined, so default bounds are used of 150') + else: + self.thresholds = thresholds def ivim_fit(self, signals, bvalues, **kwargs): """Perform the IVIM fit diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index 3a405141..49b055b0 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -6,6 +6,8 @@ import os import time +from exceptiongroup import catch + from src.wrappers.OsipiBase import OsipiBase from utilities.data_simulation.GenerateData import GenerateData #run using python -m pytest from the root folder @@ -70,7 +72,7 @@ def test_ivim_fit_saved(name, bvals, data, algorithm, xfail, kwargs, tolerances, signal = signal_helper(data["data"]) tolerances = tolerances_helper(tolerances, data) start_time = time.time() # Record the start time - fit = OsipiBase(algorithm=algorithm, **kwargs) + fit = OsipiBase(algorithm=algorithm, **kwargs) fit_result = fit.osipi_fit(signal, bvals) elapsed_time = time.time() - start_time # Calculate elapsed time def to_list_if_needed(value): @@ -113,7 +115,7 @@ def test_default_bounds_and_initial_guesses(algorithm): fit = OsipiBase(algorithm=algorithm) #assert fit.bounds is not None, f"For {algorithm}, there is no default fit boundary" #assert fit.initial_guess is not None, f"For {algorithm}, there is no default fit initial guess" - if fit.bounds is not None: + if fit.use_bounds: assert 0 <= fit.bounds[0][0] <= 0.003, f"For {algorithm}, the default lower bound of D {fit.bounds[0][0]} is unrealistic" assert 0 <= fit.bounds[1][0] <= 0.01, f"For {algorithm}, the default upper bound of D {fit.bounds[1][0]} is unrealistic" assert 0 <= fit.bounds[0][1] <= 1, f"For {algorithm}, the default lower bound of f {fit.bounds[0][1]} is unrealistic" @@ -123,7 +125,7 @@ def test_default_bounds_and_initial_guesses(algorithm): assert 0 <= fit.bounds[0][3] <= 1, f"For {algorithm}, the default lower bound of S {fit.bounds[0][3]} is unrealistic; note data is normaized" assert 1 <= fit.bounds[1][3] <= 1000, f"For {algorithm}, the default upper bound of S {fit.bounds[1][3]} is unrealistic; note data is normaized" assert fit.bounds[1][0] <= fit.bounds[0][2], f"For {algorithm}, the default upper bound of D {fit.bounds[1][0]} is higher than lower bound of D* {fit.bounds[0][2]}" - if fit.initial_guess is not None: + if fit.use_initial_guess: assert 0.0008 <= fit.initial_guess[0] <= 0.002, f"For {algorithm}, the default initial guess for D {fit.initial_guess[0]} is unrealistic" assert 0 <= fit.initial_guess[1] <= 0.5, f"For {algorithm}, the default initial guess for f {fit.initial_guess[1]} is unrealistic" assert 0.003 <= fit.initial_guess[2] <= 0.1, f"For {algorithm}, the default initial guess for Ds {fit.initial_guess[2]} is unrealistic" @@ -169,5 +171,37 @@ def test_bounds(name, bvals, data, algorithm, xfail, kwargs, tolerances, request assert bounds[0][2] <= fit_result['Dp'] <= bounds[1][2], f"Result {fit_result['Dp']} out of bounds for data: {name}" # S0 is not returned as argument... #assert bounds[0][3] <= fit_result['S0'] <= bounds[1][3], f"Result {fit_result['S0']} out of bounds for data: {name}" + '''if fit.use_initial_guess: + passD=False + passf=False + passDp=False + fit = OsipiBase(algorithm=algorithm, bounds=bounds, initial_guess=[0.0007, 0.25, 0.015, 1.2], **kwargs) + try: + fit_result = fit.osipi_fit(signal, bvals) + if fit_result['D'] == 0: + passD=True + except: + passD=True + fit = OsipiBase(algorithm=algorithm, bounds=bounds, initial_guess=[0.001, 0.19, 0.015, 1.2], **kwargs) + try: + fit_result = fit.osipi_fit(signal, bvals) + if fit_result['f'] == 0: + passf=True + except: + passf = True + try: + fit = OsipiBase(algorithm=algorithm, bounds=bounds, initial_guess=[0.001, 0.25, 0.009, 1.2], **kwargs) + if fit_result['Dp'] == 0: + passDp = True + except: + passDp = True + assert passD, f"Fit still passes when initial guess D is out of fit bounds; potentially initial guesses not respected for: {name}" + assert passf, f"Fit still passes when initial guess f is out of fit bounds; potentially initial guesses not respected for: {name}" + assert passDp, f"Fit still passes when initial guess Ds is out of fit bounds; potentially initial guesses not respected for: {name}" ''' + + + + + \ No newline at end of file From 01d51aa22730e8558d30a489bb5ce7da2b3e354a Mon Sep 17 00:00:00 2001 From: Oliver Gurney-Champion Date: Mon, 17 Mar 2025 09:20:35 +0100 Subject: [PATCH 10/10] error python 3.11 fixed(?) --- tests/IVIMmodels/unit_tests/test_ivim_fit.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/tests/IVIMmodels/unit_tests/test_ivim_fit.py b/tests/IVIMmodels/unit_tests/test_ivim_fit.py index 49b055b0..fb58d13b 100644 --- a/tests/IVIMmodels/unit_tests/test_ivim_fit.py +++ b/tests/IVIMmodels/unit_tests/test_ivim_fit.py @@ -3,13 +3,8 @@ import pytest import json import pathlib -import os import time - -from exceptiongroup import catch - from src.wrappers.OsipiBase import OsipiBase -from utilities.data_simulation.GenerateData import GenerateData #run using python -m pytest from the root folder def signal_helper(signal):