       Re: Finding local maxima in multidimensional array (efficiently)

• To: mathgroup at smc.vnet.net
• Subject: [mg114224] Re: Finding local maxima in multidimensional array (efficiently)
• From: Ray Koopman <koopman at sfu.ca>
• Date: Sat, 27 Nov 2010 03:48:05 -0500 (EST)
• References: <ico3gv\$olk\$1@smc.vnet.net>

On Nov 26, 2:52 am, Julian Francis <julian.w.fran... at gmail.com> wrote:
> Dear all,
>
> How do I find the local maxima in a 4 Dimensional array.
>
> e.g.
>
> randTable = Table[Random[], {a, 1, 40}, {b, 1, 40}, {c, 1, 40}, {d, 1,
> 40}];
>
> in this case there are 80 neighbours for each point (except at the
> edges)
>
> This is my best attempt, NB I take the highest 0.1% to reduce the
>
> In:= threshold = Quantile[Flatten[randTable], .999]
>
> Out= 0.998989
>
> In:= top = Position[randTable, x_ /; x > threshold];
>
> In:= arrayBoundsCheck[array_, index_] :=
>  If[index == {}, True,
>   First[index] > 0 && First[index] <= Length[array] &&
>    arrayBoundsCheck[array[[First[index]]], Rest[index]]]
>
> In:= maxKernel =
>   DeleteCases[
>    Flatten[Table[{h, w, y, x}, {h, -1, 1}, {w, -1, 1}, {y, -1,
>       1}, {x, -1, 1}], 3], {0, 0, 0, 0}];
>
> In:= localMax[line_, array_] :=
>  Extract[array, line] >
>   Max[Extract[array,
>     Select[Map[Function[s, s + line], maxKernel],
>      arrayBoundsCheck[array, #] &]]]
>
> In:= Timing[results = Select[top, localMax[#, randTable] &];]
>
> Out= {80.59, Null}
>
> Taking the top 0.1% as I have done is just done to try to bring the
> compute time down, but isn't necessary to the problem.
> I'm really just interested in finding a list of all the local maxima.
> If someone has a much better solution, I would love to see it.
>
> Thanks very much,
> Julian.

Here are two routines that find local maxima.

'localMaxima' is a recoded version of your algorithm.

localMaxima[data_,p_] := Module[{dims,threshold,maxKernel,localMaxQ},
dims = Dimensions@data; threshold = Quantile[Flatten@data,p];
maxKernel = Transpose@Delete[Tuples[{-1,0,1},4],41];
localMaxQ[index_] := Extract[data, index] > Max@Extract[data,
Select[Transpose[index + maxKernel],
(And@@LessEqual@@@Transpose@{{1,1,1,1},#,dims}&)]];
Select[Position[data, x_ /; x > threshold], localMaxQ]]

'localMaxima2' takes a different approach. There is no thresholding,
and the window size is variable. The second parameter, 'h',
controls the window size: j +/- h on each axis.
A local maximum must be in the center of a full window.
('localMaxima' uses h = 1 implicitly, and does not
require maxima to be in the centers of full windows.)

localMaxima2[data_,h_] := Module[{dims = Dimensions@data,
r = Range[-h,h], s = (1+(1+2h)^4)/2}, Reap[
Do[If[#[[s]] > Max@Delete[#,s], Sow@{j1,j2,j3,j4}]& @
Flatten@data[[j1+r,j2+r,j3+r,j4+r]],
{j1,1+h,dims[]-h},
{j2,1+h,dims[]-h},
{j3,1+h,dims[]-h},
{j4,1+h,dims[]-h}]][[2,1]]]

Dimensions[randata = With[{n = 40}, Table[Random[],{n},{n},{n},{n}]]]
{40,40,40,40}

Timing@Length[result = localMaxima[randata,.999]]
{18.53 Second, 2448}

Timing@Length[result1 = localMaxima2[randata,1]]
{80.35 Second, 25788}

Timing@Length[result2 = localMaxima2[randata,2]]
{84.68 Second, 2629}

Timing@Length[result3 = localMaxima2[randata,3]]
{360.38 Second, 559}

• Prev by Date: Re: Problems with Mathematica 8.0 Solve
• Next by Date: Re: Matrix Form of Quadratic Equations
• Previous by thread: Finding local maxima in multidimensional array (efficiently)
• Next by thread: Re: Finding local maxima in multidimensional array (efficiently)