diff --git a/docs/images/QC/sub-001_from-T1w_to-CITI168_regqc.png b/docs/images/QC/sub-001_from-T1w_to-CITI168_regqc.png new file mode 100755 index 00000000..26dd43db Binary files /dev/null and b/docs/images/QC/sub-001_from-T1w_to-CITI168_regqc.png differ diff --git a/docs/images/QC/sub-001_hemi-L_desc-unetf3d_dice.tsv b/docs/images/QC/sub-001_hemi-L_desc-unetf3d_dice.tsv new file mode 100755 index 00000000..99eff923 --- /dev/null +++ b/docs/images/QC/sub-001_hemi-L_desc-unetf3d_dice.tsv @@ -0,0 +1 @@ +0.8502973449669877 \ No newline at end of file diff --git a/docs/images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-dentate_desc-subfields_midthickness.surf.png b/docs/images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-dentate_desc-subfields_midthickness.surf.png new file mode 100755 index 00000000..109d7de6 Binary files /dev/null and b/docs/images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-dentate_desc-subfields_midthickness.surf.png differ diff --git a/docs/images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-hipp_desc-subfields_midthickness.surf.png b/docs/images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-hipp_desc-subfields_midthickness.surf.png new file mode 100755 index 00000000..aa6f6d5b Binary files /dev/null and b/docs/images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-hipp_desc-subfields_midthickness.surf.png differ diff --git a/docs/images/QC/sub-001_hemi-L_space-cropT1w_desc-subfields_atlas-histologyReference2023_dseg.png b/docs/images/QC/sub-001_hemi-L_space-cropT1w_desc-subfields_atlas-histologyReference2023_dseg.png new file mode 100755 index 00000000..8b69519c Binary files /dev/null and b/docs/images/QC/sub-001_hemi-L_space-cropT1w_desc-subfields_atlas-histologyReference2023_dseg.png differ diff --git a/docs/index.md b/docs/index.md index 0aaab5d3..68455213 100644 --- a/docs/index.md +++ b/docs/index.md @@ -51,6 +51,7 @@ pipeline/algorithms outputs/output_files outputs/visualization +outputs/QC ``` diff --git a/docs/outputs/QC.md b/docs/outputs/QC.md new file mode 100644 index 00000000..b7e29fe7 --- /dev/null +++ b/docs/outputs/QC.md @@ -0,0 +1,56 @@ +# Quality Control + +A complex pipeline can have multiple points of failure, so its imprtant to check the QC folder to ensure results make sense. This will contain (for example) the following files: + + hippunfold/ + └── sub-001 + └── qc + ├── sub-001_from-T1w_to-CITI168_regqc.png + ├── sub-001_hemi-L_desc-unetf3d_dice.tsv + ├── sub-001_hemi-L_space-cropT1w_desc-subfields_atlas-histologyReference2023_dseg.png + ├── sub-001_hemi-L_space-T1w_den-0p5mm_label-dentate_desc-subfields_midthickness.surf.png + └── sub-001_hemi-L_space-T1w_den-0p5mm_label-hipp_desc-subfields_midthickness.surf.png + +These files will be mirrored in the right (`hemi-R`) hemisphere, but we will skip that here for brevity. + +## Automated checks + +A very fast quality check is to simply look at `sub-001_hemi-L_desc-unetf3d_dice.tsv` and ensure that this number is greater than `0.7`. Typically this criterion can be applied to all subjects, and subjects `<0.7` can simply be automatically discarded. A more in-depth explanation for this is provided at the bottom of this page. + +To better understand some common points of failure, let's walk through the other files contained here: + +## Catch-all failures + +A good way to be sure that nothing has failed is also to look at the final results, as in the following snapshot images. + +`sub-001_hemi-L_space-cropT1w_desc-subfields_atlas-histologyReference2023_dseg.png`: + +![image](../images/QC/sub-001_hemi-L_space-cropT1w_desc-subfields_atlas-histologyReference2023_dseg.png) + +Here we see a volumetric output of subfield parcellations. If you have some hippocampal expertise then you may be able to see right away whether the results make sense. Some common errors include missed parts of the hippocampus and, sometimes, parts of collateral sulcus are included as part of the hippocapus. This is because of gross or "catastrophic" failures in the UNet tissue classification, which should be caught be the above automated check. Still, this collage of sample slices should give a good sense for whether/how the segmentation has failed. + +`sub-001_hemi-L_space-T1w_den-0p5mm_label-hipp_desc-subfields_midthickness.surf.png`: + +![image](../images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-hipp_desc-subfields_midthickness.surf.png) + +`sub-001_hemi-L_space-T1w_den-0p5mm_label-dentate_desc-subfields_midthickness.surf.png`: + +![image](../images/QC/sub-001_hemi-L_space-T1w_den-0p5mm_label-dentate_desc-subfields_midthickness.surf.png) + +In these two images we see the output hippocampal and dentate gyrus surfaces, which should roughly resemble those shown in the HippUnfold manuscript (though they may have some differences in shape and/or sulcal/gyral patterning due to inter-individual variance!). + +Sometimes the dentate gyrus surface may show a few large triangles that obscure the rest of the surface, as in the example here. This typically arises when just a few vertices are badly misplaced, and so typically this doesn't have much impact on quantitative results but it can look quite ugly in a figure. This can arise because the structure is small and contains few voxels, none of which have a laplace field close to a given value and so the xyz coordinate of a given laplace value is extrapolated badly. We will try to fix this issue in future releases. + +## Registration/cropping failures + +Though it seems simple, just getting the image roughly aligned to a template so we can then crop around the left and right hippocampi for subsequent processing (in `space-corobl`) can sometimes fail. We can examine this in the file `sub-001_from-T1w_to-CITI168_regqc.png`: + +![image](../images/QC/sub-001_from-T1w_to-CITI168_regqc.png) + +In this image, the subject tissue boundaries (e.g. grey matter, white matter, CSF, skull) are overlaid on a standard template. This is a linear alignment, so while gyri/sulci may not be well aligned, we should still see fairly good overlap that is sufficient to derive the position of the hippocampus in the image. This most often fails if the input images are in an unexpected orientation (e.g. this sometimes arises in ex-vivo scanning), and can sometimes be imporved by first running FSL's [reorient2std](https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/Fslutils). + +## Checking the automated check + +The `Dice > 0.7` rule is a good automatic check because HippUnfold typically shows high sensitivity to details in hippocampal structure, but when it fails it often does so catastrophically. These failures are very obvious, and can often be shown by this simple Dice score that compares a whole mask of the hippocampus to a whole mask of a hippocampus generated by a more conventional but less precise measure: diffeomorphic registration to a standard template. The Dice score should be good between the two methods, but will never be close to perfect (`>0.9`) since the HippUnfold method will differ as its more precise. + +This issue typically only occurs in ~1% of cases. diff --git a/docs/requirements.txt b/docs/requirements.txt index 713e48d4..32329cea 100644 --- a/docs/requirements.txt +++ b/docs/requirements.txt @@ -1,7 +1,8 @@ docutils==0.17 -sphinx-argparse +sphinx-argparse==0.4.0 sphinx_rtd_theme>=1.0.0 snakebids>=0.4.0 -myst-parser -sphinxcontrib.asciinema +myst-parser==1.0.0 +sphinxcontrib.asciinema==0.3.7 pulp<2.8.0 + diff --git a/hippunfold/config/snakebids.yml b/hippunfold/config/snakebids.yml index 4c9033c2..9342df5e 100644 --- a/hippunfold/config/snakebids.yml +++ b/hippunfold/config/snakebids.yml @@ -128,7 +128,7 @@ parse_args: --version: help: 'Print the version of HippUnfold' action: version - version: "1.4.2-pre.35" + version: "1.4.2-pre.36" --modality: @@ -365,6 +365,7 @@ singularity: autotop: 'docker://khanlab/hippunfold_deps:v0.5.0' xfm_identity: resources/etc/identity_xfm.txt +xfm_identity_itk: resources/etc/identity_xfm_itk.txt #templates enabled for template-based segmentation are here # TODO: should also perhaps include modalities avaialble, and any custom crop_native_res settings diff --git a/hippunfold/pipeline_description.json b/hippunfold/pipeline_description.json index 5b4b890c..be792a81 100644 --- a/hippunfold/pipeline_description.json +++ b/hippunfold/pipeline_description.json @@ -5,7 +5,7 @@ "GeneratedBy": [ { "Name": "hippunfold", - "Version": "1.4.2-pre.35", + "Version": "1.4.2-pre.36", "CodeURL": "https://github.com/khanlab/hippunfold", "Author": "Jordan DeKraker & Ali Khan", "AuthorEmail": "ali.khan@uwo.ca" diff --git a/hippunfold/resources/etc/identity_xfm_itk.txt b/hippunfold/resources/etc/identity_xfm_itk.txt new file mode 100644 index 00000000..fd2f22cf --- /dev/null +++ b/hippunfold/resources/etc/identity_xfm_itk.txt @@ -0,0 +1,5 @@ +#Insight Transform File V1.0 +#Transform 0 +Transform: MatrixOffsetTransformBase_double_3_3 +Parameters: 1 0 0 0 1 0 0 0 1 0 0 0 +FixedParameters: 0 0 0 diff --git a/hippunfold/workflow/rules/preproc_t2.smk b/hippunfold/workflow/rules/preproc_t2.smk index e410355c..81965af4 100644 --- a/hippunfold/workflow/rules/preproc_t2.smk +++ b/hippunfold/workflow/rules/preproc_t2.smk @@ -248,8 +248,6 @@ rule reg_t2_to_t1: def get_inputs_compose_t2_xfm_corobl(wildcards): if config["t1_reg_template"]: - # xfm0: t2 to t1 - # xfm1: t1 to corobl t2_to_t1 = ( bids( root=work, @@ -262,65 +260,23 @@ def get_inputs_compose_t2_xfm_corobl(wildcards): type_="itk" ), ) - t1_to_cor = ( - bids( - root=work, - datatype="warps", - **config["subj_wildcards"], - suffix="xfm.txt", - from_="T1w", - to="corobl", - desc="affine", - type_="itk" - ), - ) - return {"t2_to_t1": t2_to_t1, "t1_to_cor": t1_to_cor} - else: - # xfm0: t2 to template - t2_to_std = ( - bids( - root=work, - datatype="warps", - **config["subj_wildcards"], - suffix="xfm.txt", - from_="T2w", - to=config["template"], - desc="affine", - type_="itk" - ), - ) - - # xfm1: template to corobl - template_dir = Path(download_dir) / "template" / config["template"] - return {"t2_to_std": t2_to_std, "template_dir": template_dir} - - -def get_cmd_compose_t2_xfm_corobl(wildcards, input): - if config["t1_reg_template"]: - # xfm0: t2 to t1 - xfm0 = input.t2_to_t1 - # xfm1: t1 to corobl - xfm1 = input.t1_to_cor - else: - # xfm0: t2 to template - xfm0 = input.t2_to_std - # xfm1: template to corobl - xfm1 = Path(input.template_dir) / config["template_files"][config["template"]][ - "xfm_corobl" - ].format(**wildcards) - - return "c3d_affine_tool -itk {xfm0} -itk {xfm1} -mult -oitk {output}" + t2_to_t1 = Path(workflow.basedir).parent / config["xfm_identity_itk"] + to_corobl = ( + Path(download_dir) + / "template" + / config["template"] + / config["template_files"][config["template"]]["xfm_corobl"] + ) + return {"t2_to_t1": t2_to_t1, "to_corobl": to_corobl} # now have t2 to t1 xfm, compose this with t1 to corobl xfm rule compose_t2_xfm_corobl: input: unpack(get_inputs_compose_t2_xfm_corobl), - params: - cmd=get_cmd_compose_t2_xfm_corobl, output: - t2_to_cor=bids( + bids( root=work, datatype="warps", **config["subj_wildcards"], @@ -330,22 +286,12 @@ rule compose_t2_xfm_corobl: desc="affine", type_="itk" ), - log: - bids( - root="logs", - **config["subj_wildcards"], - suffix="composecorobl.txt", - from_="T2w", - to="corobl", - desc="affine", - type_="itk" - ), container: config["singularity"]["autotop"] group: "subj" shell: - "{params.cmd} > {log}" + "c3d_affine_tool -itk {input.t2_to_t1} -itk {input.to_corobl} -mult -oitk {output}" # if already have t2w in T1w space, then we don't need to use composed xfm: diff --git a/pyproject.toml b/pyproject.toml index 28a39a57..1809f783 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [tool.poetry] name = "hippunfold" -version = "1.4.2-pre.35" +version = "1.4.2-pre.36" description = "BIDS App for Hippocampal AutoTop (automated hippocampal unfolding and subfield segmentation)" authors = ["Jordan DeKraker & Ali Khan "]