[R] Fast evaluation of functions in 3D domains

[Rolf Turner]

On 23/03/15 23:44, Lluis.Hurtado at uv.es wrote:

> Dear all,
>
> I am currently working with the spatstat package with 3D samples. I am trying to
> evaluate a non analytical function over the window that encloses the sample and I
> need to know which is the fastest way of doing it.

I gather that you mean that you want to evaluate the *integral* of the 
function over the window.

What do you mean by "non analytical"?  You have, it would appear, some 
way of calculating the function at an arbitrary point (x,y,z) in 
3-space, which makes it at least in some some sense "analytical".

> The function input is a 3 coordinate position in the window (x,y,z) and a list of
> parameters (a,b,c). The output is a numerical value.
>
> n <- function(x,y,z,a,b,c)
>
> But I need to do it over the whole volume.
>
> For 2 dimensions it can be done with
>
> A <- as.im(function,window,parameters)
> norm <- integral.im(A).

There is nothing very sophisticated about the way integral.im() works. 
Just look at the code.  The operative line of code is:

    a <- with(f, sum(v, na.rm = TRUE) * xstep * ystep)

> For 3 dimensions I have tried to pass an array of a grid covering the window (like a
> quadrature scheme) and then summing up the output array, but I would like to know if
> there is any faster way of integrating the function.

It sounds, at least vaguely, as though you are imitating in 3D what 
integral.im() does in 2D.  So, given that you are doing it right, you 
are doing as well as integral.im() does.

The calculations should be pretty fast since sum() is pretty fast.  The 
bulk of the computational burden is likely to be the evaluation of the 
function at all of the voxel centres.

Note that this is a very unsophisticated method of numerical quadrature. 
  The integrand is being approximated by a step function.  The method 
seems to work well enough for the uses that we put it to in spatstat. 
Whether it works well enough for your purposes depends on what your 
purposes are.

cheers,

Rolf Turner

-- 
Rolf Turner
Technical Editor ANZJS
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276
Home phone: +64-9-480-4619