athenaepi / dmipy Goto Github PK

To allow fitting of data without b0-measurements it should be possible to add S0 as an optimized parameter.

from dmipy.distributions.distribute_models import SD2BinghamDistributed

ValueError Traceback (most recent call last)
in
----> 1 from dmipy.distributions.distribute_models import SD2BinghamDistributed
2 # bingham_dispersed_bundle = SD2BinghamDistributed(models=[stick, zeppelin])
3 # bingham_dispersed_bundle.parameter_names

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in find_and_load(name, import)

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in find_and_load_unlocked(name, import)

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in _load_unlocked(spec)

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in _load_backward_compatible(spec)

~/anaconda3/lib/python3.6/site-packages/dmipy-0.1.dev0-py3.6.egg/dmipy/distributions/distribute_models.py in
----> 1 from . import distributions
2 from collections import OrderedDict
3 from itertools import chain
4 from ..utils.spherical_convolution import sh_convolution
5 from ..utils.utils import T1_tortuosity, parameter_equality

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in find_and_load(name, import)

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in find_and_load_unlocked(name, import)

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in _load_unlocked(spec)

~/anaconda3/lib/python3.6/importlib/_bootstrap.py in _load_backward_compatible(spec)

~/anaconda3/lib/python3.6/site-packages/dmipy-0.1.dev0-py3.6.egg/dmipy/distributions/distributions.py in
37 log_bingham_normalization_splinefit = np.load(
38 join(DATA_PATH,
---> 39 "bingham_normalization_splinefit.npz"), encoding='bytes')['arr_0']
40
41 inverse_sh_matrix_kernel = {

~/anaconda3/lib/python3.6/site-packages/numpy/lib/npyio.py in getitem(self, key)
260 return format.read_array(bytes,
261 allow_pickle=self.allow_pickle,
--> 262 pickle_kwargs=self.pickle_kwargs)
263 else:
264 return self.zip.read(key)

~/anaconda3/lib/python3.6/site-packages/numpy/lib/format.py in read_array(fp, allow_pickle, pickle_kwargs)
690 # The array contained Python objects. We need to unpickle the data.
691 if not allow_pickle:
--> 692 raise ValueError("Object arrays cannot be loaded when "
693 "allow_pickle=False")
694 if pickle_kwargs is None:

ValueError: Object arrays cannot be loaded when allow_pickle=False

Dipy's visualization updated so replace all hcp examples like the one in mc-mdi.

Speaks for itself.

['S4SphereGaussianPhaseApproximation_1_diameter',
'G1Ball_1_lambda_iso',
'C1Stick_1_mu',
'C1Stick_1_lambda_par',
'partial_volume_0',
'partial_volume_1',
'partial_volume_2']

Does the partial_volume_0/1/2 represent tumor, extra-cell, and vascular respectively?

The explanation for parameter's maps followed is from example_noddi_watson.ipynb
_

"Above we see the familiar maps. Remember that partial_volume_0/1/2 represent fractions of models in the order they were given to the MultiCompartmentMicrostructureModel object. In this case, 0 represents the fraction of Ball (CSF), and the 1 represents the Stick (intra-axonal) and 2 represents the zeppelin (extra-axonal).

In the bottom-right ODI map you can see that the concentration is highest in the corpus callosum and in coherent white matter bundles. Finally, the 𝜆∥ and 𝜆iso maps are fixed to the values we gave them."

_

I think the partial_volume_0 is from G1ball, represents fraction of CSF

I have two questions:
**1. what's the parameters for partial_volume_1 and partial_volume_2 in the explaination

Does it refers to 'partial_volume_1' and 'SD1WatsonDistributed_1_partial_volume_0'?
or refers to:
'SD1WatsonDistributed_1_partial_volume_0' and 'SD1WatsonDistributed_1_SD1Watson_1_mu'

2. Is the map for Stick (intra-axonal) parameter equal to NDI (neurite density index ) map?**

Here's the watson parameters:
NODDI_mod.parameter_names

['G1Ball_1_lambda_iso',
'SD1WatsonDistributed_1_SD1Watson_1_mu',
'SD1WatsonDistributed_1_partial_volume_0',
'SD1WatsonDistributed_1_SD1Watson_1_odi',
'partial_volume_0',
'partial_volume_1']

It would be appreciated that if you could help, : )

Extension of FittedMCSH model to fit the spherical FOD using one or more parametric distributions (watson/bingam), essentially making a bridge to standard MC models. This way you can get an estimate of the dispersion in multiple crossings from e.g. MT-CSD without having to do inefficient NODDIx stuff.

