Closed maxbiostat closed 5 years ago
Hi Luiz,
Thanks for your email. I am glad you've found the package to be useful.
You are absolutely right, in that the problem is mainly of numerical instability, and can easily be avoided by replacing the negative eigenvalues with a small number. This is exactly what I have now done in the GitHub version of the package (I had the code written up, I just hadn't pushed it through yet). I will make more changes to the usability of the package over the summer break, and update the CRAN version as well. Till then, I hope this fix works for you.
Kind regards, Dootika
On Thu, May 2, 2019 at 8:19 PM Luiz Max F. Carvalho notifications@github.com<mailto:notifications@github.com> wrote:
First of all, thanks for this excellent package. What follows is not a complaint but a doubt. I'm trying to compute mESS for a particular MCMC run and I keep getting the warning message: You either need more samples or x is not full column rank, which is mcse.multi trying to tell me something is off with the covariance matrix. When I compute the eigen values, I get something like this:
eigen(cc, only.values = TRUE) $values [1] 1.831821e+03 3.489175e+02 6.867367e+01 4.431231e+01 1.894525e+01 1.736171e+01 3.760838e-02 2.132381e-02 1.203988e-02 [10] 2.598515e-03 2.441643e-03 6.132745e-04 5.422644e-04 4.971468e-04 4.429832e-04 4.231605e-04 4.144702e-04 3.477708e-04 [19] 2.920540e-04 2.866941e-04 2.245915e-08 2.896050e-09 5.610856e-18 2.122652e-18 3.862830e-19 1.115808e-19 -9.999266e-19 [28] -3.026382e-18 -7.703938e-18
So you can see that the eigen values are technically negative, but very close to zero. My question is then this: would it be OK if mcse.multi just replaced those values with zeros and carried on? I'm aware that zero/and or eigen values are a sign of perfect correlations in the matrix, but this is not the case here, it's just a matter of estimating accuracy, I guess. I have tried creating a much larger mcmc run by re-sampling the matrix rows (after burn-in/warmup) but I still get the same thing. Also tried removing a bunch of columns (parameters) but it mostly doesn't work.
— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/statvats/mcmcse/issues/1, or mute the threadhttps://github.com/notifications/unsubscribe-auth/ABSZGSBRDXETSAPVRS4UOE3PTM5D5ANCNFSM4HKFBWEA.
-- Dootika Vats NSF Postdoctoral Fellow Department of Statistics University of Warwick http://warwick.ac.uk/dvats
Thanks. I'm replying on GitHub and what I'm about to post may be a bit unwieldy on email.
When I do:
$cov
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
[1,] 4.508690e+02 1839.43792585 2.590367e+02 -5.410359e-03 -6.058960e-03 -6.041934e+01 -7.971453e+01 6.977809e+00 4.384873e-02
[2,] 1.839438e+03 8393.20239013 1.201866e+03 -2.804821e-02 -2.811205e-02 -2.603808e+02 -2.817874e+02 3.567029e+01 8.906421e-02
[3,] 2.590367e+02 1201.86611571 2.438743e+02 -3.931628e-03 -3.962773e-03 -5.763222e+01 1.561851e+00 1.014230e+01 1.106768e-02
[4,] -5.410359e-03 -0.02804821 -3.931628e-03 1.490894e-07 1.127215e-07 4.756293e-04 2.009420e-04 -3.439321e-04 6.202363e-07
[5,] -6.058960e-03 -0.02811205 -3.962773e-03 1.127215e-07 1.069559e-07 7.674596e-04 7.510838e-04 -6.185671e-05 1.063666e-07
[6,] -6.041934e+01 -260.38080078 -5.763222e+01 4.756293e-04 7.674596e-04 2.802018e+02 3.211920e+01 1.698255e+01 -1.189702e-02
[7,] -7.971453e+01 -281.78735376 1.561851e+00 2.009420e-04 7.510838e-04 3.211920e+01 1.245128e+03 3.050797e+01 -7.429194e-02
[8,] 6.977809e+00 35.67028976 1.014230e+01 -3.439321e-04 -6.185671e-05 1.698255e+01 3.050797e+01 1.028821e+02 4.686248e-02
[9,] 4.384873e-02 0.08906421 1.106768e-02 6.202363e-07 1.063666e-07 -1.189702e-02 -7.429194e-02 4.686248e-02 1.640945e-03
[10,] 1.373399e-01 0.68191961 1.078437e-01 -2.722611e-06 -2.429850e-06 -4.089632e-02 1.433887e-01 4.133326e-03 -3.005101e-04
[11,] -2.253640e-01 -1.02812349 -1.496875e-01 4.896488e-06 4.052750e-06 1.324366e-01 -7.009156e-02 -3.997834e-02 -9.210270e-04
[12,] 4.417537e-02 0.25713968 3.077613e-02 -2.794113e-06 -1.729267e-06 -7.964326e-02 9.948416e-04 -1.101746e-02 -4.092132e-04
[13,] 9.865805e-02 0.46075242 9.201851e-02 -2.023048e-06 -1.498290e-06 -3.625517e-02 -1.506605e-02 -3.782062e-02 4.213200e-05
[14,] -1.665354e-01 -0.76373146 -1.002798e-01 3.136272e-06 3.763151e-06 3.722900e-02 5.250734e-02 5.559508e-02 2.758718e-05
[15,] 1.896294e-02 0.09666447 3.928307e-02 -6.904775e-07 -1.011327e-06 1.727573e-02 -1.936268e-01 5.827012e-03 9.690377e-05
[16,] 4.891444e-02 0.20631456 -3.102179e-02 -4.227462e-07 -1.253534e-06 -1.824956e-02 1.561855e-01 -2.360147e-02 -1.666229e-04
[17,] 7.993203e-02 0.31685138 6.455440e-02 -1.205433e-06 -8.761334e-07 -4.818307e-02 6.458796e-02 4.213532e-02 3.402873e-05
[18,] 1.244104e-01 0.42303223 6.489774e-02 3.156440e-07 -8.095224e-07 4.154222e-03 2.032742e-01 -2.870947e-02 2.530136e-04
[19,] -1.278612e-01 -0.38375027 -7.787572e-02 3.661545e-07 1.272861e-06 5.261163e-02 -1.159130e-01 -6.311637e-02 -1.095936e-04
[20,] -7.648128e-02 -0.35613334 -5.157641e-02 5.236341e-07 4.127951e-07 -8.582787e-03 -1.519491e-01 4.969051e-02 -1.774487e-04
[21,] -1.171872e-01 -0.13122804 4.337744e-02 -3.804692e-07 -2.387482e-07 -4.091649e-01 4.492675e-01 -1.438156e-01 -2.349281e-04
[22,] -3.387417e-02 1.19200922 4.476041e-01 -8.851648e-06 -7.423540e-06 -4.217971e-02 8.383751e-01 2.777596e-01 1.729724e-03
[23,] -1.034520e+00 -3.95886575 -1.088741e-02 1.862063e-05 1.207094e-05 -3.740257e-01 1.346322e+00 9.764847e-01 3.063859e-04
[24,] 1.779478e-01 0.60640952 4.358449e-02 -2.969944e-06 -2.765881e-06 -1.677077e-02 -1.752024e-01 3.906886e-02 3.143540e-04
[25,] 1.033370e-01 0.41593169 5.533594e-03 -2.559259e-06 -2.125947e-06 -4.284474e-02 -3.907553e-02 -4.274739e-02 -2.418951e-04
[26,] -2.812848e-01 -1.02234121 -4.911808e-02 5.529203e-06 4.891828e-06 5.961551e-02 2.142779e-01 3.678531e-03 -7.245888e-05
[27,] 5.338433e-01 1.81922855 1.307535e-01 -8.909832e-06 -8.297642e-06 -5.031231e-02 -5.256072e-01 1.172066e-01 9.430621e-04
[28,] 3.100111e-01 1.24779508 1.660078e-02 -7.677777e-06 -6.377841e-06 -1.285342e-01 -1.172266e-01 -1.282422e-01 -7.256854e-04
[29,] -8.438543e-01 -3.06702363 -1.473543e-01 1.658761e-05 1.467548e-05 1.788465e-01 6.428338e-01 1.103559e-02 -2.173766e-04
[,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18]
[1,] 1.373399e-01 -2.253640e-01 4.417537e-02 9.865805e-02 -1.665354e-01 1.896294e-02 4.891444e-02 7.993203e-02 1.244104e-01
[2,] 6.819196e-01 -1.028123e+00 2.571397e-01 4.607524e-01 -7.637315e-01 9.666447e-02 2.063146e-01 3.168514e-01 4.230322e-01
[3,] 1.078437e-01 -1.496875e-01 3.077613e-02 9.201851e-02 -1.002798e-01 3.928307e-02 -3.102179e-02 6.455440e-02 6.489774e-02
[4,] -2.722611e-06 4.896488e-06 -2.794113e-06 -2.023048e-06 3.136272e-06 -6.904775e-07 -4.227462e-07 -1.205433e-06 3.156440e-07
[5,] -2.429850e-06 4.052750e-06 -1.729267e-06 -1.498290e-06 3.763151e-06 -1.011327e-06 -1.253534e-06 -8.761334e-07 -8.095224e-07
[6,] -4.089632e-02 1.324366e-01 -7.964326e-02 -3.625517e-02 3.722900e-02 1.727573e-02 -1.824956e-02 -4.818307e-02 4.154222e-03
[7,] 1.433887e-01 -7.009156e-02 9.948416e-04 -1.506605e-02 5.250734e-02 -1.936268e-01 1.561855e-01 6.458796e-02 2.032742e-01
[8,] 4.133326e-03 -3.997834e-02 -1.101746e-02 -3.782062e-02 5.559508e-02 5.827012e-03 -2.360147e-02 4.213532e-02 -2.870947e-02
[9,] -3.005101e-04 -9.210270e-04 -4.092132e-04 4.213200e-05 2.758718e-05 9.690377e-05 -1.666229e-04 3.402873e-05 2.530136e-04
[10,] 1.050178e-03 -3.445706e-04 -3.985687e-04 -1.005398e-04 8.149446e-05 4.100269e-05 -2.195737e-05 -5.589291e-05 2.474203e-04
[11,] -3.445706e-04 1.803853e-03 -5.270515e-04 -1.889325e-04 1.326876e-04 -1.755335e-04 2.317784e-04 1.450493e-04 -3.256021e-04
[12,] -3.985687e-04 -5.270515e-04 1.343180e-03 2.473403e-04 -2.417693e-04 3.762708e-05 -4.319812e-05 -1.231851e-04 -1.748318e-04
[13,] -1.005398e-04 -1.889325e-04 2.473403e-04 1.413405e-03 -2.919275e-04 -2.440725e-04 -8.686220e-04 -3.900474e-05 -1.546670e-04
[14,] 8.149446e-05 1.326876e-04 -2.417693e-04 -2.919275e-04 1.188409e-03 -4.051392e-04 -4.839557e-04 2.332907e-05 -1.170877e-04
[15,] 4.100269e-05 -1.755335e-04 3.762708e-05 -2.440725e-04 -4.051392e-04 1.307508e-03 -6.501708e-04 -7.952558e-05 -2.422816e-05
[16,] -2.195737e-05 2.317784e-04 -4.319812e-05 -8.686220e-04 -4.839557e-04 -6.501708e-04 2.015264e-03 9.520125e-05 2.959828e-04
[17,] -5.589291e-05 1.450493e-04 -1.231851e-04 -3.900474e-05 2.332907e-05 -7.952558e-05 9.520125e-05 8.814971e-04 -2.489529e-04
[18,] 2.474203e-04 -3.256021e-04 -1.748318e-04 -1.546670e-04 -1.170877e-04 -2.422816e-05 2.959828e-04 -2.489529e-04 8.256280e-04
[19,] -1.889365e-04 3.057019e-04 -7.171783e-06 1.109094e-04 4.079833e-05 1.009384e-04 -2.526462e-04 -3.195247e-04 -2.892870e-04
[20,] -2.590832e-06 -1.251491e-04 3.051887e-04 8.276226e-05 5.296030e-05 2.815325e-06 -1.385379e-04 -3.075418e-04 -2.822571e-04
[21,] 2.567673e-04 2.111862e-04 -2.330253e-04 2.035211e-04 -1.997493e-05 -2.350390e-04 5.149287e-05 -3.520934e-04 3.898528e-04
[22,] 4.034841e-04 -1.883254e-03 -2.499540e-04 7.716808e-04 -4.697250e-04 1.012733e-03 -1.314689e-03 -1.642789e-03 -1.046851e-04
[23,] 1.065124e-03 -4.070370e-04 -9.644725e-04 -1.023165e-03 -2.893475e-04 5.775349e-04 7.349774e-04 -5.769161e-04 1.738323e-03
[24,] -4.503362e-04 -1.891402e-04 3.251224e-04 -1.147575e-04 4.305817e-05 -4.225157e-04 4.942151e-04 2.090487e-04 8.016047e-05
[25,] 2.333969e-04 -4.737678e-05 5.587498e-05 -1.012431e-04 9.219809e-05 8.904560e-06 1.404118e-07 6.195715e-05 1.418050e-04
[26,] 2.169393e-04 2.365169e-04 -3.809973e-04 2.160006e-04 -1.352563e-04 4.136112e-04 -4.943555e-04 -2.710059e-04 -2.219655e-04
[27,] -1.351009e-03 -5.674205e-04 9.753671e-04 -3.442726e-04 1.291745e-04 -1.267547e-03 1.482645e-03 6.271462e-04 2.404814e-04
[28,] 7.001908e-04 -1.421303e-04 1.676250e-04 -3.037292e-04 2.765943e-04 2.671368e-05 4.212354e-07 1.858714e-04 4.254151e-04
[29,] 6.508179e-04 7.095508e-04 -1.142992e-03 6.480018e-04 -4.057688e-04 1.240833e-03 -1.483067e-03 -8.130177e-04 -6.658965e-04
[,19] [,20] [,21] [,22] [,23] [,24] [,25] [,26] [,27]
[1,] -1.278612e-01 -7.648128e-02 -1.171872e-01 -3.387417e-02 -1.034520e+00 1.779478e-01 1.033370e-01 -2.812848e-01 5.338433e-01
[2,] -3.837503e-01 -3.561333e-01 -1.312280e-01 1.192009e+00 -3.958866e+00 6.064095e-01 4.159317e-01 -1.022341e+00 1.819229e+00
[3,] -7.787572e-02 -5.157641e-02 4.337744e-02 4.476041e-01 -1.088741e-02 4.358449e-02 5.533594e-03 -4.911808e-02 1.307535e-01
[4,] 3.661545e-07 5.236341e-07 -3.804692e-07 -8.851648e-06 1.862063e-05 -2.969944e-06 -2.559259e-06 5.529203e-06 -8.909832e-06
[5,] 1.272861e-06 4.127951e-07 -2.387482e-07 -7.423540e-06 1.207094e-05 -2.765881e-06 -2.125947e-06 4.891828e-06 -8.297642e-06
[6,] 5.261163e-02 -8.582787e-03 -4.091649e-01 -4.217971e-02 -3.740257e-01 -1.677077e-02 -4.284474e-02 5.961551e-02 -5.031231e-02
[7,] -1.159130e-01 -1.519491e-01 4.492675e-01 8.383751e-01 1.346322e+00 -1.752024e-01 -3.907553e-02 2.142779e-01 -5.256072e-01
[8,] -6.311637e-02 4.969051e-02 -1.438156e-01 2.777596e-01 9.764847e-01 3.906886e-02 -4.274739e-02 3.678531e-03 1.172066e-01
[9,] -1.095936e-04 -1.774487e-04 -2.349281e-04 1.729724e-03 3.063859e-04 3.143540e-04 -2.418951e-04 -7.245888e-05 9.430621e-04
[10,] -1.889365e-04 -2.590832e-06 2.567673e-04 4.034841e-04 1.065124e-03 -4.503362e-04 2.333969e-04 2.169393e-04 -1.351009e-03
[11,] 3.057019e-04 -1.251491e-04 2.111862e-04 -1.883254e-03 -4.070370e-04 -1.891402e-04 -4.737678e-05 2.365169e-04 -5.674205e-04
[12,] -7.171783e-06 3.051887e-04 -2.330253e-04 -2.499540e-04 -9.644725e-04 3.251224e-04 5.587498e-05 -3.809973e-04 9.753671e-04
[13,] 1.109094e-04 8.276226e-05 2.035211e-04 7.716808e-04 -1.023165e-03 -1.147575e-04 -1.012431e-04 2.160006e-04 -3.442726e-04
[14,] 4.079833e-05 5.296030e-05 -1.997493e-05 -4.697250e-04 -2.893475e-04 4.305817e-05 9.219809e-05 -1.352563e-04 1.291745e-04
[15,] 1.009384e-04 2.815325e-06 -2.350390e-04 1.012733e-03 5.775349e-04 -4.225157e-04 8.904560e-06 4.136112e-04 -1.267547e-03
[16,] -2.526462e-04 -1.385379e-04 5.149287e-05 -1.314689e-03 7.349774e-04 4.942151e-04 1.404118e-07 -4.943555e-04 1.482645e-03
[17,] -3.195247e-04 -3.075418e-04 -3.520934e-04 -1.642789e-03 -5.769161e-04 2.090487e-04 6.195715e-05 -2.710059e-04 6.271462e-04
[18,] -2.892870e-04 -2.822571e-04 3.898528e-04 -1.046851e-04 1.738323e-03 8.016047e-05 1.418050e-04 -2.219655e-04 2.404814e-04
[19,] 1.008608e-03 -3.935302e-04 3.595293e-04 2.245927e-04 -2.432710e-03 -5.244156e-04 -1.403359e-04 6.647515e-04 -1.573247e-03
[20,] -3.935302e-04 9.894766e-04 -3.972886e-04 1.522881e-03 1.271304e-03 2.352064e-04 -6.342627e-05 -1.717801e-04 7.056191e-04
[21,] 3.595293e-04 -3.972886e-04 9.352445e-03 3.557897e-03 1.461179e-03 -1.842573e-03 3.727714e-04 1.469802e-03 -5.527719e-03
[22,] 2.245927e-04 1.522881e-03 3.557897e-03 7.660959e-02 5.147712e-03 3.982529e-04 -1.595288e-03 1.197035e-03 1.194759e-03
[23,] -2.432710e-03 1.271304e-03 1.461179e-03 5.147712e-03 1.060842e-01 1.898953e-03 -2.555431e-03 6.564775e-04 5.696859e-03
[24,] -5.244156e-04 2.352064e-04 -1.842573e-03 3.982529e-04 1.898953e-03 3.514591e-03 -6.163964e-06 -3.486650e-03 1.047844e-02
[25,] -1.403359e-04 -6.342627e-05 3.727714e-04 -1.595288e-03 -2.555431e-03 -6.163964e-06 8.673456e-04 -8.558023e-04 -1.849189e-05
[26,] 6.647515e-04 -1.717801e-04 1.469802e-03 1.197035e-03 6.564775e-04 -3.486650e-03 -8.558023e-04 4.369519e-03 -1.045995e-02
[27,] -1.573247e-03 7.056191e-04 -5.527719e-03 1.194759e-03 5.696859e-03 1.047844e-02 -1.849189e-05 -1.045995e-02 3.163132e-02
[28,] -4.210078e-04 -1.902788e-04 1.118314e-03 -4.785865e-03 -7.666292e-03 -1.849189e-05 2.585899e-03 -2.567407e-03 -5.547568e-05
[29,] 1.994255e-03 -5.153403e-04 4.409405e-03 3.591106e-03 1.969433e-03 -1.045995e-02 -2.567407e-03 1.302736e-02 -3.137985e-02
[,28] [,29]
[1,] 3.100111e-01 -8.438543e-01
[2,] 1.247795e+00 -3.067024e+00
[3,] 1.660078e-02 -1.473543e-01
[4,] -7.677777e-06 1.658761e-05
[5,] -6.377841e-06 1.467548e-05
[6,] -1.285342e-01 1.788465e-01
[7,] -1.172266e-01 6.428338e-01
[8,] -1.282422e-01 1.103559e-02
[9,] -7.256854e-04 -2.173766e-04
[10,] 7.001908e-04 6.508179e-04
[11,] -1.421303e-04 7.095508e-04
[12,] 1.676250e-04 -1.142992e-03
[13,] -3.037292e-04 6.480018e-04
[14,] 2.765943e-04 -4.057688e-04
[15,] 2.671368e-05 1.240833e-03
[16,] 4.212354e-07 -1.483067e-03
[17,] 1.858714e-04 -8.130177e-04
[18,] 4.254151e-04 -6.658965e-04
[19,] -4.210078e-04 1.994255e-03
[20,] -1.902788e-04 -5.153403e-04
[21,] 1.118314e-03 4.409405e-03
[22,] -4.785865e-03 3.591106e-03
[23,] -7.666292e-03 1.969433e-03
[24,] -1.849189e-05 -1.045995e-02
[25,] 2.585899e-03 -2.567407e-03
[26,] -2.567407e-03 1.302736e-02
[27,] -5.547568e-05 -3.137985e-02
[28,] 7.806111e-03 -7.702221e-03
[29,] -7.702221e-03 3.932567e-02
$vol
[1] 0.00456343
$est
treeModel.rootHeight treeLength constant.popSize ucld.mean meanRate CP1.kappa
6.979236e+01 2.640154e+02 3.021028e+01 9.495392e-04 9.387344e-04 1.371482e+01
CP2.kappa CP3.kappa CP1.frequencies1 CP1.frequencies2 CP1.frequencies3 CP1.frequencies4
2.098898e+01 1.967292e+01 3.038562e-01 1.501525e-01 3.497484e-01 1.962429e-01
CP2.frequencies1 CP2.frequencies2 CP2.frequencies3 CP2.frequencies4 CP3.frequencies1 CP3.frequencies2
2.661797e-01 2.292657e-01 2.163105e-01 2.882441e-01 3.006934e-01 2.391483e-01
CP3.frequencies3 CP3.frequencies4 CP1.alpha CP2.alpha CP3.alpha CP1.nu
2.072808e-01 2.528775e-01 6.077391e-02 1.274148e-01 7.368951e-01 1.420660e-01
CP2.nu CP3.nu CP1.mu CP2.mu CP3.mu
5.605684e-02 8.018772e-01 4.261980e-01 1.681705e-01 2.405632e+00
$nsim
[1] 8002
$method
[1] "bm"
$large
[1] FALSE
$size
[1] 89
it's fine. But when I do
> multiESS(dt)
[1] NaN
It seems to me that the problem is that both var_mat
and covmat
need to be adjusted. Is that right?
Hmm, most likely the issue is then that your output matrix yields a non positive-definite var(output)
. This shouldn't usually happen, unless your MCMC output is highly correlated. If you know why it's happening, I would suggest doing the fix on your own, since the multiESS
function isn't too complicated. I wouldn't want to make that fix part of the package just yet, unless I have thought about it a bit more.
That's fair. And yes, the MCMC has some highly correlated blocks that can cause problems, although I'm not positive it should not be possible to compute mESS
. I'll write a custom version of multiESS
and use it. If I spot anything I'll drop you an email. Closing the issue for the time being. Thanks again.
First of all, thanks for this excellent package. What follows is not a complaint but a doubt. I'm trying to compute
mESS
for a particular MCMC run and I keep getting the warning message:You either need more samples or x is not full column rank
, which ismcse.multi
trying to tell me something is off with the covariance matrix. When I compute the eigen values, I get something like this:So you can see that the eigen values are technically negative, but very close to zero. My question is then this: would it be OK if
mcse.multi
just replaced those values with zeros and carried on? I'm aware that zero/and or eigen values are a sign of perfect correlations in the matrix, but this is not the case here, it's just a matter of estimating accuracy, I guess. I have tried creating a much larger mcmc run by re-sampling the matrix rows (after burn-in/warmup) but I still get the same thing. Also tried removing a bunch of columns (parameters) but it mostly doesn't work.