Noise estimation

[araikes]

Hi @robbert-harms,

I appreciate this toolbox and I am having good success using this with ex-vivo mouse brains. I do have a question to ensure that I am doing this "the right way." In the modeling, the noise standard deviation is computed within the mask. If I have previously denoised my data (e.g., with DIPY's patch2self), how does this affect MDT's modeling (particularly for NODDI)?

My general workflow is:

• Denoise (DIPY)
• Unringing (DIPY)
• Eddy correction (TORTOISE)
• Bias correction (Mrtrix)
• MDT for NODDI & Kurtosis

[robbert-harms]

Hi Araikes,

Happy to hear MDT is of use to you in your study. Your question is spot on, and an important one. Noise modeling is far too often neglected in the papers I review.

In your case I would recommend modifying the NODDI ex-vivo model to use a Gaussian noise model instead of the default Offset-Gaussian. And then either fix the noise noise standard deviation to 1 or let it be the default auto-estimation. To change the likelihood model you would have to redefine the model. To do so, please unzip this file:
NODDI.py.zip

and place it in your .mdt/<version>/components/user/composite_models folder, you can then use the NODDI_ExVivo_GaussianNoise in your analysis.

The rationale behind this is as follows. In the default NODDI model, MDT uses the OffsetGaussian noise model. This noise model works under the assumption that your diffusion data is Rician distributed (and OffsetGaussian is a well-known approximation to Rician).

Edit: the below part was part of the original answer and has some nuances. Please see messages below.
In your data, you are denoising before using it in MDT, as such the Rician assumption no longer holds as the noise is nullified). As such, there is theoretically no need for further noise modeling. Yet, maximum likelihood estimation requires you to define a likelihood/noise model and for that I recommend the Gaussian model. Gaussian, because parameter estimation using Gaussian likelihood is (theoretically) insensitive to the noise level used (there might be small numerical errors when using a very high or very low noise level).

Using Gaussian you have the freedom to choose your noise standard deviation. You could fix it to a value like 1 or sqrt(2) for consistency when using different slices of the data. The only reason why you would want to set it to the actual (small) noise level remaining is that having a true estimate of the Gaussian noise level enables you to use the BIC/AICc maps for parameter comparisons over models.

I hope that helps a bit. If there are parts unclear, please ask.

Best wishes,

Robbert

[araikes]

@robbert-harms

Thanks for this. I will give this a test and let you know.

Are any noise estimation changes needed if running kurtosis?

[robbert-harms]

Hi Araikes,

Same story, copy the model and adapt the likelihood function. You can create a new Kurtosis.py file in your mdt composite models folder and add something like:

from mdt import CompositeModelTemplate

class Kurtosis_Gaussian(CompositeModelTemplate):
    """The standard Kurtosis model with in vivo defaults."""
    model_expression = '''
        S0 * KurtosisTensor
    '''
    volume_selection = {'b': [(0, 3.0e9 + 0.1e9)]}
    likelihood_function = 'Gaussian'

Best,

Robbert

[ghammad]

Hello @robbert-harms

Sorry to jump in this discussion but as you mentioned, noise estimation is important and often neglected.

Similarly to @araikes, our standard DWI processing pipeline does include denoising and unringing, using the MRtrix3 package.
According to MRTrix documentation for the denoise function, it does not correct for "non-Gaussian noise biases present in magnitude-reconstructed MR images".
However, in this MRTrix thread, I can read the following from J. Tournier:

Yes, magnitude data has a non-Gaussian noise distribution at very low SNR, but it is nonetheless very well approximated as Gaussian with an offset (and slightly lower standard deviation), which means it can still be denoised very effectively by dwidenoise – though that will undeniably not reduce the Rician-induced bias

Therefore, my (naive) understanding is that, unless complex images are used (for which Gaussian noise is expected), magnitude images, after denoising, remain affected by a residual non-Gaussian noise.

In that case, would the strategy to use a Gaussian noise model in MDT still be advisable?

Last point; since MRTrix denoise function outputs a noise image which, according to this thread gives

estimate of the standard deviation S of additive Gaussian noise in the magnitude images

is It possible to make use a such image to constrain the noise model parameters in MDT?

Thank you in advance for your help and for making MDT available.

Best regards,

Grégory

[robbert-harms]

Hi @ghammad,

Thank you for pointing me towards the MRTrix thread. I see that I overestimated the dwidenoise function, in my understanding it was also able to cope with the Rician bias. Given the new information I would retract my previous advice and still go with the good-old OffsetGaussian noise model (I also updated my answer above).

It is indeed possible to use a Gaussian noise matrix instead of a scalar. In the Python API this is done by providing (for example) a nifti file to the input data:

input_data = mdt.load_input_data(
    'data.nii',
    'bvecs.prtcl',
    'mask.nii',
    noise_std='sigma_CP.nii'
)

Best,

Robbert

[ghammad]

@robbert-harms Thank you for the prompt reply. For completeness, I checked the DIPY package and apparently there are multiple denoise techniques available. So, I don't know whether the one used by @araikes does reduce the Rician bias or not.
For MRtrix, (using the MP-PCA technique) apparently it does not.

Best regards,

Grégory

[araikes]

I use DIPY's patch2self. I'll reach out to the devs and see how it handles Rician noise.

[ghammad]

Hello @araikes, any news about the DIPY's patch2self function?

Thank you.