dwidenoise: Add option -demean #2363
Conversation
This adds command-line option -demean, which subtracts the mean signal across the shortest matrix axis, essentially pre-allocating it as a signal component. This improves stability of rank determination when the coefficient of variance of the input data is small (e.g. fMRI data). Re-introduces code from �99d71856.
|
Base branch for merge is base of #2340, since much of the experimentation done around this was based around observations only made by utilising the enhancements included there. I would however suggest on resolving the pair of these that the command-line options be arranged into groups. |
|
Suspiciously relevant Twitter post appeared on my feed... |
Conflicts: cmd/dwidenoise.cpp
There was a problem hiding this comment.
Seems good to me. @dchristiaens, are you happy with these changes?
|
Sorry for not responding earlier. I meant to test this at some point, but I forgot. The code looks OK to me, except for the minor issue about the |
|
No worries, it's not like I've been particularly diligent of late... What issue are your referring to related to rankmap? I don't see any comment here that would match the review you mention... |
|
The other option here (from another suspiciously relevant twitter thread, though I didn't grab the link at the time) is to explicitly Z-transform all input data prior to PCA, which would supersede this proposal. |
Perhaps I messed up Github reviews... I was referring to this line, on which I had commented:
Do you mean this "double centring" you mentioned earlier? I'm not convinced this adds any value over the current approach of demeaning on the shortest matrix dimension, but I haven't done any testing on this... |
The double-centering would address the issue of selecting which matrix dimension to demean, and therefore potentially having sudden changes in behaviour if the number of volumes is just below vs. just above the kernel size. But it comes with a bit of extra computational cost because you have to iterate to do so. I was thinking instead just doing a single Z-transform along the shortest dimension; so not just demeaning as currently implemented here but also multiplying to unit variance as a kind of pre-conditioning. The link I came across was essentially just stating that it's a very common step prior to PCA. I'm not the expert here, I'm just verifying that we make an informed decision. Also, currently demeaning is done via an explicit command-line option, whereas if we were to double-centre or just Z-transform, we'd probably do it all the time and remove the command-line option. |
|
My 2¢... First off, it's important to note we're not doing a PCA analysis, at least not in the regular sense of the term. We won't be looking at the principal components themselves, we're only interested in those that look like noise. I think that means some of the considerations for PCA won't be applicable here. At least, that's my understanding, might be wrong here... So with that in mind: if demeaning is only there to help with numerical stability, simply demeaning one of the dimensions should be plenty enough to get 99% of the desired benefit. Double centering strikes me as overkill in this context, especially if it comes with additional computational overhead. As to going for a full Z transform, I think that would explicitly break the assumption of constant noise variance in the data (at least, over the patch being denoised). Again, might be appropriate for a PCA, but in our case, I expect doing this would lead to the distribution of singular values for pure noise no longer being characterised by the MP distribution... Regarding @dchristiaens's concerns with the rank when demeaning applied, I'm not sure either, but I agree that a quick simulation with random matrices should confirm one way or the other. |
|
It occurred to me that there's an additional possibility here for removing non-zero mean effects prior to decomposition. If the input data possess a DWI shell structure, then what one could choose to do is: for each voxel, rather than regressing out the mean across all volumes, instead regress out the mean per shell. (This is even more likely that a basic mean regression to have a detrimental impact on the estimator; so would need to be tested) |
By default, if the input image series contains a gradient table that appears to be arranged into shells, then the demeaning will be computed per shell per voxel; otherwise a single mean will be computed per voxel. Demeaning can also be disabled using the -demean option. The code responsible for the demeaning has been integrated with the code for phase demodulation and the variance-stabilising transform as "preconditioning". The presence of demeaning does not affect the PCA or the noise level estimation in any way; it only influences the _reported_ signal rank. Supersedes #2363.
|
Closing in favour of #3023. |
This adds command-line option
-demean, which subtracts the mean signal across the shortest matrix axis, essentially pre-allocating it as a signal component.This improves stability of rank determination when the coefficient of variance of the input data is small (e.g. fMRI data).
Re-introduces code from 99d7185.
Magnitude data denoising
Rank images:
Complex data denoising
Rank images:
(Big bright spot at focus point is issue introduced by use of sub-of-squares reconstruction; should have used the SENSE1 data for demonstration...)