Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Gtm tortuosity #84

Closed
wants to merge 3 commits into from
Closed

Conversation

matteofrigo
Copy link
Contributor

Tortuosity is affected by the signal fraction bias. Let's do something about it.

This commit adds the `s0_responses` parameter to the
`set_tortuous_parameter` function in the `DistributeModel` class.
This allows a correct estimation of the volume fraction of the IC
compartment which is not affected by the signal fraction bias.
This commit changes the way in which the S0 of every tissue is passed to
the model constructors. Now it is allowed to pass None, if the
corresponding model has been defined using the `S0_tissue_responses`
parameter. All of this happens in the `MultiCompartmentModel.fit`
function.

Other changes involve the proper definition of the tortuosity constraint
that have been addressed in the parent commits.
# and the corresponding S0 is treated as the weighted
# average of the S0 tissues modelled by the "submodels".
s0 = 0.0
for k, mm in enumerate(m.models):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Nicely done, but while this may work here at the top level to the second level - just for this project - it means that the top level now has the logic to infer what is the S0 of a lower level distributed model.

In my opinion this is not so clean, as the we should be able to ask the distributed model what is it's S0 based on its signal volume fractions. This has the added benefit of it being to (potentially) generalize to more than 2 layers - we want to have the code as modular as possible.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It makes sense.
If I understand the flow correctly, this could be done by moving the definition of S0_responses to the __init__ method, where self.S0_tissue_responses will be defined inferring the different S0s from the "submodels".
Am I right?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What I mean is that we should be able to say S0 = m.S0_response(volume_fractions) to any submodel, instead of having to manually sum up the s0 values in line 1229.

Hold on, I how is line 1229 even getting the 'partial_volume' values? mm is just a DistributedModel instance right? not a dict? I don't thing this works like this, or am I wrong?

@@ -1212,8 +1216,24 @@ def fit(self, acquisition_scheme, data,
start = time()
mt_fractions = np.empty(
np.r_[N_voxels, self.N_models], dtype=float)

S0_responses = []
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We know what will be the length of S0_responses, so it's cleaner to initiate it just with S0_responses = np.zeros(self.N_models) or something.

@@ -275,7 +275,8 @@ def set_fixed_parameter(self, parameter_name, value):

def set_tortuous_parameter(self, lambda_perp,
lambda_par,
volume_fraction_intra):
volume_fraction_intra,
S0_responses=None):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think the S0_responses should not be given as a parameter of the set_tortuous_parameter function, but should be given to the SD1_WatsonDistributedModel as a second argument, and before a property as self.S0_responses, just as in the upper MC-models.

You could, potentially, add a parameter here like use_S0_debiasing_if_possible or something to just be able to choose here.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For the moment I would opt for self.S0_tissue_responses for consistency with other functions.
Do you like set_tortuous_parameter(..., use_tissue_S0=True) for triggering it?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ah of course keep things consistent. I just meant that I didn't like the multiple places to give S0_tissue_responses, which would make it possible to do things differently in different places as the same time.

perfect on the option.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

See my comment above, we do need to find a way to nicely define which model is associated with what S0 response for intra or extra.

Perhaps we can make the signature even simpler and do def set_tortuous_parameter(model_intra, model_extra, use_tissue_S0=True), and accept the model definitions themselves as input parameters, like

stick=C1Stick()
zepp = G2Zeppelin()
watson = SD1WatsonDistributedModel(models = [stick, sepp], S0s=[1, 2])
watson.set_tortuous_parameter(model_intra=stick, model_extra=zepp, use_tissue_S0=True)

Internally we can just check for the appropriate model with isinstance and get the names for the parameters and associated S0s... what do you think?

@@ -186,6 +186,8 @@ def T1_tortuosity(lambda_par, vf_intra, vf_extra=None):
intra-axonal volume fraction [0, 1].
vf_extra : float, (optional)
extra-axonal volume fraction [0, 1].
S0_responses : list, (optional)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

indicate the length of the list

fraction_intra = vf_intra
fraction_extra = 1.0 - fraction_intra if vf_extra is None else vf_extra

if S0_responses is not None:
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

add detailed comment either here or in the doctring of the tortuosity function here - add a preliminary citation to your abstract if you want.

@@ -558,7 +558,8 @@ def _add_fixed_parameter_array(self, parameter_name, parameter_array):
def set_tortuous_parameter(self, lambda_perp_parameter_name,
lambda_par_parameter_name,
volume_fraction_intra_parameter_name,
volume_fraction_extra_parameter_name):
volume_fraction_extra_parameter_name,
S0_responses=None):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Again (see below) I don't think we should be able to give S0 responses twice for any MC-model - it should be unique in the definition.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a shortcut for knowing which of the S0s in self.S0_tissue_responses list corresponds to the intra and which one to the extra?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I didn't even thing about this when I first coded this... I especially kept the parameter names as input because I didn't want to have to hardcode the ordering or definitions of models to get the tortuosity regardless.

Actually we do need to find a solution here - what you were doing with adding the S0_responses as a parameter may be the easiest way to make this work after all... Either this or we need to define a nice way to just point to the models that are intra and extra, and then make the function just find which parameters and S0s are associated with this...

@@ -1212,8 +1216,24 @@ def fit(self, acquisition_scheme, data,
start = time()
mt_fractions = np.empty(
np.r_[N_voxels, self.N_models], dtype=float)

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You need to add this also for MC-SM models.

Copy link
Collaborator

@rutgerfick rutgerfick left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

see the individual comments.

@matteofrigo
Copy link
Contributor Author

This PR is replaced by #110 .

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

2 participants