[R] Singularity problem

[David Winsemius]

On Mar 16, 2011, at 1:32 PM, Berend Hasselman wrote:

>
> Peter Langfelder wrote:
>>
>> On Wed, Mar 16, 2011 at 8:28 AM, Feng Li <m at feng.li> wrote:
>>> Dear R,
>>>
>>> If I have remembered correctly, a square matrix is singular if and  
>>> only
>>> if
>>> its determinant is zero. I am a bit confused by the following code  
>>> error.
>>> Can someone give me a hint?
>>>
>>>> a <- matrix(c(1e20,1e2,1e3,1e3),2)
>>>> det(a)
>>> [1] 1e+23
>>>> solve(a)
>>> Error in solve.default(a) :
>>>  system is computationally singular: reciprocal condition number =  
>>> 1e-17
>>>
>>
>> You are right, a matrix is mathematically singular iff its  
>> determinant
>> is zero. However, this condition is useless in practice since in
>> practice one cares about the matrix being "computationally" singular,
>> i.e. so close to singular that it cannot be inverted using the
>> standard precision of real numbers. And that's what your matrix is
>> (and the error message you got says so).
>>
>> You can write your matrix as
>>
>> a = 1e20 * matrix (c(1, 1e-18, 1e-17, 1e-17), 2, 2)
>>
>> Compared to the first element, all of the other elements are nearly
>> zero, so the matrix is numerically nearly singular even though the
>> determinant is 1e23. A better measure of how numerically unstable the
>> inversion of a matrix is is the condition number which IIRC is
>> something like the largest eigenvalue divided by the smallest
>> eigenvalue.
>>
>
> svd(a) indicates the problem.
>
> largest singular value / smallest singular value=1e17  (condition  
> number)
> --> reciprocal condition number=1e-17
> and the standard solve can't handle that.

Actually it can if you relax the default tolerance settings:

 > solve(a, tol=1e-21)
       [,1]   [,2]
[1,]  1e-20 -1e-20
[2,] -1e-21  1e-03
 > a%*%solve(a, tol=1e-21)
     [,1] [,2]
[1,]    1    0
[2,]    0    1

(I posted a similar reply to the OP soon after his complaint hit the  
list, but the system seems to have eaten it.)

-- 

David Winsemius, MD
West Hartford, CT
>
> (pivoted) QR decomposition does help. And so does SVD.
>
> Berend
>
> --
> View this message in context: http://r.789695.n4.nabble.com/Singularity-problem-tp3382093p3382465.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at 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.

David Winsemius, MD
West Hartford, CT