Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Imagine I have a series of particles floating around in $\mathrm{R}^3$, where, at any given time, particles can pop into and out of existence. Here, for some number of time points, $(t_1, \dots, t_N)$, I generate an unordered list of the positions of each of the particles present in the chamber: $(p_{(k,1)}, p_{(k,2)}, \dots) \in l_k$ s.t. I have a set of lists of particle positions, $(l_1, \dots, l_N) \in L$, each list corresponding to each time point.

I would like to generate a set of pruned lists $(q_1, \dots, q_N) \in Q$, where each $q_i \subseteq l_i$ consists of the set of particle positions in $l_i$ that likely arose from measuring the same particle that

  • existed and was measured at $\alpha$ previous time points. If, for example, we set $\alpha = 2$ and we measure a particle at time $t_i$, we would only place that measurement in $q_i$ if a very similar measurement was made at times $t_{(i-1)}$ and $t_{(i-2)}$.

  • also existed and was measured for the next $\beta$ time points. If, for example, we set $\beta = 3$ and we measure a particle at time $t_i$, we would only place that measurement in $q_i$ if a very similar measurement was made at times $t_{(i+1)}$, $t_{(i+2)}$, and $t_{(i+3)}$.

If the particles were frozen in place, generating some $q_i \subseteq l_i$ would be a simple matter calculating a few intersections between $l_i$ (using the Intersection[] function in Mathematica) under the assumption that the probability one particle pops out of existence and another arises in exactly the same place in-between measurements is arbitrarily small. However, here the particles are moving, so we need to "guess" that some measured particle position at time $t_i$ and either times $t_{(i-\delta)}$ or $t_{(i+\delta)}$ (for small enough $\delta$) corresponds to the same particle if the Euclidean distance between the measured positions is below some threshold $\tau$. So if we set $\beta = 2$, we would need to see that the particle only moves a threshold distance $\tau$ from the current time point to the next, and then only moves a distance of at most $\tau$ in the interval between that time point to the one after (i.e., for some time point $t_k$ we require that some $\lvert p_{(k, \dots)} - p_{(k+1, \dots)}\rvert \leq \tau$ and that $\lvert p_{(k+1, \dots)} - p_{(k+2, \dots)}\rvert \leq \tau$).

Said another way:

Let $(l_1, l_2, \dots, l_N) \in L$ be an array of lists. For $\alpha = 1$ and $\beta = 1$, I'd like to remove all elements in each $l_i$ that are not within some threshold Euclidean distance $\tau$ of any element in the previous list, and any element in the next list in the array. If we set $\beta = 2$, I would like to prune away elements where there is no element, $p_j$ in list $l_{(i+1)}$ that is at most a Euclidean distance $\tau$ away, and where there is no element in list $l_{(i+2)}$ that is at most a Euclidean distance $\tau$ away from $p_j$.

Provided a set of lists $(l_1, \dots, l_N) \in L$, and values for $\alpha$, $\beta$, and $\tau$, is there an elegant and fast way to generate my pruned lists $(q_1, \dots, q_N) \in Q$ (beyond writing a bunch of Do or For loops)? I am happy to handle the boundaries of the list, i.e., $l_1$ and $l_N$, in whatever manner is most convenient.

::: Let's generate a simple example :::

Let

 l1 = {{87.02,892},{1,2},{65,23},{99,99}};
 l2 = {{666,666},{777,777},{123,432},{876,28},{39,192},{87.01,892}};
 l3 = {{87.0004,892},{8683,12},{7673,777},{7662,2120}};
 l4 = {{40,40},{50,50},{87.002,892.001}};

Let $k = 3$, $\alpha = 2$, $\beta = 1$, and $\tau = 0.01$. Looking at list $l_2$, we notice that the only element in $l_3$ that is within a Euclidean distance $\tau$ of an element in $l_2$ is $p_0 = (87.0004, 892)$, displaced to the position $p_{(-1)} = (87.01, 892)$ in $l_2$. We now need to look at $l_1$ and check whether there is an element within a Euclidean distance of $\tau$ of $p_{(-1)} = (87.01, 892)$, and indeed we find that there is some $p_{(-2)} = (87.02, 892)$.

Next, we need to check that there exists some $p_1$ in $l_4$ that is within a Euclidean distance $\tau$ of $p_0$, and indeed we find that $p_1 = (87.002,892.001)$ in $l_4$ satisfies this requirement. Looking through the rest of $l_3$, we find that no other elements satisfy the above set of requirements. Therefore, we can write $q_3 = (p_0) = ((87.0004, 892))$.

We repeat this procedure to generate a pruned list $q_i$ for each $l_i$, and handle the boundary conditions in whatever manner is easiest.

Please note that if we set $\tau = 0$, we could easily just write

