RE: [R] eigenvalues of a circulant matrix

From: Mulholland, Tom <Tom.Mulholland_at_dpi.wa.gov.au>
Date: Tue 03 May 2005 - 17:06:14 EST


Well since I know nothing about this topic I have lurked so far, but here's my two bob's worth.

Firstly I tried to make sense of Brian's initial reply. I have got no idea who Bellman is and you have not referenced (his/her) work in a way I can access the issues you refer to. So I assumed that's exactly what Brian was talking about.

Secondly.

toeplitz(1:4)

     [,1] [,2] [,3] [,4]

[1,]    1    2    3    4
[2,]    2    1    2    3
[3,]    3    2    1    2
[4,]    4    3    2    1

require(magic)
 circulant(4)

     [,1] [,2] [,3] [,4]

[1,]    1    2    3    4
[2,]    4    1    2    3
[3,]    3    4    1    2
[4,]    2    3    4    1

So they are obviously two different things. Although I think you may have implied (not stated) that the particular combination you were using resulted in both being exactly the same.

It does appear as if in this case the (X) matrix is circulant. But then I'm no expert in even such simple things.

Then I had no idea where I was going. So I tried the variations in eigen.

I ran you code

x<-scan("h:/t.txt")
y<-x[c(109:216,1:108)]
X<-toeplitz(y)

 and then

> X[is.na(X)]

numeric(0)

So I didn't get any NAs

t1 <- eigen(X)$vectors
t2 <- eigen(X,symmetric = TRUE)$vectors

> identical(t1,t2)

[1] TRUE
>

Then

t2 <- eigen(X,symmetric = TRUE,EISPACK = TRUE)$vectors
> identical(t1,t2)

[1] FALSE
>

So there'e obviously more than one way of getting the vectors. Does the second one make more sense to you?

I also noticed in the eigen help that there are references to issues such as "IEEE 754 arithmetic","(They may also differ between methods and between platforms.)" and "or Hermitian if complex". All of these are out of my competence but they do signal to me that there are issues which may relate to hardware, digital arithmetic and other things of that ilk.

I added the comment about complex because I have a vague idea that they are related to imaginary parts that you refer to.

So not coming to any conclusion that makes sense to me, and given that there are often threads about supposed inaccuracies that have answers such as the digits you see are not always what are held by the machine I set my options(digits = 22) and noticed that some of the numbers are still going at the 22 decimal place suggesting that the machine might be incapable of producing perfectly accurate results using digital arithmetic.

My other big sphere of ignorance is complex numbers.

So I tried
X<-toeplitz(complex(real = y))
t1 <- eigen(X)$vectors

> t1[1:20]

 [1]  0.068041577278880341+0i -0.068041577140546913+0i  0.068041576864811659+0i -0.068041576452430155+0i
 [5]  0.068041575907139579+0i -0.068041575231135451+0i  0.068041574435267163+0i -0.068041573525828514+0i
 [9]  0.068041572538722991+0i -0.068041571498323253+0i  0.068041570619888622+0i -0.068041570256170081+0i
[13] 0.068041568759931989+0i -0.068041566476633147+0i 0.068041563560502477+0i -0.068041560000305007+0i [17] 0.068041555538765813+0i -0.068041549792984865+0i 0.068041544123969511+0i -0.068041537810956801+0i
> t2[1:20]
 [1]  0.068041381743976906 -0.068041381743976850  0.068041381743976781 -0.068041381743976753  0.068041381743976587
 [6] -0.068041381743976725  0.068041381743976920 -0.068041381743976836  0.068041381743976892 -0.068041381743976781
[11]  0.068041381743976781 -0.068041381743977392  0.068041381743976725 -0.068041381743976753  0.068041381743976753
[16] -0.068041381743976698  0.068041381743976587 -0.068041381743976642  0.068041381743976698 -0.068041381743976490

>

Which is again different. I have no idea what I'm doing but you do seem to get slightly different answers depending upon which method you use. I do not know if one is superior to the others or where one draws the line in terms of accuracy.

Tom

> -----Original Message-----
> From: r-help-bounces@stat.math.ethz.ch
> [mailto:r-help-bounces@stat.math.ethz.ch]On Behalf Of Globe Trotter
> Sent: Tuesday, 3 May 2005 10:51 AM
> To: r-help@stat.math.ethz.ch
> Subject: Re: [R] eigenvalues of a circulant matrix
>
>
> OK, here we go:
>
> I am submitting two attachments. The first is the datafile
> called kinv used to
> create my circulant matrix, using the following commands:
>
>
> x<-scan("kinv")
> y<-x[c(109:1,0:108)]
> X=toeplitz(y)
> eigen(X)
> write(X,ncol=216,file="test.dat")
>
> reports the following columns full of NaN's: 18, 58, 194,
> 200. (Note that
> eigen(X,symmetric=T) makes no difference and I get the same as above).
>
> The second attachment contains only the eigenvectors obtained
> on calling a
> LAPACK routine directly (from C). The eigenvalues are
> essentially the same as
> that obtained using R. Here, I use the LAPACK-recommended
> double precision
> routine dspevd() routine for symmetric matrices in packed
> storage format. Note
> the absence of the NaN's....I would be happy to send my C
> programs to whoever
> is interested.
>
> I am using
>
> :~> uname -a
> Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686
> i686 i386 GNU/Linux
>
> and R.2.0.1.
>
> Many thanks and best wishes!
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Tue May 03 17:15:31 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:31:32 EST