On 11.02.2013 15:00, Camarda, Carlo Giovanni wrote:
Dear Uwe,

thanks for the response.
I knew my matrix was almost singular, but is there any way to implement the 
algorithm is solve(base) with a spam-matrix?

Not sure, I think you want another threshold for the error. But I do not know if the algorithm is robust enough. Better ask the spam maintainer who probably knows much better what is actually used here. There is a lot in the documentation, you may find the answer there as well.


As workaround, I might add a ridge penalty:

     ApP <- A+10^-1*diag.spam(nrow(A))
     solve(ApP)

but this would completely jeopardize other parts of my model.

Well, just ignoring sensible defaults and relying on numerics that are not accurate with such a condition may "jeopardize" the results - just in another way.

Best,
Uwe Ligges


Thanks again,
Giancarlo



________________________________________
From: Uwe Ligges [lig...@statistik.tu-dortmund.de]
Sent: Saturday, February 09, 2013 8:04 PM
To: Camarda, Carlo Giovanni
Cc: r-h...@stat.math.ethz.ch
Subject: Re: [R]  impossible to invert a spam-object, but possible when it's a 
matrix-object

On 05.02.2013 13:29, Camarda, Carlo Giovanni wrote:
Dear R-users,

a question concerning sparse matrices in package "spam" (spam_0.29-2).

On one hand I have a spam object (n X n) from which I cannot compute the 
inverse. On the other hand, if I convert this object in a plain matrix, I can 
find the inverse without any problem.

Specifically I get the following error message:
    Error in chol.spam(a, ...) :
    Singularity problem when calculating the Cholesky factor.

Obviously I get similar behaviour when I compute determinants of these objects.

Please see below a toy-example which I created based on my actual problem.


Different algorithms are used to compute the inverse - and you are
rather close to singularity: kappa(A) results in 59638727.

Best,
Uwe Ligges








Thanks in advance for any help,
GC


library(spam)
## creating a spam matrix A
ent <- c(2312.12324929972,-2000,1000,-2000,1000,-2000,
           5031.91011235955,-2000,-2000,1000,-1,1000,-2000,
           2049.8595505618,-2000,1000,-1,-2000,5036.89635854342,
           -2000,1000,-2000,1,-2000,-2000,8119.66292134831,-2000,
           -2000,-2000,1000,-2000,5058.83426966292,-2000,-1,1000,
           -2000,2051.85434173669,-2000,1000,1,1000,-2000,-2000,
           5043.87640449438,-2000,1,1000,-2000,1000,-2000,
           2110.68820224719,-1,1,0,-1,1,-1,1)
rr <- as.integer(c(1,6,12,18,24,29,35,41,47,52,55,57,59))
cc <- as.integer(c(1,2,3,4,7,1,2,3,5,8,10,1,2,3,6,9,11,1,
                     4,5,6,7,10,2,4,5,6,8,3,4,5,6,9,12,1,4,
                     7,8,9,11,2,5,7,8,9,12,3,6,7,8,9,2,4,10,3,7,6,8))
di <- as.integer(c(12,12))
A <-   new("spam",
             entries=ent,
             colindices=cc,
             rowpointers=rr,
             dimension=di)
## calculate the determinant
det(A)
## computes the inverse
solve(A)

## transform A in a plain matrix
A1 <- matrix(c(A), di)
## calculate the determinant
det(A1)
## computes the inverse
solve(A1)


version
                 _
platform       i686-pc-linux-gnu
arch           i686
os             linux-gnu
system         i686, linux-gnu
status
major          2
minor          14.1
year           2011
month          12
day            22
svn rev        57956
language       R
version.string R version 2.14.1 (2011-12-22)

----------
This mail has been sent through the MPI for Demographic Research.  Should you 
receive a mail that is apparently from a MPI user without this text displayed, 
then the address has most likely been faked. If you are uncertain about the 
validity of this message, please check the mail header or ask your system 
administrator for assistance.

______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.



______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Reply via email to