From 8cfba0610b03d0b8ca20b30037baef5f86b7bd23 Mon Sep 17 00:00:00 2001 From: "Christopher J. Markiewicz" Date: Sat, 26 Jan 2019 21:32:50 -0500 Subject: [PATCH 1/3] BUG: Demonstrate resampling bug --- bids/variables/tests/test_variables.py | 28 ++++++++++++++++++++++++++ 1 file changed, 28 insertions(+) diff --git a/bids/variables/tests/test_variables.py b/bids/variables/tests/test_variables.py index 339ad8553..74876489c 100644 --- a/bids/variables/tests/test_variables.py +++ b/bids/variables/tests/test_variables.py @@ -1,5 +1,6 @@ from bids.layout import BIDSLayout import pytest +import os from os.path import join from bids.tests import get_test_data_path from bids.variables import (merge_variables, DenseRunVariable, SimpleVariable, @@ -7,7 +8,9 @@ from bids.variables.entities import RunInfo import numpy as np import pandas as pd +import nibabel as nb import uuid +import json def generate_DEV(name='test', sr=20, duration=480): @@ -174,3 +177,28 @@ def test_filter_simple_variable(layout2): assert merged.filter({'nonexistent': 2}, strict=True) is None merged.filter({'acquisition': 'fullbrain'}, inplace=True) assert merged.to_df().shape == (40, 9) + + +@pytest.mark.parametrize( + "TR, nvols", + [(2.00000, 251), + (2.000001, 251)]) +def test_resampling_edge_case(tmpdir, TR, nvols): + tmpdir.chdir() + os.makedirs('sub-01/func') + with open('sub-01/func/sub-01_task-task_events.tsv', 'w') as fobj: + fobj.write('onset\tduration\tval\n1\t0.1\t1\n') + with open('sub-01/func/sub-01_task-task_bold.json', 'w') as fobj: + json.dump({'RepetitionTime': TR}, fobj) + + dataobj = np.zeros((5, 5, 5, nvols), dtype=np.int16) + affine = np.diag((2.5, 2.5, 2.5, 1)) + img = nb.Nifti1Image(dataobj, affine) + img.header.set_zooms((2.5, 2.5, 2.5, TR)) + img.to_filename('sub-01/func/sub-01_task-task_bold.nii.gz') + + layout = BIDSLayout('.', validate=False) + coll = load_variables(layout).get_collections('run')[0] + dense_var = coll.variables['val'].to_dense(coll.sampling_rate) + regressor = dense_var.resample(1.0 / TR).values + assert regressor.shape == (nvols, 1) From c44fedc2a6d610bc6c5a69561b85a3d3185b6b4f Mon Sep 17 00:00:00 2001 From: Adina Wagner Date: Tue, 29 Jan 2019 14:44:19 -0500 Subject: [PATCH 2/3] FIX: revert to length of self.Index when dimensions differ slightly, blow for larger deviations --- bids/variables/variables.py | 21 ++++++++++++++++++++- 1 file changed, 20 insertions(+), 1 deletion(-) diff --git a/bids/variables/variables.py b/bids/variables/variables.py index 0f98548bb..ce7e9fa79 100644 --- a/bids/variables/variables.py +++ b/bids/variables/variables.py @@ -448,12 +448,31 @@ def resample(self, sampling_rate, inplace=False, kind='linear'): self.index = self._build_entity_index(self.run_info, sampling_rate) x = np.arange(n) - num = int(np.ceil(n * sampling_rate / old_sr)) + # In Index we trust! + num = len(self.index) + # _build_entity_index computation above does it per run, + # which possibly provides a more stable solution and yadadada + num_computed = int(np.ceil(n * sampling_rate / old_sr)) + num_diff = abs(num - num_computed) + if num_diff > 1: + raise RuntimeError( + "Our internal assumptions about resampling are too faulty " + "to continue.") + elif num_diff: + assert num_diff == 1 # the only way to get here + import warnings + warnings.warn( + "Probably due to multiple runs we ran into a slight " + "divergence between obtained number of samples in the index (computed " + "per each run) and overall estimate down/up-sampled samples. " + "But that is ok") + 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 319f81cf3325b5fba472e8e1fa24740b1c0ee514 Mon Sep 17 00:00:00 2001 From: Adina Wagner Date: Fri, 1 Feb 2019 13:37:25 -0500 Subject: [PATCH 3/3] rely exclusively on self.index length during resampling --- bids/variables/variables.py | 18 ------------------ 1 file changed, 18 deletions(-) diff --git a/bids/variables/variables.py b/bids/variables/variables.py index ce7e9fa79..6cd06d9a1 100644 --- a/bids/variables/variables.py +++ b/bids/variables/variables.py @@ -448,25 +448,7 @@ def resample(self, sampling_rate, inplace=False, kind='linear'): self.index = self._build_entity_index(self.run_info, sampling_rate) x = np.arange(n) - # In Index we trust! num = len(self.index) - # _build_entity_index computation above does it per run, - # which possibly provides a more stable solution and yadadada - num_computed = int(np.ceil(n * sampling_rate / old_sr)) - num_diff = abs(num - num_computed) - if num_diff > 1: - raise RuntimeError( - "Our internal assumptions about resampling are too faulty " - "to continue.") - elif num_diff: - assert num_diff == 1 # the only way to get here - import warnings - warnings.warn( - "Probably due to multiple runs we ran into a slight " - "divergence between obtained number of samples in the index (computed " - "per each run) and overall estimate down/up-sampled samples. " - "But that is ok") - from scipy.interpolate import interp1d f = interp1d(x, self.values.values.ravel(), kind=kind)