The CVXPY optimizer is generally applicable but not super fast for CSD optimization.
An alternative is the classic tournier07 optimizer - especially when wrapping the iteration using Numba.

• should be able to take any single-component model as long as all parameters are fixed.
• if multiple models are given, it is only accepted if the volume fractions are fixed and the combined response function should be composed internally.
• error warning if parameters are not set.
• iteration should be numba-compatible

MT-CSD doens't necessarily enforce unity on the volume fractions when fitting the signal.
Since no penalty on non-unity fractions is imposed, this is the same as fitting the normalized signal E, enforcing unity on the volume fractions and fitting only the shape of the tissue responses, and then rescaling the volume fractions to correct for the differences in S0 tissue responses.

The current parameter linking system works but is difficult to interact with and add custom functions into.

• Parameter and optimization link descriptions should be passed as dictionaries.
• each dictionary describes the mathematical relation that is imposed, which parameter are input, which parameter are removed or added.
• Linking and optimization dictionaries should have similar setup, but optimization dictionaries need to somehow be more complicated to include optimization ranges for the new optimization parameter.
• In relation to replaying parameter links from spherical mean models to regular or sh-models, it should be able to pass a list of these dictionaries and the model will then just replay these links.

In relation to fixing parameters

• parameter should not just be removed from the parameter list, (meaning if you call the same function again it will give an error since the parameter was removed, unless the model is re-initialized). The system should be changed so that parameter are recatagorized somehow and a summary function can be called which lists which parameter are optimized / fixed / linked.

The linking functions sometimes vary whether they are given to distributed models or MC models. This should be resolved.

Astro-models are just spherical mean representations of cylinder models. This should be emphasized in the documentation.

Following the recent publication of VERDICT-AMICO, It seems that AMICO can be implemented as a generalised optimizer instead of a model-wise one.

Necessary for simple iterative process to do single shell (b0 + 1 shell) multi-tissue CSD.

https://www.researchgate.net/publication/301766619_A_novel_iterative_approach_to_reap_the_benefits_of_multi-tissue_CSD_from_just_single-shell_b0_diffusion_MRI_data

I'm trying smt_noddi use my own data

(second fit)
"Cannot estimate signal grid with voxel-dependent x0_vector."
Can i ignore this ?

General function to set ranges without having to separately set internal scales and ranges. Input should be in SI units and internally the range is divided by the max of the ranges and the scale is the max of the input ranges.

Small extension to FittedMultiCompartmentSphericalHarmonicsModel to estimate the apparent fiber density.

Raffelt, David, et al. "Apparent fibre density: a novel measure for the analysis of diffusion-weighted magnetic resonance images." Neuroimage 59.4 (2012): 3976-3994.

Multi-Tissue CSD calibration should return model-free models (no parameters), but still be able to return rotational harmonics and spherical mean for the acquisition scheme from which they were generated.

Alternative FOD optimizer taking into account non-Gaussian noise and enabling total variation spatial regularization.

• Should click into the MultiCompartmentSphericalHarmonics model, but would need to be adapted to allow for iterating over the entire image for the total variation term.
• Perhaps TV should be an option for any MC framework in Dmipy.

Canales-Rodríguez, Erick J., et al. "Spherical deconvolution of multichannel diffusion MRI data with non-Gaussian noise models and spatial regularization." PloS one 10.10 (2015): e0138910.

It's not possible to import dmipy.utils in a python 3 environment due to use of print as a statement.

Right now a fractional parameter ranges between [0-1] * some_parameter, but that means that it can also take the value 0, which is typically not the case in some_parameter.

fractional parameter should be defined to also take into account the parameter_ranges of some_parameter, i.e. value_min + [0-1] * (some_parameter - value_min), where some_parameter ranges between value_min and some value_max.

Implement equation (2) of [1] after IVIM parameter estimation.

[1] Lemke, A., Laun, F. B., Simon, D., Stieltjes, B., & Schad, L. R. (2010). An in vivo verification of the intravoxel incoherent motion effect in diffusion‐weighted imaging of the abdomen. Magnetic resonance in medicine, 64(6), 1580-1585.

With #51 in place reproduce example similar to axcaliber one with activeax.

slide 60 of :
https://www.cs.jhu.edu/~misha/Spring15/17.pdf

and explanation of anisotropic response function.

Following discussion in #67 and suggestion by @thijsdhollander.

Option to clip SDM above certain number to avoid noisy CSF response functions.

