[ENH] start bidsifying output

docs/source/changes.md

@@ -12,6 +12,8 @@ Released MONTH YEAR
 
 ### Changes
 
+- [ENH] Make output more BIDS compliant. (@Remi-Gau)
+
 ## 0.4.0
 
 Released August 2023

docs/source/outputs.md

@@ -1,7 +1,79 @@
 # Outputs
 
+The output of this app aims to follow the guideline
+of the [BIDS extension proposal 17 - Generic BIDS connectivity data schema](https://bids.neuroimaging.io/bep017).
+
+Metadata files content is described in this BIDS extension proposal.
+
+## participant level
+
 When performing `participant` level analysis, the output is a HDF5 file per participant that was passed to `--participant_label` or all subjects under `bids_dir`.
-The output file name is: `sub-<participant_id>_atlas-<atlas_name>_desc-<denoising_strategy>.h5`
 
-When performing `group` level analysis, the file will contain time series and connectomes of each subject, as well as group average connectomes. The output is a HDF5 file per participant that was passed to `--participant_label` or all subjects under `bids_dir`.
-The output file name is: `atlas-<atlas_name>_desc-<denoising_strategy>.h5`
+The output file name is:
+
+```
+sub-<participant_id>/func/sub-<participant_id>_atlas-<atlas_name>_meas-PearsonCorrelation_desc-<denoising_strategy>_relmat.h5
+```
+
+## group level
+
+When performing `group` level analysis,
+the file will contain time series and connectomes of each subject,
+as well as group average connectomes.
+
+The output is a HDF5 file per participant that was passed to `--participant_label` or all subjects under `bids_dir`.
+
+The output file name is: `group/atlas-<atlas_name>_meas-PearsonCorrelation_desc-<denoising_strategy>.h5`
+
+## Examples
+
+```
+├── dataset_description.json
+├── meas-PearsonCorrelation_desc-simple_relmat.json
+├── group
+│   └── atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+├── sub-1
+│   └── func
+│       └── sub-1_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+└── sub-2
+    └── func
+        └── sub-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+```
+
+### BIDS compliant output
+
+If the `--output_to_bids` is passed to the command line then each input file will have associated output file.
+
+The output file name is:
+
+```
+sub-<participant_id>/[ses-<ses_id>]/func/<source_filename>_atlas-<atlas_name>_meas-PearsonCorrelation_desc-<denoising_strategy>_relmat.h5
+```
+
+Example:
+
+```
+giga_connectome/data/test_data/output
+├── dataset_description.json
+├── meas-PearsonCorrelation_desc-simple_relmat.json
+├── group
+│   └── atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+├── sub-1
+│   ├── ses-timepoint1
+│   │   └── func
+│   │       ├── sub-1_ses-timepoint1_task-probabilisticclassification_run-01_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+│   │       └── sub-1_ses-timepoint1_task-probabilisticclassification_run-02_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+│   └── ses-timepoint2
+│       └── func
+│           ├── sub-1_ses-timepoint2_task-probabilisticclassification_run-01_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+│           └── sub-1_ses-timepoint2_task-probabilisticclassification_run-02_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+└── sub-2
+    ├── ses-timepoint1
+    │   └── func
+    │       ├── sub-2_ses-timepoint1_task-probabilisticclassification_run-01_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+    │       └── sub-2_ses-timepoint1_task-probabilisticclassification_run-02_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+    └── ses-timepoint2
+        └── func
+            ├── sub-2_ses-timepoint2_task-probabilisticclassification_run-01_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+            └── sub-2_ses-timepoint2_task-probabilisticclassification_run-02_space-MNI152NLin2009cAsym_res-2_atlas-Schaefer20187Networks_meas-PearsonCorrelation_desc-simple_relmat.h5
+```

docs/source/usage.md

@@ -67,3 +67,29 @@ In a `json` file, define the customised atlas. We will use the atlas above as an
 ```
 
 See examples in [`giga_connectome/data`](https://github.com/SIMEXP/giga_connectome/tree/main/giga_connectome/data).
+
+
+## Example
+
+You can run a demo of the bids app by downloading some test data.
+
+Run the following from the root of the repository.
+
+```bash
+pip install tox
+tox -e test_data
+```
+
+```bash
+giga_connectome \
+    --atlas Schaefer20187Networks \
+    --denoise-strategy simple \
+    --standardize zscore \
+    --bids-filter giga_connectome/data/test_data/bids_filter.json  \
+    --reindex-bids \
+    --output-to-bids \
+    --calculate-intranetwork-average-correlation \
+    giga_connectome/data/test_data/ds000017-fmriprep22.0.1-downsampled-nosurface \
+    giga_connectome/data/test_data/output \
+    participant
+```

giga_connectome/atlas.py

@@ -17,7 +17,7 @@
 PRESET_ATLAS = ["DiFuMo", "MIST", "Schaefer20187Networks"]
 
 
-def load_atlas_setting(atlas: Union[str, Path, dict]):
+def load_atlas_setting(atlas: Union[str, Path, dict]) -> dict:
     """Load atlas details for templateflow api to fetch.
     The setting file can be configured for atlases not included in the
     templateflow collections, but user has to organise their files to

giga_connectome/postprocess.py

@@ -18,6 +18,7 @@
 
 def run_postprocessing_dataset(
     strategy: dict,
+    atlas: str,
     resampled_atlases: List[Union[str, Path]],
     images: List[BIDSImageFile],
     group_mask: Union[str, Path],
@@ -26,6 +27,7 @@ def run_postprocessing_dataset(
     output_path: Path,
     analysis_level: str,
     calculate_average_correlation: bool = False,
+    output_to_bids: bool = False,
 ) -> None:
     """
     Generate subject and group level timeseries and connectomes.
@@ -61,6 +63,9 @@ def run_postprocessing_dataset(
     strategy : dict
         Parameters for `load_confounds_strategy` or `load_confounds`.
 
+    atlas : dict
+        Atlas settings.
+
     resampled_atlases : list of str or pathlib.Path
         Atlas niftis resampled to the common space of the dataset.
 
@@ -78,16 +83,14 @@ def run_postprocessing_dataset(
         Smoothing kernel size, passed to nilearn masker.
 
     output_path:
-        Full path to output file, named in the following format:
-        output_dir / atlas-<atlas>_desc-<strategy_name>.h5
+        Full path to output directory.
 
     analysis_level : str
         Level of analysis, either "participant" or "group".
 
     calculate_average_correlation : bool
         Whether to calculate average correlation within each parcel.
     """
-    atlas = output_path.name.split("atlas-")[-1].split("_")[0]
     atlas_maskers, connectomes = {}, {}
     for atlas_path in resampled_atlases:
         if isinstance(atlas_path, str):
@@ -103,8 +106,19 @@ def run_postprocessing_dataset(
     # transform data
     gc_log.info("Processing subject")
 
-    for img in tqdm(images):
+    # single output file per subject when output is not BIDS compliant
+    if not output_to_bids:
+        subject, _, _ = utils.parse_bids_name(images[0].path)
+        filename = utils.output_filename(
+            Path(images[0].filename).stem,
+            atlas["name"],
+            strategy["name"],
+            output_to_bids,
+        )
+        connectome_path = output_path / subject / "func" / filename
+        utils.check_path(connectome_path)
 
+    for img in tqdm(images):
         print()
         gc_log.info(f"Processing image:\n{img.filename}")
 
@@ -114,8 +128,25 @@ def run_postprocessing_dataset(
         )
         # parse file name
         subject, session, specifier = utils.parse_bids_name(img.path)
+
+        # one output file per input file when output is BIDS compliant
+        if output_to_bids:
+            connectome_path = output_path / subject
+            if session:
+                connectome_path = connectome_path / session
+            filename = utils.output_filename(
+                Path(img.filename).stem,
+                atlas["name"],
+                strategy["name"],
+                output_to_bids,
+            )
+            connectome_path = connectome_path / "func" / filename
+            utils.check_path(connectome_path)
+
         for desc, masker in atlas_maskers.items():
-            attribute_name = f"{subject}_{specifier}_atlas-{atlas}_desc-{desc}"
+            attribute_name = (
+                f"{subject}_{specifier}_atlas-{atlas['name']}_desc-{desc}"
+            )
             if not denoised_img:
                 time_series_atlas, correlation_matrix = None, None
 
@@ -137,8 +168,8 @@ def run_postprocessing_dataset(
             connectomes[desc].append(correlation_matrix)
 
             # dump to h5
-            flag = _set_file_flag(output_path)
-            with h5py.File(output_path, flag) as f:
+            flag = _set_file_flag(connectome_path)
+            with h5py.File(connectome_path, flag) as f:
                 group = _fetch_h5_group(f, subject, session)
                 timeseries_dset = group.create_dataset(
                     f"{attribute_name}_timeseries", data=time_series_atlas
@@ -150,22 +181,36 @@ def run_postprocessing_dataset(
                     f"{attribute_name}_connectome", data=correlation_matrix
                 )
 
-        gc_log.info(f"Saved to:\n{output_path}")
+        gc_log.info(f"Saved to:\n{connectome_path}")
 
     if analysis_level == "group":
+        connectome_path = (
+            output_path
+            / "group"
+            / utils.output_filename(
+                "",
+                atlas["name"],
+                strategy["name"],
+                output_to_bids,
+            )
+        )
+        utils.check_path(connectome_path)
 
         gc_log.info("Create group connectome")
+        gc_log.info(connectome_path)
 
         for desc in connectomes:
             average_connectome = np.mean(
                 np.array(connectomes[desc]), axis=0
             ).astype(np.float32)
-            with h5py.File(output_path, "a") as f:
+            with h5py.File(connectome_path, "a") as f:
                 f.create_dataset(
-                    f"atlas-{atlas}_desc-{desc}_connectome",
+                    f"atlas-{atlas['name']}_desc-{desc}_connectome",
                     data=average_connectome,
                 )
 
+        gc_log.info(f"Saved to:\n{connectome_path}")
+
 
 def _set_file_flag(output_path: Path) -> str:
     """Find out if new file needs to be created."""

giga_connectome/run.py

@@ -104,6 +104,14 @@ def global_parser() -> argparse.ArgumentParser:
         "pipeline (option A). The default is False.",
         action="store_true",
     )
+    parser.add_argument(
+        "--output-to-bids",
+        help="If requested the output will be more BIDS compliant. "
+        "This means that there will be one output file, for each input file."
+        "Otherwise all subject output will be saved into a single file."
+        "The default is False.",
+        action="store_true",
+    )
     parser.add_argument(
         "--verbosity",
         help="""