Skip to content

Issue 918 bregman functional #1267

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

Merged
merged 14 commits into from
Feb 12, 2018

Conversation

aringh
Copy link
Member

@aringh aringh commented Jan 14, 2018

Finally found some time to implement a Bregman distance functional, see #918.

  • Add Bregman distance functional.
  • Add it to the general test suit for functionals.
  • Add more specific tests, checking e.g. the use of a separate subgradient operator
  • Add example of using Bregman functional, e.g. an example of Bregman-TV

@pep8speaks
Copy link

pep8speaks commented Jan 14, 2018

Checking updated PR...

No PEP8 issues.

Comment last updated on February 12, 2018 at 13:45 Hours UTC

@aringh
Copy link
Member Author

aringh commented Jan 14, 2018

This is good to go for review. An example can come later, it doesn't even have to be in this PR.

In #918, @ozanoktem suggested having something like functional.bregman. Is this something we want to pursue?

Copy link
Member

@adler-j adler-j left a comment

Choose a reason for hiding this comment

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

Exceptional PR, easy to read and to the point with all required components.

I basically nit-picked on some style issues and some small implementation stuff

in a point :math:`x` is given by

.. math::
D_f(x) = f(x) - f(y) - \langle \partial f(y), x - y \\rangle.
Copy link
Member

Choose a reason for hiding this comment

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

Missing y on L.H.S.?

optional argument `subgradient_op` is not given, the functional
needs to implement `functional.gradient`.
point : element of ``functional.domain``
The point from which to define the Bregman distance
Copy link
Member

Choose a reason for hiding this comment

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

Missing .

--------
Example of initializing the Bregman distance functional:

>>> space = odl.uniform_discr(0, 2, 14)
Copy link
Member

Choose a reason for hiding this comment

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

14 is rather arbitrary, prefer 10 or so

Copy link
Member Author

Choose a reason for hiding this comment

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

It works with 8, 9, 11, 12, 13, and 14 points... but with 10 points apparently Bregman_dist(space.zero()) returns 1.9999999999999998 and therefore fails... is there a way to make "almost equal" in doc-tests?

>>> point = space.one()
>>> Bregman_dist = odl.solvers.BregmanDistance(l2_squared, point)

This is gives the shifted L2 norm squared ||x - 1||:
Copy link
Member

Choose a reason for hiding this comment

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

missing ^2

''.format(functional))
else:
# Check that given subgradient is an operator that maps from the
# domain of the functional to the domain of the functional
Copy link
Member

Choose a reason for hiding this comment

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

to itself

