# Calculating the median difference between elements with a (particular pair) of consecutive integer indices residing in the same sublist

NOTE: I formally made a serious mistake in the first example provided. pointLists[[1]] had an extra element, and we should have a guarantee that Length[pointLists[[i]]] == Length[indexLists[[i]]]. I apologize for wasting people's time with this.

I have a list pointLists of sublists of 2D coordinates, where an N = 3 example looks like this:

pointLists =
{{
{131.048, 243.364}, {131.046, 243.321}, {131.037, 243.357},
{130.931, 243.391}, {130.909, 243.377}}, {{164.911, 244.039}, {164.929, 243.942},
{164.74, 244.083}}, {{98.1685, 239.618}, {98.1913, 239.6}, {98.1528, 239.623}
}};


I also have another list indexLists (where indexLists[[i]] corresponds to pointLists[[i]]) which assigns an integer value $\geq 1$ to each coordinate in each sublist. An example of indexLists for the above example looks like this:

indexLists = {{1, 2, 3, 5, 6}, {1, 2, 4}, {1, 2, 3}};


What I'd like to calculate here, as quickly as possible, is a list differenceList where a position k in differenceList corresponds to the mean or median difference between coordinates in the same sublist of pointLists where the first coordinate and the second coordinate have values of k and (k-1) in indexLists. Sometimes, however, there will be no examples in any of pointLists satisfying this criterion. In this case I'd like to simply set differenceList[[k]] = "NULL" (or really anything distinct and addressable).

Here's another way of thinking about the question:

Imagine, for example, that we have a bunch of cells growing in different colonies/plates, and their color changes over time. But we can't always measure the change for each colony at every time point, just whenever we can.

Each time we measure a color, we add it on the end of a sublist corresponding to the colony in pointLists and we time stamp the addition by adding the time to a position at the same index as pointLists in indexLists. We then ask... for some time point (e.g. $T = 5$), what's the median change from the previous timepoint (e.g. $T = 4$)?

And now imagine that I unfortunately have no control over this data formatting. I just have a structure like pointLists and indexLists to work with from measurement device.

Descriptive Example -

Since the description provided is a little wordy, let's plug through the above example in a step-wise fashion (this is optimized for easy reading, and should in no way constrain any algorithm that can achieve the same result):

(k = 1) We set differenceList[[1]] = {0,0} to handle this special case.

(k = 2) For k=2 we notice that pointLists[[1]], pointLists[[2]], and pointLists[[3]] have coordinates with indexLists values of $k = 2$ and $1$. We compute the differences as, for the first sublist: {131.046, 243.321} - {131.048, 243.364} = {-0.002, -0.043}, for the second sublist: {164.929, 243.942} - {164.911, 244.039} = {0.018, -0.097}, and for the third sublist: {98.1913, 239.6} - {98.1685, 239.618} = {0.0228, -0.018}. We then set: differenceList[[2]] = Median[{{-0.002, -0.043}, {0.018, -0.097}, {0.0228, -0.018}}] = {0.018, -0.043}.

(k = 3) For k=3 we notice that pointLists[[1]] and pointLists[[3]] have coordinates with indexLists values of $k = 3$ and $2$ (and that pointLists[[2]] does NOT have any coordinates satisfying this criterion). We compute the differences as - for the first sublist: {131.037, 243.357} - {131.046, 243.321} = {-0.009, 0.036}, and for the third sublist: {98.1528, 239.623} - {98.1913, 239.6} = {-0.0385, 0.023}. We then set: differenceList[[3]] = Median[{{-0.009, 0.036}, {-0.0385, 0.023}}] = {-0.02375, 0.0295}.

(k = 4) For k=4 we notice that none of the sublists have coordinates with an associated index value of $4$ AND ALSO a coordinate with an associated index value of $3$ in indexLists, so we set differenceList[[4]] = "NULL".

