diff --git a/docs/builtin.rst b/docs/builtin.rst index 921d2a66a..3c6e0b758 100644 --- a/docs/builtin.rst +++ b/docs/builtin.rst @@ -182,6 +182,16 @@ Available - Calculate (f)ALFF and aggregate using spheres placed on coordinates - Done - 0.0.1 + * - :class:`junifer.markers.EdgeCentricFCParcels` + - Calculate edge-centric functional connectivity over parcellation, as found in + `Jo et al. (2021) `_ + - Done + - 0.0.2 + * - :class:`junifer.markers.EdgeCentricFCSpheres` + - Calculate edge-centric functional connectivity over spheres placed on coordinates, + as found in `Jo et al. (2021) `_ + - Done + - 0.0.2 Planned ~~~~~~~ @@ -199,10 +209,6 @@ Planned * - Permutation entropy, Range entropy, Multiscale entropy and Hurst exponent - Calculate Permutation entropy, Range entropy, Multiscale entropy and Hurst exponent - :gh:`61` - * - EdgeCentricFC - - Calculate edge-centric functional connectivity - - :gh:`64` - Parcellations ------------- diff --git a/docs/changes/latest.inc b/docs/changes/latest.inc index 0fa9aa7d2..84b6ae143 100644 --- a/docs/changes/latest.inc +++ b/docs/changes/latest.inc @@ -32,6 +32,9 @@ Enhancements - Organize functional connectivity markers in ``junifer.markers.functional_connectivity`` (:gh:`107` by `Synchon Mandal`_). +- Add :class:`junifer.markers.EdgeCentricFCParcels` and :class:`junifer.markers.EdgeCentricFCSpheres` + (:gh:`64` by `Leonard Sasse`_). + Bugs ~~~~ diff --git a/junifer/markers/__init__.py b/junifer/markers/__init__.py index de37c5a26..badbd2f30 100644 --- a/junifer/markers/__init__.py +++ b/junifer/markers/__init__.py @@ -14,6 +14,8 @@ from .functional_connectivity import ( FunctionalConnectivityParcels, FunctionalConnectivitySpheres, CrossParcellationFC, + EdgeCentricFCParcels, + EdgeCentricFCSpheres, ) from .reho import ReHoParcels, ReHoSpheres from .falff import ( diff --git a/junifer/markers/functional_connectivity/__init__.py b/junifer/markers/functional_connectivity/__init__.py index 8763a661d..e1eb9fab5 100644 --- a/junifer/markers/functional_connectivity/__init__.py +++ b/junifer/markers/functional_connectivity/__init__.py @@ -6,3 +6,5 @@ from .functional_connectivity_parcels import FunctionalConnectivityParcels from .functional_connectivity_spheres import FunctionalConnectivitySpheres from .crossparcellation_functional_connectivity import CrossParcellationFC +from .edge_functional_connectivity_parcels import EdgeCentricFCParcels +from .edge_functional_connectivity_spheres import EdgeCentricFCSpheres diff --git a/junifer/markers/functional_connectivity/edge_functional_connectivity_parcels.py b/junifer/markers/functional_connectivity/edge_functional_connectivity_parcels.py new file mode 100644 index 000000000..eb02f56c0 --- /dev/null +++ b/junifer/markers/functional_connectivity/edge_functional_connectivity_parcels.py @@ -0,0 +1,90 @@ +"""Provide class for edge-centric functional connectivity using parcels.""" + +# Authors: Leonard Sasse +# Synchon Mandal +# License: AGPL + +from typing import Any, Dict, List, Optional, Union + +from ...api.decorators import register_marker +from ..parcel_aggregation import ParcelAggregation +from ..utils import _ets +from .functional_connectivity_base import FunctionalConnectivityBase + + +@register_marker +class EdgeCentricFCParcels(FunctionalConnectivityBase): + """Class for edge-centric FC using parcellations. + + Parameters + ---------- + parcellation : str or list of str + The name(s) of the parcellation(s). Check valid options by calling + :func:`junifer.data.parcellations.list_parcellations`. + agg_method : str, optional + The method to perform aggregation of BOLD time series. + Check valid options in :func:`junifer.stats.get_aggfunc_by_name` + (default "mean"). + agg_method_params : dict, optional + Parameters to pass to the aggregation function. Check valid options in + :func:`junifer.stats.get_aggfunc_by_name` (default None). + cor_method : str, optional + The method to perform correlation. Check valid options in + :class:`nilearn.connectome.ConnectivityMeasure` + (default "covariance"). + cor_method_params : dict, optional + Parameters to pass to the correlation function. Check valid options in + :class:`nilearn.connectome.ConnectivityMeasure` (default None). + mask : str, optional + The name of the mask to apply to regions before extracting signals. + Check valid options by calling :func:`junifer.data.masks.list_masks` + (default None). + name : str, optional + The name of the marker. If None, will use the class name (default + None). + + References + ---------- + .. [1] Jo et al. (2021) + Subject identification using + edge-centric functional connectivity + doi: https://doi.org/10.1016/j.neuroimage.2021.118204 + + """ + + def __init__( + self, + parcellation: Union[str, List[str]], + agg_method: str = "mean", + agg_method_params: Optional[Dict] = None, + cor_method: str = "covariance", + cor_method_params: Optional[Dict] = None, + mask: Optional[str] = None, + name: Optional[str] = None, + ) -> None: + self.parcellation = parcellation + super().__init__( + agg_method=agg_method, + agg_method_params=agg_method_params, + cor_method=cor_method, + cor_method_params=cor_method_params, + mask=mask, + name=name, + ) + + def aggregate(self, input: Dict[str, Any]) -> Dict: + """Perform parcel aggregation.""" + parcel_aggregation = ParcelAggregation( + parcellation=self.parcellation, + method=self.agg_method, + method_params=self.agg_method_params, + mask=self.mask, + on="BOLD", + ) + + bold_aggregated = parcel_aggregation.compute(input) + ets, edge_names = _ets( + bold_aggregated["data"], bold_aggregated["columns"] + ) + + return dict(data=ets, columns=edge_names) diff --git a/junifer/markers/functional_connectivity/edge_functional_connectivity_spheres.py b/junifer/markers/functional_connectivity/edge_functional_connectivity_spheres.py new file mode 100644 index 000000000..6fcabffbe --- /dev/null +++ b/junifer/markers/functional_connectivity/edge_functional_connectivity_spheres.py @@ -0,0 +1,97 @@ +"""Provide class for edge-centric functional connectivity using spheres.""" + +# Authors: Leonard Sasse +# Synchon Mandal +# License: AGPL + +from typing import Any, Dict, Optional + +from ...api.decorators import register_marker +from ..sphere_aggregation import SphereAggregation +from ..utils import _ets, raise_error +from .functional_connectivity_base import FunctionalConnectivityBase + + +@register_marker +class EdgeCentricFCSpheres(FunctionalConnectivityBase): + """Class for edge-centric FC using coordinates (spheres). + + Parameters + ---------- + coords : str + The name of the coordinates list to use. See + :func:`junifer.data.coordinates.list_coordinates` for options. + radius : float, optional + The radius of the sphere in mm. If None, the signal will be extracted + from a single voxel. See :class:`nilearn.maskers.NiftiSpheresMasker` + for more information (default None). + agg_method : str, optional + The aggregation method to use. + See :func:`junifer.stats.get_aggfunc_by_name` for more information + (default None). + agg_method_params : dict, optional + The parameters to pass to the aggregation method (default None). + cor_method : str, optional + The method to perform correlation using. Check valid options in + :class:`nilearn.connectome.ConnectivityMeasure` (default "covariance"). + cor_method_params : dict, optional + Parameters to pass to the correlation function. Check valid options in + :class:`nilearn.connectome.ConnectivityMeasure` (default None). + mask : str, optional + The name of the mask to apply to regions before extracting signals. + Check valid options by calling :func:`junifer.data.masks.list_masks` + (default None). + name : str, optional + The name of the marker. By default, it will use + KIND_EdgeCentricFCSpheres where KIND is the kind of data it + was applied to (default None). + + References + ---------- + .. [1] Jo et al. (2021) + Subject identification using + edge-centric functional connectivity + doi: https://doi.org/10.1016/j.neuroimage.2021.118204 + + """ + + def __init__( + self, + coords: str, + radius: Optional[float] = None, + agg_method: str = "mean", + agg_method_params: Optional[Dict] = None, + cor_method: str = "covariance", + cor_method_params: Optional[Dict] = None, + mask: Optional[str] = None, + name: Optional[str] = None, + ) -> None: + self.coords = coords + self.radius = radius + if radius is None or radius <= 0: + raise_error(f"radius should be > 0: provided {radius}") + super().__init__( + agg_method=agg_method, + agg_method_params=agg_method_params, + cor_method=cor_method, + cor_method_params=cor_method_params, + mask=mask, + name=name, + ) + + def aggregate(self, input: Dict[str, Any]) -> Dict: + """Perform sphere aggregation.""" + sphere_aggregation = SphereAggregation( + coords=self.coords, + radius=self.radius, + method=self.agg_method, + method_params=self.agg_method_params, + mask=self.mask, + on="BOLD", + ) + bold_aggregated = sphere_aggregation.compute(input) + ets, edge_names = _ets( + bold_aggregated["data"], bold_aggregated["columns"] + ) + + return dict(data=ets, columns=edge_names) diff --git a/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py new file mode 100644 index 000000000..5257980a5 --- /dev/null +++ b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py @@ -0,0 +1,60 @@ +"""Provide tests for edge-centric functional connectivity using parcels.""" + +# Authors: Leonard Sasse +# Synchon Mandal +# License: AGPL + +from pathlib import Path + +from nilearn import datasets, image + +from junifer.markers.functional_connectivity import EdgeCentricFCParcels +from junifer.storage import SQLiteFeatureStorage + + +def test_EdgeCentricFCParcels(tmp_path: Path) -> None: + """Test EdgeCentricFCParcels. + + Parameters + ---------- + tmp_path : pathlib.Path + The path to the test directory. + + """ + # get a dataset + ni_data = datasets.fetch_spm_auditory(subject_id="sub001") + fmri_img = image.concat_imgs(ni_data.func) # type: ignore + + # Check empirical correlation method parameters + efc = EdgeCentricFCParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + cor_method_params={"empirical": True}, + ) + all_out = efc.fit_transform({"BOLD": {"data": fmri_img, "meta": {}}}) + + out = all_out["BOLD"] + + # for 16 ROIs we should get (16 * (16 -1) / 2) edges in the ETS + n_edges = int(16 * (16 - 1) / 2) + assert "data" in out + assert "row_names" in out + assert "col_names" in out + assert out["data"].shape[0] == n_edges + assert out["data"].shape[1] == n_edges + assert len(set(out["row_names"])) == n_edges + assert len(set(out["col_names"])) == n_edges + + # check correct output + assert efc.get_output_type("BOLD") == "matrix" + + uri = tmp_path / "test_fc_parcellation.sqlite" + # Single storage, must be the uri + storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") + meta = {"element": {"subject": "test"}, "dependencies": {"numpy"}} + input = {"BOLD": {"data": fmri_img, "meta": meta}} + all_out = efc.fit_transform(input, storage=storage) + + features = storage.list_features() + assert any( + x["name"] == "BOLD_EdgeCentricFCParcels" for x in features.values() + ) diff --git a/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py new file mode 100644 index 000000000..6b9ad17f6 --- /dev/null +++ b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py @@ -0,0 +1,73 @@ +"""Provide tests for edge-centric functional connectivity using spheres.""" + +# Authors: Leonard Sasse +# Synchon Mandal +# License: AGPL + +from pathlib import Path + +from nilearn import datasets, image + +from junifer.markers.functional_connectivity import EdgeCentricFCSpheres +from junifer.storage import SQLiteFeatureStorage + + +def test_EdgeCentricFCSpheres(tmp_path: Path) -> None: + """Test EdgeCentricFCSpheres. + + Parameters + ---------- + tmp_path : pathlib.Path + The path to the test directory. + + """ + # get a dataset + ni_data = datasets.fetch_spm_auditory(subject_id="sub001") + fmri_img = image.concat_imgs(ni_data.func) # type: ignore + + efc = EdgeCentricFCSpheres( + coords="DMNBuckner", radius=5.0, cor_method="correlation" + ) + all_out = efc.fit_transform({"BOLD": {"data": fmri_img, "meta": {}}}) + + out = all_out["BOLD"] + + # There are six DMNBuckner coordinates, so + # for 6 ROIs we should get (6 * (6 -1) / 2) edges in the ETS + n_edges = int(6 * (6 - 1) / 2) + assert "data" in out + assert "row_names" in out + assert "col_names" in out + assert out["data"].shape[0] == n_edges + assert out["data"].shape[1] == n_edges + assert len(set(out["row_names"])) == n_edges + assert len(set(out["col_names"])) == n_edges + + # check correct output + assert efc.get_output_type("BOLD") == "matrix" + + # Check empirical correlation method parameters + efc = EdgeCentricFCSpheres( + coords="DMNBuckner", + radius=5.0, + cor_method="correlation", + cor_method_params={"empirical": True}, + ) + + meta = { + "element": {"subject": "sub001"}, + "dependencies": {"nilearn"}, + } + all_out = efc.fit_transform({"BOLD": {"data": fmri_img, "meta": meta}}) + + uri = tmp_path / "test_fc_parcellation.sqlite" + # Single storage, must be the uri + storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") + meta = {"element": {"subject": "test"}, "dependencies": {"numpy"}} + input = {"BOLD": {"data": fmri_img, "meta": meta}} + all_out = efc.fit_transform(input, storage=storage) + + features = storage.list_features() + assert any( + x["name"] == "BOLD_EdgeCentricFCSpheres" for x in features.values() + ) diff --git a/junifer/markers/tests/test_marker_utils.py b/junifer/markers/tests/test_marker_utils.py index 67fc9b7cb..bc4b63d3f 100644 --- a/junifer/markers/tests/test_marker_utils.py +++ b/junifer/markers/tests/test_marker_utils.py @@ -7,6 +7,7 @@ # License: AGPL import numpy as np +import pytest from junifer.markers.utils import _ets @@ -25,5 +26,24 @@ def test_ets() -> None: n_time, n_rois = bold_ts.shape n_edges = int(n_rois * (n_rois - 1) / 2) + # test without labels edge_ts = _ets(bold_ts) assert edge_ts.shape == (n_time, n_edges) + + # test with labels + roi_labels = [f"Label_{x}" for x in range(n_rois)] + edge_ts, edge_labels = _ets(bold_ts, roi_labels) + assert edge_ts.shape == (n_time, n_edges) + assert len(edge_labels) == n_edges + + +def test_ets_incorrect_label_list() -> None: + """Test edge-wise timeseries computing function with incorrect labels.""" + bold_ts = np.arange(30).reshape(10, 3) + # one label is missing, the bold time series suggests three ROIs + roi_labels = ["label 1", "label 2"] + + with pytest.raises( + ValueError, match="List of roi names does not correspond" + ): + _ets(bold_ts, roi_labels) diff --git a/junifer/markers/utils.py b/junifer/markers/utils.py index 4c009b5c7..270c582d4 100644 --- a/junifer/markers/utils.py +++ b/junifer/markers/utils.py @@ -7,7 +7,7 @@ # Federico Raimondo # License: AGPL -from typing import Any, Callable, Dict, Type, Union +from typing import Any, Callable, Dict, List, Type, Union import numpy as np import pandas as pd @@ -55,7 +55,9 @@ def singleton(cls: Type) -> Type: return get_instance -def _ets(bold_ts: np.ndarray) -> np.ndarray: +def _ets( + bold_ts: np.ndarray, roi_names: Union[None, List[str]] = None +) -> np.ndarray: """Compute the edge-wise time series based on BOLD time series. Take a timeseries of brain areas, and calculate timeseries for each @@ -66,12 +68,21 @@ def _ets(bold_ts: np.ndarray) -> np.ndarray: ---------- bold_ts : np.ndarray BOLD time series (time x ROIs) + roi_names : List[str] or None + List containing the names of the ROIs. + Order of the ROI names should correspond to order of the columns + in bold_ts. If None (default), only the edge-wise time series are + returned, without corresponding edge labels. Returns ------- - np.ndarray + ets : np.ndarray edge-wise time series, i.e. estimate of functional connectivity at each time point. + edge_names : List[str] + List of edge names corresponding to columns in the + edge-wise time series. This is only returned if the roi_names + are specified. References ---------- @@ -88,7 +99,21 @@ def _ets(bold_ts: np.ndarray) -> np.ndarray: # indices of unique edges (lower triangle) u, v = np.tril_indices(n_roi, k=-1) # Compute the ETS - return timeseries[:, u] * timeseries[:, v] + ets = timeseries[:, u] * timeseries[:, v] + # Obtain the corresponding edge labels if specified else return + if roi_names is None: + return ets + else: + if len(roi_names) != n_roi: + raise_error( + "List of roi names does not correspond " + "to the number of ROIs in the timeseries!" + ) + roi_names = np.array(roi_names) + edge_names = [ + "~".join([x, y]) for x, y in zip(roi_names[u], roi_names[v]) + ] + return ets, list(edge_names) def _correlate_dataframes(