q3 = Intersection[l1,l2,l3,l4]

To see this, try

q3 = Intersection[Round[l1], Round[l2], Round[l3], Round[l4]]

As expected, the result is

{{87, 892}}

If we want to grab the non-rounded value(s) in $l_3$ corresponding to this result, we can use the script

q3RealData = Flatten[Table[Nearest[l3, q3[[k]]], {k, 1, Length[q3]}], 1];

Which returns

{{87.0004,892}}

Dealing with $\tau > 0$ in the manner described above in an efficient way (rather than perhaps extending the rounding trick above) is the tricky part of this question.

share|improve this question
    
If you have control over the generation of the initial array, it seems to me that it would be a lot easier to number the particles rather than repeatedly have to guess which is which. Instead of the $l_i$ one could have $(n_i, r_i)$ pairs with $n_i$ being the list of particle numbers and $r_i$ the list of positions of these particles. –  Oleksandr R. Jul 16 '13 at 5:34
    
@OleksandrR. The problem would be very simple if that were possible. Unfortunately, I have no way of ab initio guessing the particle numbers (particles spontaneously appear and disappear between time points). –  PerpetuumImmobile Jul 16 '13 at 5:55
    
What process produces the initial array? The fact that the individual particles are not conserved is not necessarily a problem. All you would need is for them to be distinguishable, the easiest way to achieve that seemingly being to number them. –  Oleksandr R. Jul 16 '13 at 7:11
    
@OleksandrR. The initial array is produced by sampling the state of the physical system at different time intervals. If I number the particles in one array, I would then need to go to the second array, check for my Euclidean distance constraint, and number the particles accordingly. I suppose a solution to my problem could be to scan through the arrays and number the particles in this fashion, then to run MemberQ checks to generate the $q_i$. –  PerpetuumImmobile Jul 16 '13 at 7:18
    
I see. So, by "sampling the state of the physical system", you mean that this is experimental data? Sorry, I had assumed it was a simulation. –  Oleksandr R. Jul 16 '13 at 7:24

2 Answers 2

Would this be helpful (where t is your threshold):

fun[t_] := 
 Intersection[l1, l2, l3, l4, 
  SameTest -> (EuclideanDistance[#1, #2] < t &)]
share|improve this answer

Here is a possible solution. Please, consider that it needs more checks, with some longer sequence of lists.

Here is the input I considered, changing your examples so to have at least a couple of points in the pruned list.

l1 = {{87.02, 892}, {1, 2}, {65, 23}, {99, 99}, {777.02, 777.02}};
l2 = {{666, 666}, {777.01, 777.01}, {123, 432}, {876, 28}, {39, 192}, {87.01, 892}};
l3 = {{87, 892}, {8683, 12}, {777, 777}, {7662, 2120}};
l4 = {{777.01, 777.01}, {40, 40}, {50, 50}, {87.002, 892.001}};

First of all I define a distance function to pick all elments up from a list according to the sample and tolerance

(* DistanceTest gives a list of True/False according to which elements in the list 
   variable are close to sample with a tolerance *)
DistanceTest[sample_, tolerance_,  list_] := 
Map[Apply[And, 
          {IntervalMemberQ[Interval[sample[[1]] + {-tolerance, tolerance}], #[[1]]], 
           IntervalMemberQ[Interval[sample[[2]] + {-tolerance, tolerance}], #[[2]]]}]&,
    list]

then I created the general function

PruneLists[listsPositions_List, a_, b_, t_] :=
Module[{q, k, currentlist, sample},
  Table[(* q is the current pruned list *)
        q = Partition[listsPositions[[k]], 1];
        If[k > 1,
           q = Partition[Map[First,Last[Map[
                 Function[{currentlist},
                           If[q =!= {},
                           q = Map[Function[sample, 
         Join[sample, Pick[currentlist, DistanceTest[Last[sample], t, currentlist]]]], q];
         q = DeleteCases[q, {{_, _}}]]],
         Reverse[Take[listsPositions, {Max[k-a, 1], k-1}]]]]], 1]];

         If[k < Length[listsPositions],
            q = Map[First,Last[Map[
                  Function[{currentlist},
                            If[q =!= {},
                               q = Map[Function[sample, 
         Join[sample, Pick[currentlist, DistanceTest[Last[sample], t, currentlist]]]], q];
         q = DeleteCases[q, {{_, _}}]]],
    Take[listsPositions, {k + 1, Min[k + b, Length[listsPositions]]}]]]]],
{k, 1, Length[listsPositions]}]]

this is an example:

PruneLists[{l1, l2, l3, l4}, 2, 1, 0.01]

If you discover any error, let me know. I'm sure there are other ways to do the same things in a better and more elegant way, but this is what I get today ;-)

share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.