(k = 5) For k=5 we notice that, while pointLists[[1]] has a coordinate with an index value in indexLists of $5$, and pointLists[[2]] has a coordinate with an index value in indexLists of $4$, two coordinates with indices $4$ and $5$ not present in the same sublist, and there are no other examples where a sublist has both values, so we set differenceList[[5]] = "NULL".

(k = 6) For k=6 we notice that pointLists[[1]] has coordinates with indices $k = 6$ and $k = 5$ in pointLists, and that this is a lone example, so we set differenceList[[6]] = {130.909, 243.377} - {130.931, 243.391} = {-0.022, -0.014}.

Done. And our result is:

differenceList =
{{0,0}, {0.018, -0.043}, {-0.02375, 0.0295}, "NULL", "NULL", {-0.022, -0.014}}


Is there a clever way to achieve the same result as above for very large instances of pointLists (in terms of both the number and length of component sublists)? If it matters, I'd like to optimize for speed and not necessarily memory usage.

Also, since I'd like to be able to rerun the algorithm after updating pointLists, is there some way to quickly calculate the differenceList value at specific positions (i.e. for specific values of $k$) without recalculating the entire list?

Let's do one more quick example where we calculate differenceList[[5]] for new versions of pointLists and indexLists.

pointLists = {{{122.534, 191.041}, {122.735, 191.023}, {122.692, 190.999}, {122.742, 191.051}}, {{16.2591, 182.679}, {16.2429, 182.681}, {16.3327, 182.726}, {16.4354, 182.733}, {16.4144, 182.737}, {16.3727, 182.759}, {16.3141, 182.783}, {16.3692, 182.77}, {16.4068, 182.745}, {16.3346, 182.712}}};

indexLists = {{1, 3, 4, 5}, {1, 2, 3, 4, 5, 6, 7, 8, 9, 10}};


Here's a description of how to to calculate differenceList[[5]] in the new example:

From pointLists[[1]] we calculate: {122.742, 191.051} - {122.692, 190.999} = {0.05, 0.052}. From pointLists[[2]] we calculate: {16.4144, 182.737} - {16.4354, 182.733} = {-0.021, 0.004}. And finally, we take the median of those values: differenceList[[5]] = Median[{{0.05, 0.052},{-0.021, 0.004}}] = {0.0145, 0.028}.

Here, {122.742, 191.051} & {16.4144, 182.737} have indices of 5 at the same position in the associated indexLists. And points {122.692, 190.999} & {16.4354, 182.733} have indices of four at the same position in indexLists. We calculate: {122.742, 191.051} - {122.692, 190.999} since these values are in the same sublist of pointLists (pointLists[[1]]), and {16.4144, 182.737} - {16.4354, 182.733} are both in pointLists[[2]].

So that's how we do k = 5 in the new example.

If {122.742, 191.051} OR {16.4144, 182.737} OR both are missing, we'd have differenceList[[5]] = {-0.021, 0.004} (and differenceList[[5]] = {0.05, 0.052} if {122.692, 190.999} OR {16.4354, 182.733} OR both are missing). If three of those coordinates are missing differenceList[[5]] = "NULL".

-
@PinguinDirk That's precisely right. You need coordinates with both k = 4 and k = 3 values simultaneously in at least one sublist in order to not have a NULL value in the differenceList at k = 4. –  HStoley Sep 16 '13 at 7:38
@PinguinDirk Sorry, I think I'm not understanding... pointLists[[2]] has three items, the last of which has index 4? So the point corresponding to the index k = 4 is {164.74, 244.083}. –  HStoley Sep 16 '13 at 7:46
@PinguinDirk We're only guaranteed that there will only be one point in each sublist with a particular index, that pointLists will be ordered according to the element indices in indexLists, and no guarantee that there won't be arbitrarily large gaps for the indices of consecutive points in pointLists. –  HStoley Sep 16 '13 at 7:49
@PinguinDirk I'm sorry about that :(... I've tried my best to explain with the example, but please let me know how I can improve things. –  HStoley Sep 16 '13 at 8:01
–  Pinguin Dirk Sep 16 '13 at 8:13

