Calculus.hessian result oddities

[Boylan, Ross]

The hessian produced by Calculus.hessian has some properties that seem odd to me. First, the resulting matrix is not symmetric (though it's really close). Second, it is not positive definite even though the eigenvalues are all positive. This might partly be the result of the first problem.

Comments?

julia> @time h6 = Calculus.hessian(RB.mylike, sol6[2])
206.988247 seconds (2.43 G allocations: 83.830 GB, 15.48% gc time)
10x10 Array{Float64,2}:
....
julia> issym(h6)
false
julia> isposdef(h6)
false
julia> eigvals(h6)
10-element Array{Float64,1}:
3.27548e5
41407.1
1873.44
809.515
385.469
336.118
269.82
64.4228
117.009
122.845
julia> [det(h6[1:i, 1:i]) for i in 1:10] # all submatrix determinants are positive
10-element Array{Any,1}:
583.356
2.93666e5
3.76334e7
4.26776e9
3.17529e12
6.53405e14
1.12477e17
8.54011e21
3.75658e26
6.65866e29

julia> showall(h6)
[583.356378500377 -215.9519377704782 42.707753088841706 -12.046905619427504 100.26418659571075 30.907721733245413 5.1999446295105445 -3522.0797654066855 4410.401341876791 29.123944464570087
-215.95193777047822 583.3512842422355 -12.0619061856546 42.69860798380086 -82.01778794681003 -20.00373356894301 -8.643894790795585 1484.6195212852567 -3316.905502014349 -2.7622075180025965
42.7077530888417 -12.061906185654598 131.30467049554812 -44.14799847519552 16.748039679197223 -10.98960509547607 43.02830286169418 234.68906012694512 -1258.6622801992237 9.86994636336174
-12.046905619427504 42.69860798380086 -44.147998475195514 131.35743705036847 -17.498616283805696 4.724998443749094 -31.481387268916205 -95.390157895702 209.5245372890838 12.903636963625035
100.26418659571075 -82.01778794681003 16.748039679197223 -17.498616283805696 766.6651018177663 -20.766430168130857 -2.164074090378522 610.206909770971 -6927.858446549729 -200.4409796490271
30.907721733245413 -20.00373356894301 -10.98960509547607 4.724998443749093 -20.766430168130857 209.8283757716871 -60.222027409043555 280.3888402431943 -1638.62573982097 -158.60599205247166
5.199944629510544 -8.643894790795585 43.02830286169418 -31.4813872689162 -2.164074090378524 -60.222027409043555 204.34068072875016 -3.631354965653126 314.50383342119943 30.12788292183477
-3522.079765406686 1484.6195212852565 234.6890601269451 -95.39015789570199 610.2069097709712 280.3888402431943 -3.631354965653125 102596.98982032725 -117630.58331133435 267.23733374838395
4410.401341876791 -3316.9055020143487 -1258.6622801992235 209.5245372890838 -6927.858446549729 -1638.62573982097 314.5038334211995 -117630.58331133435 265600.50878023007 5508.687364627178
29.123944464570084 -2.762207518002596 9.86994636336174 12.903636963625033 -200.4409796490271 -158.6059920524717 30.12788292183477 267.23733374838395 5508.687364627178 2125.8954716959997]

The largest absolute difference between symmetrical entries is < 5e-13.

Ross Boylan

[Kristoffer Carlsson]

Doesn't Calculu's Hessian just do finite difference? There will naturally be floating point errors.

isposdef returns false if the matrix is not hermitian, under the hood, potrf from Lapack is used which work with hermitian matrices.

[Boylan, Ross]

I believe Calculus.hessian does finite differences.  I had not expected it to computer both h[i, j] and h[j, i], since they should be the same.  That is, AFAIK they are the same mathematically; numerically some differences dependent on the order of operations seem inevitable.

The double calculation has the potential advantage of alerting you if the numerical issues are significant, but returning a matrix that lacks the expected symmetry seems a pretty high price.
Ross

---

From: julia...@googlegroups.com [julia...@googlegroups.com] on behalf of Kristoffer Carlsson [kcarl...@gmail.com]
Sent: Tuesday, May 31, 2016 1:26 AM
To: julia-users
Subject: [julia-users] Re: Calculus.hessian result oddities

[Patrick Kofod Mogensen]

File an issue :)