# Why is there no PositionFunction in Mathematica?

Too often I have seen the programs of inexperienced users greatly slowed by using Position in an iterative fashion, when far faster would have been to compute a look-up table for positions beforehand.

Mathematica provides such functionality for Nearest and Interpolation "out of the box" with the syntax of function[data] along with dedicated functions NearestFunction and InterpolatingFunction.

There is no equivalent with Position[data] producing a PositionFunction object. I observe that Position appears to be used more often than either of the other functions so this seems like a regrettable omission. Two complications I can think of are:

• the levelspec of Position

• the handling of patterns

1. What other complications are there for creating a PositionFunction?

2. How best can such a function be implemented?

3. How can the utility and performance of the function be maximized?

4. All limitations considered would such a function be valuable?
Why is there no such functionality in Mathematica?

-
I belong to the team of Position[]-users... :-( –  Rod May 21 '13 at 13:23
@RodLm I'm not sure how to read that; often Position is well used too, and I certainly use it a lot myself (not implying that makes it right). Rather, certain uses of it are quite slow. I hope to bring attention to the issue and perhaps develop a useable general function. –  Mr.Wizard May 21 '13 at 13:26
I find it funny that I really only started to use Position after participating here. I think this site has warped my sensibilities. :D –  rcollyer May 21 '13 at 13:56
@ruebenko That is not the case I'm thinking of. The usual misuse is something like: a = RandomInteger[{1, 1000}, {5000}]; Table[{i, Position[a, i]}, {i, 1, 1000}] // Timing // First when magnitudes faster would be Sort[{#[[1, 1]], #[[All, 2]]} & /@ GatherBy[Transpose[{a, Range@Length@a}], First]] // Timing // First –  Mr.Wizard May 22 '13 at 7:56
There does occur to me a 'reason' for why there is no PositionFunction, and that is that Nearest and Interpolation are "subject matter" functions whereas Position is more a language thing. E.g. hypothetically, Nearest could be constructed from a variety of data structures because all that matters is the mathematical distance between one unit and another, whereas Position is entirely about the specific list structure you are using. Not saying I'd agree with this, but there is this conceptual difference which may have something to do with the absence of a PositionFunction. –  amr May 23 '13 at 1:52

This is perhaps a place to start:

position[expr_, level_: 1] :=
With[{positionData =
SortBy[
#[[1, 1]] -> #[[All, 2]] & /@
GatherBy[Extract[expr, #, Verbatim] -> # & /@ Position[expr, _, level], First],
Min[Length /@ #[[2]]] &
] // Dispatch},
Replace[#, positionData] &
]


The second argument controls the depth of the indexing. An example:

f = position[x^2 + y^2 + q_^r_, 3];

f[2]

{{1, 2}, {2, 2}}


Due to the use of Verbatim patterns are matched literally, which deviates from normal Position behavior:

f[_]

{{3, 1, 2}, {3, 2, 2}}

-
Hm.. Are you aware that setting x=y breaks the positionFunction? i.e. f[y] -> {1,1} . Verbatim does not have HoldAll. Maybe use HoldPattern (as well as Verbatim?). I guess it would be even better if you could shield the symbols from updates. Compare, with your definitions, {Drop[Trace[f, TraceOriginal -> True], 2], Update[y]; Drop[Trace[f, TraceOriginal -> True], 2]} –  Jacob Akkerboom May 21 '13 at 16:36
@Jacob Good point; I should add that, or change to DownValues on a HoldAll function, etc. I only intended this to be a rough start. –  Mr.Wizard May 21 '13 at 18:05
Maybe I should have added that I respect that it is rough start. It seems like a good start though :) –  Jacob Akkerboom May 21 '13 at 18:08
It turns out HoldComplete is not the only thing we need. Mathematica is really thorough in storing which subexpressions contain which symbols, it seems. Very interesting :) –  Jacob Akkerboom May 21 '13 at 18:14
@rm-rf Don't miss the videos. –  Mr.Wizard May 22 '13 at 5:54

Temporary message: First of all, sorry for the big mess, especially the 20 line long paragraph. I am tired :)

Below, on second thought the Verbatim option is no good yet, as the rules should be gathered differently in the case we do not use Verbatim, but I guess the code does not hinder the proper implementation of this. Allow me to show you this work in progress.

ClearAll[positionFunction]
positionFunction[expr_, level_: 1,
OptionsPattern[{"attributes" -> None, "verbatim" -> True}]
] :=
Module[{myHold},
SetAttributes[myHold, HoldAll];
With[
{
attrs =
If[#1 === None, Unevaluated[Sequence[]], #1] &@
OptionValue["attributes"]
,
positions = Position[Unevaluated[expr], _]
},

With[
{
symbolSower = (If[
MatchQ[#1, myHold[_Symbol]] && Context @@ #1 =!= "System",
Sow[#1]]; #1) &
}
,
With[
{
reapage =
Reap[(symbolSower[#1[[1, 1]]] -> #1[[All, 2]] &) /@

GatherBy[
Partition[
Riffle[
Extract[
Unevaluated[expr],
positions,
myHold
],
positions
]
,
2
],
First
]
]
}
,
With[
{
heldRules =
If[
OptionValue["verbatim"]
,
ReplaceAll[#, myHold -> Verbatim ] &
,
DeleteCases[#, myHold, {4, 4}, Heads -> True] &
]@(HoldComplete @@ {reapage[[1]]})
,
symbolsHeld = Hold[Evaluate[reapage[[2, 1]]]][[All, All, 1]]
}
,
With[
{
completelyHeldDispatch =
ReleaseHold[
Hold[Block][symbolsHeld,
Hold[HoldComplete @@ {Dispatch @@ heldRules}]]]
}
,
Function[Null,
Replace @@ {Unevaluated[#1],
Unevaluated @@ completelyHeldDispatch}, attrs]
]
]
]
]
]
];


The dispatch tables in this version should never break. However, if we use Verbatim, somehow no DispatchTables will be made. This may be because Dispatch just does not work well with Verbatim. Example:

Dispatch[{Verbatim[List] -> {{0}}, Verbatim[c] -> {{1}},
Verbatim[d] -> {{2}}, Verbatim[e] -> {{3}}, Verbatim[f] -> {{4}},
Verbatim[g] -> {{5}}, Verbatim[{c, d, e, f, g}] -> {{}}}]


-> {Verbatim[List] -> {{0}}, Verbatim[c] -> {{1}}, Verbatim[d] -> {{2}}, Verbatim[e] -> {{3}}, Verbatim[f] -> {{4}}, Verbatim[g] -> {{5}}, Verbatim[{c, d, e, f, g}] -> {{}}}

I would expected an expression that displayed as Dispatch[..., -DispatchTables-] here.

This is (probably) not because there are to few arguments in the list of rules. You make a dispatch table from a list of 4 rules.

However, it is still nice that the lists of rules in the case of Verbatim do not break. If you wonder why they might break in the first place (and what I mean precisely), see the section "about HoldAllComplete" and this Q&A.

The Unevaluated in the second argument of Replace is necessary to prevent this breaking. I made another wrapper Unevaluated in the first argument of Replace so that we will be able to use positionFunctionRealisation[Unevaluated[arg]]. I have made it optional whether or not the resulting positionFunction will have HoldAll (in fact we can give the resulting positionFunction any attribute). The default value is None, which means that the positionFunction will have no arguments. I have renamed position positionFunction as it generates positionFunctions in the same sense that Function generates functions. I have improved positionFunction to work with Held expressions, as well as made it possible for positionFunction to be used with Unevaluated to make a positionFunction from an expression that would change due to evaluation. I have made it so that we only use Extract once, instead of mapping some function comprising of Extract of the positions. It still bothers me that we use both Extract and Position on our expression expr, even though they are similar. But hey, what can you do?

Examples:

Clear[ffff11, x, y]
x = 6;
ffff11 = positionFunction[Unevaluated[x^2 + y^2 + q_^r_], 3];

ffff11[Unevaluated[x]]


-> {{1, 1}}

We can change x, but still dispatch table works

Clear[gggg11, x, y]
x = y;
gggg11 = positionFunction[Unevaluated[x^2 + y^2 + q_^r_], 3,
"attributes" -> HoldAll];

gggg11[x]
gggg11[Evaluate[x]]


-> {{1, 1}}

-> {{2, 1}}

gggg11[_]


-> {{3, 1, 2}, {3, 2, 2}}

Now let

hhhh11 = positionFunction[x^2 + y^2 + q_^r_, Infinity,
"attributes" -> HoldAll, "verbatim" -> False]


To show that the dispatch table does not break, we have:

x = 8;
Cases[
Trace[
hhhh11[x]
,
TraceOriginal -> True],
Dispatch | Rule, Infinity
] == {}


(-> True)

The point here is that if a dispatch table will reevaluate, you will have both a case of {Rule} and a case of {Dispatch} in your trace. We find no such cases, so the dispatch table is not reevaluated. However, a Dispatch table is certainly present in this case. This can be seen by

x = 9;
Cases[Trace[hhhh11[x], TraceOriginal -> True], Dispatch, Infinity,


-> 14

See my Q&A here for more about breaking dispatch tables.

Old stuff

WARNING: everything beyond this point is outdated. The use of HoldComplete below is not so nice. However, I will keep it here for reference and as a bad example of a fix or whatever.

To fix undesirable aspect(s) of Mr.Wizards solution, we might alter the code as follows

ClearAll[position3]
position3[expr_, level_: 1] :=
With[
{
positionData =
SortBy[
#[[1, 1]] -> #[[All, 2]] & /@
GatherBy[
Extract[expr, #, Composition[HoldPattern, Verbatim]] ->
# & /@
Position[expr, _, level]
,
First
],
Min[Length /@ #[[2]]] &
] // Dispatch
}
,
ReleaseHold@HoldComplete[Replace[#, positionData]] &
]


Where I have just added HoldPattern and a HoldComplete. First of all, note that Verbatim is still necessary and it still works. We have for

Clear[g, x, y]
g = position3[x^2 + y^2 + q_^r_, 3];


that

g[_]


-> {{3, 1, 2}, {3, 2, 2}}

But we now also have

x = y;
g[y]


-> {{2, 1}}

Whereas with Mr.Wizards f we would have

f[y]
x=y;
f[y]


-> {{2,1}}

-> {{1,1}}

WARNING: this section is also outdated. Maybe you can learn something here about breaking dispatch tables though, who knows :).

We have

{Drop[Trace[g, TraceOriginal -> True], 2], Update[y]; Drop[Trace[g, TraceOriginal -> True], 2]}


-> {{}, {}}

(no matter when we last updated), whereas with Mr.Wizards f we would have

f; (*this is just to make sure the first element of the output below is an empty list*)
{Drop[Trace[f, TraceOriginal -> True], 2], Update[y];
Drop[Trace[f, TraceOriginal -> True], 2]}

(*output*)
{{}, {{HoldForm[Function]},
HoldForm[Replace[#1, {Verbatim[Plus] -> {{0}}, Verbatim[x^2] -> {{1}},
Verbatim[y^2] -> {{2}}, Verbatim[(q_)^(r_)] -> {{3}},
Verbatim[2] -> {{1, 2}, {2, 2}}, Verbatim[Power] ->
{{1, 0}, {2, 0}, {3, 0}}, Verbatim[x] -> {{1, 1}},
Verbatim[y] -> {{2, 1}}, Verbatim[q_] -> {{3, 1}},
Verbatim[r_] -> {{3, 2}}, Verbatim[Pattern] -> {{3, 1, 0}, {3, 2, 0}},
Verbatim[q] -> {{3, 1, 1}}, Verbatim[r] -> {{3, 2, 1}},
Verbatim[_] -> {{3, 1, 2}, {3, 2, 2}}}] & ]}}


So we have successfully shielded (most of the) expression (function) to which g is assigned from updates. So we do not have to reconstruct this expression every time a variable that occurs in this expression (x, y, q or r) changes. The function shielding against such tracking of symbols is HoldComplete. The reason it is not dangerous to shield from updates in this case is that the evaluation of our function should never depend on changes of symbols.

Actually on closer look this shielding does not seem to be so important in this case. The trace I seem to have underestimated Mathematica. The trace generated is not long at all, so our HoldComplete construction probably did not save much work. My guess is that the Attribute HoldAll of Function helps here (may be a wild guess). Anyway, to show what I mean, I expected something like the following to happen. Set for the moment

ClearAll[positionData]
positionData[expr_, level_: 1] :=
SortBy[
#[[1, 1]] -> #[[All, 2]] & /@
GatherBy[
Extract[expr, #, Composition[HoldPattern, Verbatim]] ->
# & /@
Position[expr, _, level]
,
First
],
Min[Length /@ #[[2]]] &
] // Dispatch
Clear[hhhh, x, y]
hhhh = positionData[x^2 + y^2 + q_^r_, 3];


For this we have

hhhh;
Drop[Trace[hhhh, TraceOriginal -> True], 2]
Update[y];
Drop[Trace[hhhh, TraceOriginal -> True], 2]//Length


-> {}

-> 16

So that here the rebuilding is more significant.

Edit It turns out such rebuilding will happen in all cases, when the function generated by position is evaluated with arguments. If we really want to block the tracking of symbols we will have to do it somewhere low in the expression tree rather than high up. /Edit

Anyway this might be something to be wary of.

-
Come on guys, at least leave a comment saying that we would rather have Replace[Unevaluated[..], Unevaluated[..]] in the resulting function than something where Replace still has to be applied. Be a good sport! :P. Ill fix soon ;). Another matter to adress is why Dispatch does not work with Verbatim in the example... –  Jacob Akkerboom May 24 '13 at 11:40
Ah yes, and Inner[List, extracts, positions, Hold] instead of Partition[Riffle[extracts, positions],2] –  Jacob Akkerboom May 24 '13 at 16:19

This is a similar approach to Mr. Wizard's I suppose, but using the function's DownValues instead of a dispatch table to store the rules.

One major difference is that this code won't work with data containing pattern objects (I guess this might be fixable with Verbatim).

The expression is traversed using MapIndexed, for each part visited the position is Sowed, with the part's value as the tag. The actual downvalues are set afterwards using the third argument of Reap (the first time I've ever used it I think).

makePositionFunction[f_Symbol, data_, level_: {-1}] := Block[{},
ClearAll[f];
Reap[
MapIndexed[Sow[#2, #1] &, data, level, Heads -> True],
_, (f[#] = #2) &];
f[other_] := Position[data, other, level]]


Example:

data = RandomInteger[1000, {3000, 20}];

makePositionFunction[pos, data];

First @ Timing[test1 = Table[pos[i], {i, 1000}]]
(* 0. *)

First @ Timing[test2 = Table[Position[data, i], {i, 1000}]]
(* 3.532 *)

test1 == test2
(* True *)


The "PositionFunction" defaults to using plain old Position for any search patterns which have not been precomputed:

pos[n_ /; n < 20] // Length
(* 1225 *)

-
+1 for the post. However, I considered MapIndexed and went with Position and Extract for one reason: held expressions (or subexpressions). To make this robust I think you'll need something additional, e.g. MapIndexed onto an explicitly held expression and then a Cases to pull out those parts. You'll also have to consider how to keep parts extracted from a held expression from evaluating. –  Mr.Wizard May 22 '13 at 14:38
@Mr.Wizard, good point, I didn't consider held expressions. I suppose one could, with sufficient code gymnastics, work with held expressions in MapIndexed - but then why bother when you've already shown how to do it neatly with Position and Extract. Great question, by the way. –  Simon Woods May 22 '13 at 16:02
@Mr.Wizard I find it very frustrating that there is no function like ScanIndexed. That would work perfectly with held expressions. I have spent a lot of time trying to find an alternative for held expressions to the construct with (Inner, my addition) Position and Extract (that I suppose will soon appear in my answer), to try do something with Scan and/or MapIndexed. Maybe we can discuss things soon? –  Jacob Akkerboom May 29 '13 at 15:27

I think the complexities of things like PatternTest will obstruct any kind of data structure for searching for general pattern matches. I think some regularity, either with the pattern as in Mr Wizard's answer or with the expression as in Leonid Shifrin's remark, will be needed in order to beat using Position. Jacob Akkerboom has pointed out problems with what to do if variables in expression change after the position function is computed.

That said, here's a start along the lines of Mr. Wizard's answer, using literal matches (literally, as you will see ;-). When he mentioned NearestFunction in the question, I thought, I wonder if there is a way to measure the "pattern-distance" between expressions. I don't think it's likely in general, but you can in the case of exact matches. Recall that Nearest uses EditDistance when the data are strings.

ClearAll[position2, positionFunction];
positionFunction[nf_NearestFunction, expr_][x_] :=
With[{arg = ToString[FullForm @ x]},
With[{pos = nf[arg]},
If[pos =!= {} && ToString[FullForm @ Extract[expr, Prepend[First @ pos, 1]]] === arg,
pos, {}]]];
Format[positionFunction[nf_NearestFunction, expr_]] := positionFunction["<>", Short[expr]];
position2[expr_, level_: 1] := positionFunction[
Nearest[ToString[FullForm @ Extract[expr, #]] -> # & /@ Position[expr, _, level]],
HoldForm[expr]]


Mr.Wizard's example:

pf = position2[x^2 + y^2 + q_^r_, 3]

positionFunction["<>", x^2 + y^2 + q_^r_]

pf[2]

{{1, 2}, {2, 2}}

pf[_]

{{3, 1, 2}, {3, 2, 2}}


Another example:

gf = position2[Plot3D[x^2 + y^2, {x, -2, 2}, {y, -2, 2}], Infinity]


gf[Polygon]

{{1, 2, 1, 1, 2, 1, 1, 0}, {1, 2, 1, 1, 2, 1, 2, 0},
{1, 2, 1, 1, 2, 1, 3, 0}, {1, 2, 1, 1, 2, 1, 4, 0}}
`
-
I love how your positionFunction looks. Also, you code is so short compared to mine ;). –  Jacob Akkerboom May 23 '13 at 21:37
@JacobAkkerboom Thanks. I think to do this right, the code might have to be long, though. –  Michael E2 May 23 '13 at 23:03