Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

OpenBUGS fails to converge on model that converges in WinBUGS. Precision limit?

As the title of this post says, when I try to run code and data that work fine in WinBUGS from R using BRugsFit (with coda=T), I get these errors:

Error in glm.fit(x = structure(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,  : 
  NA/NaN/Inf in foreign function call (arg 1)
In addition: Warning messages:
1: glm.fit: algorithm did not converge 
2: glm.fit: algorithm did not converge 
3: glm.fit: algorithm did not converge 
4: glm.fit: algorithm did not converge 
5: step size truncated due to divergence 

When I do tail() on the coda object, I get the same numbers over and over. On the other hand, when I run WinBUGS, save the coda, and load it into R, I get some stochastic variation like I would expect, and no warnings about convergence.

Here is my model (it uses the 'ones trick' to find the posteriors for the parameters of a Logistic-Makeham distribution).

model {
    for(i in 1:n){
            ones[i]<-1;
            # here I pre-calculate two quantities that occur several times 
            # in the PDF, to save a little processing time
            expbx[i] <- exp(b*x[i]); expbx1[i]<- 1/(1+sd*(expbx[i]-1));
            # below is the actual PDF
            p[i]<-(c+a*expbx[i]*expbx1[i])*exp(-c*x[i])*pow(expbx1[i],1/s);
            # the ones trick
            ones[i]~dbern(p[i]);
    }
b~dunif(0,1); d~dexp(1); c~dexp(1); s.raw~dflat();
    # a (lambda) parametrized this way because it comes out more accurate
    # s forced to be > 0
    a<-b*d; s<-abs(s.raw);
    # NOT a standard deviation, just s times d, to save processing time
    sd<-s*d;
    # I save all the parameters I'm interested in to one vector, for convenient
    # viewing of them all in the same window.
    out[1]<-a; out[2]<-b; out[3]<-c; out[4]<-s; out[5]<-d;
}

Here is a typical example of my data:

 list(n= 148 ,x=c( 1246,1175,1048,1169,1043,802,543,719,1296,817,1122,542,839,443,1536,700,834,232,596,1096,1118,957,974,1031,1149,1044,1108,
519,677,569,952,1243,970,1736,1262,1026,979,1543,1029,761,533,540,511,1150,1589,1169,1029,1248,1572,638,731,525,968,1467,1596,1077,712,1603,1
203,991,37,1775,893,993,913,1487,1186,1381,977,1247,857,786,615,733,1206,1059,1508,569,1205,754,886,1099,843,599,780,1294,1603,1242,883,1320,
507,1097,1193,218,1362,1181,1118,453,1291,972,787,1061,1097,1100,1117,1174,596,1305,1071,940,919,999,1209,1043,1161,1016,1025,750,423,732,
1389,907,1184,1275,944,1209,1073,1098,1348,976,817,557,211,961,880,1039,1287,736,1400,1757,1355,977,198,689,853,1251,767,768 ))

...and typical inits (I use 4 chains, thinning 20, burnin 2000, 20000 iterations)

list( d=0.001738,b=0.0009672,c=0.002451,s.raw=0.001511 )
list( d=0.006217,b=0.005596,c=0.00777,s.raw=0.007019 )
list( d=1.504E-05,b=4.825E-06,c=2.172E-07,s.raw=6.104E-05 )
list( d=0.3011,b=0.03552,c=0.1274,s.raw=0.2549 )

Does OpenBUGS simply round off to fewer significant digits than does WinBUGS, and if so, perhaps there is a setting I can change to make it stop doing that?

like image 554
bokov Avatar asked Jul 07 '12 23:07

bokov


1 Answers

My tentative answer to this seems to be some combination of...

  • Setting the format argument for bugsInits() and bugsData() commands to fg.
  • Parametrizing the prior distribution such that if the parameter is very small (in the negative double digits on the log scale) the reciprocal (or some other appropriate transformation) is being sampled.
  • Just using large thinning intervals (in my case, 80) and a LOT of iterations. OpenBUGS does not currently support over-relaxation, and that's that.
  • If some of the variables are categorical, do not attempt to include them in the same summary as the continuous variables.

To the responder who suggested turning off over-relaxation: the problem is that I can't turn it on, and without it, the iterations take forever. But that does seem to be the only choice at this time. I wish the WinBUGS manual was more specific about what is meant by using this feature with caution. Oh well, I guess I'll eventually need to read the paper they cite. But, since this is not even available in OpenBUGS it's semi-moot at the moment.

I'll leave my question open for a while, in case someone has a better or more detailed answer than what I've been able to come up with.

like image 73
bokov Avatar answered Oct 12 '22 23:10

bokov