From 8288587faf26e659ccd14fb6ef4f58425f6a6ba4 Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Mon, 4 Aug 2025 18:21:59 +0200 Subject: [PATCH 1/8] Actually allow SMR files to be processed --- phys2bids/phys2bids.py | 6 ++++++ phys2bids/utils.py | 2 +- 2 files changed, 7 insertions(+), 1 deletion(-) diff --git a/phys2bids/phys2bids.py b/phys2bids/phys2bids.py index 5161792bf..00d581027 100644 --- a/phys2bids/phys2bids.py +++ b/phys2bids/phys2bids.py @@ -260,6 +260,12 @@ def phys2bids( from phys2bids.io import load_gep phys_in = load_gep(infile) + elif ftype == "smr": + from phys2bids.io import load_smr + + phys_in = load_smr(infile, chtrig) + else: + raise RuntimeError(f'Unsupported extension "{ftype}"') LGR.info("Checking that units of measure are BIDS compatible") for index, unit in enumerate(phys_in.units): diff --git a/phys2bids/utils.py b/phys2bids/utils.py index 05949afa9..33a8ee8b7 100644 --- a/phys2bids/utils.py +++ b/phys2bids/utils.py @@ -9,7 +9,7 @@ LGR = logging.getLogger(__name__) -SUPPORTED_FTYPES = ("acq", "txt", "mat", "gep") +SUPPORTED_FTYPES = ("acq", "txt", "mat", "gep", "smr") def check_input_ext(filename, ext): From a37440392e9afb6e2f5d49f32b1b8f0568f8cace Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Mon, 4 Aug 2025 18:22:42 +0200 Subject: [PATCH 2/8] Update workflows --- .github/workflows/auto-release.yml | 4 ++-- .github/workflows/python-publish.yml | 9 +++++---- 2 files changed, 7 insertions(+), 6 deletions(-) diff --git a/.github/workflows/auto-release.yml b/.github/workflows/auto-release.yml index 4a36dfb68..6c9f5c158 100644 --- a/.github/workflows/auto-release.yml +++ b/.github/workflows/auto-release.yml @@ -10,7 +10,7 @@ on: jobs: auto-release: - runs-on: ubuntu-20.04 + runs-on: ubuntu-22.04 # Set skip ci to avoid loops if: "!contains(github.event.head_commit.message, 'ci skip') && !contains(github.event.head_commit.message, 'skip ci')" # Set bash as default shell for jobs @@ -19,7 +19,7 @@ jobs: shell: bash steps: - name: Checkout source - uses: actions/checkout@v3 + uses: actions/checkout@v4 with: # Fetch all history for all branches and tags fetch-depth: 0 diff --git a/.github/workflows/python-publish.yml b/.github/workflows/python-publish.yml index bd3d0d25d..03e58b58f 100644 --- a/.github/workflows/python-publish.yml +++ b/.github/workflows/python-publish.yml @@ -10,15 +10,15 @@ on: jobs: deploy: - runs-on: ubuntu-20.04 + runs-on: ubuntu-22.04 steps: - name: Checkout source - uses: actions/checkout@v3 + uses: actions/checkout@v4 - name: Set up Python - uses: actions/setup-python@v4 + uses: actions/setup-python@v5 with: - python-version: '3.10' + python-version: 3.11 - name: Install dependencies run: | python -m pip install --upgrade pip @@ -29,4 +29,5 @@ jobs: TWINE_PASSWORD: ${{ secrets.PYPI_PASSWORD }} run: | python setup.py sdist bdist_wheel + twine check dist/* twine upload dist/* From bf0c136b3fcfa1f9f2cf4b7d372574db3ba52847 Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Mon, 4 Aug 2025 18:23:10 +0200 Subject: [PATCH 3/8] Fix some string issues --- phys2bids/physio_obj.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/phys2bids/physio_obj.py b/phys2bids/physio_obj.py index 131cb7ce9..0eb88de74 100644 --- a/phys2bids/physio_obj.py +++ b/phys2bids/physio_obj.py @@ -465,7 +465,7 @@ def delete_at_index(self, idx): del self.units[idx] if self.trigger_idx == idx: - LGR.warning("Removing trigger channel - are you sure you are doing" "the right thing?") + LGR.warning("Removing trigger channel - are you sure you are doing the right thing?") self.trigger_idx = 0 def check_trigger_amount(self, thr=None, num_timepoints_expected=0, tr=0): @@ -527,7 +527,7 @@ def check_trigger_amount(self, thr=None, num_timepoints_expected=0, tr=0): num_timepoints_found = np.count_nonzero(np.ediff1d(timepoints.astype(np.int8)) > 0) if flag == 1: LGR.info( - f"The number of timepoints according to the std_thr method " + f"The number of timepoints according to the average threshold method " f"is {num_timepoints_found}. The computed threshold is {thr:.4f}" ) else: @@ -550,7 +550,8 @@ def check_trigger_amount(self, thr=None, num_timepoints_expected=0, tr=0): elif num_timepoints_found < num_timepoints_expected: timepoints_missing = num_timepoints_expected - num_timepoints_found - LGR.warning(f"Found {timepoints_missing} timepoints" " less than expected!") + LGR.warning(f"Found {timepoints_missing} timepoints less than expected!") + # if tr checks for any truthy value, including not 0 nor empty things. if tr: LGR.warning( "Correcting time offset, assuming missing " @@ -559,14 +560,14 @@ def check_trigger_amount(self, thr=None, num_timepoints_expected=0, tr=0): ) time_offset -= timepoints_missing * tr else: - LGR.warning("Can't correct time offset - you should " "specify the TR") + LGR.warning("Can't correct time offset - you should specify the TR") else: LGR.info("Found just the right amount of timepoints!") else: LGR.warning( - "The necessary options to find the amount of timepoints " "were not provided." + "The necessary options to find the amount of timepoints were not provided." ) self.thr = thr self.time_offset = time_offset From a2ba8e9fcf4ca6f5ea7de4d9424be79d7edb9863 Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Mon, 4 Aug 2025 18:23:29 +0200 Subject: [PATCH 4/8] Update setup to newer physiopy --- setup.cfg | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) diff --git a/setup.cfg b/setup.cfg index f75ceb262..771d7d1b7 100644 --- a/setup.cfg +++ b/setup.cfg @@ -2,18 +2,14 @@ name = phys2bids url = https://github.com/physiopy/phys2bids download_url = https://github.com/physiopy/phys2bids -author = phys2bids developers -maintainer = Stefano Moia +author = The Physiopy Community +maintainer = The Physiopy Community maintainer_email = physiopy.community@gmail.com classifiers = Development Status :: 5 - Production/Stable Intended Audience :: Science/Research License :: OSI Approved :: Apache Software License - Programming Language :: Python :: 3.6 - Programming Language :: Python :: 3.7 - Programming Language :: Python :: 3.8 - Programming Language :: Python :: 3.9 - Programming Language :: Python :: 3.10 + Programming Language :: Python :: 3 license = Apache-2.0 description = Python library to convert physiological data files into BIDS format long_description = file:README.md @@ -64,7 +60,7 @@ test = coverage %(interfaces)s %(style)s -all = +dev = %(doc)s %(duecredit)s %(interfaces)s From 412f642f15a5a6ac7df04b90704fe188d9b91c63 Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Tue, 5 Aug 2025 18:37:50 +0200 Subject: [PATCH 5/8] Add a brand new cluster estimation to find groups of triggers, i.e. takes or runs. --- phys2bids/cli/run.py | 25 +++++++++++ phys2bids/phys2bids.py | 9 +++- phys2bids/slice4phys.py | 91 +++++++++++++++++++++++++++++++++++++++++ 3 files changed, 124 insertions(+), 1 deletion(-) diff --git a/phys2bids/cli/run.py b/phys2bids/cli/run.py index a7770171e..574bb91a0 100644 --- a/phys2bids/cli/run.py +++ b/phys2bids/cli/run.py @@ -134,6 +134,31 @@ def _get_parser(): "or just one TR if it is consistent throughout the session.", default=None, ) + optional.add_argument( + "-esttakes", + "--estimate_takes", + dest="estimate_takes", + action="store_true", + help="Run automatic algorithm to estimate clusters of triggers, i.e. the " + "'takes' or 'runs' of fMRI. Useful when sequences were stopped and restarted, " + "or when you don't know how many triggers or trs you have in each take. " + "This might work 95% of the time. Default is False.", + default=False, + ) + optional.add_argument( + "-ci", + "--confidence-interval", + dest="ci", + # Here always as float, later it will check if the float is an integer instead. + type=float, + help="The Confidence Interval (CI) to use in the estimation of the trigger clusters. " + "The cluster algorithm considers triggers with duration (in samples) within this " + "CI as part of the same group. If CI is an integer, it will consider that amount " + "of triggers. If CI is a float and < 1, it will consider that percentage of the " + "trigger duration. CI cannot be a float > 1. Default is 1. Change to .25 if " + "there is a CMRR DWI sequence or if you are recording sub-triggers.", + default=1, + ) optional.add_argument( "-thr", "--threshold", diff --git a/phys2bids/phys2bids.py b/phys2bids/phys2bids.py index 00d581027..09615b33e 100644 --- a/phys2bids/phys2bids.py +++ b/phys2bids/phys2bids.py @@ -38,7 +38,7 @@ from phys2bids import _version, bids, utils, viz from phys2bids.cli.run import _get_parser from phys2bids.physio_obj import BlueprintOutput -from phys2bids.slice4phys import slice4phys +from phys2bids.slice4phys import estimate_ntp_and_tr, slice4phys from . import __version__ from .due import Doi, due @@ -141,6 +141,8 @@ def phys2bids( chsel=None, num_timepoints_expected=None, tr=None, + estimate_takes=False, + ci=1, thr=None, pad=9, ch_name=[], @@ -304,6 +306,11 @@ def phys2bids( LGR.info("Renaming channels with given names") phys_in.rename_channels(ch_name) + # If requested, run the automatic detection of timepoints and groups + + if estimate_takes: + num_timepoints_expected, tr = estimate_ntp_and_tr(phys_in, thr=None, ci=1) + # Checking acquisition type via user's input if tr is not None and num_timepoints_expected is not None: # Multi-run acquisition type section diff --git a/phys2bids/slice4phys.py b/phys2bids/slice4phys.py index e4e491a57..1fad45999 100644 --- a/phys2bids/slice4phys.py +++ b/phys2bids/slice4phys.py @@ -8,6 +8,97 @@ LGR = logging.getLogger(__name__) +def estimate_ntp_and_tr(phys_in, thr=None, ci=1): + """ + Find groups of trigger in a spiky signal like the trigger channel signal. + """ + LGR.info('Running automatic clustering of triggers to find timepoints and tr of each "take"') + trigger = phys_in.timeseries[phys_in.trigger_idx] + + thr = np.mean(trigger) if thr is None else thr + timepoints = trigger > thr + spikes = np.flatnonzero(np.ediff1d(timepoints.astype(np.int8)) > 0) + interspike_interval = np.diff(spikes) + unique_isi, counts = np.unique(interspike_interval, return_counts=True) + + # The following line is for python < 3.12. From 3.12, ci.is_integer() is enough. + if isinstance(ci, int) or isinstance(ci, float) and ci.is_integer(): + upper_ci_isi = unique_isi + ci + elif isinstance(ci, float) and ci < 1: + upper_ci_isi = unique_isi * (1 + ci) + elif isinstance(ci, float) and ci > 1: + raise ValueError("Confidence intervals above 1 are not supported.") + else: + raise ValueError("Confidence intervals must be either integers or floats.") + + # Loop through the uniques ISI and group them within the specified CI. + # Also compute the average TR of the group. + isi_groups = {} + average_tr = {} + k = 0 + current_group = [unique_isi[0]] + + for n, i in enumerate(range(1, len(unique_isi))): + if unique_isi[i] <= upper_ci_isi[n]: + current_group.append(unique_isi[i]) + else: + isi_groups[k] = current_group + average_tr[k] = np.mean(current_group) / phys_in.freq[0] + k += 1 + current_group = [unique_isi[i]] + + isi_groups[k] = current_group + average_tr[k] = np.mean(current_group) / phys_in.freq[0] + + # Invert the isi_group into value per group + group_by_isi = {isi: group for group, isis in isi_groups.items() for isi in isis} + + # Use the found groups to find the number of timepoints and assign the right TR + estimated_ntp = [] + estimated_tr = [] + + i = 0 + while i < interspike_interval.size - 1: + current_group = group_by_isi.get(interspike_interval[i]) + for n in range(i + 1, interspike_interval.size): + if current_group != group_by_isi.get(interspike_interval[n]): + break + # Repeat one last time outside of for loop + estimated_ntp += [n - i] + estimated_tr += [average_tr[current_group]] + i = n + + if len(estimated_ntp) < 1: + raise Exception("This should not happen. Something went very wrong.") + # The algorithm found n groups, the last of which has two timepoints less due to + # diff computations. Each real group of n>1 triggers counts one trigger less but is + # followed by a "fake" group of 1 trigger that is actually the interval to the next + # group. That does not hold if there is a real group of 1 trigger. + # Loop through the estiamtions to fix all that. + ntp = [] + tr = [] + i = 0 + + while i < len(estimated_ntp): + if estimated_ntp[i] == 1: + ntp.append(estimated_ntp[i]) + tr.append(estimated_tr[i]) + i += 1 + elif i + 1 < len(estimated_ntp): + ntp.append(estimated_ntp[i] + estimated_ntp[i + 1]) + tr.append(estimated_tr[i]) + i += 2 + else: + ntp.append(estimated_ntp[i] + 2) + tr.append(estimated_tr[i]) + i += 1 + + LGR.info( + f"The automatic clustering found {len(ntp)} groups of triggers long: {ntp} with respective TR: {tr}" + ) + return ntp, tr + + def find_takes(phys_in, ntp_list, tr_list, thr=None, padding=9): """ Find takes slicing index. From ee1822477ba1546fb0b2919cf522489a8f82ed2b Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Tue, 5 Aug 2025 18:53:40 +0200 Subject: [PATCH 6/8] numpy docstrings. --- phys2bids/cli/run.py | 8 ++++---- phys2bids/slice4phys.py | 31 ++++++++++++++++++++++++++++++- 2 files changed, 34 insertions(+), 5 deletions(-) diff --git a/phys2bids/cli/run.py b/phys2bids/cli/run.py index 574bb91a0..4045dec80 100644 --- a/phys2bids/cli/run.py +++ b/phys2bids/cli/run.py @@ -153,10 +153,10 @@ def _get_parser(): type=float, help="The Confidence Interval (CI) to use in the estimation of the trigger clusters. " "The cluster algorithm considers triggers with duration (in samples) within this " - "CI as part of the same group. If CI is an integer, it will consider that amount " - "of triggers. If CI is a float and < 1, it will consider that percentage of the " - "trigger duration. CI cannot be a float > 1. Default is 1. Change to .25 if " - "there is a CMRR DWI sequence or if you are recording sub-triggers.", + "CI as part of the same group, thus the same. If CI is an integer, it will consider " + "that amount of triggers. If CI is a float and < 1, it will consider that " + "percentage of the trigger duration. CI cannot be a float > 1. Default is 1. " + "Change to .25 if there is a CMRR DWI sequence or when recording sub-triggers.", default=1, ) optional.add_argument( diff --git a/phys2bids/slice4phys.py b/phys2bids/slice4phys.py index 1fad45999..1619e5ba5 100644 --- a/phys2bids/slice4phys.py +++ b/phys2bids/slice4phys.py @@ -11,6 +11,35 @@ def estimate_ntp_and_tr(phys_in, thr=None, ci=1): """ Find groups of trigger in a spiky signal like the trigger channel signal. + + Parameters + ---------- + phys_in : BlueprintInput object + A BlueprintInput object containing a physiological acquisition + thr : None, optional + The threshold for automatic spike detection. Default is to use the average of + the signal. + ci : int or float, optional + Confidence Interval (CI) to use in the estimation of the trigger clusters. The + cluster algorithm considers triggers with duration (in samples) within this CI + as part of the same group, thus the same. If CI is an integer, it will consider + that amount of triggers. If CI is a float and < 1, it will consider that + percentage of the trigger duration. CI cannot be a float > 1. Default is 1. + Change to .25 if there is a CMRR DWI sequence or when recording sub-triggers. + + Returns + ------- + ntp + The list of number of timepoints found for each take. + tr + The list of corresponding TR, computed as average samples per group / frequency. + + Raises + ------ + Exception + If it doesn't find at least a group. + ValueError + If CI is a float above 1 or if it is not an int or a float. """ LGR.info('Running automatic clustering of triggers to find timepoints and tr of each "take"') trigger = phys_in.timeseries[phys_in.trigger_idx] @@ -27,7 +56,7 @@ def estimate_ntp_and_tr(phys_in, thr=None, ci=1): elif isinstance(ci, float) and ci < 1: upper_ci_isi = unique_isi * (1 + ci) elif isinstance(ci, float) and ci > 1: - raise ValueError("Confidence intervals above 1 are not supported.") + raise ValueError("Confidence intervals percentages above 1 are not supported.") else: raise ValueError("Confidence intervals must be either integers or floats.") From 6324c1cba569158d142da4e8cd860e506442402b Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Mon, 11 Aug 2025 14:21:41 +0200 Subject: [PATCH 7/8] Fix issue with help --- phys2bids/cli/run.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/phys2bids/cli/run.py b/phys2bids/cli/run.py index 4045dec80..4b34f5b4a 100644 --- a/phys2bids/cli/run.py +++ b/phys2bids/cli/run.py @@ -142,7 +142,7 @@ def _get_parser(): help="Run automatic algorithm to estimate clusters of triggers, i.e. the " "'takes' or 'runs' of fMRI. Useful when sequences were stopped and restarted, " "or when you don't know how many triggers or trs you have in each take. " - "This might work 95% of the time. Default is False.", + "This might work 95%% of the time. Default is False.", default=False, ) optional.add_argument( From 3e5bea319e45b5d854fa54b6379c58fa7f34f8ba Mon Sep 17 00:00:00 2001 From: Stefano Moia Date: Mon, 11 Aug 2025 15:03:03 +0200 Subject: [PATCH 8/8] Fix more explanation of CI and also add a tiny comment to stop puzzling about why it works. --- phys2bids/cli/run.py | 7 ++++--- phys2bids/slice4phys.py | 3 ++- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/phys2bids/cli/run.py b/phys2bids/cli/run.py index 4b34f5b4a..aef761ff1 100644 --- a/phys2bids/cli/run.py +++ b/phys2bids/cli/run.py @@ -154,9 +154,10 @@ def _get_parser(): help="The Confidence Interval (CI) to use in the estimation of the trigger clusters. " "The cluster algorithm considers triggers with duration (in samples) within this " "CI as part of the same group, thus the same. If CI is an integer, it will consider " - "that amount of triggers. If CI is a float and < 1, it will consider that " - "percentage of the trigger duration. CI cannot be a float > 1. Default is 1. " - "Change to .25 if there is a CMRR DWI sequence or when recording sub-triggers.", + "that amount of samples around the distance. If CI is a float and < 1, it will " + "consider that percentage of the trigger duration. CI cannot be a float > 1. " + "Default is 1. Change to .25 if there is a CMRR DWI sequence or when recording " + "sub-triggers.", default=1, ) optional.add_argument( diff --git a/phys2bids/slice4phys.py b/phys2bids/slice4phys.py index 1619e5ba5..bec7171fc 100644 --- a/phys2bids/slice4phys.py +++ b/phys2bids/slice4phys.py @@ -23,7 +23,7 @@ def estimate_ntp_and_tr(phys_in, thr=None, ci=1): Confidence Interval (CI) to use in the estimation of the trigger clusters. The cluster algorithm considers triggers with duration (in samples) within this CI as part of the same group, thus the same. If CI is an integer, it will consider - that amount of triggers. If CI is a float and < 1, it will consider that + that amount of samples. If CI is a float and < 1, it will consider that percentage of the trigger duration. CI cannot be a float > 1. Default is 1. Change to .25 if there is a CMRR DWI sequence or when recording sub-triggers. @@ -67,6 +67,7 @@ def estimate_ntp_and_tr(phys_in, thr=None, ci=1): k = 0 current_group = [unique_isi[0]] + # np.unique returns sorted elements → unique_isi[0] == min(unique_isi), so THIS WORKS. for n, i in enumerate(range(1, len(unique_isi))): if unique_isi[i] <= upper_ci_isi[n]: current_group.append(unique_isi[i])