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

BF: Keep up with changes in scipy 0.16 #707

Merged
merged 7 commits into from Sep 14, 2015

Conversation

Projects
None yet
3 participants
@arokem
Member

arokem commented Sep 1, 2015

Recent changes in scipy use a different heuristic to find the epsilon
parameter used for RBF interpolation, causing for us issue #688. This
should fix it, by using the old heuristic.

@arokem

This comment has been minimized.

Member

arokem commented Sep 1, 2015

Note - the first commit (arokem@7c4e670) are the substantial changes. The rest is just PEP8 fixes that I did while I was in this file already.

@arokem

This comment has been minimized.

Member

arokem commented Sep 12, 2015

Oh dear. Hey @stefanv - any ideas here? Looks pretty bad.

@arokem

This comment has been minimized.

Member

arokem commented Sep 12, 2015

Just to be clear on this one - these tests would have failed with previous versions of scipy as well (this PR basically just patches this up to the previous behavior, but adds more tests). It seems that RBF is just inaccurate in many use-cases, which we weren't thinking of before.

@Garyfallidis

This comment has been minimized.

Member

Garyfallidis commented Sep 12, 2015

Oy! This is strange. Does that mean that all the functionality with interpolating on spheres is now uncertain?

@arokem

This comment has been minimized.

Member

arokem commented Sep 12, 2015

Well - it was always uncertain. It's just that we now know for certain just
how uncertain it was...

On Sat, Sep 12, 2015 at 4:39 PM, Eleftherios Garyfallidis <
notifications@github.com> wrote:

Oy! This is strange. Does that mean that all the functionality with
interpolating on sphere is now uncertain?


Reply to this email directly or view it on GitHub
#707 (comment).

@Garyfallidis

This comment has been minimized.

Member

Garyfallidis commented Sep 12, 2015

Okay so what is your suggestion. How to move forward with this?

@arokem

This comment has been minimized.

Member

arokem commented Sep 12, 2015

A few options (in order of preference, I think):

  • Spend some time understanding RBF interpolation, and fix it up (might
    need to fix the scipy implementation, see links to discussions there - IIUC
    that epsilon parameter is pretty much made up there as well...).
  • Put a big fat warning on RBF for users to check the results of
    interpolation when they use this, and relax test requirements.
  • Drop RBF interpolation altogether.

I can spend some time on option 1, but not before this Thursday. Even then,
progress might be slow. What do you think?

On Sat, Sep 12, 2015 at 4:44 PM, Eleftherios Garyfallidis <
notifications@github.com> wrote:

Okay so what is your suggestion. How to move forward with this?


Reply to this email directly or view it on GitHub
#707 (comment).

@Garyfallidis

This comment has been minimized.

Member

Garyfallidis commented Sep 13, 2015

Hi @arokem,
here is some feedback from my investigations.

The problem at least at my scipy 0.14 seems to be only with norm="euclidean_norm"
(norm="angle" works fine) which can be resolved using a small smoothing (here using 0.1).

Check my code snippet here

def test_interp_rbf():
    def data_func(s, a, b):
        return a * np.cos(s.theta) + b * np.sin(s.phi)

    from dipy.core.sphere import Sphere, interp_rbf
    from dipy.core.subdivide_octahedron import create_unit_hemisphere
    import numpy as np

    s0 = create_unit_hemisphere(5)
    s1 = create_unit_hemisphere(6)

    from dipy.viz import fvtk
    s0 = Sphere(xyz = np.vstack((s0.vertices, -s0.vertices)))
    s1 = Sphere(xyz = np.vstack((s1.vertices, -s1.vertices)))


    for a, b in zip([1, 2, 0.5], [1, 0.5, 2]):
        ren = fvtk.ren()
        print(a, b)
        data = data_func(s0, a, b)
        expected = data_func(s1, a, b)
        interp_data_en = interp_rbf(data, s0, s1, smooth=.1, norm="euclidean_norm")
        interp_data_a = interp_rbf(data, s0, s1, smooth=.1, norm="angle")

        sf = np.zeros((3, s1.vertices.shape[0]))
        sf[0] = expected
        sf[1] = interp_data_en
        sf[2] = interp_data_a
        fvtk.add(ren, fvtk.sphere_funcs(sf, sphere=s1, norm=False, radial_scale=False))
        nt.assert_(np.mean(np.abs(interp_data_en - expected)) < 0.1)
        nt.assert_(np.mean(np.abs(interp_data_a - expected)) < 0.1)

        fvtk.show(ren)


if __name__ == "__main__":
    # nt.run_module_suite()
    test_interp_rbf()

image

Now if you change to the default smooth (zero) you get some weird result looking

image

