🔀 Merge branch 'mri_robust_template' into fix/longitudinal-configuration-subject_id

Author: shnizzedy

CHANGELOG.md

@@ -24,13 +24,16 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Required positional parameter "wf" in input and output of `ingress_pipeconfig_paths` function, where a node to reorient templates is added to the `wf`.
 - Required positional parameter "orientation" to `resolve_resolution`.
 - Optional positional argument "cfg" to `create_lesion_preproc`.
+- `mri_robust_template` for longitudinal template generation.
+- `max_iter` parameter for longitudinal template generation.
 
 ### Changed
 
 - Moved `pygraphviz` from requirements to `graphviz` optional dependencies group.
 - Automatically tag untagged `subject_id` and `unique_id` as `!!str` when loading data config files.
 - Made orientation configurable (was hard-coded as "RPI").
 - Disabled variant image builds.
+- Made `mri_robust_template` default implementation for longitudinal template generation.
 
 ### Fixed
 

CPAC/longitudinal/__init__.py

@@ -18,6 +18,7 @@
 
 from CPAC.utils.docs import DOCS_URL_PREFIX
 
+assert isinstance(__doc__, str)
 __doc__ += f"""
 
 See {DOCS_URL_PREFIX}/user/longitudinal

CPAC/longitudinal/preproc.py

@@ -19,7 +19,9 @@
 
 from collections import Counter
 from multiprocessing.dummy import Pool as ThreadPool
+from multiprocessing.pool import Pool
 import os
+from typing import Literal, Optional
 
 import numpy as np
 import nibabel as nib
@@ -131,27 +133,23 @@ def norm_transformation(input_mat):
 
 
 def template_convergence(
-    mat_file, mat_type="matrix", convergence_threshold=np.finfo(np.float64).eps
-):
+    mat_file: str,
+    mat_type: Literal["matrix", "ITK"] = "matrix",
+    convergence_threshold: float | np.float64 = np.finfo(np.float64).eps,
+) -> bool:
     """Check that the deistance between matrices is smaller than the threshold.
 
     Calculate the distance between transformation matrix with a matrix of no transformation.
 
     Parameters
     ----------
-    mat_file : str
+    mat_file
         path to an fsl flirt matrix
-    mat_type : str
-        'matrix'(default), 'ITK'
+    mat_type
         The type of matrix used to represent the transformations
-    convergence_threshold : float
-        (numpy.finfo(np.float64).eps (default)) threshold for the convergence
+    convergence_threshold
         The threshold is how different from no transformation is the
         transformation matrix.
-
-    Returns
-    -------
-    bool
     """
     if mat_type == "matrix":
         translation, oth_transform = read_mat(mat_file)
@@ -347,50 +345,51 @@ def flirt_node(in_img, output_img, output_mat):
 
 
 def template_creation_flirt(
-    input_brain_list,
-    input_skull_list,
-    init_reg=None,
-    avg_method="median",
-    dof=12,
-    interp="trilinear",
-    cost="corratio",
-    mat_type="matrix",
-    convergence_threshold=-1,
-    thread_pool=2,
-    unique_id_list=None,
-):
+    input_brain_list: list[str],
+    input_skull_list: list[str],
+    init_reg: Optional[list[pe.Node]] = None,
+    avg_method: Literal["median", "mean", "std"] = "median",
+    dof: Literal[12, 9, 7, 6] = 12,
+    interp: Literal["trilinear", "nearestneighbour", "sinc", "spline"] = "trilinear",
+    cost: Literal[
+        "corratio", "mutualinfo", "normmi", "normcorr", "leastsq", "labeldiff", "bbr"
+    ] = "corratio",
+    mat_type: Literal["matrix", "ITK"] = "matrix",
+    convergence_threshold: float | np.float64 = -1,
+    max_iter: int = 5,
+    thread_pool: int | Pool = 2,
+    unique_id_list: Optional[list[str]] = None,
+) -> tuple[str, str, list[str], list[str], list[str]]:
     """Create a temporary template from a list of images.
 
     Parameters
     ----------