e.g. MC-MDI cannot be fitted with MIX because mix assumes volume fractions are linear parameters. The algorithms needs to have a fall-back to just using differential evolution on the model if this is the case, instead of the error it has now.

Currently, the Watson and Bingham distributions are automatically parameterized using the ODI and beta_fraction parameters, but these are in fact optimization parameters just to facilitate optimization of the actual \kappa and \beta parameters.

• Using these optimization parameters (with their respective optimization ranges) should be optional.

Likely due to some kind of overflow error some cylinder models don't give meaninful signal attenuation values when giving huge diameters. There should be a catch for this.

When a gradient table is converted from dmipy to dipy format and delta is None or Delta is None, the devoted function gives the following error.

In [8]: gtab_dmipy2dipy(gtab_dmipy)                                                 
None
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-8-89848d3de32f> in <module>
----> 1 gtab_dmipy2dipy(scheme)

/home/user/dmipy/dmipy/core/acquisition_scheme.py in gtab_dmipy2dipy(dmipy_gradient_table)
    907     elif len(np.unique(delta)) == 1:
--> 908         delta = delta[0]
    909 
    910     if len(np.unique(Delta)) > 1:

TypeError: 'NoneType' object is not subscriptable

and this happens also at line 911 for Delta.

A simple hack that fixes this problem is to change the checks in lines 907 and 910 as follows:

elif delta is not None and len(np.unique(delta)) == 1:
...
elif Delta is not None and len(np.unique(Delta)) == 1:

I'll put up a PR with this in the close future.

Using tissue_response_models to other frameworks outside CSD should also include the option to fit the tissue response of the model to the signal.

Can be seen as extension of #19

The median_otsu method from dipy takes as input the variable vol_idx which must be an iterable.
In the three_tissue_response_dhollander16 function this parameter is passed as a scalar.

I'll fix it and push a PR, but we should be aware of the fact that this could be happening somewhere else. It may be due to the recent dipy update (v1.0) but again, I'm not sure and this must be checked.

Currently it produces an error in python 3.

dmipy is a great tool to me, but why no dki model here? :)

delta [ms] |Delta[ms]

@rutgerfick I think that the dependencies are too lax. We need to add everything that the tool needs to import to run. Boto is just needed for the examples so I added it but is there something missing? I think, for instance, that you are calling matplotlib from within the modules.

Current implementation in dmipy/dmipy/utils/spherical_mean.py is very slow and ugly with shell-wise for-loops. In spherical mean model fitting it is often half the estimation time.

• precalculate and pre-invert multishell spherical mean observation matrix A_sm_inv in acquisition scheme
• multi-shell spherical mean estimation should be nothing more than sm = np.dot(A_sm_inv, data)

verdict_fit = verdict_mod.fit(acq_scheme, data, solver='mix', use_parallel_processing=True

/home/omnisky/anaconda3/lib/python3.6/site-packages/numpy/core/fromnumeric.py:3118: RuntimeWarning: Mean of empty slice.
out=out, **kwargs)
/home/omnisky/anaconda3/lib/python3.6/site-packages/numpy/core/_methods.py:78: RuntimeWarning: invalid value encountered in true_divide
ret, rcount, out=ret, casting='unsafe', subok=False)

With the recent update of dipy, many functions will not work.
There are two possible solutions:

• Require an old version of dipy (suboptimal)
• Do a complete revision of the parts of dmipy where dipy is used (will take longer than expected)

What do we do @rutgerfick ?

I'm attempting to set up an acquisition scheme for CSD for some data from a collaborator, and I'm having trouble finding little delta and big Delta, despite having the DICOM header. Where should I start looking for these values?

In trying to set up my analysis, one reference has beentutorial_setting_up_acquisition_scheme.ipynb. That notebook sets up the acquisition scheme for HCP data and says that
the delta and Delta times we know from the HCP documentation in seconds.
However, I'm also unable to find those values in the HCP docs (for instance, under "Diffusion Session Scan Protocol" in https://humanconnectome.org/storage/app/media/documentation/s1200/HCP_S1200_Release_Appendix_I.pdf.

The models inside a bundle model are assumed to be aligned - so they should not have individual orientations.

The input orientations should therefore be linked to a bundle_mu parameter of the bundle itself.

Currently, Dmipy can modularly generate most models proposed for PGSE-based microstructure imaging.

However, the DIAMOND model by Sherrer et al. [1] is a bit particular in that it uses a matrix-variate Gamma distribution - the implementation of which does not look trivial at all.