Some other important points to improving the tests. I think those tests should use higher sampled spheres. We were using create_unit_hemisphere(2) and (3) which provide only 9 and 33 points respectively (too discrete). Also you will need to make sure that using hemispheres is a good idea for this type of interpolation. I am afraid that the points in the edges might create problems. So, I would use full spheres to be on the safe side.

Finally, I noticed some another minor issue. In lines 550 and 552 we use is for comparing strings which is generally not recommended. I would change that to ==.

So, in summary it seems that the problem may not be so serious and it can be resolved with a tiny bit of smoothing and less discrete spheres. So, if you can update the test I think we can move forward. At the same time it would be nice to have in our mind any new developments from the scipy devs. And give them feedback as we learn more from the behaviour of these functions.

How does that sound?

@arokem

This comment has been minimized.

Member

arokem commented Sep 14, 2015

So what do you think? Should we enforce smoothing per default? Might be a
safer way to go(?) And drop the Euclidean distance? In a way, it's not
surprising that doesn't work as well as angle - Euclidean is just an
approximation of angle in this case.

I've made the following changes:

  • Smoothing is set to 0.1 per default.
  • norm is set to angle per default. If a user sets it to
    euclidean_norm a warning is issued. The euclidean will be completely
    removed as an option in a subsequent release cycle. I've added a test to
    see that this warning is raised.
  • I've set the test to use spheres, instead of hemispheres (good call), and
    these have more points than before.

On Sun, Sep 13, 2015 at 4:20 PM, Eleftherios Garyfallidis <
notifications@github.com> wrote:

Hi @arokem https://github.com/arokem,
here is some feedback from my investigations.

The problem at least at my scipy 0.14 seems to be only with
norm="euclidean_norm"
(norm="angle" works fine) which can be resolved using a small smoothing
(here using 0.1).

Check my code snippet here

def test_interp_rbf():
def data_func(s, a, b):
return a * np.cos(s.theta) + b * np.sin(s.phi)

from dipy.core.sphere import Sphere, interp_rbf
from dipy.core.subdivide_octahedron import create_unit_hemisphere
import numpy as np

s0 = create_unit_hemisphere(5)
s1 = create_unit_hemisphere(6)

from dipy.viz import fvtk
s0 = Sphere(xyz = np.vstack((s0.vertices, -s0.vertices)))
s1 = Sphere(xyz = np.vstack((s1.vertices, -s1.vertices)))


for a, b in zip([1, 2, 0.5], [1, 0.5, 2]):
    ren = fvtk.ren()
    print(a, b)
    data = data_func(s0, a, b)
    expected = data_func(s1, a, b)
    interp_data_en = interp_rbf(data, s0, s1, smooth=.1, norm="euclidean_norm")
    interp_data_a = interp_rbf(data, s0, s1, smooth=.1, norm="angle")

    sf = np.zeros((3, s1.vertices.shape[0]))
    sf[0] = expected
    sf[1] = interp_data_en
    sf[2] = interp_data_a
    fvtk.add(ren, fvtk.sphere_funcs(sf, sphere=s1, norm=False, radial_scale=False))
    nt.assert_(np.mean(np.abs(interp_data_en - expected)) < 0.1)
    nt.assert_(np.mean(np.abs(interp_data_a - expected)) < 0.1)

    fvtk.show(ren)

if name == "main":
# nt.run_module_suite()
test_interp_rbf()

[image: image]
https://cloud.githubusercontent.com/assets/134276/9839458/23f7d4f8-5a49-11e5-8b7f-72da9618fd78.png

Now if you change to the default smooth (zero) you get some weird result
looking

[image: image]
https://cloud.githubusercontent.com/assets/134276/9839503/84ae56c2-5a4a-11e5-8615-eefafd8dd7ed.png

Some other important points to improving the tests. I think those tests
should use higher sampled spheres. We were using create_unit_hemisphere(2)
and (3) which provide only 9 and 33 points respectively. Also you will need
to make sure that using hemispheres is a good idea for this type of
interpolation. I am afraid that the points in the edges might create
problems. So, I would use full spheres to be on the safe side.

Finally, I noticed some another minor issue. In lines 550 and 552 we use
is for comparing strings which is generally not recommended. I would
change that to ==.


Reply to this email directly or view it on GitHub
#707 (comment).

@arokem

This comment has been minimized.

Member

arokem commented Sep 14, 2015

Thanks for all the feedback, by they way!

On Sun, Sep 13, 2015 at 7:20 PM, Ariel Rokem arokem@gmail.com wrote:

So what do you think? Should we enforce smoothing per default? Might be a
safer way to go(?) And drop the Euclidean distance? In a way, it's not
surprising that doesn't work as well as angle - Euclidean is just an
approximation of angle in this case.

