[pymvpa] Q about faster GNB searchlight script: quick way to compute sphere voxels? Plus some Matlab code, in case it helps

[Rajeev Raizada]

Dear MVPA python + matlab folks (I'm cross-posting to both lists),

I was chatting with Yarik and Michael the other day
about their ultra-fast Gaussian Naive Bayes searchlight script.
I think it's a nice idea. Thanks for making it!
(The script I'm talking about is here in the code tree:
pymvpa/mvpa/measures/gnbsearchlight.py )

The script gains its speed by exploiting a nice trick
that I hadn't thought of:
because GNB doesn't use any info about particular
combinations of voxels, you can calculate the means and s.d.s
of all the voxels at once at the beginning of each cross-validation rep,
and then reuse those precomputed values for any searchlight sphere.

I'm still far too lacking in fluency in Python to be able to really
read through and understand individual PyMVPA scripts,
so I reimplemented the GNB part in Matlab:
http://www.nmr.mgh.harvard.edu/~raj/Matlab/raj_searchlight_gnb.m

That script is nice and fast, but it relies on the coords of the sphere voxels
having been pre-computed in advance.
That still takes a fair chunk of CPU time:
http://www.nmr.mgh.harvard.edu/~raj/Matlab/premake_sphere_coords_list.m

As far as I can tell from the gnbsearchlight.py script,
it assumes that the sphere radius is one voxel.
One-voxel-radius sphere-coords can be computed on the fly very quickly,
simply by listing the six +/-1 x,y,z voxels attached to the sphere center.

But I was wondering if you guys have a quick trick
for computing the voxel coords of searchlight spheres which
are of larger radius. Once you get past radius 1 or 2,
the lookup-table approach starts to get pretty impractical.

My best attempt so far is that script
http://www.nmr.mgh.harvard.edu/~raj/Matlab/premake_sphere_coords_list.m
It's not very good though. It still takes a long time.
At least you only need to run it once, but it's still slow.
Any advice on faster tricks greatly appreciated!

Raj