It would be great if someone could help with DIAMOND's implementation to make Dmipy's collection truly complete!

[1] Scherrer, Benoit, et al. "Characterizing brain tissue by assessment of the distribution of anisotropic microstructural environments in diffusion‐compartment imaging (DIAMOND)." Magnetic resonance in medicine 76.3 (2016): 963-977.

fvtk visualization in dipy has been deprecated.

• create function in utils for specifically vizualizing the examples that uses the new visualization tools.
• takes the fod field and the background as input and just prints the image, so the clutter in the examples is removed.

Every framework should have a general catch for when the optimization errors-out for whatever reason, and then store that voxel in an "error map" or something. Much easier to inspect weird voxels that way instead of having to debug the entire data for the one error voxel to occur again.

The dmipy.utils.estimate_spherical_mean_shell function prints the value of the spherical mean every time it is called. This information may not be interesting in every single case.

I propose to eliminate this print.

seems to be a boundary case.

Sphere implementation for Callaghan model does not pass tests yet.

The returned variables listed in the API documentation are not in the correct order.

Listed:

    Returns
    -------
    S0_wm : float,
        white matter S0 response value.
    TR2_wm_model : Dmipy TR2AnisotropicTissueResponseModel,
        ModelFree representation of white matter response.
    S0_gm : float,
        grey matter S0 response value.
    TR1_gm_model : Dmipy TR1IsotropicTissueResponseModel,
        ModelFree representation of grey matter response.
    S0_csf : float,
        csf S0 response value.
    TR1_csf_model : Dmipy TR1IsotropicTissueResponseModel,
        ModelFree representation of csf response.
    three_tissue_selection: array of size (x, y, z, 3),
        RGB mask of selected voxels used for white/grey matter and CSD.

Actual:

return ([S0_wm, S0_gm, S0_csf],
            [TR2_wm_model, TR1_gm_model, TR1_csf_model],
            three_tissue_selection)

Hi again @rutgerfick and others,

This is not so much an "issue" as just a small announcement, that might be of use for some of your developments here! So as the title implies, SS3T-CSD is now also "publicly" released from my end. Related, the 2019 version of the response function calibration (for either multi-shell as well as single-shell data) now also comes with it. The latter improves on the 2016 version mostly on the single-fibre white matter front, and while at it, delivers quite an improvement in speed as a bonus. Just a few sources, should it come in handy in any possible way:

• Release is a fork of mrtrix, called "MRtrix3Tissue"; all info on this website: https://3tissue.github.io
• Some release info (because this incorporates all mrtrix developments since RC3 as well, so about 1.5 years...): https://3tissue.github.io/doc/install/release.html (mostly if you need to know about changed command names)
• SS3T-CSD information (this will mostly interest you, I reckon): https://3tissue.github.io/doc/ss3t-csd.html . Note this now also includes a weighting mechanism for b=0 images. This is (part of the) key to controlling the "speed" of the algorithm somewhat. I noticed at some point your notebooks, basically finding what I was looking at every since the initial abstract, i.e. that the algorithm is best "stopped" at a certain point. So this version deals somewhat with that aspect. Note that the goal is also not to reproduce MSMT-CSD results at all times anymore: there's cases where SS3T-CSD fits the high b-value signal itself better if other shells don't come into play. This can lead to interesting findings ( https://www.biorxiv.org/content/10.1101/629873v1 ).
• New response function calibration is mentioned on that page too, but the essential information sits in the abstract ( https://www.researchgate.net/publication/331165168_Improved_white_matter_response_function_estimation_for_3-tissue_constrained_spherical_deconvolution ) and the ISMRM 2019 talk for some added visual intuition ( https://youtu.be/7yPSFgLt8CA ).
• Both response function calibration and SS3T-CSD are open source Python scripts, if you need any details on how things are exactly done in practice; though these scripts rely heavily on calling mrtrix commands (so might not be all that pure Python-readable, if that makes sense).

Note the SS3T-CSD implementation, while being exactly the mechanism of the 2016 abstract, does switch things a bit around in terms of what "delineates" an iteration (it starts with an "initialisation", which is essentially half an iteration, and then starts counting the next iteration as "iteration 1"). And then there's the added weighting of course as well.

I'll probably drop a bit of a wider announcement on Twitter or something within the next 24 hours; but I thought you might want to be among the first to have access to all info regardless. ;-)

Cheers,
Thijs

Le Bihan, Denis, et al. "Separation of diffusion and perfusion in intravoxel incoherent motion MR imaging." Radiology 168.2 (1988): 497-505.