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.

A list like

l = {0, 1, 2, 3, 4, 5, 7, 9, 12, 13, 18, 19}

may have subsequences of constant increase, $a_{n+1} = a_n + k$. For example:

{0, 1, 2, 3, 4, 5}
{0, 2, 4}
{1, 3, 5, 7, 9}
{1, 7, 13, 19}

How can these sublist be found with Mathematica?

I am particularly curious whether there are any reasonably fast functional solutions.

share|improve this question
    
I'm a bit confused; your example list is already ordered, so any elements drawn in order will also be ordered. Could you revise your question to use a list that is not ordered, and give applicable examples? –  Mr.Wizard Aug 20 '12 at 9:49
1  
@Mr.Wizard I think the OP wants to identify sequences a(n+1)=a(n)+k –  belisarius Aug 20 '12 at 10:23
4  
@Verde I wish that had been made more clear. –  Mr.Wizard Aug 20 '12 at 10:29

2 Answers 2

up vote 11 down vote accepted

Idea of the solution

The following attempt is based on every sequence represented as

{elnum, last, step}

where elnum is a number of elements in the sequence, last is the last element, and step is the increment. Then, the proposed algorithm is essentially recursive: when you add a new element of the original list to a list of already processed ones, you first check all already accumulated subsequences and update them if it can serve as their last number, and then add new pairs of any preceding element of the list smaller than a given one, and a given one, as starting new sequences. The result you get is a set of sequences expressed in this form, from which we filter out those sequences containing only two numbers. The complexity of this algorithm should be (according to my estimate) roughly O(n^3), where n is the size of the original list.

Top-level recursive implementation

Here is a tail-recursive implementation of this idea:

ClearAll[incrseqs];
incrseqs[l_List] := incrseqs[{}, {}, l];
incrseqs[accum_List, seqs_List, {}] :=
  Cases[seqs, {elnum_, _, _} /; elnum > 1];
incrseqs[accum_List, seqs_List, left : {new_, ___}] :=
  incrseqs[
     Append[accum, new],
     Block[{replaced},
      replaced[_] = False;
      Join[
        Replace[seqs,
          {elnum_, last_, step_} /; last + step == new :>
             (replaced[last] = True; {elnum + 1, new, step}),
          {1}],
        Cases[accum, 
           num_ /; ! replaced[num] && num < new :> {1, new, new - num}]
      ]
    ],
    Rest[left]
  ];

For example:

l = {0, 1, 2, 3, 4, 5, 7, 9, 12, 13, 18, 19};
incrseqs[l]

(* {{5,5,1},{2,4,2},{4,9,2},{2,7,3},{3,13,4},{3,19,6},{2,12,5},{2,18,9},{2,19,7}} *)

We will need to reconstruct the sequences of original list elements from these sequence specifications. Here are the top-level functions for that.This one will reconstruct the original sub-sequences (I could have implemented it somewhat more efficiently at the expenseof clarity):

Clear[reconstructSubsequnces];
reconstructSubsequnces[seqspecs_List] :=
  Replace[
     seqspecs, 
     {elnum_, last_, step_} :> (last - step*Range[elnum, 0, -1]),
     {1}
  ];

and the final one to bring it all together:

ClearAll[getLinearSubsequences];
getLinearSubsequences[l_List, subsFun_] :=
   reconstructSubsequnces[Union[subsFun[l]]];

Here is how it can be used:

getLinearSubsequences[l, incrseqs]

(* 
   {{0,2,4},{1,4,7},{2,7,12},{0,9,18},{5,12,19},{1,5,9,13},{1,7,13,19},
      {1,3,5,7,9},{0,1,2,3,4,5}}
*)

where I chose the syntax where I pass a name of the subsequence-finding function explicitly, for reasons to be clarified below.

Compiled version

While the following will be implementing an iterative version of the same algorithm and thus will have the same asymptotic complexity, we will use Compile to get a speed-up of a factor 20-30. So, here is the ugly and boring compiled code:

