9
$\begingroup$

I am looping through the 3 million elements of a time-series (v1, o1) and corresponding values (v2, o2).

  • v1 contains a uniformly increasing timestamp (e.g. 1 to 3,000,000 in increments of 1)
  • o1 contains a non-uniformly increasing timestamp (includes gaps, meaning o1 is a shorter list)
  • v2 contains zeroes
  • o2 contains real-type values

The purpose is to assign the o2 values to v2 for each v1 timestamp. In case of a gap, a previous o2 value will be used.

Could this code be optimized?

Sample input

 v1 = {1, 2, 3, 4, 5, 6}
 o1 = {2, 3, 5}
 v2 = {0, 0, 0, 0, 0, 0}
 o2 = {0.3, 0.5, 0.1}

Sample output

 v2 = {0, 0.3, 0.5, 0.5, 0.1, 0.1}

Code

v1 = Range[1, 3000000, 1];
v2 = ConstantArray[0., 3000000];
o1 = Delete[v1, {#} & /@ RandomInteger[3000000, 500000]];
o2 = RandomReal[{0, 1}, Length[o1]];

Length[o1]
Length[v1]

y = 1;
i = 1;                              (* Revision: Set i to 1 *)
  Monitor[
    While[i<Length[v1],             (* Revision: substitute Do for While *)
      If[SameQ[v1[[i]], o1[[y]]],
        v2[[i]] = o2[[y]];
        i++;                     (* Revision: i++ *)
        y++,
        If[v1[[i]] < o1[[y]],
          v2[[i]] = o2[[y - 1]];
          i++,                   (* Revision: i++ *)
          y++]]
    ], {i, y}] // Timing

Big thank-you to all the participants emphasizing the importance of clear problem formulation, as well as providing with inputs, and expected outputs.

$\endgroup$
  • 6
    $\begingroup$ Can you please provide some example data or a procedure to generate simulated data. $\endgroup$ – s0rce Sep 2 '15 at 2:23
  • 2
    $\begingroup$ Don't do it that way - mutating a step at a time is... slow. Make a zero list of correct length. Populate the indices at o1 with o1. Fold Max over that. Use result to pull values. All are one-step operations, should handle millions of entries in a fraction of a second... $\endgroup$ – ciao Sep 2 '15 at 4:26
  • $\begingroup$ @ciao Are you assuming o2 is increasing? $\endgroup$ – Dr. belisarius Sep 2 '15 at 5:26
  • $\begingroup$ @belisarius: you mean o1? Then yes - it's time steps, so unless we're back to the future... if you meant o2 - not sure what you mean, it's just data values from my read. $\endgroup$ – ciao Sep 2 '15 at 5:28
  • $\begingroup$ @ciao I've read o2. Silly me $\endgroup$ – Dr. belisarius Sep 2 '15 at 5:32
11
$\begingroup$

Now that OP has accepted an answer, implicitly making the results from it correct, here's a much faster method than the fastest there (as in 3x to nearly order of magnitude faster, depending on amount of missing points):

v2 = c2 = ConstantArray[0, Length@v1];
v2[[o1]] = o2;
c2[[o1]] = o1;
c2[[Range[1, First@o1]]] = Range[1, First@o1];
v2 = v2[[With[{xx = c2}, Compile[{}, FoldList[Max, xx]]][]]];

Had this flash when nodding off, should be much faster still:

c2 = ConstantArray[0, Length@v1];
c2[[o1]] = 1;
If[c2[[1]] == 0, c2[[1]] = 1; o2 = Prepend[o2, 0]];
c2 = Accumulate@c2;
v2 = o2[[c2]];

Take a look at SquareOne's neat post of various methods and timings...

$\endgroup$
  • $\begingroup$ Faster indeed ! ;) Very instructive answer (for me at least) ! (As requested, see timings in my post). $\endgroup$ – SquareOne Sep 5 '15 at 17:00
  • $\begingroup$ @SquareOne: Glad you like it, give update a whirl - simpler, and much faster again in my quick tests... $\endgroup$ – ciao Sep 5 '15 at 22:31
  • $\begingroup$ @Matvei Kruglyak :appreciate the accept, glad you find answer useful. I'd humbly ask however that you return accept to SquareOne - their answer is most educational, shows effort, and they were kind enough to mention mine so if a reader has performance as top priority, they'll see mine referenced. I mark my answers as wiki so I get no points, I'd rather a deserving new poster like SquareOne get them... $\endgroup$ – ciao Sep 6 '15 at 3:23
  • $\begingroup$ Added the new faster solution to my timings. Nice ! Also I had to correct your previous timings because I forgot to include the beginning part of your code, sorry ;( However, it is only a little bit slower. $\endgroup$ – SquareOne Sep 6 '15 at 10:14
  • $\begingroup$ I did actually ask OP to accept your answer as performance is clearly the goal. You can mention my post for the detailed timings and comparison with other approaches. $\endgroup$ – SquareOne Sep 6 '15 at 10:18
9
$\begingroup$

Here are 3 other approaches.

I'll first test their output with very basic input:

n = 10;
o1 = {3, 7, 9}; o2 = {2.2, 3.3, 4.4}

It is not clear what the OP exactly needs at the extremities of the new list: for example here do the values for time steps 1 and 2 have to be 0 or 2.2 ? Same question for time step 10 ? The different approaches below show some differences regarding to this, but they could be easily adapted.

TimeSeries

Normal@Quiet@
  TimeSeriesResample[TimeSeries[o2, {o1}], {1, n, 1}, 
   ResamplingMethod -> {"Interpolation", InterpolationOrder -> 0}]

{{1, 2.2}, {2, 2.2}, {3, 2.2}, {4, 2.2}, {5, 2.2}, {6, 2.2}, {7, 3.3}, {8, 3.3}, {9, 4.4}, {10, 4.4}}

FoldList & SparseArray

(*v2 = FoldList[If[#2 == 0, #1, #2] &, Normal@SparseArray[o1 -> o2, n]]*)
v2 = FoldList[Max, Normal@SparseArray[o1 -> o2, n]]

{0, 0, 2.2, 2.2, 2.2, 2.2, 3.3, 3.3, 4.4, 4.4}

To compare to the previous result with the time series form:

Thread@{Range@n, v2}

{{1, 0}, {2, 0}, {3, 2.2}, {4, 2.2}, {5, 2.2}, {6, 2.2}, {7, 3.3}, {8, 3.3}, {9, 4.4}, {10, 4.4}}

ConstantArray & Differences

v2 = Join @@ MapThread[ConstantArray, Join @@@ {{{0}, o2}, {{First@o1 - 1}, 
 Differences@o1, {n - Last@o1 + 1}}}]

{0, 0, 2.2, 2.2, 2.2, 2.2, 3.3, 3.3, 4.4, 4.4}

Thread@{Range@n, v2}

{{1, 0}, {2, 0}, {3, 2.2}, {4, 2.2}, {5, 2.2}, {6, 2.2}, {7, 3.3}, {8,3.3}, {9, 4.4}, {10, 4.4}}

Performance

(Comparison of all the user's answers)

Here are some timings on my macbook for the following parameters:

SeedRandom[409];
n = 3000000;
n1 = 2500000; (*and n1=500000*)
v1 = Range[n];
o1 = Union@RandomInteger[{1, n}, n1];
o2 = Reverse[o1/100 // N];

(Timings are in seconds):

enter image description here

Ciao has provided the fastest answer ! It's a very instructive answer as it shows how fast can be working with Part and list of indices ([[{i1,i2,...}]]) and to always deal as much as possible with integers (FoldList[Max,...]).

Also, as commented below, the built-in TimeSeriesResample function is 40.4/0.3 =about 130x slower than Ciao's ...

$\endgroup$
  • $\begingroup$ I don't think any of these is are the fastest possible way, but seeing as OP code is broken (or question is underspecified), moot at this point... $\endgroup$ – ciao Sep 4 '15 at 22:15
  • $\begingroup$ Innovative use of Join, MapThread, ConstantArray and Differences. I get it to be about an order of magnitude faster than the original loop method. $\endgroup$ – Jack LaVigne Sep 4 '15 at 23:13
  • $\begingroup$ If you have time, please add my answer to timings - curious of results on your hardware, it is much faster on mine... $\endgroup$ – ciao Sep 5 '15 at 5:43
  • $\begingroup$ BTW, if you do re-test, might want to do it reflective of OP example - you have the opposite: OP keeps most of the data, you are dropping most of it - should change n1 to about 2500000 to mimic OP.... I'd venture the order of performance of the different solutions might get juggled if the OP example is followed $\endgroup$ – ciao Sep 5 '15 at 6:22
  • 2
    $\begingroup$ The interesting question is why TimeSeries is so slow. Is this a bug? Does WRI care? $\endgroup$ – Rolf Mertig Sep 6 '15 at 7:30
5
$\begingroup$

Expanding @ciao's comment above:

n = 3 10^6;
n1 = 10000;
v1 = Range[n];
o1 = Union@RandomInteger[{1, n}, n1];
o2 = Reverse[ o1/1000 // N]; (*Anything*)

getV2[v1_, o1_, o2_] := Module[{f, v2},
  v2 = 0 v1;
  v2[[o1]] = o1;
  (f[#1] = #2) & @@@ Transpose[{o1, o2}];
  f[0] = o2[[1]];
  v2 = f /@ FoldList[Max, v2]
  ]

getV2[v1, o1, o2]; // Timing
$\endgroup$
  • 1
    $\begingroup$ The beauty of direct references... +1 o/c! $\endgroup$ – ciao Sep 2 '15 at 5:47
  • $\begingroup$ I don't think this produces what the OP is asking for. I think he is asking for a zero-order interpolation. $\endgroup$ – m_goldberg Sep 2 '15 at 14:34
  • $\begingroup$ @m_goldberg And I don´t think what you do think :) Who knows! $\endgroup$ – Dr. belisarius Sep 2 '15 at 19:14
  • 1
    $\begingroup$ @m_goldberg This answer substitues the previous valid non-zero value for zero values in v2 (a direct substitution of o2 values is used). I applied Matvei Kruglyak code to a short version of example data and it matched belisarius's result. I believe it answers the question correctly. $\endgroup$ – Jack LaVigne Sep 2 '15 at 19:51
3
$\begingroup$

My interpretation of this question is that the OP is asking for a zero-order interpolation of his data. There for I propose the following solution.

First a proof of concept using about 50 points based on belisarius' data model.

SeedRandom[42];
  Module[{interp, n, n1, o1, o2},
    n = 50;
    n1 = 15; 
    o1 = {1, Sequence @@ Union@RandomInteger[{2, n - 1}, n - n1], n}; 
    o2 = Reverse[o1/1000 // N];
    interp = Interpolation[Transpose[{o1, o2}], InterpolationOrder -> 0];
    v2 = Table[interp[i], {i, 1, n}]];

ListPlot[v2, Filling -> Bottom]

plot

Here is timing with 3*^6 points with about 1000 gaps punched into it.

SeedRandom[42];
Module[{interp, n, n1, o1, o2},
  n = 3*^6;
  n1 = 1*^3; 
  o1 = {1, Sequence @@ Union@RandomInteger[{2, n - 1}, n - n1], n}; 
  o2 = Reverse[o1/1000 // N];
  Timing[
    interp = Interpolation[Transpose[{o1, o2}], InterpolationOrder -> 0]; 
    Table[interp[i], {i, 1, n}];]]
{11.7651, Null}

This generate the required list about 8500 times faster than the OP's estimate of 30 elements per second.

$\endgroup$
  • $\begingroup$ When I make a table of a short example o1 with gaps and associated o2, it becomes apparent that rather than getting the previous value it uses the next value. Here there is a gap between 3 and 8 Table[{i, Interpolation[{{1, 1}, {2, 5}, {3, 7}, {8, 2}, {9, 3}, {10, 1}}, InterpolationOrder -> 0][i]}, {i, 1, 10}] and the result {{1, 5}, {2, 5}, {3, 7}, {4, 2}, {5, 2}, {6, 2}, {7, 2}, {8, 2}, {9, 3}, {10, 1}} shows that 4 thru 8 are filled with 2's. Kruglyak's question asked for a previous value (check his code on a short sample). Can your answer be modified to provide that? $\endgroup$ – Jack LaVigne Sep 2 '15 at 19:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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