Ok this makes a lot of sense, thank you very much Ilai! Cheers Guillermo On Fri, Feb 24, 2012 at 12:16 AM, ilai <ke...@math.montana.edu> wrote:
> On Thu, Feb 23, 2012 at 8:32 PM, Adan Jordan-Garza > <ajordangarza2...@my.fit.edu> wrote: > > Hello Ilai, > > thank you very much for your response, > > can I bother you a little further? > > What do you mean my model is not identifiable? > > You should read up on model identifiably or consult your local > statistician for a complete explanation, especially if you are to be > doing more analysis. An incomplete explanation just regarding your > bugs code: you are estimating "a" (an overall mean) and the three > levels b[1:3]. That's too much, you need to remove "a". Think what you > are asking of bugs to estimate (I'll make up some numbers): for bugs, > the solution > a = 10 > b[1] = 0 > b[2] = -2 > b[3] = -5 > is no different than > a = 8 > b[1] = 2 > b[2] = 0 > b[3] = -3 > or > a = 5 > b[1] = 5 > b[2] = 3 > b[3] = 0 > > This is why you're getting large intervals - even with convergence > etc. the MCMC is simply skipping between these options with no way to > "identify" which is it since it knows the means for depths 1:3, but > nothing about "a". a could be 100 and all others -99,-95 etc. > Makes sense? > > It is, for e.g., the reason why summary(glm(y~Depth)) gave you > intercept (in R the first level) and estimates for the DIFFERENCES of > all other factor levels from the first. It was not just to make it > harder on you, it's because the same idea holds for the frequentist > approach. > > It runs ok (10000 > > iteractions) no autocorrelation, and it does converge but the credible > > intervals are too big. Yes Depth has 3 levels and I coded them as 1, 2, > and > > 3. In all the examples I have seen for winBugs they use this type of > coding > > for categorical predictors. > > Like I said, winBugs may have some feature where it will build the > design matrix internally, so you can pass categorical variables as > column vectors (like R) but I seriously doubt it can guess your > variable coded 1,2,3 is categorical. If you can reference an example > like that from winBugs examples (not some random code from the > Internet) I'd really appreciate it. > What I think is happening is you're fitting > y1 = 1b[1] + e1 > y2 = 3b[3] + e2 > ... > y232 = ... + e232 > > But why should b[3] be 3 times b[1] ? that makes no sense for > categories. And the estimates are completely off (as you've noted in > your post) not just large variance. > > I have read about the matrix of dummy variables > > but I am not sure if I am supposed to code that explicitly in winBugs or > how > > to do it. I already got the dummy matrix from R. > > There are R2winbugs and other packages to let you run winbugs from R > or pass the data. I don't use them so I don't know. At the end what > you want is Depth to be a 3 column matrix, either the same as glm to > get estimate of Depth1,Depth2-1,Depth3-1 : > 1 0 0 > 1 0 0 > 1 1 0 > 1 1 0 > 1 0 1 > 1 0 1 > > Or you could have an estimate for every depth: > 1 0 0 > 1 0 0 > 0 1 0 > 0 1 0 > 0 0 1 > 0 0 1 > > And your code: > model > { > for (i in 1:232) > { > Recruitment[i]~dpois(lambda[i]) > log(lambda[i])<- b[1]*Depth[,1]+b[2]*Depth[,2]+b[3]*Depth[,3] > } > b[1]~dnorm(0,0.01) > b[2]~dnorm(0,0.01) > b[3]~dnorm(0,0.01) > } > > Now compare the result to glm. > > > from R. Do you recommend JAGS over winBugs? > > I don't use winBugs because I'm not on a Win OS. I don't think JAGS > has an advantage for simple models like this. > > Good luck, and please cc the r-help list in your answer so the thread > is complete. > > Elai. > > > > Thank you very much for your time, I really appreciate it. > > > > This is my whole code in winBugs: > > > > model > > { > > for (i in 1:232) > > { > > Recruitment[i]~dpois(lambda[i]) > > log(lambda[i])<-a+b[Depth[i]]*Depth[i] > > } > > a~dnorm(0,0.01) > > b[1]~dnorm(0,0.01) > > b[2]~dnorm(0,0.01) > > b[3]~dnorm(0,0.01) > > } > > list(a=0, b=c(0,0,0)) > > > > Depth[] Recruitment[] > > > > 1 0 > > > > 1 0 > > > > 1 1 > > > > 1 0 > > > > 1 1 > > > > 1 1 > > > > 1 1 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 1 > > > > 1 1 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 8 > > > > 1 0 > > > > 1 0 > > > > 1 2 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 1 > > > > 1 2 > > > > 1 1 > > > > 1 0 > > > > 1 1 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 0 > > > > 1 8 > > > > 2 2 > > > > 2 0 > > > > 2 0 > > > > 2 4 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 1 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 1 > > > > 2 1 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 2 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 1 > > > > 2 2 > > > > 2 1 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 1 > > > > 2 1 > > > > 2 0 > > > > 2 1 > > > > 2 1 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 1 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 0 > > > > 2 1 > > > > 2 1 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 1 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 1 > > > > 3 5 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 1 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 2 > > > > 3 0 > > > > 3 1 > > > > 3 1 > > > > 3 1 > > > > 3 0 > > > > 3 2 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 1 > > > > 3 0 > > > > 3 1 > > > > 3 1 > > > > 3 0 > > > > 3 1 > > > > 3 0 > > > > 3 1 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 2 > > > > 3 3 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 1 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 0 > > > > 3 3 > > > > 3 0 > > > > 3 0 > > > > 3 1 > > > > 3 2 > > > > 3 4 > > > > 3 0 > > > > 3 0 > > > > 3 3 > > > > 3 1 > > > > 3 0 > > > > 3 0 > > END > > > > > > > > > > On Thu, Feb 23, 2012 at 1:19 PM, ilai <ke...@math.montana.edu> wrote: > >> > >> Adan, > >> How many levels does Depth have? my wild guess: 3 and your bugs model > >> is not identifiable. > >> > >> Second, I think you may have a critical error in the way you formatted > >> the data for the bugs model. From your code it looks like you are just > >> using the factor Depth and not a design matrix of dummy variables. > >> > >> I may be wrong with respect to WinBugs (I use JAGS), but if Depth is > >> denoted as, for e.g., "low","med","high" wouldn't your multiply > >> operation "...*Depth[i] " on line 5 fail ? > >> > >> More likely Depth is denoted "1","2","3" and WinBugs thinks it's > >> numerical. Well, in that case clearly coefficients for this model > >> don't make any sense (you'll only need one b for the slope). You can > >> use model.matrix(~Depth) to get the proper format for your data. > >> > >> Hope this solves it. Next time, knowing n.chains n.iter and if they > >> achieved convergence (with different starting values) can help sort > >> through these sort of issues. > >> > >> Cheers > >> > >> Elai > >> > >> On Thu, Feb 23, 2012 at 9:57 AM, Adan Jordan-Garza > >> <ajordangarza2...@my.fit.edu> wrote: > >> > Hi, > >> > I am running a model with count data and one categorical predictor > >> > (simple > >> > model for me to understand it fully), I did in R a glm like this: > >> > glm(Recruitment~Depth, family=poisson). I get the coefficientes and > >> > confidence intervals and all is ok. But then I want to do the same > model > >> > with Bayesian stats, here is my code: > >> > > >> > model > >> > { for (i in 1:232) > >> > { > >> > Recruitment[i]~dpois(lambda[i]) > >> > log(lambda[i])<-a+b[Depth[i]]*Depth[i] > >> > } > >> > a~dnorm(0,0.000001) > >> > b[1]~dnorm(0,0.000001) > >> > b[2]~dnorm(0,0.000001) > >> > b[3]~dnorm(0,0.0000001) > >> > } > >> > list(a=0, b=c(0,0,0)) > >> > > >> > I have two problems: 1) the resulting credible intervals for the > >> > coefficients (a, b1, b2 and b3) are HUGE don t make any reasonable > >> > sense; > >> > 2) Using OpenBugs and Winbugs I get different results, > >> > > >> > if anyone can help me I appreciate a lot your time, > >> > > >> > thanks > >> > > >> > Guillermo > >> > > >> > [[alternative HTML version deleted]] > >> > > >> > ______________________________________________ > >> > R-help@r-project.org mailing list > >> > https://stat.ethz.ch/mailman/listinfo/r-help > >> > PLEASE do read the posting guide > >> > http://www.R-project.org/posting-guide.html<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 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.