Integrated approach results

[Nigel/essence]

I thought I would start a new thread as previous ones have got a bit long and off topic.

I have an algorithm which integrates NUTS, DREAM, and BFGS.

I present some initial results for a 64 dimensional highly anisotropic normal function.

In all cases, except where noted, I use about 100,000 function/derivative calculations. That is about the limit which is realistcally possible on my 'real' problem.

I use a burn-in of 10% of the total steps.

I calculate the mean and variance for the correlated direction. The analytical solution is mean 0.0, variance 10.0.

The results are very promising, but this is not surprising as the Hessian is constant. I will do some similar tests for Rosenbrock, where initial tests have shown a less clear picture, but I also want to find something for which there is an analytical solution to the cdf.....

This will be about the limit of my reporting results for this normal function - I will not be doing any more tuning or analysis of any more diagnostics, as I have other priorities.

The results for 10 repeats using different methods is:

NUTS 1 chain

mean 0.6167839377810526 var 5.508288511201127
mean 0.64074412865635 var 11.07387882387111
mean 0.5786005354458039 var 6.148639408380258
mean -0.553956160334462 var 10.271694324726974
mean -0.16574704176567948 var 10.291984898452851
mean -0.060295048642573346 var 8.135330174540236
mean 0.4206766620046686 var 8.351722833561766
mean -0.17667388800409556 var 10.844319803519346
mean 0.49097034819048385 var 8.163976158049639
mean 0.26417087200287426 var 11.565764687975726

NUTS 8 chains

mean 0.4614363575853935 var 6.971332591994299
mean 0.9822252978746535 var 11.810216948945522
mean -0.11363361170981276 var 8.440949718773384
mean -0.61272042063882 var 6.102013042650599
mean 0.7894698582375683 var 9.983850258988891
mean -0.07726062935466374 var 9.144181870153783
mean 0.19658631945875066 var 13.312626944908295
mean 0.36325926914230855 var 9.300629588537122
mean -0.07778737748101591 var 10.247546740479583
mean -0.7081519752101677 var 12.226224758092478

DREAM 16 chains

mean 1.4076046022483477 var 3.4720014716004544
mean 0.43813970357130333 var 4.97361589402284
mean -0.5185127761980877 var 13.362818781371592
etc.

DREAM 64 chains

mean 0.20177980329374878 var 24.828242114779634
mean -0.04612499550161708 var 10.821707647927411
mean -0.5109623173955199 var 9.907169552479512
mean 0.14950364967701685 var 16.800493640956425
etc.

NUTS + DREAM + BFGS 8 chains

mean -0.029567726201169306 var 10.010209977347648
mean -0.04594260971671796 var 10.491858555279794
mean -0.08142132600787867 var 10.12482862462074
mean -0.10936468607870985 var 10.096719326961662
mean -0.09224625982895168 var 9.736906762100327
mean 0.095397111290218 var 9.909153675080598
mean -0.10611103958593468 var 10.25922187711361
mean -0.022282443304893616 var 9.157979424033764
mean 0.003361049690066385 var 9.483483797914294
mean -0.08266872503596824 var 9.743713768346252

[Nigel/essence]

Here are some results for Rosenbrock. 64 dimensions.

I do the integration:

INTEGRATE[R exp(-R)]

i.e. my log likelihood is rosenbrock, and I am integration a function rosenbrock.

I don't think this has any analytical solution - if anybody can suggest an integration of the form:

INTEGRATE[f exp(-R)]

as a test case, I would be very grateful...

My results are below. As I say, much less clear situation. The amount of burn-in is clearly important, and the use of BFGS is less productive. This is easy to understand in such a complex multi modal situation. i do not use a changing Hessian during leapfrogs a la Riemann paper.

Both NUTS and NUTS+DREAM+BFGS settle on a mean of about 34 - 40, but a wider range in variance of about 40 - 80. Of course, there is a danger that all the methods get stuck in a local minimum and report a variance which is too small.