I've made the following changes:

  • Smoothing is set to 0.1 per default.
  • norm is set to angle per default. If a user sets it to
    euclidean_norm a warning is issued. The euclidean will be completely
    removed as an option in a subsequent release cycle. I've added a test to
    see that this warning is raised.
  • I've set the test to use spheres, instead of hemispheres (good call),
    and these have more points than before.

On Sun, Sep 13, 2015 at 4:20 PM, Eleftherios Garyfallidis <
notifications@github.com> wrote:

Hi @arokem https://github.com/arokem,
here is some feedback from my investigations.

The problem at least at my scipy 0.14 seems to be only with
norm="euclidean_norm"
(norm="angle" works fine) which can be resolved using a small smoothing
(here using 0.1).

Check my code snippet here

def test_interp_rbf():
def data_func(s, a, b):
return a * np.cos(s.theta) + b * np.sin(s.phi)

from dipy.core.sphere import Sphere, interp_rbf
from dipy.core.subdivide_octahedron import create_unit_hemisphere
import numpy as np

s0 = create_unit_hemisphere(5)
s1 = create_unit_hemisphere(6)

from dipy.viz import fvtk
s0 = Sphere(xyz = np.vstack((s0.vertices, -s0.vertices)))
s1 = Sphere(xyz = np.vstack((s1.vertices, -s1.vertices)))


for a, b in zip([1, 2, 0.5], [1, 0.5, 2]):
    ren = fvtk.ren()
    print(a, b)
    data = data_func(s0, a, b)
    expected = data_func(s1, a, b)
    interp_data_en = interp_rbf(data, s0, s1, smooth=.1, norm="euclidean_norm")
    interp_data_a = interp_rbf(data, s0, s1, smooth=.1, norm="angle")

    sf = np.zeros((3, s1.vertices.shape[0]))
    sf[0] = expected
    sf[1] = interp_data_en
    sf[2] = interp_data_a
    fvtk.add(ren, fvtk.sphere_funcs(sf, sphere=s1, norm=False, radial_scale=False))
    nt.assert_(np.mean(np.abs(interp_data_en - expected)) < 0.1)
    nt.assert_(np.mean(np.abs(interp_data_a - expected)) < 0.1)

    fvtk.show(ren)

if name == "main":
# nt.run_module_suite()
test_interp_rbf()

[image: image]
https://cloud.githubusercontent.com/assets/134276/9839458/23f7d4f8-5a49-11e5-8b7f-72da9618fd78.png

Now if you change to the default smooth (zero) you get some weird result
looking

[image: image]
https://cloud.githubusercontent.com/assets/134276/9839503/84ae56c2-5a4a-11e5-8615-eefafd8dd7ed.png

Some other important points to improving the tests. I think those tests
should use higher sampled spheres. We were using create_unit_hemisphere(2)
and (3) which provide only 9 and 33 points respectively. Also you will need
to make sure that using hemispheres is a good idea for this type of
interpolation. I am afraid that the points in the edges might create
problems. So, I would use full spheres to be on the safe side.

Finally, I noticed some another minor issue. In lines 550 and 552 we use
is for comparing strings which is generally not recommended. I would
change that to ==.


Reply to this email directly or view it on GitHub
#707 (comment).

arokem added some commits Sep 1, 2015

BF: Keep up with changes in scipy 0.16
Recent changes in scipy use a different heuristic to find the epsilon
parameter used for RBF interpolation, causing for us issue #688. This
should fix it, by using the old heuristic.
PEP8
While I'm here, some PEP8 fixes on this file.
BF: Don't use Euclidean norm as distance function.
Instead, warn that it will be deprecated in future versions.

Additional changes in testing.

@arokem arokem force-pushed the arokem:fix-688 branch from 94e02d6 to d332cfb Sep 14, 2015

@Garyfallidis

This comment has been minimized.

Member

Garyfallidis commented Sep 14, 2015

Great! All good here.
I think this is good as it is and it doesn't affect any other parts of the code.

But put in your TODO list (after Thursday) to update the api_changes.rst file by saying that the default parameter has changed and some norms will be deprecated. You should keep an issue open for scipy 0.16 if necessary.

I will go ahead and merge this as it is now. That will make it easier to know which of the other PRs are failing for a good reason and not because of RBF.

Oh and best of luck with your deadline. Thank you for looking into this.

Garyfallidis added a commit that referenced this pull request Sep 14, 2015

Merge pull request #707 from arokem/fix-688
BF: Keep up with changes in scipy 0.16

@Garyfallidis Garyfallidis merged commit d5f7540 into nipy:master Sep 14, 2015

1 check passed

continuous-integration/travis-ci/pr The Travis CI build passed
Details
@stefanv

This comment has been minimized.

Contributor

stefanv commented Sep 14, 2015

I guess it's not too surprising that the Euclidean norm performed badly. This looks like a reasonable solution!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment