On Feb 18, 2015, at 1:13 PM, C W wrote: > Thanks Thierry for the pointer, that's explains the problem. > > Is there anything I can do about the matrix instability or numerical > inaccuracy?
There are matrix methods in the Rmpfr package that support increased precision, but it is implemented with S4 methods. You would probably need to retool the numDeriv functions to use the mpfrMatrix-class. -- david. > > Mike > > On Wed, Feb 18, 2015 at 11:57 AM, Thierry Onkelinx <thierry.onkel...@inbo.be >> wrote: > >> Have a look at FAQ 7.31 >> >> ir. Thierry Onkelinx >> Instituut voor natuur- en bosonderzoek / Research Institute for Nature and >> >> 2015-02-18 17:44 GMT+01:00 C W <tmrs...@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.hu...@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@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@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@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. David Winsemius Alameda, CA, USA ______________________________________________ R-help@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.