My starting points are normally distributed, mean zero variance one.

In my 'real' problem I would start off with much better starting values, using different methods. I am sure that problem with a long burn-in being needed is because of poor starting points, but I am trying to expose rather than hide problems.

I report the mean and variance, and the best value found for the Rosenbrock.

Here are the results:

NUTS 8 chains about 1M evals 10% burn-in

mean -38.32194964546049 var 180.85784031453113 max -12.745943852056271
mean -36.96379605126499 var 98.04004027249516 max -12.974918658789422
mean -38.12864613158856 var 147.44671112004062 max -13.837908616082387
mean -38.159644536804755 var 153.72261230270965 max -13.275055140742301

NUTS 8 chains about 3.1M evals 10% burn-in

mean -34.16639714810331 var 40.20492239264718 max -12.381343005832425
mean -34.10459925366785 var 45.11544218967106 max -11.661992656079981

DREAM 8 chains about 1M evals 10% burn-in

mean -118.83615260058633 var 913.6830693057509 max -75.83356719103655

NUTS + DREAM + BFGS 8 chains about 1M evals 10% burn-in

mean -49.39746340870195 var 335.17351106501746 max -10.995533790280497
mean -57.65654459582138 var 395.7876069033441 max -15.362306214992628
mean -59.93414340908645 var 406.5943820013581 max -12.757985650451177
mean -67.26903667503059 var 363.94361228173875 max -16.254710101753794
mean -52.84015964344487 var 381.06917514480347 max -12.772920133367546

NUTS + DREAM + BFGS 8 chains about 0.8M evals 20% burn-in

mean -59.14301337295896 var 355.0887918356966 max -14.646649593603888

NUTS + DREAM + BFGS 8 chains about 1.6M evals 20% burn-in

mean -44.180026781591174 var 269.94025655645714 max -13.443530270269894

NUTS + DREAM + BFGS 8 chains about 2M evals 20% burn-in

mean -38.27556925476736 var 120.19687543307032 max -12.24281558468894

NUTS + DREAM + BFGS 8 chains about 3M evals 20% burn-in

mean -36.057110754520345 var 80.47409109338219 max -12.310630009278363

NUTS + DREAM + BFGS 8 chains about 1.2M evals 50% burn-in

mean -40.094223644415976 var 155.88316817778218 max -14.08746181547723

NUTS + DREAM + BFGS 32 chains about 1.3M evals 20% burn-in

mean -85.0723409382521 var 119.52936467289965 max -24.227365226518508

NUTS + DREAM + BFGS 32 chains about 3M evals 20% burn-in

mean -69.79074789137864 var 203.54625482131542 max -14.625246406735375

NUTS + DREAM + BFGS 8 chains about 3.6M evals 20% burn-in

mean -34.41540527828443 var 45.187602930163166 max -9.746282477730288

NUTS + DREAM + BFGS 8 chains about 8.5M evals 20% burn-in

mean -35.147930108756675 var 52.46126388490143 max -12.15626247216249

[Bob Carpenter]

Did you say how you were combining these techniques
in earlier messages? Is BFGS being used to init?

- Bob

[Nigel/essence]

No, I haven't publicly stated how i am doing it, but have given some hints.

I am using BFGS throughout, not just burn-in, and I think I preserve detailed balance but am happy to submit to experts.

At present any more detailed discussion on how I do it I would prefer to keep private by email rather than on this public group.

[Nigel/essence]

Some more Rosenbrock results. DREAM with parallel tempering (PT) and a basic HMC with PT.

HMC with PT seems to confirm the best of previous results, which is encouraging, as PT should help multi modes.

DREAM, with or without PT, gives a consistently higher mean.

Who knows which is 'correct'!

Although DREAM, by itself, does not do very well compared with NUTS on the anisotropic normal case.

But that doesn't mean that DREAM does not do better in a multi modal situation.

You will see in these results classic cases where repeatability does NOT show convergence.

I haven't yet implemented NUTS + PT, and probably won't do so.

