Ok this makes a lot of sense, thank you very much Ilai!
Cheers
Guillermo
On Fri, Feb 24, 2012 at 12:16 AM, ilai <keren@math.montana.edu> wrote:
> On Thu, Feb 23, 2012 at 8:32 PM, Adan Jordan-Garza
> <ajordangarza2009@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 <keren@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
> >> <ajordangarza2009@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]]