I'm sure this isn't the fastest way, but it's a start:

Median /@
((SplitBy[ Sort@Flatten[{Position[indexLists, #], Position[indexLists, # - 1]}, 1], First] /.
{{x_Integer, y_Integer}} .. ->  Sequence[] /.
{{a_Integer, b_Integer}, {a_Integer, d_Integer}} :>
(pointLists[[a, indexLists[[a, d]]]] - pointLists[[a, indexLists[[a, b]]]])) & /@
Range@Max@Flatten@indexLists) /. Median[{}] -> Null


{Null, {0.018, -0.043}, {-0.02375, 0.0295}, Null, Null, {-0.022, -0.014}}

You can set the first element to {0, 0} if you wish.

-
I'm getting some errors with this scheme... –  HStoley Sep 16 '13 at 6:46
Try: pointLists = {{{122.534, 191.041}, {122.735, 191.023}, {122.692, 190.999}, {122.742, 191.051}}, {{16.2591, 182.679}, {16.2429, 182.681}, {16.3327, 182.726}, {16.4354, 182.733}, {16.4144, 182.737}, {16.3727, 182.759}, {16.3141, 182.783}, {16.3692, 182.77}, {16.4068, 182.745}, {16.3346, 182.712}}}; –  HStoley Sep 16 '13 at 6:46
And indexLists = {{1, 3, 4, 5}, {1, 2, 3, 4, 5, 6, 7, 8, 9, 10}}; –  HStoley Sep 16 '13 at 6:46
I see a: "Part 5 of {{122.534,191.041},{122.735,191.023},{122.692,190.999},{122.742,191.051}} does not exist." error, and some other such messages. –  HStoley Sep 16 '13 at 6:48
it works for me –  Pinguin Dirk Sep 16 '13 at 8:04

You could insert zeros in for the missing data so that it makes a rectangular packed array. The some of the operations can be made efficient. If most of the data is not missing, this might not be very wasteful. The trade-off is the use of packed arrays.

In my solution, I'm not sure I'm constructing the padded arrays in the most efficient way. Almost half the time is spent setting these up in the first timing example below. Around 40% is taken up by Median. The rest of it is pretty efficient by comparison.

medDiff[pointLists_, indexLists_] :=
paddedPts = ConstantArray[0., {Length @ pointLists, Max[indexLists], 2}];
paddedIdx = ConstantArray[0, {Length @ indexLists, Max[indexLists]}];
Do[
{i, Length@indexLists}];

picklist = Most@# Rest@# &@ Transpose @ paddedIdx;
Median /@ Pick[Transpose @ Differences[paddedPts, {0, 1}], picklist, 1]
]


Example:

medDiff[pointLists, indexLists]
(* {{0.018, -0.043}, {-0.02375, 0.0295}, Median[{}], Median[{}], {-0.022, -0.014}} *)


One can use Median[{}] to mark the null results or apply a replacement rule such as % /. _Median -> Null. One can also Prepend {0, 0} if desired.

Timing

SeedRandom[1];
coloniesN = 100;
dataN = 10^4;
len = RandomInteger[{coloniesN - 10, coloniesN}, dataN]; (* up to 10 missed meas. *)
points = Table[RandomReal[{100, 110}, {len[[i]], 2}], {i, dataN}];
indices = Table[Sort@RandomSample[Range@coloniesN, len[[i]]], {i, dataN}];

medDiff[points, indices] /. _Median -> Null // AbsoluteTiming // First
(* 0.339532 *)


We calculate all the differences and Pick out the valid ones. How much time is wasted depends on the density of the valid measurements. The array paddedIdx really works like a binary mask indicating which elements in paddedPts are valid. Offsetting this by one place, I compute a mask of which differences are valid my multiplying Most and Rest of each row to store in picklist, which passed to Pick.

It is about 80 times faster than Pinguin Dirk's solution:

doStuff[points, indices] // AbsoluteTiming // First
(* 29.595675 *)


If we increase the number of missed measurements, medDiff is only 45 times faster:

SeedRandom[1];
coloniesN = 100;
dataN = 10^4;
len = RandomInteger[{coloniesN / 10, coloniesN}, dataN]; (* up to 90 missed meas. *)
points = Table[RandomReal[{100, 110}, {len[[i]], 2}], {i, dataN}];
indices = Table[Sort@RandomSample[Range@coloniesN, len[[i]]], {i, dataN}];

medDiff[points, indices] /. _Median -> Null // AbsoluteTiming // First
doStuff[points, indices] // AbsoluteTiming // First

(*  0.229826
10.418164 *)

-
Very nice answer - I shall rethink my answer again, when I find time. I wrote my idea based on a slightly different question, thus it is by no means efficient. But I'd love to try n see if I can get at least a bit closer to your solution, with regard to timings. –  Pinguin Dirk Sep 17 '13 at 6:53
@PinguinDirk Thanks. There's certainly some waste in my solution, so theoretically a faster one might be possible. –  Michael E2 Sep 17 '13 at 13:12

I'm also sure that this, too, isn't the fastest way, but I hope it'll help me/us understand the problem a bit better: (adjusted to the edited version)

doStuff[points_, indices_] :=
Module[{pos},
pos = Pick[#1, #2, 1] & @@@ Transpose@{Rest /@ indices, Differences /@ indices};
(Median /@
Apply[#1 - #2 &,
Map[Extract[points, #] &,
(Position[pos, #] /. {x_Integer, y_Integer} :>
{{x, Position[indices[[x]], #][[1, 1]]},
{x, Position[indices[[x]], #][[1, 1]] - 1}} & /@
Range[2, Max@pos]), {2}], {2}] /. Median[{}] -> "NULL") ~Prepend~ {0, 0}]


So:

doStuff[pointLists, indexLists]


{{0, 0}, {0.018, -0.043}, {-0.02375, 0.0295}, "NULL", "NULL", {-0.022, -0.014}}

How do I do it? In pos, I get the positions of relevant numbers (i.e. where we have an index that is only 1 bigger than the previous). Using that, I extract the relevant values of the pointsLists, subtract, median.

Note on edit: please be aware that I wrote the code on the first version and then adapted to the edit. Thus, it is definitely suboptimal. But before adjusting, I want to make sure I get the output right.

-
Thanks @MichaelE2, fixed –  Pinguin Dirk Sep 17 '13 at 6:02
f[pointLists_, indexLists_, k_ /; k > 1] :=
Module[{thr, col, cln},
col = If[
And[MemberQ[#[[2]], k],
MemberQ[#[[2]],
k - 1]], #[[1, First@Position[#[[2]], k]]] - #[[1,
First@Position[#[[2]], k - 1]]], {}] & /@ thr;
cln = col /. {} -> Sequence[] ;
If[cln == {}, "Null", Median[cln]]
]
f[pointLists_, indexLists_, 1] := {0, 0}


Using:

pointLists = {{{131.048, 243.364}, {131.046, 243.321}, {131.037,
243.357}, {130.931, 243.391}, {130.909, 243.377}}, {{164.911,
244.039}, {164.929, 243.942}, {164.74, 244.083}}, {{98.1685,
239.618}, {98.1913, 239.6}, {98.1528, 239.623}}};
indexLists = {{1, 2, 3, 5, 6}, {1, 2, 4}, {1, 2, 3}};


then

Table[f[pointLists, indexLists, j], {j, 6}]


yields:

{{-0.02375, 0.0295}}, "Null", "Null", {{-0.022, -0.014}}}

-