1
$\begingroup$

Wolfram isn't helping me much so I am curious if there are other programs out there.

I don't know what degree it is, but I have a series of numbers and I'd like to determine the linear recurrence function/coefficients so I can find any value I want via matrix exponentiation.

Are there such programs for interpolating, and how would I do it?

  • 0
    Mathematica does this2012-11-18
  • 0
    I cannot afford Mathematica and I've already used the free trial in the past. I think it uses the same FindLinearRecurrence() function that Wolfram online does, which timed out on me2012-11-18
  • 0
    I believe Excel has some functions you might find useful.2012-11-18
  • 0
    @Manos such as what? Solver? Sometimes Solver can't find an answer even when I know one exists (to which it'll find the solution if I bootstrap the numbers first)2012-11-18

3 Answers 3

3

You want the book The Book of Numbers by Conway and Guy. A description of the algorithm is at ALGORITHM but is not typeset.

  • 0
    Thanks for the response but this technically doesn't accomplish what I am after. I know what I have is some form of linear recurrence, just not how to extract the coefficients I can use in matrix exponentiation2012-11-18
  • 0
    @user49873, this gives you an upper bound on the number of coefficients needed. Once you know that, call it $n,$ and you have roughly $2n$ elements of the sequence, you can solve for the vector $c_1, \cdots, c_n$ with $n$ linear equations. How would you solve for the coefficients for a degree 2 (similar to Fibonacci) if you had some 5 consecutive entries, and you knew how to solve $Ac = b$ for a given fixed square matrix $A$ and column vector $b$? For degree 2 is is just a 2 by 2 system.2012-11-18
  • 0
    @user49873, Look, solve the degree 2 recurrence $1,1,5,29,169,985, \ldots.$ You have coefficients $c,d.$ Then you have $c + d = 5, \; \; c + 5 d = 29.$2012-11-18
  • 0
    so for the Fibonacci sequence, is this the correct way to use that method? http://i.imgur.com/hRmxa.png2012-11-18
  • 0
    @user49873, back from shopping. Yes, your treatment of the Fibonacci sequence is correct. Do try my little problem with $1,1,5,29,169,$ it's easy and illustrates the linear algebra part of finding the explicit coefficients. You will want to see the book, as it appears Conway does not give quite enough detail of how to deal with premature occurrences of $0$ on the web page. Of course, if you have a large number of consecutive sequence values, you can just go farther out to a place with no unfortunate early $0$'s.2012-11-18
  • 0
    I did. However, I also tried this method with my current sequence and the numbers got bigger instead of smaller for some reason.2012-11-19
  • 0
    It seems to me, that the Conway-method is in principle Gaussian elimination on the matrix M in the Pari/GP-script in my answer below. Perhaps the two methods can be translated into each other...2012-11-19
  • 0
    @GottfriedHelms, it could be. Conway indicates a second method on that web page as well. I think the OP wants the actual coefficients, which is just some linear algebra once you accept a bound on the degree of the recurrence.2012-11-19
  • 0
    When I first had seen Conways method (somewhere in the 70ies or 80ies in the german version of the Scientific American, I think) I was much impressed, but did never check it thoroughly and so never really used it. I'd like it, if it could be expressed using the matrix-ansatz in such a simple translation. Perhaps tomoorow I'll take a deeper look into it...2012-11-19
2

It can easily be done using for instance the freely available program Pari/GP. Here is a short excerpt of an example session.
We define the linear recurrence $x_{k+1} = B x_k + A x_{k-1} \quad $
For an example we set $[A,B]=[1,4] \quad $, then produce testdata beginning with $[x_1,x_2]=[1,2] \quad $. This is the vector of example data $[9, 38, 161, 682, 2889, 12238] \quad $ for which we want to find the recurrence-coefficients A and B (and possibly more if we do not know the order of recurrence beforehand).

