Skip to content

TEST: Allow SVD columns to flip sign #3405

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 2 commits into from
Nov 6, 2021
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
108 changes: 68 additions & 40 deletions nipype/algorithms/tests/test_CompCor.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,41 @@
from ..confounds import CompCor, TCompCor, ACompCor


def close_up_to_column_sign(a, b, rtol=1e-05, atol=1e-08, equal_nan=False):
"""SVD can produce sign flips on a per-column basis."""
kwargs = dict(rtol=rtol, atol=atol, equal_nan=equal_nan)
if np.allclose(a, b, **kwargs):
return True

ret = True
for acol, bcol in zip(a.T, b.T):
ret &= np.allclose(acol, bcol, **kwargs) or np.allclose(acol, -bcol, **kwargs)
if not ret:
break

return ret


@pytest.mark.parametrize(
"a, b, close",
[
([[0.1, 0.2], [0.3, 0.4]], [[-0.1, 0.2], [-0.3, 0.4]], True),
([[0.1, 0.2], [0.3, 0.4]], [[-0.1, 0.2], [0.3, -0.4]], False),
],
)
def test_close_up_to_column_sign(a, b, close):
a = np.asanyarray(a)
b = np.asanyarray(b)
assert close_up_to_column_sign(a, b) == close
# Sign flips of all columns never changes result
assert close_up_to_column_sign(a, -b) == close
assert close_up_to_column_sign(-a, b) == close
assert close_up_to_column_sign(-a, -b) == close
# Trivial case
assert close_up_to_column_sign(a, a)
assert close_up_to_column_sign(b, b)


class TestCompCor:
"""Note: Tests currently do a poor job of testing functionality"""

Expand Down Expand Up @@ -42,11 +77,11 @@ def setup_class(self, tmpdir):

def test_compcor(self):
expected_components = [
["-0.1989607212", "-0.5753813646"],
["0.5692369697", "0.5674945949"],
["-0.6662573243", "0.4675843432"],
["0.4206466244", "-0.3361270124"],
["-0.1246655485", "-0.1235705610"],
[-0.1989607212, -0.5753813646],
[0.5692369697, 0.5674945949],
[-0.6662573243, 0.4675843432],
[0.4206466244, -0.3361270124],
[-0.1246655485, -0.1235705610],
]

self.run_cc(
Expand All @@ -73,11 +108,11 @@ def test_compcor(self):

def test_compcor_variance_threshold_and_metadata(self):
expected_components = [
["-0.2027150345", "-0.4954813834"],
["0.2565929051", "0.7866217875"],
["-0.3550986008", "-0.0089784905"],
["0.7512786244", "-0.3599828482"],
["-0.4500578942", "0.0778209345"],
[-0.2027150345, -0.4954813834],
[0.2565929051, 0.7866217875],
[-0.3550986008, -0.0089784905],
[0.7512786244, -0.3599828482],
[-0.4500578942, 0.0778209345],
]
expected_metadata = {
"component": "CompCor00",
Expand Down Expand Up @@ -111,11 +146,11 @@ def test_tcompcor(self):
self.run_cc(
ccinterface,
[
["-0.1114536190", "-0.4632908609"],
["0.4566907310", "0.6983205193"],
["-0.7132557407", "0.1340170559"],
["0.5022537643", "-0.5098322262"],
["-0.1342351356", "0.1407855119"],
[-0.1114536190, -0.4632908609],
[0.4566907310, 0.6983205193],
[-0.7132557407, 0.1340170559],
[0.5022537643, -0.5098322262],
[-0.1342351356, 0.1407855119],
],
"tCompCor",
)
Expand All @@ -138,11 +173,11 @@ def test_compcor_no_regress_poly(self):
pre_filter=False,
),
[
["0.4451946442", "-0.7683311482"],
["-0.4285129505", "-0.0926034137"],
["0.5721540256", "0.5608764842"],
["-0.5367548139", "0.0059943226"],
["-0.0520809054", "0.2940637551"],
[0.4451946442, -0.7683311482],
[-0.4285129505, -0.0926034137],
[0.5721540256, 0.5608764842],
[-0.5367548139, 0.0059943226],
[-0.0520809054, 0.2940637551],
],
)

Expand Down Expand Up @@ -225,27 +260,20 @@ def run_cc(
assert os.path.getsize(expected_file) > 0

with open(ccresult.outputs.components_file, "r") as components_file:
if expected_n_components is None:
expected_n_components = min(
ccinterface.inputs.num_components, self.fake_data.shape[3]
)
header = components_file.readline().rstrip().split("\t")
components_data = np.loadtxt(components_file, delimiter="\t")

if expected_n_components is None:
expected_n_components = min(
ccinterface.inputs.num_components, self.fake_data.shape[3]
)

assert header == [
f"{expected_header}{i:02d}" for i in range(expected_n_components)
]

components_data = [line.rstrip().split("\t") for line in components_file]

# the first item will be '#', we can throw it out
header = components_data.pop(0)
expected_header = [
expected_header + "{:02d}".format(i)
for i in range(expected_n_components)
]
for i, heading in enumerate(header):
assert expected_header[i] in heading

num_got_timepoints = len(components_data)
assert num_got_timepoints == self.fake_data.shape[3]
for index, timepoint in enumerate(components_data):
assert len(timepoint) == expected_n_components
assert timepoint[:2] == expected_components[index]
assert components_data.shape == (self.fake_data.shape[3], expected_n_components)
assert close_up_to_column_sign(components_data[:, :2], expected_components)

if ccinterface.inputs.save_metadata:
expected_metadata_file = ccinterface._list_outputs()["metadata_file"]
Expand Down