From 9af82a703cafb8c275b6e6054e534a07d9657e34 Mon Sep 17 00:00:00 2001 From: Taylor Salo Date: Tue, 23 Jul 2024 10:00:58 -0400 Subject: [PATCH] Plot CIFTI maps on subject-specific surfaces when available (#1208) --- docs/outputs.rst | 4 + xcp_d/data/reports-spec-func.yml | 74 +++++++++++- xcp_d/data/reports-spec.yml | 74 +++++++++++- xcp_d/interfaces/plotting.py | 191 +++++++++++++++++++++---------- xcp_d/workflows/anatomical.py | 21 ++++ xcp_d/workflows/base.py | 8 ++ xcp_d/workflows/cifti.py | 13 +++ xcp_d/workflows/connectivity.py | 60 ++++++++-- xcp_d/workflows/restingstate.py | 95 +++++++++++---- 9 files changed, 435 insertions(+), 105 deletions(-) diff --git a/docs/outputs.rst b/docs/outputs.rst index 78ac50126..034d178b8 100644 --- a/docs/outputs.rst +++ b/docs/outputs.rst @@ -125,6 +125,10 @@ Surface mesh files If the ``--warp-surfaces-native2std`` option is selected, and reconstructed surfaces are available in the preprocessed dataset, then these surfaces will be warped to fsLR space at 32k density. +The resulting mesh files will reflect the subject's morphology with the same geometry and density +as fsLR-32k surfaces, which may be useful for visualizing fsLR-space derivatives on a subject's +brain. + .. code-block:: xcp_d/ diff --git a/xcp_d/data/reports-spec-func.yml b/xcp_d/data/reports-spec-func.yml index 26908ec5a..641903eb8 100755 --- a/xcp_d/data/reports-spec-func.yml +++ b/xcp_d/data/reports-spec-func.yml @@ -30,7 +30,22 @@ sections: datatype: figures desc: coverage suffix: bold - caption: Parcellation coverage. + caption: | + Parcellation coverage. + subtitle: Coverage + - bids: + datatype: figures + desc: coverageParcellatedStandard + suffix: bold + caption: | + Parcellation coverage, overlaid on the fsLR template surface. + subtitle: Coverage + - bids: + datatype: figures + desc: coverageParcellatedSubject + suffix: bold + caption: | + Parcellation coverage, overlaid on the subject's surface warped to the fsLR template. subtitle: Coverage - bids: {datatype: figures, desc: connectivityplot, suffix: bold} caption: This plot shows heatmaps from ROI-to-ROI correlations from one to four atlases. @@ -44,9 +59,18 @@ sections: subtitle: ALFF - bids: datatype: figures - desc: alffSurfacePlot + desc: alffSurfaceStandard suffix: bold - caption: ALFF, or amplitude of low frequency fluctuations. + caption: | + ALFF, or amplitude of low frequency fluctuations, overlaid on the fsLR template surface. + subtitle: ALFF + - bids: + datatype: figures + desc: alffSurfaceSubject + suffix: bold + caption: | + ALFF, or amplitude of low frequency fluctuations, overlaid on the subject's surface + warped to the fsLR template. subtitle: ALFF - bids: datatype: figures @@ -54,6 +78,21 @@ sections: suffix: bold caption: Parcellated ALFF. subtitle: ALFF + - bids: + datatype: figures + desc: alffParcellatedStandard + suffix: bold + caption: | + ALFF, or amplitude of low frequency fluctuations, overlaid on the fsLR template surface. + subtitle: ALFF + - bids: + datatype: figures + desc: alffParcellatedSubject + suffix: bold + caption: | + ALFF, or amplitude of low frequency fluctuations, overlaid on the subject's surface + warped to the fsLR template. + subtitle: ALFF - bids: datatype: figures desc: rehoVolumetricPlot @@ -63,9 +102,19 @@ sections: subtitle: ReHo - bids: datatype: figures - desc: rehoSurfacePlot + desc: rehoSurfaceStandard suffix: bold - caption: ReHo, or regional homogeneity. + caption: | + ReHo, or regional homogeneity. + Overlaid on the fsLR midthickness surface. + subtitle: ReHo + - bids: + datatype: figures + desc: rehoSurfaceSubject + suffix: bold + caption: | + ReHo, or regional homogeneity. + Overlaid on the participant's midthickness surface mapped onto fsLR space. subtitle: ReHo - bids: datatype: figures @@ -73,6 +122,21 @@ sections: suffix: bold caption: Parcellated ReHo. subtitle: ReHo + - bids: + datatype: figures + desc: rehoParcellatedStandard + suffix: bold + caption: | + ReHo, or regional homogeneity, overlaid on the fsLR template surface. + subtitle: ReHo + - bids: + datatype: figures + desc: rehoParcellatedSubject + suffix: bold + caption: | + ReHo, or regional homogeneity, overlaid on the subject's surface + warped to the fsLR template. + subtitle: ReHo - name: About reportlets: - bids: {datatype: figures, desc: about, suffix: bold} diff --git a/xcp_d/data/reports-spec.yml b/xcp_d/data/reports-spec.yml index 6a05dad30..70a09d1a3 100755 --- a/xcp_d/data/reports-spec.yml +++ b/xcp_d/data/reports-spec.yml @@ -40,6 +40,20 @@ sections: suffix: bold caption: Parcellation coverage. subtitle: Coverage + - bids: + datatype: figures + desc: coverageParcellatedStandard + suffix: bold + caption: | + Parcellation coverage, overlaid on the fsLR template surface. + subtitle: Coverage + - bids: + datatype: figures + desc: coverageParcellatedSubject + suffix: bold + caption: | + Parcellation coverage, overlaid on the subject's surface warped to the fsLR template. + subtitle: Coverage - bids: {datatype: figures, desc: connectivityplot, suffix: bold} caption: This plot shows heatmaps from ROI-to-ROI correlations from one to four atlases. subtitle: Correlation Heatmaps from A Subset of Atlases @@ -48,13 +62,23 @@ sections: desc: alffVolumetricPlot suffix: bold caption: ALFF, or amplitude of low frequency fluctuations. - Overlaid on T1W image with same entities as the original image. + Overlaid on T1w image with same entities as the original image. subtitle: ALFF - bids: datatype: figures - desc: alffSurfacePlot + desc: alffSurfaceStandard suffix: bold - caption: ALFF, or amplitude of low frequency fluctuations. + caption: | + ALFF, or amplitude of low frequency fluctuations. + Overlaid on the fsLR midthickness surface. + subtitle: ALFF + - bids: + datatype: figures + desc: alffSurfaceSubject + suffix: bold + caption: | + ALFF, or amplitude of low frequency fluctuations. + Overlaid on the participant's midthickness surface mapped onto fsLR space. subtitle: ALFF - bids: datatype: figures @@ -62,6 +86,21 @@ sections: suffix: bold caption: Parcellated ALFF. subtitle: ALFF + - bids: + datatype: figures + desc: alffParcellatedStandard + suffix: bold + caption: | + ALFF, or amplitude of low frequency fluctuations, overlaid on the fsLR template surface. + subtitle: ALFF + - bids: + datatype: figures + desc: alffParcellatedSubject + suffix: bold + caption: | + ALFF, or amplitude of low frequency fluctuations, overlaid on the subject's surface + warped to the fsLR template. + subtitle: ALFF - bids: datatype: figures desc: rehoVolumetricPlot @@ -71,9 +110,19 @@ sections: subtitle: ReHo - bids: datatype: figures - desc: rehoSurfacePlot + desc: rehoSurfaceStandard suffix: bold - caption: ReHo, or regional homogeneity. + caption: | + ReHo, or regional homogeneity. + Overlaid on the fsLR midthickness surface. + subtitle: ReHo + - bids: + datatype: figures + desc: rehoSurfaceSubject + suffix: bold + caption: | + ReHo, or regional homogeneity. + Overlaid on the participant's midthickness surface mapped onto fsLR space. subtitle: ReHo - bids: datatype: figures @@ -81,6 +130,21 @@ sections: suffix: bold caption: Parcellated ReHo. subtitle: ReHo + - bids: + datatype: figures + desc: rehoParcellatedStandard + suffix: bold + caption: | + ReHo, or regional homogeneity, overlaid on the fsLR template surface. + subtitle: ReHo + - bids: + datatype: figures + desc: rehoParcellatedSubject + suffix: bold + caption: | + ReHo, or regional homogeneity, overlaid on the subject's surface + warped to the fsLR template. + subtitle: ReHo - name: About reportlets: - bids: {datatype: figures, desc: about, suffix: bold} diff --git a/xcp_d/interfaces/plotting.py b/xcp_d/interfaces/plotting.py index e9c364f77..7abf8b8ce 100644 --- a/xcp_d/interfaces/plotting.py +++ b/xcp_d/interfaces/plotting.py @@ -905,10 +905,27 @@ class _PlotCiftiParcellationInputSpec(BaseInterfaceInputSpec): usedefault=True, desc="Maximum value for the colormap.", ) + base_desc = traits.Str( + mandatory=False, + default_value="", + usedefault=True, + desc="Base description for the output file.", + ) + lh_underlay = File( + exists=True, + mandatory=False, + desc="Left hemisphere underlay.", + ) + rh_underlay = File( + exists=True, + mandatory=False, + desc="Right hemisphere underlay.", + ) class _PlotCiftiParcellationOutputSpec(TraitedSpec): out_file = File(exists=True, desc="Output file.") + desc = traits.Str(desc="Description of the output file.") class PlotCiftiParcellation(SimpleInterface): @@ -921,24 +938,30 @@ def _run_interface(self, runtime): assert len(self.inputs.in_files) == len(self.inputs.labels) assert len(self.inputs.cortical_atlases) > 0 - rh = str( - get_template( - template="fsLR", - hemi="R", - density="32k", - suffix="midthickness", - extension=".surf.gii", + if not (isdefined(self.inputs.lh_underlay) and isdefined(self.inputs.rh_underlay)): + self._results["desc"] = f"{self.inputs.base_desc}ParcellatedStandard" + rh = str( + get_template( + template="fsLR", + hemi="R", + density="32k", + suffix="midthickness", + extension=".surf.gii", + ) ) - ) - lh = str( - get_template( - template="fsLR", - hemi="L", - density="32k", - suffix="midthickness", - extension=".surf.gii", + lh = str( + get_template( + template="fsLR", + hemi="L", + density="32k", + suffix="midthickness", + extension=".surf.gii", + ) ) - ) + else: + self._results["desc"] = f"{self.inputs.base_desc}ParcellatedSubject" + rh = self.inputs.rh_underlay + lh = self.inputs.lh_underlay # Create Figure and GridSpec. # One subplot for each file. Each file will then have four subplots, arranged in a square. @@ -973,13 +996,17 @@ def _run_interface(self, runtime): subplot.set_axis_off() vmin, vmax = self.inputs.vmin, self.inputs.vmax + threshold = 0.01 if vmin == vmax: + threshold = None + # Define vmin and vmax based on all of the files vmin, vmax = np.inf, -np.inf for cortical_file in cortical_files: img_data = nb.load(cortical_file).get_fdata() vmin = np.min([np.nanmin(img_data), vmin]) vmax = np.max([np.nanmax(img_data), vmax]) + vmin = 0 for i_file in range(n_files): subplot = subplots[i_file] @@ -1011,6 +1038,7 @@ def _run_interface(self, runtime): plot_surf_stat_map( lh, lh_surf_data, + threshold=threshold, vmin=vmin, vmax=vmax, hemi="left", @@ -1024,6 +1052,7 @@ def _run_interface(self, runtime): plot_surf_stat_map( rh, rh_surf_data, + threshold=threshold, vmin=vmin, vmax=vmax, hemi="right", @@ -1037,6 +1066,7 @@ def _run_interface(self, runtime): plot_surf_stat_map( lh, lh_surf_data, + threshold=threshold, vmin=vmin, vmax=vmax, hemi="left", @@ -1050,6 +1080,7 @@ def _run_interface(self, runtime): plot_surf_stat_map( rh, rh_surf_data, + threshold=threshold, vmin=vmin, vmax=vmax, hemi="right", @@ -1095,15 +1126,27 @@ class _PlotDenseCiftiInputSpec(BaseInterfaceInputSpec): mandatory=True, desc="CIFTI file to plot.", ) - name_source = File( - exists=False, - mandatory=True, - desc="File to use as the name source. Unused but retained for compatibility.", + base_desc = traits.Str( + mandatory=False, + default_value="", + usedefault=True, + desc="Base description for the output file.", + ) + lh_underlay = File( + exists=True, + mandatory=False, + desc="Left hemisphere underlay.", + ) + rh_underlay = File( + exists=True, + mandatory=False, + desc="Right hemisphere underlay.", ) class _PlotDenseCiftiOutputSpec(TraitedSpec): out_file = File(exists=True, desc="Output file.") + desc = traits.Str(desc="Description of the output file.") class PlotDenseCifti(SimpleInterface): @@ -1113,30 +1156,52 @@ class PlotDenseCifti(SimpleInterface): output_spec = _PlotDenseCiftiOutputSpec def _run_interface(self, runtime): - rh = str( - get_template( - template="fsLR", - hemi="R", - density="32k", - suffix="midthickness", - extension=".surf.gii", + if not (isdefined(self.inputs.lh_underlay) and isdefined(self.inputs.rh_underlay)): + self._results["desc"] = f"{self.inputs.base_desc}SurfaceStandard" + rh = str( + get_template( + template="fsLR", + hemi="R", + density="32k", + suffix="midthickness", + extension=".surf.gii", + ) ) - ) - lh = str( - get_template( - template="fsLR", - hemi="L", - density="32k", - suffix="midthickness", - extension=".surf.gii", + lh = str( + get_template( + template="fsLR", + hemi="L", + density="32k", + suffix="midthickness", + extension=".surf.gii", + ) ) - ) + else: + self._results["desc"] = f"{self.inputs.base_desc}SurfaceSubject" + rh = self.inputs.rh_underlay + lh = self.inputs.lh_underlay cifti = nb.load(self.inputs.in_file) cifti_data = cifti.get_fdata() cifti_axes = [cifti.header.get_axis(i) for i in range(cifti.ndim)] - fig, axes = plt.subplots(figsize=(4, 4), ncols=2, nrows=2, subplot_kw={"projection": "3d"}) + # Create Figure and GridSpec. + fig = plt.figure(constrained_layout=False) + fig.set_size_inches(6.5, 6) + # Add an additional column for the colorbar + gs = GridSpec(1, 2, figure=fig, width_ratios=[1, 0.05]) + subplot_gridspec = gs[0, 0] + subplot = fig.add_subplot(subplot_gridspec) + colorbar_gridspec = gs[0, 1] + + subplot.set_axis_off() + + # Create 4 Axes (2 rows, 2 columns) from the subplot + gs_inner = GridSpecFromSubplotSpec(2, 2, subplot_spec=subplot_gridspec) + inner_subplots = [ + fig.add_subplot(gs_inner[i, j], projection="3d") for i in range(2) for j in range(2) + ] + lh_surf_data = surf_data_from_cifti( cifti_data, cifti_axes[1], @@ -1148,13 +1213,8 @@ def _run_interface(self, runtime): "CIFTI_STRUCTURE_CORTEX_RIGHT", ) - vmax = np.max([np.max(lh_surf_data), np.max(rh_surf_data)]) - vmin = np.min([np.min(lh_surf_data), np.min(rh_surf_data)]) - - axes[0, 0].set_rasterized(True) - axes[0, 1].set_rasterized(True) - axes[1, 0].set_rasterized(True) - axes[1, 1].set_rasterized(True) + vmax = np.nanmax([np.nanmax(lh_surf_data), np.nanmax(rh_surf_data)]) + vmin = np.nanmin([np.nanmin(lh_surf_data), np.nanmin(rh_surf_data)]) plot_surf_stat_map( lh, @@ -1164,32 +1224,35 @@ def _run_interface(self, runtime): hemi="left", view="lateral", engine="matplotlib", + cmap="cool", colorbar=False, - axes=axes[0, 0], + axes=inner_subplots[0], figure=fig, ) plot_surf_stat_map( - lh, - lh_surf_data, + rh, + rh_surf_data, vmin=vmin, vmax=vmax, - hemi="left", - view="medial", + hemi="right", + view="lateral", engine="matplotlib", + cmap="cool", colorbar=False, - axes=axes[1, 0], + axes=inner_subplots[1], figure=fig, ) plot_surf_stat_map( - rh, - rh_surf_data, + lh, + lh_surf_data, vmin=vmin, vmax=vmax, - hemi="right", - view="lateral", + hemi="left", + view="medial", engine="matplotlib", + cmap="cool", colorbar=False, - axes=axes[0, 1], + axes=inner_subplots[2], figure=fig, ) plot_surf_stat_map( @@ -1200,12 +1263,24 @@ def _run_interface(self, runtime): hemi="right", view="medial", engine="matplotlib", + cmap="cool", colorbar=False, - axes=axes[1, 1], + axes=inner_subplots[3], figure=fig, ) - axes[0, 0].set_title("Left Hemisphere", fontsize=10) - axes[0, 1].set_title("Right Hemisphere", fontsize=10) + + inner_subplots[0].set_title("Left Hemisphere", fontsize=10) + inner_subplots[1].set_title("Right Hemisphere", fontsize=10) + + for ax in inner_subplots: + ax.set_rasterized(True) + + # Create a ScalarMappable with the "cool" colormap and the specified vmin and vmax + sm = ScalarMappable(cmap="cool", norm=Normalize(vmin=vmin, vmax=vmax)) + + colorbar_ax = fig.add_subplot(colorbar_gridspec) + # Add a colorbar to colorbar_ax using the ScalarMappable + fig.colorbar(sm, cax=colorbar_ax) self._results["out_file"] = fname_presuffix( self.inputs.in_file, diff --git a/xcp_d/workflows/anatomical.py b/xcp_d/workflows/anatomical.py index 77f3c1c59..0f0cc5275 100644 --- a/xcp_d/workflows/anatomical.py +++ b/xcp_d/workflows/anatomical.py @@ -305,6 +305,8 @@ def init_postprocess_surfaces_wf( fsnative space, this workflow will warp them to fsLR space. If process-surfaces is enabled and the mesh files are already in fsLR space, they will be copied to the output directory. + These fsLR-space mesh files retain the subject's morphology, + and are thus useful for visualizing fsLR-space statistical derivatives on the subject's brain. As long as process-surfaces is enabled and mesh files (in either space) are available, HCP-style midthickness, inflated, and very-inflated surfaces will be generated from them. @@ -394,6 +396,17 @@ def init_postprocess_surfaces_wf( ), name="inputnode", ) + outputnode = pe.Node( + niu.IdentityInterface( + fields=[ + "lh_midthickness", + "rh_midthickness", + ], + ), + name="outputnode", + ) + workflow.add_nodes([outputnode]) # outputnode may not be used + workflow.__desc__ = "" if abcc_qc and mesh_available: @@ -458,6 +471,8 @@ def init_postprocess_surfaces_wf( workflow.connect([ (inputnode, hcp_surface_wfs["lh"], [("lh_pial_surf", "inputnode.name_source")]), (inputnode, hcp_surface_wfs["rh"], [("rh_pial_surf", "inputnode.name_source")]), + (hcp_surface_wfs["lh"], outputnode, [("outputnode.midthickness", "lh_midthickness")]), + (hcp_surface_wfs["rh"], outputnode, [("outputnode.midthickness", "rh_midthickness")]), ]) # fmt:skip if mesh_available and standard_space_mesh: @@ -761,6 +776,11 @@ def init_generate_hcp_surfaces_wf(name="generate_hcp_surfaces_wf"): name="inputnode", ) + outputnode = pe.Node( + niu.IdentityInterface(fields=["midthickness"]), + name="outputnode", + ) + generate_midthickness = pe.Node( SurfaceAverage(), name="generate_midthickness", @@ -772,6 +792,7 @@ def init_generate_hcp_surfaces_wf(name="generate_hcp_surfaces_wf"): ("pial_surf", "surface_in1"), ("wm_surf", "surface_in2"), ]), + (generate_midthickness, outputnode, [("out_file", "midthickness")]), ]) # fmt:skip ds_midthickness = pe.Node( diff --git a/xcp_d/workflows/base.py b/xcp_d/workflows/base.py index bdd2d9f9c..e55055f0a 100644 --- a/xcp_d/workflows/base.py +++ b/xcp_d/workflows/base.py @@ -481,6 +481,14 @@ def init_single_subject_wf(subject_id: str): ]), ]) # fmt:skip + if config.workflow.process_surfaces or (config.workflow.abcc_qc and mesh_available): + workflow.connect([ + (postprocess_surfaces_wf, postprocess_bold_wf, [ + ("outputnode.lh_midthickness", "inputnode.lh_midthickness"), + ("outputnode.rh_midthickness", "inputnode.rh_midthickness"), + ]), + ]) # fmt:skip + if config.execution.atlases: workflow.connect([ (load_atlases_wf, postprocess_bold_wf, [ diff --git a/xcp_d/workflows/cifti.py b/xcp_d/workflows/cifti.py index 6e87a217a..fdac8186b 100644 --- a/xcp_d/workflows/cifti.py +++ b/xcp_d/workflows/cifti.py @@ -154,6 +154,9 @@ def init_postprocess_cifti_wf( "atlases", "atlas_files", "atlas_labels_files", + # for plotting, if the anatomical workflow was used + "lh_midthickness", + "rh_midthickness", ], ), name="inputnode", @@ -280,6 +283,10 @@ def init_postprocess_cifti_wf( alff_wf = init_alff_wf(name_source=bold_file, TR=TR, mem_gb=mem_gbx) workflow.connect([ + (inputnode, alff_wf, [ + ("lh_midthickness", "inputnode.lh_midthickness"), + ("rh_midthickness", "inputnode.rh_midthickness"), + ]), (prepare_confounds_wf, alff_wf, [ ("outputnode.temporal_mask", "inputnode.temporal_mask"), ]), @@ -291,6 +298,10 @@ def init_postprocess_cifti_wf( reho_wf = init_reho_cifti_wf(name_source=bold_file, mem_gb=mem_gbx) workflow.connect([ + (inputnode, reho_wf, [ + ("lh_midthickness", "inputnode.lh_midthickness"), + ("rh_midthickness", "inputnode.rh_midthickness"), + ]), (denoise_bold_wf, reho_wf, [ ("outputnode.censored_denoised_bold", "inputnode.denoised_bold"), ]), @@ -373,6 +384,8 @@ def init_postprocess_cifti_wf( ("atlases", "inputnode.atlases"), ("atlas_files", "inputnode.atlas_files"), ("atlas_labels_files", "inputnode.atlas_labels_files"), + ("lh_midthickness", "inputnode.lh_midthickness"), + ("rh_midthickness", "inputnode.rh_midthickness"), ]), (prepare_confounds_wf, connectivity_wf, [ ("outputnode.temporal_mask", "inputnode.temporal_mask"), diff --git a/xcp_d/workflows/connectivity.py b/xcp_d/workflows/connectivity.py index a1007c5b6..1cf912f31 100644 --- a/xcp_d/workflows/connectivity.py +++ b/xcp_d/workflows/connectivity.py @@ -522,6 +522,10 @@ def init_functional_connectivity_nifti_wf(mem_gb, name="connectivity_wf"): def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivity_wf"): """Extract CIFTI time series. + This will parcellate the CIFTI file using the selected atlases and compute functional + connectivity between all regions for the selected atlases. + It will also parcellate ReHo and ALFF maps if they are provided. + Workflow Graph .. workflow:: :graph2use: orig @@ -610,6 +614,9 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit "atlases", "atlas_files", "atlas_labels_files", + # for plotting, if the anatomical workflow is enabled + "lh_midthickness", + "rh_midthickness", ], ), name="inputnode", @@ -655,19 +662,27 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit cortical_atlases = select_atlases(atlases=config.execution.atlases, subset="cortical") if cortical_atlases: plot_coverage = pe.Node( - PlotCiftiParcellation(cortical_atlases=cortical_atlases, vmin=0, vmax=1), + PlotCiftiParcellation( + base_desc="coverage", + cortical_atlases=cortical_atlases, + vmin=0, + vmax=1, + ), name="plot_coverage", mem_gb=mem_gb["resampled"], ) workflow.connect([ - (inputnode, plot_coverage, [("atlases", "labels")]), + (inputnode, plot_coverage, [ + ("atlases", "labels"), + ("lh_midthickness", "lh_underlay"), + ("rh_midthickness", "rh_underlay"), + ]), (parcellate_bold_wf, plot_coverage, [("outputnode.coverage_cifti", "in_files")]), ]) # fmt:skip ds_plot_coverage = pe.Node( DerivativesDataSink( base_directory=output_dir, - desc="coverage", datatype="figures", ), name="ds_plot_coverage", @@ -675,7 +690,10 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit ) workflow.connect([ (inputnode, ds_plot_coverage, [("name_source", "source_file")]), - (plot_coverage, ds_plot_coverage, [("out_file", "in_file")]), + (plot_coverage, ds_plot_coverage, [ + ("out_file", "in_file"), + ("desc", "desc"), + ]), ]) # fmt:skip # Reduce the CIFTI before calculating correlations @@ -829,12 +847,19 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit if cortical_atlases: plot_parcellated_reho = pe.Node( - PlotCiftiParcellation(cortical_atlases=cortical_atlases), + PlotCiftiParcellation( + base_desc="reho", + cortical_atlases=cortical_atlases, + ), name="plot_parcellated_reho", mem_gb=mem_gb["resampled"], ) workflow.connect([ - (inputnode, plot_parcellated_reho, [("atlases", "labels")]), + (inputnode, plot_parcellated_reho, [ + ("atlases", "labels"), + ("lh_midthickness", "lh_underlay"), + ("rh_midthickness", "rh_underlay"), + ]), (parcellate_reho_wf, plot_parcellated_reho, [ ("outputnode.parcellated_cifti", "in_files"), ]), @@ -843,7 +868,6 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit ds_plot_reho = pe.Node( DerivativesDataSink( base_directory=output_dir, - desc="rehoParcellated", datatype="figures", ), name="ds_plot_reho", @@ -851,7 +875,10 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit ) workflow.connect([ (inputnode, ds_plot_reho, [("name_source", "source_file")]), - (plot_parcellated_reho, ds_plot_reho, [("out_file", "in_file")]), + (plot_parcellated_reho, ds_plot_reho, [ + ("desc", "desc"), + ("out_file", "in_file"), + ]), ]) # fmt:skip if bandpass_filter: @@ -875,12 +902,19 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit if cortical_atlases: plot_parcellated_alff = pe.Node( - PlotCiftiParcellation(cortical_atlases=cortical_atlases), + PlotCiftiParcellation( + base_desc="alff", + cortical_atlases=cortical_atlases, + ), name="plot_parcellated_alff", mem_gb=mem_gb["resampled"], ) workflow.connect([ - (inputnode, plot_parcellated_alff, [("atlases", "labels")]), + (inputnode, plot_parcellated_alff, [ + ("atlases", "labels"), + ("lh_midthickness", "lh_underlay"), + ("rh_midthickness", "rh_underlay"), + ]), (parcellate_alff_wf, plot_parcellated_alff, [ ("outputnode.parcellated_cifti", "in_files"), ]), @@ -889,7 +923,6 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit ds_plot_alff = pe.Node( DerivativesDataSink( base_directory=output_dir, - desc="alffParcellated", datatype="figures", ), name="ds_plot_alff", @@ -897,7 +930,10 @@ def init_functional_connectivity_cifti_wf(mem_gb, exact_scans, name="connectivit ) workflow.connect([ (inputnode, ds_plot_alff, [("name_source", "source_file")]), - (plot_parcellated_alff, ds_plot_alff, [("out_file", "in_file")]), + (plot_parcellated_alff, ds_plot_alff, [ + ("out_file", "in_file"), + ("desc", "desc"), + ]), ]) # fmt:skip return workflow diff --git a/xcp_d/workflows/restingstate.py b/xcp_d/workflows/restingstate.py index 0fb0fae24..5e6db8ebb 100644 --- a/xcp_d/workflows/restingstate.py +++ b/xcp_d/workflows/restingstate.py @@ -83,6 +83,11 @@ def init_alff_wf( using a Lomb-Scargle periodogram :footcite:p:`lomb1976least,scargle1982studies,townsend2010fast,taylorlomb`. + This workflow will also generate a plot of the ALFF map. + For CIFTI data, the plot will be overlaid on the midthickness surface- + either the subject's surface warped to fsLR space (when the anatomical workflow is enabled) + or the fsLR 32k midthickness surface template. + References ---------- .. footbibliography:: @@ -117,11 +122,20 @@ def init_alff_wf( """ inputnode = pe.Node( - niu.IdentityInterface(fields=["denoised_bold", "bold_mask", "temporal_mask"]), + niu.IdentityInterface( + fields=[ + "denoised_bold", + "bold_mask", + "temporal_mask", + # only used for CIFTI data if the anatomical workflow is enabled + "lh_midthickness", + "rh_midthickness", + ], + ), name="inputnode", ) outputnode = pe.Node( - niu.IdentityInterface(fields=["alff", "smoothed_alff", "alffplot"]), + niu.IdentityInterface(fields=["alff", "smoothed_alff"]), name="outputnode", ) @@ -132,22 +146,50 @@ def init_alff_wf( name="alff_compt", n_procs=omp_nthreads, ) - - plot_interface = PlotDenseCifti if (file_format == "cifti") else PlotNifti - alff_plot = pe.Node( - plot_interface(name_source=name_source), - name="alff_plot", - ) workflow.connect([ (inputnode, alff_compt, [ ("denoised_bold", "in_file"), ("bold_mask", "mask"), ("temporal_mask", "temporal_mask"), ]), - (alff_compt, alff_plot, [("alff", "in_file")]), (alff_compt, outputnode, [("alff", "alff")]) ]) # fmt:skip + # Plot the ALFF map + ds_alff_plot = pe.Node( + DerivativesDataSink( + base_directory=output_dir, + source_file=name_source, + datatype="figures", + ), + name="ds_alff_plot", + run_without_submitting=False, + ) + + if file_format == "cifti": + alff_plot = pe.Node( + PlotDenseCifti(base_desc="alff"), + name="alff_plot", + ) + workflow.connect([ + (inputnode, alff_plot, [ + ("lh_midthickness", "lh_underlay"), + ("rh_midthickness", "rh_underlay"), + ]), + (alff_plot, ds_alff_plot, [("desc", "desc")]), + ]) # fmt:skip + else: + alff_plot = pe.Node( + PlotNifti(name_source=name_source), + name="alff_plot", + ) + ds_alff_plot.inputs.desc = "alffVolumetricPlot" + + workflow.connect([ + (alff_compt, alff_plot, [("alff", "in_file")]), + (alff_plot, ds_alff_plot, [("out_file", "in_file")]), + ]) # fmt:skip + if smoothing: # If we want to smooth if file_format == "nifti": workflow.__desc__ = workflow.__desc__ + ( @@ -206,18 +248,6 @@ def init_alff_wf( (fix_cifti_intent, outputnode, [("out_file", "smoothed_alff")]), ]) # fmt:skip - ds_alff_plot = pe.Node( - DerivativesDataSink( - base_directory=output_dir, - source_file=name_source, - desc="alffSurfacePlot" if (file_format == "cifti") else "alffVolumetricPlot", - datatype="figures", - ), - name="ds_alff_plot", - run_without_submitting=False, - ) - workflow.connect([(alff_plot, ds_alff_plot, [("out_file", "in_file")])]) - return workflow @@ -263,6 +293,13 @@ def init_reho_cifti_wf( ------- reho ReHo in a CIFTI file. + + Notes + ----- + This workflow will also generate a plot of the ReHo map. + The plot will be overlaid on the midthickness surface- + either the subject's surface warped to fsLR space (when the anatomical workflow is enabled) + or the fsLR 32k midthickness surface template. """ workflow = Workflow(name=name) workflow.__desc__ = """ @@ -279,7 +316,7 @@ def init_reho_cifti_wf( omp_nthreads = config.nipype.omp_nthreads inputnode = pe.Node( - niu.IdentityInterface(fields=["denoised_bold"]), + niu.IdentityInterface(fields=["denoised_bold", "lh_midthickness", "rh_midthickness"]), name="inputnode", ) outputnode = pe.Node( @@ -335,15 +372,20 @@ def init_reho_cifti_wf( n_procs=omp_nthreads, ) reho_plot = pe.Node( - PlotDenseCifti(name_source=name_source), + PlotDenseCifti(base_desc="reho"), name="reho_cifti_plot", ) + workflow.connect([ + (inputnode, reho_plot, [ + ("lh_midthickness", "lh_underlay"), + ("rh_midthickness", "rh_underlay"), + ]), + ]) # fmt:skip ds_reho_plot = pe.Node( DerivativesDataSink( base_directory=output_dir, source_file=name_source, - desc="rehoSurfacePlot", datatype="figures", ), name="ds_reho_plot", @@ -364,7 +406,10 @@ def init_reho_cifti_wf( (subcortical_reho, merge_cifti, [("out_file", "volume_all")]), (merge_cifti, outputnode, [("out_file", "reho")]), (merge_cifti, reho_plot, [("out_file", "in_file")]), - (reho_plot, ds_reho_plot, [("out_file", "in_file")]), + (reho_plot, ds_reho_plot, [ + ("out_file", "in_file"), + ("desc", "desc"), + ]), ]) # fmt:skip return workflow