DREAM + PT 32 chains 1.6M evals 20% burn-in

mean -89.03613681593912 var 20.37232148860403 max -78.26125190907436

DREAM + PT 8 chains 1.6M evals 20% burn-in

mean -87.29606387170365 var 36.54471184315753 max -71.01381090163459

DREAM + PT 8 chains 3.2M evals 20% burn-in

mean -80.37555180252886 var 37.44132872934097 max -59.90315831549844

HMC + PT 8 chains 2.2M evals 20% burn-in

mean -33.900531259141786 var 43.46642338297693 max -12.398650999438178
mean -33.72633201427938 var 41.235161335690236 max -12.338493763754895

HMC + PT 8 chains 5.4M evals 20% burn-in

mean -34.207301920478926 var 46.092148080574816 max -11.447937206861475
mean -34.29046547649496 var 45.52222161821214 max -11.771148503270533

[Nigel/essence]

Going back to the anisotropic normal, here are some results but without BFGS.

By comparing it is clear BFGS is helping a lot.

However, in the Rosenbrock case, BFGS does not seem to help at all. This may be just be a symptom of the multi modal situation.

NUTS + DREAM + NO BFGS 8 chains

mean -0.4527932902703773 var 10.75539596413109
mean -0.3816234498994248 var 10.160075542755328
mean 0.1372849012659407 var 9.034381922699243
mean 0.5268646590066911 var 8.382138840809764
mean -0.21538359081460606 var 9.301996360294126
mean 0.15122484873784195 var 9.187865185808029
mean 0.9613678534996046 var 13.392641695868104
mean -0.3801296569336813 var 12.353874961597565

[Matt Hoffman]

Is your BFGS method based on Zhang and Sutton's algorithm? It's not
clear to me that there's any good (non-Riemannian) alternative to that
approach that preserves detailed balance.

If you are doing something along those lines, then these results make
perfect sense. The Zhang and Sutton method is great, but (like DREAM,
and unlike RMHMC) it can only use global curvature information. The
Rosenbrock function's curvature varies wildly depending on your
location, so it's not surprising that there's no single Hessian that
gives huge performance gains. The multivariate normal, on the other
hand, has constant curvature, so this ceases to be an issue.

Matt

> --
>
>

[Nigel/essence]

No, I'm not doing anything like Zhang and Sutton.

I'm pretty sure I'm maintaining detailed balance.

I am using both NUTS and DREAM in a single scheme.

But your analysis is correct - I am using a 'global' Hessian, rather than changing it along leapfrogs, so there are obvious limitations for Rosenbrock. It would be different if all the chains were converging to the same mode.

Interestingly enough, when I try to change it along leapfrogs, or when calculating an initial epsilon, it very quickly starts exploding and showing numerical problems with very long steps. So it isn't just detailed balance which is thrown, the whole thing becomes unstable.

The choices on the changes in theta and changes in gradient, for input to updating the Hessian, are a particular set I have chosen, and there may be alternatives which maintain balance and are better for improving the NUTS aspects of the integrated algorithm - i.e. are a better 'local' approximation.

Whether my 'real' problem is as difficult as Rosenbrock I will see in the next day or two.

Sorry for teasing - if you want to know more, contact me by email.

At least my results so far show that there is a class of problems for which my integrated approach seems to perform significantly better than NUTS.

[Nigel/essence]

I've done some work making the BFGS more stable. The results are better, about the same as pure NUTS on the Rosenbrock problem (and much better than NUTS on the anisotropic normal problem).

NUTS + DREAM + BFGS 8 chains about 1M evals 10% burn-in

mean -36.43522539621117 var 133.67800542315763 max -13.003913686446266
mean -36.67449996668811 var 120.83460224912051 max -13.606218581461587
mean -41.04687058887709 var 198.23511056363145 max -13.504380553927058
mean -38.7698179246861 var 161.4843531420904 max -11.995581919827114
mean -39.9360899907532 var 183.14895200882984 max -12.164018884187646