-    input_brain_list : list of str
+    input_brain_list
         list of brain images paths
-    input_skull_list : list of str
+    input_skull_list
         list of skull images paths
-    init_reg : list of Node
+    init_reg
         (default None so no initial registration performed)
         the output of the function register_img_list with another reference
         Reuter et al. 2012 (NeuroImage) section "Improved template estimation"
         doi:10.1016/j.neuroimage.2012.02.084 uses a ramdomly
         selected image from the input dataset
-    avg_method : str
-        function names from numpy library such as 'median', 'mean', 'std' ...
-    dof : integer (int of long)
-        number of transform degrees of freedom (FLIRT) (12 by default)
-    interp : str
-        ('trilinear' (default) or 'nearestneighbour' or 'sinc' or 'spline')
+    avg_method
+        function names from numpy library
+    dof
+        number of transform degrees of freedom (FLIRT)
+    interp
         final interpolation method used in reslicing
-    cost : str
-        ('mutualinfo' or 'corratio' (default) or 'normcorr' or 'normmi' or
-         'leastsq' or 'labeldiff' or 'bbr')
+    cost
         cost function
-    mat_type : str
-        'matrix'(default), 'ITK'
+    mat_type
         The type of matrix used to represent the transformations
-    convergence_threshold : float
+    convergence_threshold
         (numpy.finfo(np.float64).eps (default)) threshold for the convergence
         The threshold is how different from no transformation is the
         transformation matrix.
+    max_iter
+        Maximum number of iterations if transformation does not converge
     thread_pool : int or multiprocessing.dummy.Pool
         (default 2) number of threads. You can also provide a Pool so the
         node will be added to it to be run.
@@ -496,7 +495,14 @@ def template_creation_flirt(
     and the loop stops when this temporary template is close enough (with a transformation
     distance smaller than the threshold) to all the images of the precedent iteration.
     """
-    while not converged:
+    iterator = 1
+    iteration = 0
+    if max_iter == -1:
+        # make iteration < max_iter always True
+        iterator = 0
+        iteration = -2
+    while not converged and iteration < max_iter:
+        iteration += iterator
         temporary_brain_template, temporary_skull_template = create_temporary_template(
             input_brain_list=output_brain_list,
             input_skull_list=output_skull_list,
@@ -628,6 +634,7 @@ def subject_specific_template(
                     "cost",
                     "mat_type",
                     "convergence_threshold",
+                    "max_iter",
                     "thread_pool",
                     "unique_id_list",
                 ],

CPAC/longitudinal/robust_template.py

@@ -0,0 +1,18 @@
+# -*- coding: utf-8 -*-
+# Copyright (C) 2024  C-PAC Developers
+
+# This file is part of C-PAC.
+
+# C-PAC is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Lesser General Public License as published by the
+# Free Software Foundation, either version 3 of the License, or (at your
+# option) any later version.
+
+# C-PAC is distributed in the hope that it will be useful, but WITHOUT
+# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+
+# You should have received a copy of the GNU Lesser General Public
+# License along with C-PAC. If not, see <https://www.gnu.org/licenses/>.
+"""Create longitudinal template using ``mri_robust_template``."""

CPAC/longitudinal/wf/__init__.py

@@ -0,0 +1,25 @@
+# Copyright (C) 2024  C-PAC Developers
+
+# This file is part of C-PAC.
+
+# C-PAC is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Lesser General Public License as published by the
+# Free Software Foundation, either version 3 of the License, or (at your
+# option) any later version.
+
+# C-PAC is distributed in the hope that it will be useful, but WITHOUT
+# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+
+# You should have received a copy of the GNU Lesser General Public
+# License along with C-PAC. If not, see <https://www.gnu.org/licenses/>.
+"""Workflows for longitudinal preprocessing."""
+
+from CPAC.utils.docs import DOCS_URL_PREFIX
+
+assert isinstance(__doc__, str)
+__doc__ += f"""
+
+See {DOCS_URL_PREFIX}/user/longitudinal
+"""  # noqa: A001