EHN detact different affine matrix from input (#19)

* EHN detact different affine matrix from input

* add gitattribute

* ADD codecov.yml

Author: htwangtw

.gitattributes

@@ -0,0 +1 @@
+.git_archival.txt export-subst

README.md

@@ -1,5 +1,6 @@
 [![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://opensource.org/licenses/MIT)
-[![codecov](https://codecov.io/gh/SIMEXP/giga_connectome/branch/main/graph/badge.svg?token=P4EGV7NKZ8)](https://codecov.io/gh/SIMEXP/giga_connectome)
+[![codecov](https://codecov.io/github/SIMEXP/giga_connectome/branch/main/graph/badge.svg?token=TYE4UURNTQ)](https://codecov.io/github/SIMEXP/giga_auto_qc)
+
 # giga_connectome
 
 BIDS App to generate connectome from fMRIPrep outputs.

codecov.yml

@@ -0,0 +1,26 @@
+---
+# similar to scikit-learn .codecov.yml
+# coverage:
+#   status:
+#     project:
+#       default:
+#         # Commits pushed to main should not make the overall
+#         # project coverage decrease by more than 1%:
+#         target: auto
+#         threshold: 1%
+#     patch:
+#       default:
+#         # Be tolerant on slight code coverage diff on PRs to limit
+#         # noisy red coverage status on github PRs.
+#         target: auto
+#         threshold: 1%
+comment:                  # this is a top-level key
+  layout: reach, diff, flags, files
+  behavior: default
+  require_changes: false  # if true: only post the comment if coverage changes
+  require_base: no        # [yes :: must have a base report to post]
+  require_head: yes       # [yes :: must have a head report to post]
+
+ignore:
+- '*/tests/'    # ignore folders related to testing
+- '*/data/'

giga_connectome/mask.py

@@ -1,15 +1,21 @@
 import os
 import re
 import json
-from typing import Optional, Union, List
+from typing import Optional, Union, List, Tuple
 
 from pathlib import Path
 from tqdm import tqdm
 import nibabel as nib
-
 from nilearn.masking import compute_multi_epi_mask
-from nilearn.image import resample_to_img, new_img_like, get_data, math_img
+from nilearn.image import (
+    resample_to_img,
+    new_img_like,
+    get_data,
+    math_img,
+    load_img,
+)
 from nibabel import Nifti1Image
+import numpy as np
 from scipy.ndimage import binary_closing
 
 from pkg_resources import resource_filename
@@ -74,6 +80,7 @@ def generate_group_mask(
     template: str = "MNI152NLin2009cAsym",
     templateflow_dir: Optional[Path] = None,
     n_iter: int = 2,
+    verbose: int = 1,
 ) -> Nifti1Image:
     """
     Generate a group EPI grey matter mask, and overlaid with a MNI grey
@@ -98,6 +105,9 @@ def generate_group_mask(
         Number of repetitions of dilation and erosion steps performed in
         scipy.ndimage.binary_closing function.
 
+    verbose :
+        Level of verbosity.
+
     Keyword Arguments
     -----------------
     Used to filter the cirret
@@ -109,7 +119,12 @@ def generate_group_mask(
     nibabel.nifti1.Nifti1Image
         EPI (grey matter) mask for the current group of subjects.
     """
-    # TODO: subject native space grey matter mask???
+    if verbose > 1:
+        print(f"Found {len(imgs)} masks")
+    if exclude := _check_mask_affine(imgs, verbose):
+        imgs, __annotations__ = _get_consistent_masks(imgs, exclude)
+        if verbose > 1:
+            print(f"Remaining: {len(imgs)} masks")
 
     # templateflow environment setting to get around network issue
     if templateflow_dir and templateflow_dir.exists():
@@ -258,3 +273,107 @@ def _load_config(atlas: Union[str, Path, dict]) -> dict:
     else:
         raise ValueError(f"Invalid input: {atlas}")
     return atlas_config
+
+
+def _get_consistent_masks(
+    mask_imgs: List[Union[Path, str, Nifti1Image]], exclude: List[int]
+) -> Tuple[List[int], List[str]]:
+    """Create a list of masks that has the same affine.
+
+    Parameters
+    ----------
+
+    mask_imgs :
+        The original list of functional masks
+
+    exclude :
+        List of index to exclude.
+
+    Returns
+    -------
+    List of str
+        Functional masks with the same affine.
+
+    List of str
+        Identidiers of scans with a different affine.
+    """
+    weird_mask_identifiers = []
+    odd_masks = np.array(mask_imgs)[np.array(exclude)]
+    odd_masks = odd_masks.tolist()
+    for odd_file in odd_masks:
+        identifier = Path(odd_file).name.split("_space")[0]
+        weird_mask_identifiers.append(identifier)
+    cleaned_func_masks = set(mask_imgs) - set(odd_masks)
+    cleaned_func_masks = list(cleaned_func_masks)
+    return cleaned_func_masks, weird_mask_identifiers
+
+
+def _check_mask_affine(
+    mask_imgs: List[Union[Path, str, Nifti1Image]], verbose: int = 1
+) -> Union[list, None]:
+    """Given a list of input mask images, show the most common affine matrix
+    and subjects with different values.
+
+    Parameters
+    ----------
+    mask_imgs : :obj:`list` of Niimg-like objects
+        See :ref:`extracting_data`.
+        3D or 4D EPI image with same affine.
+
+    verbose :
+        Level of verbosity.
+
+    Returns
+    -------
+
+    List or None
+        Index of masks with odd affine matrix. Return None when all masks have
+        the same affine matrix.
+    """
+    # save all header and affine info in hashable type...
+    header_info = {"affine": []}
+    key_to_header = {}
+    for this_mask in mask_imgs:
+        img = load_img(this_mask)
+        affine = img.affine
+        affine_hashable = str(affine)
+        header_info["affine"].append(affine_hashable)
+        if affine_hashable not in key_to_header:
+            key_to_header[affine_hashable] = affine
+
+    if isinstance(mask_imgs[0], Nifti1Image):
+        mask_imgs = np.arange(len(mask_imgs))
+    else:
+        mask_imgs = np.array(mask_imgs)
+    # get most common values
+    common_affine = max(
+        set(header_info["affine"]), key=header_info["affine"].count
+    )
+    if verbose > 0:
+        print(
+            f"We found {len(set(header_info['affine']))} unique affine "
+            f"matrices. The most common one is "
+            f"{key_to_header[common_affine]}"
+        )
+    odd_balls = set(header_info["affine"]) - {common_affine}
+    if not odd_balls:
+        return None
+
+    exclude = []
+    for ob in odd_balls:
+        ob_index = [
+            i for i, aff in enumerate(header_info["affine"]) if aff == ob
+        ]
+        if verbose > 1:
+            print(
+                "The following subjects has a different affine matrix "
+                f"({key_to_header[ob]}) comparing to the most common value: "
+                f"{mask_imgs[ob_index]}."
+            )
+        exclude += ob_index
+    if verbose > 0:
+        print(
+            f"{len(exclude)} out of {len(mask_imgs)} has "
+            "different affine matrix. Ignore when creating group mask."
+        )
+    return sorted(exclude)

giga_connectome/tests/test_mask.py

@@ -1,22 +1,55 @@
 import pytest
 import numpy as np
-from giga_connectome.mask import generate_group_mask
+from giga_connectome import mask
 from nilearn import datasets
+from nibabel import Nifti1Image
 
 
 def test_generate_group_mask():
     """Generate group epi grey matter mask and resample atlas."""
     data = datasets.fetch_development_fmri(n_subjects=3)
     imgs = data.func
 
-    group_epi_mask = generate_group_mask(imgs)
+    group_epi_mask = mask.generate_group_mask(imgs)
     # match the post processing details: https://osf.io/wjtyq
     assert group_epi_mask.shape == (50, 59, 50)
-    diff_tpl = generate_group_mask(imgs, template="MNI152NLin2009aAsym")
+    diff_tpl = mask.generate_group_mask(imgs, template="MNI152NLin2009aAsym")
     assert diff_tpl.shape == (50, 59, 50)
 
     # test bad inputs
     with pytest.raises(
         ValueError, match="TemplateFlow does not supply template blah"
     ):
-        generate_group_mask(imgs, template="blah")
+        mask.generate_group_mask(imgs, template="blah")
+
+
+def test_check_mask_affine():
+    """Check odd affine detection."""
+
+    img_base = np.zeros([5, 5, 6])
+    processed_vol = img_base.copy()
+    processed_vol[2:4, 2:4, 2:4] += 1
+    processed = Nifti1Image(processed_vol, np.eye(4))
+    weird = Nifti1Image(processed_vol, np.eye(4) * np.array([1, 1, 1.5, 1]).T)
+    weird2 = Nifti1Image(processed_vol, np.eye(4) * np.array([1, 1, 1.6, 1]).T)
+    exclude = mask._check_mask_affine(
+        [processed, processed, processed, processed, weird, weird, weird2],
+        verbose=2,
+    )
+    assert len(exclude) == 3
+    assert exclude == [4, 5, 6]
+
+
+def test_get_consistent_masks():
+    """Check odd affine detection."""
+    mask_imgs = [
+        f"sub-{i + 1:2d}_task-rest_space-MNI_desc-brain_mask.nii.gz"
+        for i in range(10)
+    ]
+    exclude = [1, 2, 5]
+    (
+        cleaned_func_masks,
+        weird_mask_identifiers,
+    ) = mask._get_consistent_masks(mask_imgs, exclude)
+    assert len(cleaned_func_masks) == 7
+    assert len(weird_mask_identifiers) == 3

giga_connectome/tests/test_metadata.py

@@ -10,7 +10,7 @@
 
 
 @pytest.mark.skipif(
-    IN_GITHUB_ACTIONS, reason="Test doesn't work in Github Actions."
+    IN_GITHUB_ACTIONS, reason="Test data is not set up on Github Actions."
 )
 def test_get_metadata():
     """Check the function can load from different fMRIPrep dataset."""

giga_connectome/workflow.py

@@ -25,7 +25,7 @@ def workflow(args):
     output_dir = args.output_dir
     working_dir = args.work_dir
     analysis_level = args.analysis_level
-    
+
     subjects = utils.get_subject_lists(args.participant_label, bids_dir)
     strategy = get_denoise_strategy_parameters(args.denoise_strategy)
     atlas = load_atlas_setting(args.atlas)
@@ -43,7 +43,7 @@ def workflow(args):
     # https://github.com/nipreps/fmriprep/blob/689ad26811cfb18771fdb8d7dc208fe24d27e65c/fmriprep/cli/parser.py#L72
     fmriprep_bids_layout = BIDSLayout(
         root=bids_dir,
-        database_path=working_dir,
+        database_path=bids_dir,
         validate=False,
         derivatives=True,
     )