Given for instance 6 numbers which are consecutive iterations, we could restate the problem as a 3x3-matrix equation: $$ \small \begin{array} {} & \cdot & \begin{bmatrix} A \\ B \\ C \end{bmatrix} \\ \begin{bmatrix} 9 & 38 & 161 \\ 38&161&682 \\ 161&682&2889 \end{bmatrix} & = & \begin{bmatrix} 682 \\ 2889 \\ 12238 \end{bmatrix}\end{array}$$ If the lhs-matrix has smaller rank than 3, then we need only two coefficients $[A,B] \quad $ instead of three $[A,B,C] \quad $.

The Pari/GP code to solve the problem:

\\pari/gp  
\\ leading ? indicates a command to Pari/GP, 
\\ leading % and a number indicates the (numbered) answer 

\\ we introduce the given vector of numbers
? v= [9, 38, 161, 682, 2889, 12238]
%141 = [9, 38, 161, 682, 2889, 12238]

\\ of six consecutive values we can generate a 3x4-matrix
? M=matrix(3,4,r,c,v[r+c-1])
%142 = 
[9 38 161 682]
[38 161 682 2889]
[161 682 2889 12238]

\\ the rank of the matrix gives the order of the linear recurrence
? rk=matrank(M)
 %143 = 2

\\ so we'll have two coefficients and need a 2x2 and a 2x1 matrix to solve
? MX=matrix(rk,rk,r,c,M[r,c])
 %144 = 
[9 38]
[38 161]

? MY=matrix(rk,1,r,c,M[r,rk+1])
 %145 = 
[161]
[682]

\\ the following command solves the problem:
? AB = MX^-1 * MY
 %146 = 
[1]
[4]

The resulting matrix contains the coefficients A and B.

As a complete function this could be written:

{solvelinrec(v)=local(rs,cs,M,rk,MX,MY,AB);
  rs=#v \ 2; cs = rs+1 ;
  M=matrix(rs,cs,r,c,v[r+c-1]);
  rk=matrank(M);
  MX=matrix(rk,rk,r,c,M[r,c]);
  MY=vectorv(rk,r,M[r,rk+1]);
  AB = matsolve(MX,MY); \\ MX^-1 * MY ;
 return(AB); }

and then

? solvelinrec([9, 38, 161, 682, 2889, 12238])
%173 = [1, 4]~

or, receiving a order-three-recurrent sequence:

 ? solvelinrec([4, 5, 13, 11, 46, 11, 184, -94, 841, -1033, 4303, -7594])
 %178 = [1, 4, -1]~

thus in this case we have $$x_k = -1 x_{k-1} + 4 x_{k-2} + 1 x_{k-3}$$

0

There are more powerful tools in Pari/GP, for example the LLL algorithm. Here is a small program that eats a vector of integers and spits out the rational generating function or 0 if it can't find it. Note that the coefficients of the denominator polynom determine the coefficients of the linear recurrence:

ggf(v)=local(l,m,p,q,B);l=length(v);B=floor(l/2);if(B<3,return(0));m=matrix(B,B,x,y,v[x-y+B+1]);q=qflll(m,4)[1];if(length(q)==0,return(0));p=sum(k=1,B,x^(k-1)*q[k,1]);q=Pol(Pol(vector(l,n,v[l-n+1]))*p+O(x^(B+1)));if(polcoeff(p,0)<0,q=-q;p=-p);q=q/p;p=Ser(q+O(x^(l+1)));for(m=1,l,if(polcoeff(p,m-1)!=v[m],return(0)));q

Example:

? ggf([1,3,7,15,31,63,127,255])
1/(2*x^2 - 3*x + 1)

So, $a_{n}=3a_{n-1}-2a_{n-2},$ starting $a_0=1, a_2=3$.

Note also that my program needs $2d$ values to determine the coefficients of a recurrence of order $d$. This is just to have a check on the result and could perhaps be optimized if really necessary.

Have fun.