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

Add positivity constraint on the propagator #338

Merged
merged 5 commits into from Apr 11, 2014

Conversation

Projects
None yet
4 participants
@maurozucchelli
Contributor

maurozucchelli commented Mar 13, 2014

This pull request add the possibility to enforce positiveness in the EAP in the set of points defined by a grid of a certain radius. The number of points depend on the size of the grid. This option improves the results but makes the fitting a lot slower than without it, so it is better to chose the size of the grid wisely.

G = None
h = None
else:
G = self.cache_get('shore_matrix_G', key=(self.pos_grid, self.pos_radius))

This comment has been minimized.

@demianw

demianw Mar 14, 2014

Contributor

I would give the cached matrix a more declarative name than shore_matrix_G. Something like shore_matrix_positive_constraint

This comment has been minimized.

@Garyfallidis
pos_grid : int,
Grid that define the points of the EAP in which we want to enforce
positiveness, if None no constrain is imposed. The parameter
constrain_e0 must be set to True.

This comment has been minimized.

@demianw

demianw Mar 14, 2014

Contributor

Maybe this warning "The parameter constrain_e0" must be set to true could be checked in the code and through a ValueError exception in the case that pos_grid > 0 and constrain_e0 is False

Personally, I rather have a constrain_positivity=False and the pos_grid and pos_radius with default values that make sense.

Then modify the code below such that constrain_e0 and constrain_positivity can be used separately.

This comment has been minimized.

@maurozucchelli

maurozucchelli Mar 17, 2014

Contributor

I also thought of adding another flag for the positiveness but I think that too many parameters make the code less simple to use and to understand.
We can do it of course but in the case I have to copy all your fitting down below adding only the positivity which is quite redundant... I will give it a try and push it if I see that is better!

@demianw

This comment has been minimized.

Contributor

demianw commented Mar 14, 2014

Nice job, I just made a couple comments

@@ -108,6 +110,13 @@ def __init__(self,
square root of the b-value.
constrain_e0 : bool,
Constrain the optimization such that E(0) = 1.
pos_grid : int,
Grid that define the points of the EAP in which we want to enforce
positiveness, if None no constrain is imposed. The parameter

This comment has been minimized.

@arokem

arokem Mar 14, 2014

Member

Typos => "positivity", "no constraint"

positiveness, if None no constrain is imposed. The parameter
constrain_e0 must be set to True.
pos_radius : float,
Radius of the grid of the EAP in which enforce positiveness in

This comment has been minimized.

@arokem

arokem Mar 14, 2014

Member

Typo => "positivity"

v,t = create_rspace(self.pos_grid, self.pos_radius)
lg = int(np.floor(self.pos_grid**3 / 2))
psi = shore_matrix_pdf(self.radial_order, self.zeta, t[:lg])#2456])

This comment has been minimized.

@arokem

arokem Mar 14, 2014

Member

Please remove the "#2456" etc. We don't need that here, right?

@Garyfallidis Garyfallidis changed the title from Add positivity constrain on the propagator to Add positivity constraint on the propagator Mar 17, 2014

@@ -108,6 +110,13 @@ def __init__(self,
square root of the b-value.
constrain_e0 : bool,
Constrain the optimization such that E(0) = 1.

This comment has been minimized.

@Garyfallidis

Garyfallidis Mar 17, 2014

Member

typo: constraint

psi = shore_matrix_pdf(self.radial_order, self.zeta, t[:lg])
G = cvxopt.matrix(-1*psi)
self.cache_set('shore_matrix_G', (self.pos_grid, self.pos_radius), G)
h = cvxopt.matrix((1e-10)*np.ones((lg)),(lg,1))

This comment has been minimized.

@Garyfallidis

Garyfallidis Mar 17, 2014

Member

pep8: add spaces between operators.

@Garyfallidis

This comment has been minimized.

Member

Garyfallidis commented Mar 17, 2014

@maurozucchelli when you say it makes it slower can you give us some timings? How much slower and perhaps how much better the result is? Is there anything that can be done to make it faster can the psi matrix be smaller or something? Or can some more things be cached in advance?

@maurozucchelli

This comment has been minimized.

Contributor

maurozucchelli commented Mar 18, 2014

@Garyfallidis I made some tests:
Creating the matrix takes 1 second for a grid of size 11 and up to 30 seconds for a grid of size 35, but the matrix is then cached so this operation is made only for the first voxel than takes 5e-06 seconds to take the matrix from the cache for the other voxels.
The problem is the optimization, here a list of times:
Gridsize | Solution time (s)
11 | 0.07 per voxel
17 | 0.26 per voxel
21 | 0.45 per voxel
35 | 1.70 per voxel
For a full brain a gridsize of 11 makes the fitting last 2:30 hours which is ok I think.
The positive constraint makes avoid the presence of negative peaks on the ODF, the fitting is more robust to noise and all the maps looks better.

maurozucchelli added some commits Mar 18, 2014

@Garyfallidis

This comment has been minimized.

Member

Garyfallidis commented Apr 11, 2014

Thx!

Garyfallidis added a commit that referenced this pull request Apr 11, 2014

Merge pull request #338 from maurozucchelli/positive_constrain
Add positivity constraint on the propagator

@Garyfallidis Garyfallidis merged commit d5e6426 into nipy:master Apr 11, 2014

1 check passed

default The Travis CI build passed
Details
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment