Skip to content

Commit

Permalink
Update the QC metrics (PennLINC#958)
Browse files Browse the repository at this point in the history
  • Loading branch information
tsalo authored Oct 11, 2023
1 parent e0f0848 commit 67be955
Show file tree
Hide file tree
Showing 3 changed files with 80 additions and 51 deletions.
71 changes: 48 additions & 23 deletions xcp_d/interfaces/plotting.py
Original file line number Diff line number Diff line change
Expand Up @@ -277,7 +277,7 @@ class QCPlots(SimpleInterface):
output_spec = _QCPlotsOutputSpec

def _run_interface(self, runtime):
# Load confound matrix and load motion with motion filtering
# Load confound matrix and load motion without motion filtering
confounds_df = pd.read_table(self.inputs.fmriprep_confounds_file)
preproc_motion_df = load_motion(
confounds_df.copy(),
Expand All @@ -297,6 +297,7 @@ def _run_interface(self, runtime):
censoring_df = pd.read_table(self.inputs.temporal_mask)
tmask_arr = censoring_df["framewise_displacement"].values
num_censored_volumes = int(tmask_arr.sum())
num_retained_volumes = int((tmask_arr == 0).sum())

# Apply temporal mask to interpolated/full data
rmsd_censored = rmsd[tmask_arr == 0]
Expand Down Expand Up @@ -382,40 +383,57 @@ def _run_interface(self, runtime):

# Calculate QC measures
mean_fd = np.mean(preproc_fd_timeseries)
mean_rms = np.nanmean(rmsd_censored) # first value can be NaN if no dummy scans
mean_fd_post_censoring = np.mean(postproc_fd_timeseries)
mean_relative_rms = np.nanmean(rmsd_censored) # first value can be NaN if no dummy scans
mean_dvars_before_processing = np.mean(dvars_before_processing)
mean_dvars_after_processing = np.mean(dvars_after_processing)
motionDVCorrInit = np.corrcoef(preproc_fd_timeseries, dvars_before_processing)[0][1]
motionDVCorrFinal = np.corrcoef(postproc_fd_timeseries, dvars_after_processing)[0][1]
fd_dvars_correlation_initial = np.corrcoef(preproc_fd_timeseries, dvars_before_processing)[
0, 1
]
fd_dvars_correlation_final = np.corrcoef(postproc_fd_timeseries, dvars_after_processing)[
0, 1
]
rmsd_max_value = np.nanmax(rmsd_censored)

# A summary of all the values
qc_values_dict.update(
{
"meanFD": [mean_fd],
"relMeansRMSMotion": [mean_rms],
"relMaxRMSMotion": [rmsd_max_value],
"meanDVInit": [mean_dvars_before_processing],
"meanDVFinal": [mean_dvars_after_processing],
"mean_fd": [mean_fd],
"mean_fd_post_censoring": [mean_fd_post_censoring],
"mean_relative_rms": [mean_relative_rms],
"max_relative_rms": [rmsd_max_value],
"mean_dvars_initial": [mean_dvars_before_processing],
"mean_dvars_final": [mean_dvars_after_processing],
"num_dummy_volumes": [dummy_scans],
"num_censored_volumes": [num_censored_volumes],
"nVolsRemoved": [dummy_scans],
"motionDVCorrInit": [motionDVCorrInit],
"motionDVCorrFinal": [motionDVCorrFinal],
"num_retained_volumes": [num_retained_volumes],
"fd_dvars_correlation_initial": [fd_dvars_correlation_initial],
"fd_dvars_correlation_final": [fd_dvars_correlation_final],
}
)

qc_metadata = {
"meanFD": {
"mean_fd": {
"LongName": "Mean Framewise Displacement",
"Description": (
"Average framewise displacement without any motion parameter filtering. "
"This value includes high-motion outliers, but not dummy volumes. "
"FD is calculated according to the Power definition."
),
"Units": "mm",
"Units": "mm / volume",
"Term URL": "https://doi.org/10.1016/j.neuroimage.2011.10.018",
},
"mean_fd_post_censoring": {
"LongName": "Mean Framewise Displacement After Censoring",
"Description": (
"Average framewise displacement without any motion parameter filtering. "
"This value does not include high-motion outliers or dummy volumes. "
"FD is calculated according to the Power definition."
),
"Units": "mm / volume",
"Term URL": "https://doi.org/10.1016/j.neuroimage.2011.10.018",
},
"relMeansRMSMotion": {
"mean_relative_rms": {
"LongName": "Mean Relative Root Mean Squared",
"Description": (
"Average relative root mean squared calculated from motion parameters, "
Expand All @@ -424,7 +442,7 @@ def _run_interface(self, runtime):
),
"Units": "arbitrary",
},
"relMaxRMSMotion": {
"max_relative_rms": {
"LongName": "Maximum Relative Root Mean Squared",
"Description": (
"Maximum relative root mean squared calculated from motion parameters, "
Expand All @@ -433,44 +451,51 @@ def _run_interface(self, runtime):
),
"Units": "arbitrary",
},
"meanDVInit": {
"mean_dvars_initial": {
"LongName": "Mean DVARS Before Postprocessing",
"Description": (
"Average DVARS (temporal derivative of root mean squared variance over "
"voxels) calculated from the preprocessed BOLD file, after dummy scan removal."
),
"TermURL": "https://doi.org/10.1016/j.neuroimage.2011.02.073",
},
"meanDVFinal": {
"mean_dvars_final": {
"LongName": "Mean DVARS After Postprocessing",
"Description": (
"Average DVARS (temporal derivative of root mean squared variance over "
"voxels) calculated from the denoised BOLD file."
),
"TermURL": "https://doi.org/10.1016/j.neuroimage.2011.02.073",
},
"num_dummy_volumes": {
"LongName": "Number of Dummy Volumes",
"Description": (
"The number of non-steady state volumes removed from the time series by XCP-D."
),
},
"num_censored_volumes": {
"LongName": "Number of Censored Volumes",
"Description": (
"The number of high-motion outlier volumes censored by XCP-D. "
"This does not include dummy volumes."
),
},
"nVolsRemoved": {
"LongName": "Number of Dummy Volumes",
"num_retained_volumes": {
"LongName": "Number of Retained Volumes",
"Description": (
"The number of non-steady state volumes removed from the time series by XCP-D."
"The number of volumes retained in the denoised dataset. "
"This does not include dummy volumes or high-motion outliers."
),
},
"motionDVCorrInit": {
"fd_dvars_correlation_initial": {
"LongName": "FD-DVARS Correlation Before Postprocessing",
"Description": (
"The Pearson correlation coefficient between framewise displacement and DVARS "
"(temporal derivative of root mean squared variance over voxels), "
"after removal of dummy volumes, but before removal of high-motion outliers."
),
},
"motionDVCorrFinal": {
"fd_dvars_correlation_final": {
"LongName": "FD-DVARS Correlation After Postprocessing",
"Description": (
"The Pearson correlation coefficient between framewise displacement and DVARS "
Expand Down
36 changes: 20 additions & 16 deletions xcp_d/interfaces/report.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,12 +30,13 @@
\t\t<ul class="elem-desc">
\t\t\t<li>BOLD volume space: {space}</li>
\t\t\t<li>Repetition Time (TR): {TR:.03g}</li>
\t\t\t<li>Mean Framewise Displacement: {meanFD}</li>
\t\t\t<li>Mean Relative RMS Motion: {meanRMS}</li>
\t\t\t<li>Max Relative RMS Motion: {maxRMS}</li>
\t\t\t<li>Mean Framewise Displacement: {mean_fd}</li>
\t\t\t<li>Mean Relative RMS Motion: {mean_relative_rms}</li>
\t\t\t<li>Max Relative RMS Motion: {max_relative_rms}</li>
\t\t\t<li>DVARS Before and After Processing : {dvars_before_after}</li>
\t\t\t<li>Correlation between DVARS and FD Before and After Processing : {corrfddv}</li>
\t\t\t<li>Number of Volumes Censored : {volcensored}</li>
\t\t\t<li>Correlation between DVARS and FD Before and After Processing :
{fd_dvars_correlation}</li>
\t\t\t<li>Number of Volumes Censored : {num_vols_censored}</li>
\t\t</ul>
"""

Expand Down Expand Up @@ -145,25 +146,28 @@ class FunctionalSummary(SummaryInterface):
def _generate_segment(self):
space = get_entity(self.inputs.bold_file, "space")
qcfile = pd.read_csv(self.inputs.qc_file)
meanFD = str(round(qcfile["meanFD"][0], 4))
meanRMS = str(round(qcfile["relMeansRMSMotion"][0], 4))
maxRMS = str(round(qcfile["relMaxRMSMotion"][0], 4))
dvars = f"{round(qcfile['meanDVInit'][0], 4)}, {round(qcfile['meanDVFinal'][0], 4)}"
mean_fd = str(round(qcfile["mean_fd"][0], 4))
mean_relative_rms = str(round(qcfile["mean_relative_rms"][0], 4))
max_relative_rms = str(round(qcfile["max_relative_rms"][0], 4))
dvars = (
f"{round(qcfile['mean_dvars_initial'][0], 4)}, "
f"{round(qcfile['mean_dvars_final'][0], 4)}"
)
fd_dvars_correlation = (
f"{round(qcfile['motionDVCorrInit'][0], 4)}, "
f"{round(qcfile['motionDVCorrFinal'][0], 4)}"
f"{round(qcfile['fd_dvars_correlation_initial'][0], 4)}, "
f"{round(qcfile['fd_dvars_correlation_final'][0], 4)}"
)
num_vols_censored = str(round(qcfile["num_censored_volumes"][0], 4))

return QC_TEMPLATE.format(
space=space,
TR=self.inputs.TR,
meanFD=meanFD,
meanRMS=meanRMS,
maxRMS=maxRMS,
mean_fd=mean_fd,
mean_relative_rms=mean_relative_rms,
max_relative_rms=max_relative_rms,
dvars_before_after=dvars,
corrfddv=fd_dvars_correlation,
volcensored=num_vols_censored,
fd_dvars_correlation=fd_dvars_correlation,
num_vols_censored=num_vols_censored,
)


Expand Down
24 changes: 12 additions & 12 deletions xcp_d/utils/qcmetrics.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,15 +46,15 @@ def compute_registration_qc(bold2t1w_mask, anat_brainmask, bold2template_mask, t
template_mask_arr = nb.load(template_mask).get_fdata()

reg_qc = {
"coregDice": [dice(bold2t1w_mask_arr, t1w_mask_arr)],
"coregPearson": [pearson(bold2t1w_mask_arr, t1w_mask_arr)],
"coregCoverage": [overlap(bold2t1w_mask_arr, t1w_mask_arr)],
"normDice": [dice(bold2template_mask_arr, template_mask_arr)],
"normPearson": [pearson(bold2template_mask_arr, template_mask_arr)],
"normCoverage": [overlap(bold2template_mask_arr, template_mask_arr)],
"coreg_dice": [dice(bold2t1w_mask_arr, t1w_mask_arr)],
"coreg_correlation": [pearson(bold2t1w_mask_arr, t1w_mask_arr)],
"coreg_overlap": [overlap(bold2t1w_mask_arr, t1w_mask_arr)],
"norm_dice": [dice(bold2template_mask_arr, template_mask_arr)],
"norm_correlation": [pearson(bold2template_mask_arr, template_mask_arr)],
"norm_overlap": [overlap(bold2template_mask_arr, template_mask_arr)],
}
qc_metadata = {
"coregDice": {
"coreg_dice": {
"LongName": "Coregistration Sørensen-Dice Coefficient",
"Description": (
"The Sørensen-Dice coefficient calculated between the binary brain masks from the "
Expand All @@ -64,7 +64,7 @@ def compute_registration_qc(bold2t1w_mask, anat_brainmask, bold2template_mask, t
),
"Term URL": "https://en.wikipedia.org/wiki/S%C3%B8rensen%E2%80%93Dice_coefficient",
},
"coregPearson": {
"coreg_correlation": {
"LongName": "Coregistration Pearson Correlation",
"Description": (
"The Pearson correlation coefficient calculated between the binary brain masks "
Expand All @@ -74,7 +74,7 @@ def compute_registration_qc(bold2t1w_mask, anat_brainmask, bold2template_mask, t
),
"Term URL": "https://en.wikipedia.org/wiki/Pearson_correlation_coefficient",
},
"coregCoverage": {
"coreg_overlap": {
"LongName": "Coregistration Coverage Metric",
"Description": (
"The Szymkiewicz-Simpson overlap coefficient calculated between the binary brain "
Expand All @@ -83,7 +83,7 @@ def compute_registration_qc(bold2t1w_mask, anat_brainmask, bold2template_mask, t
),
"Term URL": "https://en.wikipedia.org/wiki/Overlap_coefficient",
},
"normDice": {
"norm_dice": {
"LongName": "Normalization Sørensen-Dice Coefficient",
"Description": (
"The Sørensen-Dice coefficient calculated between the binary brain masks from the "
Expand All @@ -93,7 +93,7 @@ def compute_registration_qc(bold2t1w_mask, anat_brainmask, bold2template_mask, t
),
"Term URL": "https://en.wikipedia.org/wiki/S%C3%B8rensen%E2%80%93Dice_coefficient",
},
"normPearson": {
"norm_correlation": {
"LongName": "Normalization Pearson Correlation",
"Description": (
"The Pearson correlation coefficient calculated between the binary brain masks "
Expand All @@ -103,7 +103,7 @@ def compute_registration_qc(bold2t1w_mask, anat_brainmask, bold2template_mask, t
),
"Term URL": "https://en.wikipedia.org/wiki/Pearson_correlation_coefficient",
},
"normCoverage": {
"norm_overlap": {
"LongName": "Normalization Overlap Coefficient",
"Description": (
"The Szymkiewicz-Simpson overlap coefficient calculated between the binary brain "
Expand Down

0 comments on commit 67be955

Please sign in to comment.