[R] Numerical stability of eigenvalue and hessian matrix in R

[C W]

Thanks Thierry for the pointer, that's explains the problem.

Is there anything I can do about the matrix instability or numerical
inaccuracy?

Mike

On Wed, Feb 18, 2015 at 11:57 AM, Thierry Onkelinx <thierry.onkelinx at inbo.be
> wrote:

> Have a look at FAQ 7.31
>
> ir. Thierry Onkelinx
> Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
> Forest
> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
> Kliniekstraat 25
> 1070 Anderlecht
> Belgium
>
> To call in the statistician after the experiment is done may be no more
> than asking him to perform a post-mortem examination: he may be able to say
> what the experiment died of. ~ Sir Ronald Aylmer Fisher
> The plural of anecdote is not data. ~ Roger Brinner
> The combination of some data and an aching desire for an answer does not
> ensure that a reasonable answer can be extracted from a given body of data.
> ~ John Tukey
>
> 2015-02-18 17:44 GMT+01:00 C W <tmrsg11 at gmail.com>:
>
>> Hi Ben and JS,
>>
>> Thanks for the reply.
>>
>> I tried using: hessian(func = h_x, x, method = "complex"), it gives zero,
>> that's good.
>>
>> # R code
>>
>> > hess.h <- hessian(func = h_x, x, method = "complex")
>>
>> > mat <- h_x(x)*hess.h - grad(h_x, x) %o% grad(h_x, x)
>>
>> > mat
>>
>>         [,1]    [,2]     [,3]    [,4]
>>
>> [1,] 2060602       0        0       0
>>
>> [2,]       0 2060602        0       0
>>
>> [3,]       0       0 -4039596 -816080
>>
>> [4,]       0       0  -816080 4039596
>>
>>
>> But later I do,
>>
>> > eigen(mat)
>>
>> $values
>>
>> [1] -4121204  4121204  2060602  2060602
>>
>> $vectors
>>
>>             [,1]        [,2] [,3] [,4]
>>
>> [1,]  0.00000000  0.00000000    1    0
>>
>> [2,]  0.00000000  0.00000000    0    1
>>
>> [3,] -0.99503719  0.09950372    0    0
>>
>> [4,] -0.09950372 -0.99503719    0    0
>>
>>
>> And here is the problem,
>>
>> > eigen(mat)$values[2] == 4121204
>>
>> [1] FALSE
>>
>> > eigen(mat)$values[2] - 4121204
>>
>> [1] -5.494803e-08
>>
>> Why doesn't the second value equal to 412104?  How do I overcome that?
>>
>> Thanks,
>>
>> Mike
>>
>> On Wed, Feb 18, 2015 at 9:13 AM, JS Huang <js.huang at protective.com>
>> wrote:
>>
>> > Hi,
>> >
>> >   Since all entries in your hessian matrix and grad vector are
>> integers, I
>> > suggest you execute the following for mat assignment.
>> >
>> > > mat <- round(h_x(x),digits=0)*round(hess.h,digits=0) - round(grad(h_x,
>> > > x),digits=0) %o% round(grad(h_x, x),digits=0)
>> > > mat
>> >          [,1]     [,2]     [,3]     [,4]
>> > [1,]        0        0        0 -4080400
>> > [2,]        0  7920000        0 -1600000
>> > [3,]        0        0 12160400        0
>> > [4,] -4080400 -1600000        0 -7920000
>> >
>> >
>> >
>> > --
>> > View this message in context:
>> >
>> http://r.789695.n4.nabble.com/Numerical-stability-of-eigenvalue-and-hessian-matrix-in-R-tp4703443p4703456.html
>> > Sent from the R help mailing list archive at Nabble.com.
>> >
>> > ______________________________________________
>> > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> > 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.
>> >
>>
>>         [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> 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.
>>
>
>

	[[alternative HTML version deleted]]