[R] eigenvalues of a circulant matrix

Globe Trotter itsme_410 at yahoo.com
Mon May 2 14:57:56 CEST 2005


Dear Professor Ripley:

Lets do this professionally, shall we?

--- Prof Brian Ripley <ripley at stats.ox.ac.uk> wrote:
> On Sun, 1 May 2005, someone who didn't give his name wrote:
> 
> > It is my understanding that the eigenvectors of a circulant matrix are 
> > given as follows:
> >
> > 1,omega,omega^2,....,omega^{p-1}
> >
> > where the matrix has dimension given by p x p and omega is one of p complex
> > roots of unity. (See Bellman for an excellent discussion on this).
> 
> What is the relevance of this?  Also, your reference is useless to us, 
> which is important as this all hinges on your definitions.

Bellman is an excellent book on the topic and that is what I was alluding to,
that you can calculate the eigendecomposition by hand without any costly
computations, actually.

> 
> > The matrix created by the attached row and obtained using the following 
> > commands indicates no imaginary parts for the eigenvectors. It appears 
> > that the real values are close, but not exactly so, and there is no 
> > imaginary part whatsoever.
> >
> > x<-scan("kinv.dat")       #length(x) = 216
> > y<-x[c(109:216,1:108)]
> > X<-toeplitz(y)
> > eigen(X)$vectors
> 
> We don't have "kinv.dat", but X is not circulant as usually defined.

Sorry about kinv.dat -- in the e-mail that came back to me, it read "kinv",
btw.

I am unclear why you say that "X is not circulant as usually defined" -- do you
think you could clarify? It is true I use a Toeplitz matrix to set this up, but
how does that matter? The end result in this case is still a circulant matrix
that is symmetric, is it not? I would like to know why my result is not
circulant here.

> > Note that the eigenvectors are correct, and they are indeed real, 
> > because X is symmetric.
> >
> > Is this a bug in R? Any insight if not, please!
> 
> Well, first R calls LAPACK or EISPACK, so it would be a bug in one of 
> those.  But in so far as I understand you, X is a real symmetric matrix, 
> and those have real eigenvalues and eigenvectors.

Yes, I know that R calls LAPACK (which now contains EISPACK, btw). But I also
know that LAPACK contains complex eigendecomposition routines in addition to
double precision ones and it would need to be used if there is reason to
believe that the result is complex valued. (In particular ZGESDD would do it.)

The eigendecomposition of a matrix is unique. Whatever you think of Bellman,
the book does show how the eigenvectors of a circulant matrix are given by the
complex roots of unity as given above. We have therefore exhibited an
eigendecomposition without actually going through major computations (which is
good, because statistical computing is best when you use it sparingly). 
Why then does the result differ from that in R, and why by so much? (After all,
the eigendecomposition is unique, or is that only fpr real matrices?)

> I think you are confused about the meaning of Toeplitz and circulant.

Unclear, but would like to hear about your views on the actual differences in
this specific example.

> Compare
> 
> http://mathworld.wolfram.com/CirculantMatrix.html
> http://mathworld.wolfram.com/ToeplitzMatrix.html
> 
> and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix.

In my case, my matrix is symmetric and the result is a circulant matrix.

> There is a very regretable tendency here for people to assume their 
> lack of understanding is `a bug in R'.

True, but bugs in software are not exactly rare. Though R does have very few
bugs and which is why I recommend the software to every Tom.

Besides I asked a question here because I was confused....

Many thanks and best wishes!

> -- 
> Brian D. Ripley,                  ripley at stats.ox.ac.uk
> Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
> University of Oxford,             Tel:  +44 1865 272861 (self)
> 1 South Parks Road,                     +44 1865 272866 (PA)
> Oxford OX1 3TG, UK                Fax:  +44 1865 272595
>




More information about the R-help mailing list