def __repr__(self):
'''Return ``repr(self)``.'''
return '{}({!r}, {!r}, {!r}, {!r})'.format(self.__class__.__name__,
self.domain,
Copy link
Member

Choose a reason for hiding this comment

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

This argument is not in init

Copy link
Member

Choose a reason for hiding this comment

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

It's more maintainable to use the helper functions for this:

posargs = [functional, point]
subgrad_op_default = getattr(self.functional, 'gradient', None)
optargs = [('subgradient_op', self.subgradient_op, subgrad_op_default)]
inner_str = signature_string(posargs, optargs, sep=',\n')
return '{}(\n{}\n)'.format(self.__class__.__name__, indent(inner_str))

The utilities signature_string and indent are in odl.util.
(BTW, I'm working on making this MUCH better).

def test_bregman_functional_no_gradient(space):
"""Test that the Bregman distance functional fails if the underlying
functional does not have a gradient and no subgradient operator is
given."""
Copy link
Member

Choose a reason for hiding this comment

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

Move """ to newline

Copy link
Member

Choose a reason for hiding this comment

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

Also first line should be "sufficient" and fit in one line, the rest should elaborate


def test_bregman_functional_l2_squared(space, sigma):
"""Test for the Bregman distance functional, using l2 norm squared as
underlying functional."""
Copy link
Member

Choose a reason for hiding this comment

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

""" on newline

assert all_almost_equal(bregman_dist(x), expected_func(x))

# Gradient evaluation
assert all_almost_equal(bregman_dist(x), expected_func(x))
Copy link
Member

Choose a reason for hiding this comment

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

Copy-paste-error here? No gradient


l2_sq = odl.solvers.L2NormSquared(space)
point = noise_element(space)
subgrad_op = odl.ScalingOperator(space, 2.0)
Copy link
Member

Choose a reason for hiding this comment

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

Why does this need to be provided?

Copy link
Member Author

Choose a reason for hiding this comment

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

It does not need to be provided. It is a relic since I wrote this test before the other test, and wanted to check that it in fact return the expected gradient also when a subdifferential operator is given. But given the test above that now tests this, I will remove it.

@aringh aringh force-pushed the issue-918__bregman_functional branch from 89de84d to a2e2f75 Compare January 27, 2018 13:10
@aringh
Copy link
Member Author

aringh commented Jan 27, 2018

Question: do we want to be able to update the point of the Bregman distance? In iterative reconstructions using the Bregman distance, as far as I know one normally updates this point between iterations. This could of course be achieved by simply creating a new functional each time one needs to update the point, but it should also be possible to do something like

    @point.setter
    def point(self, new_point):
        self.__point = new_point
        self.__subgrad_eval = self.subgradient_op(self.__point)
        self.__constant = (-self.__functional(self.__point) +
                           self.__subgrad_eval.inner(self.__point))
        self.__bregman_dist = FunctionalQuadraticPerturb(
            self.__functional, linear_term=-self.__subgrad_eval,
            constant_coeff=self.__constant)

Is this to much, or do we want to go for that?

Apart from this, I think that the other issues have been addressed in the last commits.

@ozanoktem
Copy link
Contributor

Updating the point of the Bregman distance will as @aringh mentioned occur when using it in iterative reconstructions. I don't know what overhead one generates by creating a new functional at each iterate, but I think it is sensible to be able to update the point.

@aringh
Copy link
Member Author

aringh commented Jan 27, 2018

I think that the overhead is approximately the same between this way and if the user creates a completely new Bregman distance functional. At least with this code since I create a new FunctionalQuadraticPerturb when the point is updated. For me it is more a question of user-friendliness vs. to much built-in intelligence. I think it looks nice from a user-persective, but I would like to hear what @kohr-h and @adler-j has to say about it.

A "counter example", maybe advocating that we should not do this, is that we for example create a new proximal operator each time the step size sigma is updated.

EDIT: if you like the idea, here is an implementation which I can just merge into the current branch.

Copy link
Member

@kohr-h kohr-h left a comment

Choose a reason for hiding this comment

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

The comments are by and large of cosmetic nature.

Regarding the question about point.setter: I'm skeptical to mutating a Functional without a clear need. I definitely prefer making a new functional each time.

That brings me to a larger point: I'm not convinced that the Bregman functional should be a class. There are two main issues leading me to this conclusion:

  1. The initialization of the Bregman functional can be quite expensive since it involves evaluation of the original functional and its gradient (or the subgradient operator).
  2. The operator more or less only stores the quadratic perturbation and caches some stuff.

Regarding the first point, I'd favor an implementation that initializes the two performance-critical caches to None and sets them during the first _call.

However I think we should also consider the alternative of providing a function instead of a class, similar to what we do with proximal. A baseline implementation would be this one, and users (or concrete functionals) could choose to override it with something optimized.

I'm also in favor of adding a bregman method (or property) to Functional due to its generic nature and usefulness.

Opinions?

:math:`D_f(\cdot, y)` in a point :math:`x` is given by

.. math::
D_f(x, y) = f(x) - f(y) - \langle \partial f(y), x - y \\rangle.
Copy link
Member

Choose a reason for hiding this comment

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

I think you need \\langle

Copy link
Member

Choose a reason for hiding this comment

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

See my comment on #1243: If you prefix the docstring with an r (for "raw string") you can write normal LaTeX without using double-backslashes all the time.

[Bur2016] Burger, M. *Bregman Distances in Inverse Problems and Partial
Differential Equation*. In: Advances in Mathematical Modeling, Optimization
and Optimal Control, 2016. p. 3-33. Arxiv version online at
https://arxiv.org/abs/1505.05191
Copy link
Member

Choose a reason for hiding this comment

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

Suggestion: Turn this "Arxiv version ..." sentence into a clickable link in the text above, like

For mathematical details, see `[Bur2016]<https://arxiv.org/abs/1505.05191>`_ .

and leave the rest as-is. That way, you have the clickable link in the text and the full citation in the "References".

https://arxiv.org/abs/1505.05191
"""

def __init__(self, functional, point, subgradient_op=None):
Copy link
Member

@kohr-h kohr-h Jan 27, 2018

Choose a reason for hiding this comment

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

I would pledge for the shorter subgrad instead of subgradient_op. The docstring makes clear that it's an Operator.

The point from which to define the Bregman distance.
subgradient_op : `Operator`, optional
The operator that takes an element in `functional.domain` and
returns a subgradient of the functional in that point.
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps mention that it is a subgradient, i.e., an arbitrary element from the subdifferential.

Parameters
----------
functional : `Functional`
The functional on which to base the Bregman distance. If the
Copy link
Member

Choose a reason for hiding this comment

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

We tend to not use "The" in the parameter list, only in the "Returns" because there things are more narrowly specified.


@property
def functional(self):
"""The underlying functional for the Bregman distance."""
Copy link
Member

Choose a reason for hiding this comment

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

underlying, same below


def __init__(self, func, quadratic_coeff=0, linear_term=None):
def __init__(self, func, quadratic_coeff=0, linear_term=None,
constant_coeff=0):
Copy link
Member

Choose a reason for hiding this comment

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

constant or constant_term

@@ -1025,9 +1039,14 @@ def convex_conj(self):
the convex conjugate of :math:`f`:

.. math::
(f(x) + <y, x>)^* (x) = f^*(x - y).
(f(x) + <y, x>)^* (x^*) = f^*(x^* - y).
Copy link
Member

Choose a reason for hiding this comment

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

(f + \\langle y, \cdot \\rangle)^*

conjugate of :math:`f + c` is by definition

.. math::
(f(x) + c)^* (x^*) = f^*(x^*) - c.
Copy link
Member

Choose a reason for hiding this comment

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

f(x) -> f

(f(x) + <y, x>)^* (x) = f^*(x - y).
(f(x) + <y, x>)^* (x^*) = f^*(x^* - y).

For reference on the identity used, see [KP2015]. Moreover, the convex
Copy link
Member

Choose a reason for hiding this comment

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

A clickable link would be cool.

@ozanoktem
Copy link
Contributor

After thinking some on how ODL should handle notion of Bregman distance, and despite the arguments provided earlier by @adler-j, I think the best option would be to introduce it as a property to the functional. An advantage is that the Bregman distance respects functional calculus much like the derivative. Ideally, if a developer does not explicitly write a routine for computing the Bregman distance for a given functional, then ODL should automatically resort to the definition for computing it.

@aringh
Copy link
Member Author

aringh commented Jan 31, 2018

Then I will just move the BregmanDistance class to functional.py, not expose it in __all__, and add a bregman-property to the Functional class (that works much like proximal). Sounds good?

Also: should I leave the option of supplying a subgradient?

@kohr-h
Copy link
Member

kohr-h commented Jan 31, 2018

Also: should I leave the option of supplying a subgradient?

That's a very good question. One way to handle that would be to do it in a similar way as in derivative, where you supply the point and get back an operator for that point. For such a construction, it would be easy to add an option for a custom subgradient.
If we do it the proximal way, i.e., as a property that returns a factory function, it's also possible if that function accepts an optional subgrad argument.

I actually prefer the more lightweight solution (the derivative way), since I think in hindsight that the proximal factory stuff is too complicated for what it does. We could actually have something like in mathematics: proximal is the prox operator for step size 1.0, and the user (or optimization method) simply does tau * func.proximal. We know how to construct that guy as a proximal_scaling.

Sounds good?

To me it does, but maybe we should have more votes on this. I'm thinking of @adler-j and @mehrhardt in particular.

@adler-j
Copy link
Member

adler-j commented Jan 31, 2018

That's a very good question. One way to handle that would be to do it in a similar way as in derivative, where you supply the point and get back an operator for that point. For such a construction, it would be easy to add an option for a custom subgradient.
If we do it the proximal way, i.e., as a property that returns a factory function, it's also possible if that function accepts an optional subgrad argument.

Lets go for the functional.bregman implementation, people seem to be in favor of it and it looks more clean.

I actually prefer the more lightweight solution (the derivative way), since I think in hindsight that the proximal factory stuff is too complicated for what it does. We could actually have something like in mathematics: proximal is the prox operator for step size 1.0, and the user (or optimization method) simply does tau * func.proximal. We know how to construct that guy as a proximal_scaling.

I actually tried implementing this at some point, wasting quite a bit of time, and it turns out that it is non-trivial to do it as you say. Basically this line:

return self.functional.proximal(sigma * self.scalar)

is what screws with us. In this one, FunctionalLeftScalarMult.proximal calls Functional.proximal but with a scaled argument. Basically one cannot compute the proximal of a * F in closed form given only a and the proximal of F.

To me it does, but maybe we should have more votes on this. I'm thinking of @adler-j and @mehrhardt in particular.

Go ahead

@kohr-h
Copy link
Member

kohr-h commented Jan 31, 2018

Basically one cannot compute the proximal of a * F in closed form given only a and the proximal of F.

Ah crap, you need the prox of the right-scaled functional which you cannot derive from the prox with scalar 1. Alright, never mind.

@aringh
Copy link
Member Author

aringh commented Feb 2, 2018

I have updated the PR. Good to go for a new review round :-)

I also tried to add an example for Bregman-TV in a separate branch here. If we want such an example, I would appreciate if someone could take a look at it. Maybe @ozanoktem that knows how it is supposed to work?

Copy link
Member

@adler-j adler-j left a comment

Choose a reason for hiding this comment

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

Beautiful pull request, just some minor comments w.r.t. style and we'll go!

----------
point : element of ``functional.domain``
Point from which to define the Bregman distance.
subgrad : `Operator`, optional
Copy link
Member

Choose a reason for hiding this comment

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

Am i too rough if i propose we call this gradient? Throughout all of ODL we use gradient while we really mean sub-gradient. Even the default hints at this

Copy link
Member Author

Choose a reason for hiding this comment

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

We can go for gradient, no problem with me :-)

Copy link
Member Author

Choose a reason for hiding this comment

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

Actually, trying to go through and making this change it turns out that we would probably like to distinguish it from the gradient. Because this is the (sub)gradient of the "underlying" functional. Meanwhile, the BregmanDistance is going to have a gradient of its own (which is a translated version of what is now called subgrad, see line 1505), and we don't want to confusing these two.

Any other naming suggestions? gradient_op? Or stick with subgrad?

Copy link
Member

Choose a reason for hiding this comment

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

I'm not in favor of changing subgrad to gradient anyway, it's less precise IMO. To make clear that it's a property of the underlying functional, we could call it func_subgrad or so.

"""Initialize a new instance.

Parameters
----------
func : `Functional`
Function corresponding to ``f``.
quadratic_coeff : ``domain.field`` element, optional
Coefficient of the quadratic term.
Coefficient of the quadratic term. Default: zero.
Copy link
Member

Choose a reason for hiding this comment

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

zero -> 0

"zero" may be interpreted in so many ways, e.g. zero vector, zero real, etc. 0 is a more clear to me

constant = func.domain.field.element(constant)
if constant.imag != 0:
raise ValueError(
"Complex-valued constant coefficient is not supported.")
Copy link
Member

Choose a reason for hiding this comment

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

Add ` around "constant" to indicate it is a parameter

self.functional,
self.quadratic_coeff,
self.linear_term)
return '{}({!r}, {!r}, {!r}, {!r})'.format(self.__class__.__name__,
Copy link
Member

Choose a reason for hiding this comment

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

Perhaps @kohr-h could give some hints on how to use the utility here

Copy link
Member

Choose a reason for hiding this comment

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

Nah, don't bother for now. I'm refurbishing this stuff anyway in #1276

@@ -1268,6 +1333,186 @@ def __str__(self):
return '{}.convex_conj'.format(self.convex_conj)


class BregmanDistance(Functional):
Copy link
Member

Choose a reason for hiding this comment

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

Wouldn't "BregmanDivergence" be a more correct term?

Copy link
Member

Choose a reason for hiding this comment

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

I had the same thought initially, but the main reference (Burger et al (?)) calls it Bregman distance. Divergence leans more on the statistical side I assume, so IMO let's go with the terms the inverse problems and imaging people use.

self.__subgrad = subgrad

self.__initialized = False
self.__subgrad_eval = None
Copy link
Member

Choose a reason for hiding this comment

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

I don't think there is any reason to set these to None?


def _initialize(self):
"""Helper function to initialize the functional."""

Copy link
Member

Choose a reason for hiding this comment

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

Remove extra blankline

Copy link
Member

Choose a reason for hiding this comment

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

Also there is a very nice pattern that you can use to save some lines. Rename this to initialize_if_needed and then move the if not self.initialized line inside here. That way you can save some duplication.

def __repr__(self):
'''Return ``repr(self)``.'''
posargs = [self.functional, self.point]
subgrad_op_default = getattr(self.functional, 'gradient', None)
Copy link
Member

Choose a reason for hiding this comment

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

Functional always implements gradient (since it is an Functional), so the getattr here is redundant

Copy link
Member Author

Choose a reason for hiding this comment

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

I do not think I fully understand how signature_string works...

But even the current implementation can get into trouble. If I uses, e.g., an IndicatorBox which I supplied with an appropriate subgradient operator in the creation of the BregmanDistance as follows

>>> import odl
>>> space = odl.uniform_discr(0,1,10)
>>> indic_func = odl.solvers.IndicatorBox(space, 0.2, 0.8)
>>> point = 0.5 * space.one()
>>> subgrad = odl.ZeroOperator(space)
>>> breg = odl.solvers.BregmanDistance(indic_func, point, subgrad)
>>> breg.__repr__()

I get the following output

Traceback (most recent call last):

  File "<ipython-input-1-dcb69fc015ba>", line 7, in <module>
    breg.__repr__()

  File "[...]/odl/odl/solvers/functional/functional.py", line 1503, in __repr__
    subgrad_op_default = getattr(self.functional, 'gradient', None)

  File "[...]/odl/odl/solvers/functional/functional.py", line 97, in gradient
    ''.format(self))

NotImplementedError: no gradient implemented for functional IndicatorBox(uniform_discr(0.0, 1.0, 10), 0.2, 0.8)

However self.func_subgrad will, by the way that __init__ works, always have a value. Either it is supplied separately or it is set to be the gradient. Thus one could just change to

        optargs = [('func_subgrad', self.func_subgrad, None)]

where None will never be used, right?

Copy link
Member

Choose a reason for hiding this comment

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

This solution should be fine. It will make repr always print out the func_subgrad even if the user didn't specify it, but the logic is complex enough that you may want to skip it.
If you don't want this behavior, you probably need to do the gettattr thing to get the default and use that instead of None in optargs.

ind_func = odl.solvers.IndicatorNonnegativity(space)
point = noise_element(space)

# Indicator function has no gradient, hence one cannot create a bregman
Copy link
Member

Choose a reason for hiding this comment

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

Frankly we should just add gradient=0 for the indicator function. It makes sense, no?

This would be for later though

Copy link
Member Author

Choose a reason for hiding this comment

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

The indicator function on a set has gradient=0 on that set. But outside the set, where it takes value inf, I think that the subgradient is empty. So we could return EmptySet in these cases, but I don't know if that would make anyone happier or if it would just confuse people...?

Copy link
Member

Choose a reason for hiding this comment

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

I think what you're talking about @aringh is the subdifferential, not subgradient. The subgradient is not a set, it's an (arbitrary) element from the subdifferential. So rather than returning an empty set, the subgrad operator should in principle fail if we're outside the set of the indicator function, since logically, we're trying to return an element from the empty set.

I don't know if that trouble is worth going through. Is there a case where we would need this?

Copy link
Member

Choose a reason for hiding this comment

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

Frankly we should just add gradient=0 for the indicator function. It makes sense, no?

I think this is mathematically wrong since 0 is not a subgradient in points outside the set of the indicator function. So IMO we shouldn't do this.

# Indicator function has no gradient, hence one cannot create a bregman
# distance functional
with pytest.raises(NotImplementedError):
odl.solvers.BregmanDistance(ind_func, point)
Copy link
Member

Choose a reason for hiding this comment

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

We need to do tests for both functiona.bregman and the pure functional odl.solvers.BregmanDistance, at least one that verifies that it points to the later I'd say.

Copy link
Member Author

Choose a reason for hiding this comment

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

My idea, when keeping these test after implementing the functional.bregman, was that they still test the BregmanDistance class. Not for all functionals, but some kind of test for the class itself.

I also added a test for functional.bregman to the suite of tests for functionals, here. My idea was that this should cover all available functionals and see that there .bregman works as expected. By that we also cover the case if one does a special optimized implementation of the Bregman distance for some functionals later (e.g., the squared L2 norm should just be a translated functional).

Is that good enough, or would you like to add some other tests or change some of the current ones?

@adler-j
Copy link
Member

adler-j commented Feb 2, 2018

W.r.t. the bregman TV example, I think it looks quite good but I'd make a separate PR after this one for it. Overall, if it works I'd like it in!


@property
def subgrad(self):
"""The subgradient operator for the Bregman distance."""
Copy link
Member Author

Choose a reason for hiding this comment

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

Related to a comment on the (sub)gradient of the underlying functional and the gradient of the Bregman distance: from this description it is not obvious what is meant. It is suppose to be the gradient of the underlying functional, Suggestion: change to The ... used in the Bregman distance for all three above. Ok?

Copy link
Member

Choose a reason for hiding this comment

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

Maybe "used to define the Bregman distance"? Generally a good suggestion.

@aringh
Copy link
Member Author

aringh commented Feb 7, 2018

The question about tests is not addressed, but for that I would need some more exact pointers on what needs to be changed. Otherwise, good to go in?

@adler-j
Copy link
Member

adler-j commented Feb 7, 2018

The question about tests is not addressed, but for that I would need some more exact pointers on what needs to be changed. Otherwise, good to go in?

Excuse my bad memory, but what pointers would you like here?

@aringh
Copy link
Member Author

aringh commented Feb 7, 2018

You made a comment about the tests just above

We need to do tests for both functiona.bregman and the pure functional odl.solvers.BregmanDistance, at least one that verifies that it points to the later I'd say.

Summarizing my answer from above: we have two tests for BregmanDistance in default_functionals_test.py, and for all functionals in the test suite in functional_test.py we test the .bregman functionality. I just wander if that is sufficient, or if you want changes in the test code? If so, what kind of changes?

@aringh aringh added the merge? label Feb 8, 2018
@aringh
Copy link
Member Author

aringh commented Feb 12, 2018

I will poke this one: is it good to go in? Other changes?

Also: @adler-j or @kohr-h will have to explicitly merge it.

@adler-j
Copy link
Member

adler-j commented Feb 12, 2018

Could use a rebase on master, other than that I'll have a quick look!

Copy link
Member

@adler-j adler-j left a comment

Choose a reason for hiding this comment

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

All good. Rebase with master and then we go

@kohr-h
Copy link
Member

kohr-h commented Feb 12, 2018

It fails to import build_main from sphinx. I checked, it's a Sphinx 1.7 feature, so travis-sphinx has an implicit dependency on that.

@aringh
Copy link
Member Author

aringh commented Feb 12, 2018

Here is the output from the failing doc-build:

[...]
ufuncs                    : generated odl.util.ufuncs.rst
utility                   : generated odl.util.utility.rst
vectorization             : generated odl.util.vectorization.rst
/home/travis/build/odlgroup/odl
Traceback (most recent call last):
  File "/home/travis/miniconda/envs/testenv/bin/travis-sphinx", line 11, in <module>
    sys.exit(main())
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/click/core.py", line 722, in __call__
    return self.main(*args, **kwargs)
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/click/core.py", line 697, in main
    rv = self.invoke(ctx)
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/click/core.py", line 1066, in invoke
    return _process_result(sub_ctx.command.invoke(sub_ctx))
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/click/core.py", line 895, in invoke
    return ctx.invoke(self.callback, **ctx.params)
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/click/core.py", line 535, in invoke
    return callback(*args, **kwargs)
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/click/decorators.py", line 17, in new_func
    return f(get_current_context(), *args, **kwargs)
  File "/home/travis/miniconda/envs/testenv/lib/python3.6/site-packages/travis_sphinx/build.py", line 39, in build
    if sphinx.build_main(args + [source, outdir]):
AttributeError: module 'sphinx' has no attribute 'build_main'

@kohr-h
Copy link
Member

kohr-h commented Feb 12, 2018

Weird enough, Travis does install sphinx 1.7.0 so that should be fine (still a bummer on the side of travis-sphinx).
Can you try excluding newer versions of travis-sphinx? Using travis-sphinx<2.0 in the pip install command. While editing the Travis config, can you also remove the conda version restriction? The current version should be fine.

@kohr-h
Copy link
Member

kohr-h commented Feb 12, 2018

Ah no, it's the other way around. The build_main is gone from sphinx 1.7 on, so just use sphinx<1.7 in the pip install.

@kohr-h
Copy link
Member

kohr-h commented Feb 12, 2018

See here as well

@kohr-h
Copy link
Member

kohr-h commented Feb 12, 2018

Ah shoot, you need to put it in quotes since the < sign is interpreted by the shell: "sphinx<1.7"

@kohr-h
Copy link
Member

kohr-h commented Feb 12, 2018

Looks like this one will work, except that it may complain about too much console output. Get ready for re-inserting the 2>/dev/null. And take out the conda restriction when you make the final version please.

@aringh
Copy link
Member Author

aringh commented Feb 12, 2018

If this works I will take out the conda restriction and squash all the commits into one. But I just wanted to do one thing at a time, to see what works and what doesn't ;-)

@aringh
Copy link
Member Author

aringh commented Feb 12, 2018

It passed! I will update as asked by @kohr-h and clean up last part of the commit history

Add restriction sphinx<1.7 in .travis.yml. since
build_main is gone from sphinx 1.7.

Also remove conda restriction since fixed.
@aringh aringh force-pushed the issue-918__bregman_functional branch from c8224de to eb9e7d4 Compare February 12, 2018 13:45
@aringh
Copy link
Member Author

aringh commented Feb 12, 2018

Could @kohr-h or @adler-j merge this? :-)

@adler-j adler-j merged commit a0e3a19 into odlgroup:master Feb 12, 2018
@adler-j
Copy link
Member

adler-j commented Feb 12, 2018

Done! Well done

@ozanoktem
Copy link
Contributor

ozanoktem commented Feb 21, 2018

It seems judging from the checkmarks that we still lack an example of a reconstruction algorithm that uses Bregma distances.

In https://link.springer.com/chapter/10.1007/978-3-319-01712-9_2 there is a formulation for the Bregman-)FB-EM-TV algorithm as a nested two step iteration strategy for solving the TV regularized Poisson likelihood estimation problem in eq. (20). The numerical realization of the TV correction half step can be formulated in various ways, see eqs. (29), (31), (45) and (46). A uniform numerical framework for these variants is given in eq. (100) along with table 1 on p. 121, which contains many state-of-the-art schemes as special cases. The (alternating) split Bregman algorithm for solving TV correction half step is given on eps (104a)-(104e). Perhaps we could implement that in ODL and use it to generate an example?

@adler-j
Copy link
Member

adler-j commented Feb 21, 2018

Imo that would be far too complicated for the standard examples. A more simple example is more urgently needed.

With that said, adding that example as well would be nice I guess!

@aringh
Copy link
Member Author

aringh commented Feb 21, 2018

As mentioned higher up in this conversation I have a Bregman-TV example for tomography lying around. But I am not to convinced by the actual performance of the algorithm. I will create a new PR, to make it easier to find, and if @ozanoktem would like to have a look at it that would be great :-)

@aringh
Copy link
Member Author

aringh commented Feb 21, 2018

See PR #1305

----------
point : element of ``functional.domain``
Point from which to define the Bregman distance.
subgrad : `Operator`, optional
Copy link
Contributor

Choose a reason for hiding this comment

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

When I define the Bregman distance, I normally would like to give a specific subgradient and not a subgradient operator that would compute one.

functional :math:`D_f(\cdot, y)` in a point :math:`x` is given by

.. math::
D_f(x, y) = f(x) - f(y) - \\langle \partial f(y), x - y \\rangle.
Copy link
Contributor

Choose a reason for hiding this comment

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

The notation here is non-standard. Usually, one would pick a single subgradient p from the subdifferential \partial f(y) and define D_f^p(x, y) = f(x) - f(y) - \langle p, x - y \rangle.

@adler-j
Copy link
Member

adler-j commented Mar 14, 2018

Thank you for the comments @mehrhardt, I agree with both of those.

@kohr-h
Copy link
Member

kohr-h commented Mar 14, 2018

Without looking into it, I agree with the general idea that we should try to follow standard notation and make an implementation in accordance to that.

I made a new issue #1312, please continue the discussion there.

@aringh aringh mentioned this pull request May 9, 2018
@aringh aringh deleted the issue-918__bregman_functional branch May 11, 2018 14:25
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.

6 participants