From 90b861d7a4683e65a50f22e11c8071daec3322d3 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Mon, 4 Feb 2019 17:00:12 -0500 Subject: [PATCH 01/13] ENH: Calculate TA and store in run_info --- bids/analysis/analysis.py | 13 +++++++++++-- bids/variables/entities.py | 8 +++++--- bids/variables/io.py | 26 +++++++++++++++++++++++--- 3 files changed, 39 insertions(+), 8 deletions(-) diff --git a/bids/analysis/analysis.py b/bids/analysis/analysis.py index fd3c28660..814d45f1b 100644 --- a/bids/analysis/analysis.py +++ b/bids/analysis/analysis.py @@ -388,13 +388,22 @@ def get_design_matrix(self, names=None, format='long', mode='both', if sampling_rate == 'TR': trs = {var.run_info[0].tr for var in self.collection.variables.values()} + tas = {var.run_info[0].ta for var in self.collection.variables.values()} if not trs: raise ValueError("Repetition time unavailable; specify sampling_rate " "explicitly") elif len(trs) > 1: raise ValueError("Non-unique Repetition times found ({!r}); specify " - "sampling_rate explicitly") - sampling_rate = 1. / trs.pop() + "sampling_rate explicitly".format(trs)) + TR = trs.pop() + if not tas: + warnings.warn("Acquisition time unavailable; assuming TA = TR") + tas = {TR} + elif len(tas) > 1: + raise ValueError("Non-unique acquisition times found ({!r})".format(tas)) + + sampling_rate = 1. / TR + acquisition_time = tas.pop() elif sampling_rate == 'highest': sampling_rate = None dense_df = coll.to_df(names, format='wide', diff --git a/bids/variables/entities.py b/bids/variables/entities.py index a1060f6f3..d46da3b36 100644 --- a/bids/variables/entities.py +++ b/bids/variables/entities.py @@ -36,23 +36,25 @@ class RunNode(Node): image_file (str): The full path to the corresponding nifti image. duration (float): Duration of the run, in seconds. repetition_time (float): TR for the run. + acquisition_time (float): TA for the run. task (str): The task name for this run. ''' - def __init__(self, entities, image_file, duration, repetition_time): + def __init__(self, entities, image_file, duration, repetition_time, acquisition_time): self.image_file = image_file self.duration = duration self.repetition_time = repetition_time + self.acquisition_time = acquisition_time super(RunNode, self).__init__('run', entities) def get_info(self): return RunInfo(self.entities, self.duration, self.repetition_time, - self.image_file) + self.acquisition_time, self.image_file) # Stores key information for each Run. -RunInfo = namedtuple('RunInfo', ['entities', 'duration', 'tr', 'image']) +RunInfo = namedtuple('RunInfo', ['entities', 'duration', 'tr', 'ta', 'image']) class NodeIndex(Node): diff --git a/bids/variables/io.py b/bids/variables/io.py index d35c32094..64ef16916 100644 --- a/bids/variables/io.py +++ b/bids/variables/io.py @@ -82,6 +82,24 @@ def load_variables(layout, types=None, levels=None, skip_empty=True, return dataset +def _get_timing_info(img_md): + if 'RepetitionTime' in img_md: + tr = img_md['RepetitionTime'] + if 'DelayTime' in img_md: + ta = tr - img_md['DelayTime'] + elif 'SliceTiming' in img_md: + slicetimes = sorted(img_md['SliceTiming']) + # a, b ... z + # z = final slice onset, b - a = slice duration + ta = slicetimes[-1] + slicetimes[1] - slicetimes[0] + else: + ta = tr + elif 'VolumeTiming' in img_md: + return NotImplemented + + return tr, ta + + def _load_time_variables(layout, dataset=None, columns=None, scan_length=None, drop_na=True, events=True, physio=True, stim=True, regressors=True, skip_empty=True, **selectors): @@ -141,8 +159,9 @@ def _load_time_variables(layout, dataset=None, columns=None, scan_length=None, if 'run' in entities: entities['run'] = int(entities['run']) - tr = layout.get_metadata(img_f, suffix='bold', domains=domains, - full_search=True)['RepetitionTime'] + img_md = layout.get_metadata(img_f, suffix='bold', domains=domains, + full_search=True) + tr, ta = _get_timing_info(img_md) # Get duration of run: first try to get it directly from the image # header; if that fails, try to get NumberOfVolumes from the @@ -162,7 +181,8 @@ def _load_time_variables(layout, dataset=None, columns=None, scan_length=None, raise ValueError(msg) run = dataset.get_or_create_node('run', entities, image_file=img_f, - duration=duration, repetition_time=tr) + duration=duration, repetition_time=tr, + acquisition_time=ta) run_info = run.get_info() # Process event files From 9057a7cac3ae833fbb0a86967ca49b0c827f0c55 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Tue, 5 Feb 2019 09:58:48 -0500 Subject: [PATCH 02/13] RF: Set default TA = TR, update RunInfo calls in tests --- bids/analysis/tests/test_transformations.py | 2 +- bids/variables/entities.py | 4 ++-- bids/variables/tests/test_variables.py | 2 +- 3 files changed, 4 insertions(+), 4 deletions(-) diff --git a/bids/analysis/tests/test_transformations.py b/bids/analysis/tests/test_transformations.py index 9db5625f1..83e2c779f 100644 --- a/bids/analysis/tests/test_transformations.py +++ b/bids/analysis/tests/test_transformations.py @@ -28,7 +28,7 @@ def sparse_run_variable_with_missing_values(): 'duration': [1.2, 1.6, 0.8, 2], 'amplitude': [1, 1, np.nan, 1] }) - run_info = [RunInfo({'subject': '01'}, 20, 2, 'dummy.nii.gz')] + run_info = [RunInfo({'subject': '01'}, 20, 2, 2, 'dummy.nii.gz')] var = SparseRunVariable('var', data, run_info, 'events') return BIDSRunVariableCollection([var]) diff --git a/bids/variables/entities.py b/bids/variables/entities.py index d46da3b36..0ba13d633 100644 --- a/bids/variables/entities.py +++ b/bids/variables/entities.py @@ -40,11 +40,11 @@ class RunNode(Node): task (str): The task name for this run. ''' - def __init__(self, entities, image_file, duration, repetition_time, acquisition_time): + def __init__(self, entities, image_file, duration, repetition_time, acquisition_time=None): self.image_file = image_file self.duration = duration self.repetition_time = repetition_time - self.acquisition_time = acquisition_time + self.acquisition_time = acquisition_time or repetition_time super(RunNode, self).__init__('run', entities) def get_info(self): diff --git a/bids/variables/tests/test_variables.py b/bids/variables/tests/test_variables.py index 74876489c..8cdc4db4c 100644 --- a/bids/variables/tests/test_variables.py +++ b/bids/variables/tests/test_variables.py @@ -19,7 +19,7 @@ def generate_DEV(name='test', sr=20, duration=480): ent_names = ['task', 'run', 'session', 'subject'] entities = {e: uuid.uuid4().hex for e in ent_names} image = uuid.uuid4().hex + '.nii.gz' - run_info = RunInfo(entities, duration, 2, image) + run_info = RunInfo(entities, duration, 2, 2, image) return DenseRunVariable('test', values, run_info, 'dummy', sr) From 4921a0a5e39c3842e6b86ca37fbb5abec340d053 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Thu, 7 Feb 2019 12:39:17 -0500 Subject: [PATCH 03/13] ENH: Upsample convolution, begin calculating downsample --- bids/analysis/analysis.py | 4 +++- bids/analysis/transformations/compute.py | 9 +++++++++ bids/variables/variables.py | 18 ++++++++++++++---- 3 files changed, 26 insertions(+), 5 deletions(-) diff --git a/bids/analysis/analysis.py b/bids/analysis/analysis.py index 814d45f1b..21e226bba 100644 --- a/bids/analysis/analysis.py +++ b/bids/analysis/analysis.py @@ -408,7 +408,9 @@ def get_design_matrix(self, names=None, format='long', mode='both', sampling_rate = None dense_df = coll.to_df(names, format='wide', include_sparse=include_sparse, - sampling_rate=sampling_rate, **kwargs) + sampling_rate=sampling_rate, + integration_window=acquisition_time, + **kwargs) if dense_df is not None: dense_df = dense_df.drop(['onset', 'duration'], axis=1) diff --git a/bids/analysis/transformations/compute.py b/bids/analysis/transformations/compute.py index ae51189a9..a8379cc96 100644 --- a/bids/analysis/transformations/compute.py +++ b/bids/analysis/transformations/compute.py @@ -2,6 +2,7 @@ Transformations that primarily involve numerical computation on variables. ''' +from math import gcd import numpy as np import pandas as pd from bids.utils import listify @@ -35,6 +36,14 @@ def _transform(self, var, model='spm', derivative=False, dispersion=False, if isinstance(var, SparseRunVariable): sr = self.collection.sampling_rate + # Resolve diferences in TR and TA to milliseconds + TR = int(np.round(1000. * var.run_info.tr)) + TA = int(np.round(1000. * var.run_info.ta)) + if TA is None or TA < TR: + # Use a unit that fits an whole number of times into both + # the interscan interval (TR) and the integration window (TA) + dt = gcd(TR, TA) + sr = 1000. / dt var = var.to_dense(sr) df = var.to_df(entities=False) diff --git a/bids/variables/variables.py b/bids/variables/variables.py index 87cd5301b..59f04d6d6 100644 --- a/bids/variables/variables.py +++ b/bids/variables/variables.py @@ -419,7 +419,7 @@ def _build_entity_index(self, run_info, sampling_rate): self.timestamps = pd.concat(_timestamps, axis=0, sort=True) return pd.concat(index, axis=0, sort=True).reset_index(drop=True) - def resample(self, sampling_rate, inplace=False, kind='linear'): + def resample(self, sampling_rate, integration_window=None, inplace=False, kind='linear'): '''Resample the Variable to the specified sampling rate. Parameters @@ -436,7 +436,7 @@ def resample(self, sampling_rate, inplace=False, kind='linear'): ''' if not inplace: var = self.clone() - var.resample(sampling_rate, True, kind) + var.resample(sampling_rate, integration_window, True, kind) return var if sampling_rate == self.sampling_rate: @@ -450,6 +450,14 @@ def resample(self, sampling_rate, inplace=False, kind='linear'): x = np.arange(n) num = len(self.index) + if integration_window is not None: + # Number of timesteps in existing time series + n_frames = int(np.floor(integration_window * old_sr)) + integrator = np.zeroes((n, num), dtype=np.uint8) + + # old_onsets = old_sr * self.index. + # onset_frames = sampling_rate + from scipy.interpolate import interp1d f = interp1d(x, self.values.values.ravel(), kind=kind) x_new = np.linspace(0, n - 1, num=num) @@ -458,7 +466,8 @@ def resample(self, sampling_rate, inplace=False, kind='linear'): self.sampling_rate = sampling_rate - def to_df(self, condition=True, entities=True, timing=True, sampling_rate=None): + def to_df(self, condition=True, entities=True, timing=True, sampling_rate=None, + integration_window=None): '''Convert to a DataFrame, with columns for name and entities. Parameters @@ -474,7 +483,8 @@ def to_df(self, condition=True, entities=True, timing=True, sampling_rate=None): sampled uniformly). If False, omits them. ''' if sampling_rate not in (None, self.sampling_rate): - return self.resample(sampling_rate).to_df(condition, entities) + return self.resample(sampling_rate, + integration_window=integration_window).to_df(condition, entities) df = super(DenseRunVariable, self).to_df(condition, entities) From ff319d1b416ac7568ca5c0ba83eaffb18a89211c Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Fri, 8 Feb 2019 14:09:16 -0500 Subject: [PATCH 04/13] FIX: Pass new variable args as kwargs in split() --- bids/variables/variables.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/bids/variables/variables.py b/bids/variables/variables.py index 59f04d6d6..c8b0856d3 100644 --- a/bids/variables/variables.py +++ b/bids/variables/variables.py @@ -253,8 +253,8 @@ def split(self, grouper): subsets = [] for i, (name, g) in enumerate(data.groupby(grouper)): name = '%s.%s' % (self.name, name) - args = [name, g, getattr(self, 'run_info', None), self.source] - col = self.__class__(*args) + col = self.__class__(name=name, data=g, source=self.source, + run_info=getattr(self, 'run_info', None)) subsets.append(col) return subsets From 03cbeb052a44edc6d53f960902619feb7ea4b85d Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Fri, 8 Feb 2019 14:52:33 -0500 Subject: [PATCH 05/13] ENH: Choose sampling rate to accommodate TR, TA and target SR --- bids/analysis/transformations/compute.py | 15 +++++++++++++-- 1 file changed, 13 insertions(+), 2 deletions(-) diff --git a/bids/analysis/transformations/compute.py b/bids/analysis/transformations/compute.py index a8379cc96..741169dfa 100644 --- a/bids/analysis/transformations/compute.py +++ b/bids/analysis/transformations/compute.py @@ -37,13 +37,24 @@ def _transform(self, var, model='spm', derivative=False, dispersion=False, if isinstance(var, SparseRunVariable): sr = self.collection.sampling_rate # Resolve diferences in TR and TA to milliseconds - TR = int(np.round(1000. * var.run_info.tr)) - TA = int(np.round(1000. * var.run_info.ta)) + trs = {ri.tr for ri in var.run_info} + tas = {ri.ta for ri in var.run_info} + assert len(trs) == 1 + assert len(tas) == 1 + TR = int(np.round(1000. * trs.pop())) + TA = int(np.round(1000. * tas.pop())) + SR = int(np.round(1000. / sr)) if TA is None or TA < TR: # Use a unit that fits an whole number of times into both # the interscan interval (TR) and the integration window (TA) dt = gcd(TR, TA) + if dt > SR: + # Find the nearest whole factor of dt >= SR + # This compromises on the target sampling rate and + # one that neatly bins TR and TA + dt = dt // (dt // SR) sr = 1000. / dt + var = var.to_dense(sr) df = var.to_df(entities=False) From 42384fc3c83d656c8190d2e5f64733718e773807 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Fri, 8 Feb 2019 16:57:55 -0500 Subject: [PATCH 06/13] RF: Thread integration_window through BIDSRunVariableCollection --- bids/variables/kollekshuns.py | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/bids/variables/kollekshuns.py b/bids/variables/kollekshuns.py index c94ca712d..c278fbf09 100644 --- a/bids/variables/kollekshuns.py +++ b/bids/variables/kollekshuns.py @@ -243,7 +243,7 @@ def _all_dense(self): for v in self.variables.values()]) def resample(self, sampling_rate=None, variables=None, force_dense=False, - in_place=False, kind='linear'): + in_place=False, kind='linear', **kwargs): ''' Resample all dense variables (and optionally, sparse ones) to the specified sampling rate. @@ -276,7 +276,8 @@ def resample(self, sampling_rate=None, variables=None, force_dense=False, # None if in_place; no update needed _var = var.resample(sampling_rate, inplace=in_place, - kind=kind) + kind=kind, + **kwargs) if not in_place: _variables[name] = _var @@ -289,6 +290,7 @@ def resample(self, sampling_rate=None, variables=None, force_dense=False, def to_df(self, variables=None, format='wide', sparse=True, sampling_rate=None, include_sparse=True, include_dense=True, + integration_window=None, **kwargs): ''' Merge columns into a single pandas DataFrame. @@ -343,9 +345,10 @@ def to_df(self, variables=None, format='wide', sparse=True, sampling_rate = sampling_rate or self.sampling_rate # Make sure all variables have the same sampling rate - variables = list(self.resample(sampling_rate, variables, - force_dense=True, - in_place=False).values()) + variables = list( + self.resample(sampling_rate, variables, + force_dense=True, in_place=False, + integration_window=integration_window).values()) return super(BIDSRunVariableCollection, self).to_df(variables, format, **kwargs) From 979a4a283cf33f2082d5957bf9487311b453e791 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Fri, 8 Feb 2019 16:58:37 -0500 Subject: [PATCH 07/13] ENH: Perform integration over TA --- bids/variables/variables.py | 24 +++++++++++++++--------- 1 file changed, 15 insertions(+), 9 deletions(-) diff --git a/bids/variables/variables.py b/bids/variables/variables.py index c8b0856d3..1f492abbf 100644 --- a/bids/variables/variables.py +++ b/bids/variables/variables.py @@ -447,21 +447,27 @@ def resample(self, sampling_rate, integration_window=None, inplace=False, kind=' self.index = self._build_entity_index(self.run_info, sampling_rate) - x = np.arange(n) num = len(self.index) if integration_window is not None: # Number of timesteps in existing time series n_frames = int(np.floor(integration_window * old_sr)) - integrator = np.zeroes((n, num), dtype=np.uint8) - - # old_onsets = old_sr * self.index. - # onset_frames = sampling_rate + old_times = np.arange(n) / old_sr + new_times = np.arange(num) / sampling_rate + diffs = old_times.reshape((1, -1)) - new_times.reshape((-1, 1)) + integrator = ((diffs >= 0) & (diffs < integration_window)).astype(np.uint8) + del diffs + + old_vals = self.values.values + self.values = pd.DataFrame( + integrator.dot(old_vals) / integrator.sum(1, keepdims=True)) + else: + from scipy.interpolate import interp1d + x = np.arange(n) + f = interp1d(x, self.values.values.ravel(), kind=kind) + x_new = np.linspace(0, n - 1, num=num) + self.values = pd.DataFrame(f(x_new)) - from scipy.interpolate import interp1d - f = interp1d(x, self.values.values.ravel(), kind=kind) - x_new = np.linspace(0, n - 1, num=num) - self.values = pd.DataFrame(f(x_new)) assert len(self.values) == len(self.index) self.sampling_rate = sampling_rate From 998ce3714f0875ffe601cdf59aed6bc9e3e215d0 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Mon, 11 Feb 2019 13:49:46 -0500 Subject: [PATCH 08/13] FIX: TA within 1ms of TR is just TR --- bids/variables/io.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/bids/variables/io.py b/bids/variables/io.py index 64ef16916..8718b4db5 100644 --- a/bids/variables/io.py +++ b/bids/variables/io.py @@ -91,7 +91,10 @@ def _get_timing_info(img_md): slicetimes = sorted(img_md['SliceTiming']) # a, b ... z # z = final slice onset, b - a = slice duration - ta = slicetimes[-1] + slicetimes[1] - slicetimes[0] + ta = np.round(slicetimes[-1] + slicetimes[1] - slicetimes[0], 3) + # If the "silence" is <1ms, consider acquisition continuous + if np.abs(tr - ta) < 1e-3: + ta = tr else: ta = tr elif 'VolumeTiming' in img_md: From 10af7089ba5e8126ad70666ef6c1c0da6f667c2f Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Mon, 11 Feb 2019 14:19:01 -0500 Subject: [PATCH 09/13] RF: Use sparse matrix for integration --- bids/variables/variables.py | 17 ++++++++++------- 1 file changed, 10 insertions(+), 7 deletions(-) diff --git a/bids/variables/variables.py b/bids/variables/variables.py index 1f492abbf..83d01215d 100644 --- a/bids/variables/variables.py +++ b/bids/variables/variables.py @@ -450,17 +450,20 @@ def resample(self, sampling_rate, integration_window=None, inplace=False, kind=' num = len(self.index) if integration_window is not None: - # Number of timesteps in existing time series - n_frames = int(np.floor(integration_window * old_sr)) + from scipy.sparse import lil_matrix old_times = np.arange(n) / old_sr new_times = np.arange(num) / sampling_rate - diffs = old_times.reshape((1, -1)) - new_times.reshape((-1, 1)) - integrator = ((diffs >= 0) & (diffs < integration_window)).astype(np.uint8) - del diffs + integrator = lil_matrix((num, n), dtype=np.uint8) + count = None + for i, new_time in enumerate(new_times): + cols = (old_times >= new_time) & (old_times < new_time + integration_window) + # This could be determined analytically, but dodging off-by-one errors + if count is None: + count = np.sum(cols) + integrator[i, cols] = 1 old_vals = self.values.values - self.values = pd.DataFrame( - integrator.dot(old_vals) / integrator.sum(1, keepdims=True)) + self.values = pd.DataFrame(integrator.tocsr().dot(old_vals) / count) else: from scipy.interpolate import interp1d x = np.arange(n) From b6d82719fa8370c0a3dd7ceed90e2f7c7fef55fd Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Mon, 11 Feb 2019 15:56:11 -0500 Subject: [PATCH 10/13] RF: Drop unneeded TA is None guard --- bids/analysis/transformations/compute.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/bids/analysis/transformations/compute.py b/bids/analysis/transformations/compute.py index 741169dfa..58a712848 100644 --- a/bids/analysis/transformations/compute.py +++ b/bids/analysis/transformations/compute.py @@ -44,7 +44,7 @@ def _transform(self, var, model='spm', derivative=False, dispersion=False, TR = int(np.round(1000. * trs.pop())) TA = int(np.round(1000. * tas.pop())) SR = int(np.round(1000. / sr)) - if TA is None or TA < TR: + if TA < TR: # Use a unit that fits an whole number of times into both # the interscan interval (TR) and the integration window (TA) dt = gcd(TR, TA) From 0f0531c55fe19062a02a507b616f95853e581e9e Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Mon, 11 Feb 2019 15:58:17 -0500 Subject: [PATCH 11/13] DOC: Drop missing arguments from docstring --- bids/variables/entities.py | 2 -- 1 file changed, 2 deletions(-) diff --git a/bids/variables/entities.py b/bids/variables/entities.py index 0ba13d633..edde25be5 100644 --- a/bids/variables/entities.py +++ b/bids/variables/entities.py @@ -31,13 +31,11 @@ class RunNode(Node): ''' Represents a single Run in a BIDS project. Args: - id (int): The index of the run. entities (dict): Dictionary of entities for this Node. image_file (str): The full path to the corresponding nifti image. duration (float): Duration of the run, in seconds. repetition_time (float): TR for the run. acquisition_time (float): TA for the run. - task (str): The task name for this run. ''' def __init__(self, entities, image_file, duration, repetition_time, acquisition_time=None): From aa72c4666cc3f779e697eeba574e10873525d5be Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Tue, 12 Feb 2019 14:25:21 -0500 Subject: [PATCH 12/13] FIX: Make sure acquisition_time is defined --- bids/analysis/analysis.py | 1 + 1 file changed, 1 insertion(+) diff --git a/bids/analysis/analysis.py b/bids/analysis/analysis.py index 21e226bba..d5fcf5666 100644 --- a/bids/analysis/analysis.py +++ b/bids/analysis/analysis.py @@ -386,6 +386,7 @@ def get_design_matrix(self, names=None, format='long', mode='both', kwargs['timing'] = True kwargs['sparse'] = False + acquisition_time = None if sampling_rate == 'TR': trs = {var.run_info[0].tr for var in self.collection.variables.values()} tas = {var.run_info[0].ta for var in self.collection.variables.values()} From 258f24325462e4e069c37c69821cfd03a46fd237 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Wed, 13 Feb 2019 15:59:42 -0500 Subject: [PATCH 13/13] PY2: gcd --- bids/analysis/transformations/compute.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/bids/analysis/transformations/compute.py b/bids/analysis/transformations/compute.py index da6a9b8da..3254740d3 100644 --- a/bids/analysis/transformations/compute.py +++ b/bids/analysis/transformations/compute.py @@ -2,7 +2,6 @@ Transformations that primarily involve numerical computation on variables. ''' -from math import gcd import numpy as np import pandas as pd from bids.utils import listify @@ -10,6 +9,14 @@ from bids.analysis import hrf from bids.variables import SparseRunVariable, DenseRunVariable +try: + from math import gcd +except ImportError: + def gcd(a, b): + if b == 0: + return a + return gcd(b, a % b) + class Convolve(Transformation): """Convolve the input variable with an HRF.