cmpseqs = 
  Compile[{{l, _Integer, 1}},
    Module[{seqs = Rest@{{0, 0, 0}}, sctr = 0, replaced = Table[0, {Length[l]}],
         new = {{0, 0, 0}}, repctr = 0, lnr = {0}
      },
     Do[
       sctr = 0;
       replaced = Table[0, {Length[l]}];
       repctr = 0;
       Do[
         If[seqs[[j, 2]] + seqs[[j, 3]] == l[[i]],
           replaced[[++repctr]] = seqs[[j, 2]];
           seqs[[j, 2]] = l[[i]];
           seqs[[j, 1]]++;
         ],
         {j, Length[seqs]}
       ];
       lnr = Complement[Take[l, i], Take[replaced, repctr]];
       new = Table[{1, 0, 0}, {Length[lnr]}];
       Do[
         If[ lnr[[k]] < l[[i]],
            ++sctr;
            new[[sctr, 2]] = l[[i]];
            new[[sctr, 3]] = l[[i]] - lnr[[k]];
         ],
         {k, Length[lnr]}
       ];
       seqs = Join[seqs, Take[new, sctr]],
       {i, Length[l]}
     ];
     Select[seqs, First@# > 1 &]
    ], 
    CompilationTarget -> "C",RuntimeOptions -> "Speed"];

We call it the same way as the top-level recursive function of the previous section:

cmpseqs[l]

(* {{5,5,1},{2,4,2},{4,9,2},{2,7,3},{3,13,4},{3,19,6},{2,12,5},{2,18,9},{2,19,7}} *)

and

getLinearSubsequences[l, cmpseqs]

(* 
   {{0,2,4},{1,4,7},{2,7,12},{0,9,18},{5,12,19},{1,5,9,13},{1,7,13,19},
      {1,3,5,7,9},{0,1,2,3,4,5}}
*)

Note that cmpseqs will generally (for unordered initial list with duplicates) return a different number of sequence specification lists, and in different order than incseqs. However, taking Union of both calls should yield the same results.

Benchmarks

Here are some benchmarks on larger lists:

largerTest = RandomInteger[500,{200}];
(res1 = getLinearSubsequences[largerTest,incrseqs])//Length//AbsoluteTiming
(res2 = getLinearSubsequences[largerTest,cmpseqs])//Length//AbsoluteTiming
res1==res2

(*
  {0.7958985,415}
  {0.0214844,415}
  True
*)

largeTest = RandomInteger[2000,{500}];
(res1 = getLinearSubsequences[largeTest,incrseqs])//Length//AbsoluteTiming
(res2 = getLinearSubsequences[largeTest,cmpseqs])//Length//AbsoluteTiming
res1==res2

(*
  {10.7353516,2137}
  {0.3652343,2137}
  True
*)

We get a factor 20-30 speed increase for the compiled version, which is a rather typical factor for the speed-ups possible with Compile w.r.t. decently-written (performance-wise) top-level code. I am sure that one can tune both versions more, and squeeze extra factor of 2-3 from either, but I didn't bother.

share|improve this answer

Here is my brute-force enumeration attempt:

pickLinearSubsequence[list_?VectorQ] := Module[{ss = Subsets[list, {3, Length[list]}]}, 
  Pick[ss, MatchQ[Differences[#, 2], {0 ..}] & /@ ss]]

Examples:

pickLinearSubsequence[{0, 1, 2, 3, 4, 5, 7, 9, 12, 13, 18, 19}]
   {{0, 1, 2}, {0, 2, 4}, {0, 9, 18}, {1, 2, 3}, {1, 3, 5}, {1, 4, 7}, {1, 5, 9},
    {1, 7, 13}, {2, 3, 4}, {2, 7, 12}, {3, 4, 5}, {3, 5, 7}, {5, 7, 9}, {5, 9, 13},
    {5, 12, 19}, {7, 13, 19}, {0, 1, 2, 3}, {1, 2, 3, 4}, {1, 3, 5, 7},
    {1, 5, 9, 13}, {1, 7, 13, 19}, {2, 3, 4, 5}, {3, 5, 7, 9}, {0, 1, 2, 3, 4},
    {1, 2, 3, 4, 5}, {1, 3, 5, 7, 9}, {0, 1, 2, 3, 4, 5}}

pickLinearSubsequence[{7, 3, 2, 6, 12, 1, 4, 8, 9, 5, 10, 11}]
   {{7, 6, 5}, {7, 8, 9}, {7, 9, 11}, {3, 2, 1}, {3, 6, 9}, {3, 4, 5}, {2, 6, 10},
    {6, 8, 10}, {8, 9, 10}, {9, 10, 11}, {7, 8, 9, 10}, {8, 9, 10, 11},
    {7, 8, 9, 10, 11}}
share|improve this answer
    
Thanks JM, this also answered my question, which turned out to be a duplicate.. –  cartonn Dec 5 '12 at 15:40

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.