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

[Yaroslav Halchenko]

and as promised:
> > are of larger radius. Once you get past radius 1 or 2,
> > the lookup-table approach starts to get pretty impractical.
> hm... shouldn't be too bad... I will do some timing and report back.

Here is the output -- indeed with increasing of radius it starts to take some
considerable time... but lets not forget that radius of 5 voxels is approx 4cm
diameter within the brain for typical BOLD data -- not sure how much of
spatial specificity people could claim from such analyses ;)

 Working on dataset <Dataset: 84x14042 at float64, <sa: chunks,targets>, <fa: voxel_indices>, <a: mapper>>
 Testing sphere neighborhoods
 Lets work with radius=0 voxels
[SPH_TIME] DBG{1032.626 sec}:        Training iqe
[SPH_TIME] DBG{0.091 sec}:        Training cached qe
[SPH_TIME] DBG{0.087 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{1.596 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{1.732 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.022 sec}:        done
 Lets work with radius=1 voxels
[SPH_TIME] DBG{0.021 sec}:        Training iqe
[SPH_TIME] DBG{0.703 sec}:        Training cached qe
[SPH_TIME] DBG{0.074 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{2.426 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{2.472 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.021 sec}:        done
 Lets work with radius=3 voxels
[SPH_TIME] DBG{0.022 sec}:        Training iqe
[SPH_TIME] DBG{0.071 sec}:        Training cached qe
[SPH_TIME] DBG{0.664 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{16.338 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{16.648 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.025 sec}:        done
 Lets work with radius=5 voxels
[SPH_TIME] DBG{0.024 sec}:        Training iqe
[SPH_TIME] DBG{0.096 sec}:        Training cached qe
[SPH_TIME] DBG{0.704 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{58.577 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{59.659 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.025 sec}:        done
 Testing cube neighborhoods
 Lets work with radius=0 voxels
[SPH_TIME] DBG{0.027 sec}:        Training iqe
[SPH_TIME] DBG{0.886 sec}:        Training cached qe
[SPH_TIME] DBG{0.068 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{2.123 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{1.584 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.021 sec}:        done
 Lets work with radius=1 voxels
[SPH_TIME] DBG{0.020 sec}:        Training iqe
[SPH_TIME] DBG{0.070 sec}:        Training cached qe
[SPH_TIME] DBG{0.651 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{4.924 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{5.083 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.021 sec}:        done
 Lets work with radius=3 voxels
[SPH_TIME] DBG{0.020 sec}:        Training iqe
[SPH_TIME] DBG{0.075 sec}:        Training cached qe
[SPH_TIME] DBG{0.721 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{41.490 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{41.752 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.022 sec}:        done
 Lets work with radius=5 voxels
[SPH_TIME] DBG{0.021 sec}:        Training iqe
[SPH_TIME] DBG{0.074 sec}:        Training cached qe
[SPH_TIME] DBG{0.763 sec}:        Lets query all neighbors for all features
[SPH_TIME] DBG{231.021 sec}:        Lets query all neighbors for all features 1st time using cached qe
[SPH_TIME] DBG{244.082 sec}:        Lets query all neighbors for all features 2nd time using cached qe
[SPH_TIME] DBG{0.022 sec}:        done

for the script which does only queries for all of 14k voxels without any other
computation (e.g. searchlight):

and the actual pymvpa script is below.  I've used debug facility to provide
timing. Mention that the time in {} is for how long it took to reach given a
line from previous DBG line:

   if not 'SPH_TIME' in debug.active:
        debug.register('SPH_TIME', 'Timing of sphere calls')
        debug.active = ['SPH_TIME']
        debug.metrics = ['reltime']

    verbose(1, "Working on dataset %s" % ds)
    for dist, dist_str in ((cartesian_distance, 'sphere'),
                           (absmin_distance, 'cube')):
        verbose(1, 'Testing %s neighborhoods' % dist_str)
        for radius in [0, 1, 3, 5]:
            verbose(1, 'Lets work with radius=%i voxels' % radius)
            sphere = Sphere(radius=radius, distance_func=dist)
            iqe = IndexQueryEngine(voxel_indices=sphere)
            qe = CachedQueryEngine(IndexQueryEngine(voxel_indices=sphere))

            debug('SPH_TIME', 'Training iqe')
            iqe.train(ds)
            debug('SPH_TIME', 'Training cached qe')
            qe.train(ds)
            debug('SPH_TIME', 'Lets query all neighbors for all features')
            n_qe = [iqe[i] for i in xrange(ds.nfeatures)]

            debug('SPH_TIME', 'Lets query all neighbors for all features 1st time using cached qe')
            n_qe = [qe[i] for i in xrange(ds.nfeatures)]
            debug('SPH_TIME', 'Lets query all neighbors for all features 2nd time using cached qe')
            n_qe = [qe[i] for i in xrange(ds.nfeatures)]
            debug('SPH_TIME', 'done')


-- 
                                  .-.
=------------------------------   /v\  ----------------------------=
Keep in touch                    // \\     (yoh@|www.)onerussian.com
Yaroslav Halchenko              /(   )\               ICQ#: 60653192
                   Linux User    ^^-^^    [175555]