Hi; I am having problems inverting matrices using the function solve() For example R can not invert the following matrix [,1] [,2] [,3] [,4] [,5] [1,] 25 500 11250 275000 7.106250e+06 [2,] 500 11250 275000 7106250 1.906250e+08 [3,] 11250 275000 7106250 190625000 5.247656e+09 [4,] 275000 7106250 190625000 5247656250 1.471719e+11 [5,] 7106250 190625000 5247656250 147171875000 4.184754e+12 solve(t(xxmodel)%*%(xxmodel)) Yields the following massage: Error in solve.default(t(xxmodel) %*% (xxmodel)) : singular matrix `a' in solve The above 5X5 matrix is invertible. It has non-zero eigenvalues. Could someone explain whether there is a problem in R's solve() function. Dursun Bulutoglu
Dursun.Bulutoglu@afit.edu writes:> Hi; > > I am having problems inverting matrices using the function > solve() > > For example R can not invert the following matrix > > =20 > > [,1] [,2] [,3] [,4] > [,5] > > [1,] 25 500 11250 275000 > 7.106250e+06 > > [2,] 500 11250 275000 7106250 > 1.906250e+08 > > [3,] 11250 275000 7106250 190625000 > 5.247656e+09 > > [4,] 275000 7106250 190625000 5247656250 > 1.471719e+11 > > [5,] 7106250 190625000 5247656250 147171875000 > 4.184754e+12 > > solve(t(xxmodel)%*%(xxmodel)) > > Yields the following massage: > > Error in solve.default(t(xxmodel) %*% (xxmodel)) : singular matrix `a' > in solve > > The above 5X5 matrix is invertible. It has non-zero eigenvalues. Could > someone explain whether there is a problem in R's solve() function.Please use the r-help list unless you are sure that you have found a bug (and check up on the definition of a bug in the FAQ, and also details required when reporting). You seem to be running into a generic problem with matrix inverters, namely that they are unhappy when the columns are on widely different scales. Your diagonal elements vary by a factor of about 2e11. Since detection of linear dependency has to operate with a small tolerance to account for roundoff, this can become indistinguishable from a nonsingular matrix with a large range of eigenvalues. The typical workaround would to rescale the columns of xxmodel. -- O__ ---- Peter Dalgaard Blegdamsvej 3 c/ /'_ --- Dept. of Biostatistics 2200 Cph. N (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~~~~~~~~~ - (p.dalgaard@biostat.ku.dk) FAX: (+45) 35327907
>>>>> "PD" == Peter Dalgaard BSA <p.dalgaard@biostat.ku.dk> >>>>> on 07 Jul 2003 23:59:59 +0200 writes:PD> Dursun.Bulutoglu@afit.edu writes: >> Hi; >> >> I am having problems inverting matrices using the function >> solve() >> >> For example R can not invert the following matrix >> >> =20 >> >> [,1] [,2] [,3] [,4] >> [,5] >> >> [1,] 25 500 11250 275000 >> 7.106250e+06 >> >> [2,] 500 11250 275000 7106250 >> 1.906250e+08 >> >> [3,] 11250 275000 7106250 190625000 >> 5.247656e+09 >> >> [4,] 275000 7106250 190625000 5247656250 >> 1.471719e+11 >> >> [5,] 7106250 190625000 5247656250 147171875000 >> 4.184754e+12 >> >> solve(t(xxmodel)%*%(xxmodel)) >> >> Yields the following massage: >> >> Error in solve.default(t(xxmodel) %*% (xxmodel)) : singular matrix `a' >> in solve >> >> The above 5X5 matrix is invertible. It has non-zero eigenvalues. Could >> someone explain whether there is a problem in R's solve() function. PD> Please use the r-help list unless you are sure that you have found a PD> bug (and check up on the definition of a bug in the FAQ, and also PD> details required when reporting). yes!! PD> You seem to be running into a generic problem with matrix inverters, PD> namely that they are unhappy when the columns are on widely different PD> scales. Your diagonal elements vary by a factor of about 2e11. Since PD> detection of linear dependency has to operate with a small tolerance PD> to account for roundoff, this can become indistinguishable from a PD> nonsingular matrix with a large range of eigenvalues. PD> The typical workaround would to rescale the columns of xxmodel. 2 more things: 1) my version of R *does* invert your example matrix. Probably because newer versions of R have been using LAPACK as opposed to LINPACK. 2) If you do use a pre-LAPACK version of solve(), i.e., solve.default(); there's an argument `tol = 1e-7' which allows to set a tolerance about how singularity is determined. Martin