Performance of nearest vs linear

[pelson]

From https://groups.google.com/forum/#!searchin/scitools-iris/performance/scitools-iris/HajdRLOZZQ4/ypy1FzbxtkMJ, with code from https://gist.github.com/pelson/8256b9bab1c775348492 investigate why the performance of nearest neighbour is ~2 orders of magnitude better than linear interpolation.

[duncanwp]

Any update on this? It's a bit of a killer performance wise... Thanks!

[corinnebosley]

Having undertaken a preliminary investigation, it appears that there are several different ways of interpolating a series of points, most of which refer back to the cube.interpolate method with which you specify the interpolation scheme as an argument.

It is understandable that the linear interpolation scheme is slower than the nearest neighbour scheme because much more computation is required per sample point. This is not just a mathematical process, but seemingly quite a complex problem for cubes in terms of dimension manipulation.

There are, however, some things which can be done to optimise the linear scheme:

1. If you are looking for a set of sample points, it is always best to use a list of the points in a pre-defined interpolator, like this:

sample_lats = [0.1, 1.1, 2.1]
sample_longs = [0.1, 1.1, 2.1]
scheme = iris.analysis.Linear()
interpolator = scheme.interpolator(cube, ('latitude', 'longitude'))
interpolated_cube = interpolator((sample_lats, sample_longs))

This process will be much quicker than making a separate interpolator for every point of latitude and longitude.

2. If you can help it, don't use iris.analysis.trajectory as this is the slowest interpolation method. It does give you a 1-dimensional cube of interpolated data points (where other methods return a 2-dimensional cube containing unnecessary values), but it is currently not optimised, and it is still quicker to write your own routine to extract the required values from a cube than use trajectory:

interp_cube = linear(regular_cube, sample_points)
data_points = []
for x in xrange(len(lats)):
    data_points.append(interp_cube.data[x, x])

3. If it is possible to use a grid translation of some sort to interpolate your points, you can use cube.regrid instead.

I have timed a few different tasks using the three methods and compared the results in this notebook (in case you are still interested):

https://gist.github.com/corinnebosley/7c0370c3cde62c67386d#file-interpolate_timings-ipynb

I will continue to work on this issue and report findings or progress where appropriate. I will be proposing a small fix to optimise the trajectory method in the near future.

[duncanwp]

Thanks for taking a look at this Corinne. This backs up my hunch that the new interpolate routines are much faster than the trajectory methods.

I take your point about always trying to pass the whole list of sample points in one go, but this gets very inefficient very quickly as you go to more sample points, or larger cubes. I've put together some more tests based on your findings (using only method 1) in this notebook:
https://gist.github.com/duncanwp/b8de87520450f27a6ce4

I also did some tests on doing one sample point at a time, and for lots of sample points (100 or more) this is actually quicker.

It would be great to have an efficient trajectory like routine which allows us to ask for the interpolation of a (potentially large) number of points across a cube. We do this a lot when trying to interpolate model values onto e.g. aircraft data.

[duncanwp]

I've just updated that gist to include another method I've started using in the meantime which is to cache the interpolator and only do the interpolation using the ._points method (since we don't want a cube out at the end anyway). This helps but it's still pretty slow.

[corinnebosley]

Thanks Duncan, I'm just having a look at using the ._points method in trajectory to see if that's any faster. I'm making progress, I'll let you know how I get on.

[corinnebosley]

Hi Duncan, I have hit a lot more issues trying to use ._points in the module.

I have made a pull request for some changes which do speed up the process a bit (~25% for the linear scheme and ~62% for nearest neighbour), but not nearly as much as I think we can ultimately achieve.

Having spoken to a few people in the team it looks like we might try and put some time into this at some point in the future to make the changes that this whole area needs. At that point, I would like to come and talk to you about this in person, if that's ok with you?

[duncanwp]

Absolutely, just drop me an email and we can arrange it. I'll let you know if I'm down at the MO in the meantime. Do you have a rough idea when this might be?

Thanks

[corinnebosley]

I'd like to have a chat about this while I'm still knee-deep in it, so I think if you are around, sometime in the next week or two would be good. Otherwise, just let me know when you will be in the MO.

[duncanwp]

Great, I'm not going to be down before then so just let me know when works for you.

[corinnebosley]

How about Wednesday, you around on Wednesday?

[duncanwp]

Can you email me and we can work it out: Duncan.watson-Parris (at) physics.ox.ac.ukhttp://physics.ox.ac.uk. Thanks!

On 22 Feb 2016, at 17:14, Corinne Bosley <notifications@git.colasdn.topmailto:notifications@github.com> wrote:

How about Wednesday, you around on Wednesday?

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/1767#issuecomment-187274098.