What is the difference between these very simple normal and hierarchical normal BUGS models?

classic Classic list List threaded Threaded
5 messages Options
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

What is the difference between these very simple normal and hierarchical normal BUGS models?

Adam Ryczkowski
Consider these two basic JAGS models, each having the same input.

**Model 1 ("normal")**

    var
       ns, #number of values
       x[ns], #individual normal measurments
       tau, #precision
       mea, #fitted mean
       sd; #fitted sd
   
    model{
    for(i in 1:ns) {
       x[i] ~ dnorm(mea,tau)
    }
    mea ~ dnorm(0, 0.001)
    tau~dgamma(0.001,0.001)
    sd<-1/tau^0.5
    }


**Model 2 ("meta-normal")**

    var
       ns, #number of values
       x[ns], #individual normal measurments
       tau, #precision
       a[ns], #meta-value for x
       tau2, #superfluous(?) precision parameter
       mea, #fitted mean
       sd; #fitted sd
   
    model{
    for(i in 1:ns) {
       a[i] ~ dnorm(mea,tau)
       x[i] ~ dnorm(a[i],tau2)
    }
   
    mea ~ dnorm(0, 0.001)
    tau~dgamma(0.001,0.001)
    tau2<-100 #or whatever else, but fixed
    sd<-1/tau^0.5
    }


I've tested these models with normal variables and I confirm that they yield the same results (although the second one should be less computationally efficient). The question: is there any subtle difference between these models? Like different assumptions about data?


Note: This question is a crosspost from http://stats.stackexchange.com/questions/57095/what-is-the-difference-between-these-very-simple-normal-and-hierarchical-normal.
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: What is the difference between these very simple normal and hierarchical normal BUGS models?

John K. Kruschke
Administrator

This is one of those cases in which writing down the math of Bayes' rule will probably help. If you write down Bayes' rule for the more complicated model, you may be able to see how it simplifies, or at least relates to, the less complex model.

When you say that the two models "yield the same results" what do you mean? The central tendency of the top-level mean may be the same in the two models, but it's not clear (without the math) why other aspects should be the same, such as tau. In particular, the a
[i] values should be shrunken versions of the x[i] values, and so tau in the complex model (the precision of the a[i]) should be larger than tau in the simpler model (the precision of the x[i]).


On Thu, Apr 25, 2013 at 5:59 AM, Adam Ryczkowski [via Doing Bayesian Data Analysis] <[hidden email]> wrote:
Consider these two basic JAGS models, each having the same input.

**Model 1 ("normal")**

    var
       ns, #number of values
       x[ns], #individual normal measurments
       tau, #precision
       mea, #fitted mean
       sd; #fitted sd
   
    model{
    for(i in 1:ns) {
       x[i] ~ dnorm(mea,tau)
    }
    mea ~ dnorm(0, 0.001)
    tau~dgamma(0.001,0.001)
    sd<-1/tau^0.5
    }


**Model 2 ("meta-normal")**

    var
       ns, #number of values
       x[ns], #individual normal measurments
       tau, #precision
       a[ns], #meta-value for x
       tau2, #superfluous(?) precision parameter
       mea, #fitted mean
       sd; #fitted sd
   
    model{
    for(i in 1:ns) {
       a[i] ~ dnorm(mea,tau)
       x[i] ~ dnorm(a[i],tau2)
    }
   
    mea ~ dnorm(0, 0.001)
    tau~dgamma(0.001,0.001)
    tau2<-100 #or whatever else, but fixed
    sd<-1/tau^0.5
    }


I've tested these models with normal variables and I confirm that they yield the same results (although the second one should be less computationally efficient). The question: is there any subtle difference between these models? Like different assumptions about data?


Note: This question is a crosspost from http://stats.stackexchange.com/questions/57095/what-is-the-difference-between-these-very-simple-normal-and-hierarchical-normal.


To start a new topic under Doing Bayesian Data Analysis, email [hidden email]
To unsubscribe from Doing Bayesian Data Analysis, click here.
NAML

Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: What is the difference between these very simple normal and hierarchical normal BUGS models?

Adam Ryczkowski
Thank you for your answer!

>When you say that the two models "yield the same results" what do you mean?
I ran several combinations of tau/tau2 with x being drawn from normal distribution (mean 100, sd 15) with length 10 and 40 samples. The post-hoc distribution of mea and tau was consistent with 100 and 15.

For the same 40 data points sampled from N(100; 15), from model 1 I get [<mean> (<HDI, as defined in your book>) ): mea = 101 (95.3; 106), sd = 16.7 (13.1; 20.6)
Model 2: mea = 101 (95.2; 106), sd = 16.7 (13.2; 20.7)

> the a[i] values should be shrunken versions of the x[i] values, and so tau in the complex model (the precision of the a[i]) should be larger than tau in the simpler model (the precision of the x[i]).

Well, they aren't. Really.

But in this, 3rd model, they are:
Model 3: mea = 99,1 (93,7; 105), sd = 13,1 (8,8; 17,6)

** Model 3 **

var
   ns,
   x[ns],
   a[ns],
   tau[ns],
   tau2,
   mea,
   sd2;

model{
for(i in 1:ns) {
   a[i] ~ dnorm(mea,tau2)
   x[i] ~ dnorm(a[i],tau[i])
   tau[i]~dgamma(0.001,0.001)
}
mea ~ dnorm(0, 0.001)
tau2~dgamma(0.001,0.001)
sd2<-1/tau2^0.5
}

There is also a second variation on the subject, let's call it Model 4:

** Model 4 **
Model 4: mea = 100 (95.1; 106), sd2 = 10 (0.014; 19.2), sd = 8,69 (0,0154; 18,9)

var
   ns,
   x[ns],
   a[ns],
   tau,
   tau2,
   mea,
   sd2,
   sd;

model{
for(i in 1:ns) {
   a[i] ~ dnorm(mea,tau2)
   x[i] ~ dnorm(a[i],tau)
}

mea ~ dnorm(0, 0.001)
tau~dgamma(0.001,0.001)
tau2~dgamma(0.001,0.001)
sd<-1/tau^0.5
sd2<-1/tau2^0.5
}

I just thought, that maybe the answer to this question is immediate simple and I just overlooked something. (My intuition tells me, that adding another layer of parameter hierarchy in normal model, when each meta-parameter is just a mean of an input variable, changes very little in the model).

I think, that the only difference between the models that may be relevant may surface when input data are heteroskedastic, but I'm still not certain what exactly will happen.

The problem arose, when I was reviewing various BUGS models for network meta-analysis. In one publication found BUGS model similar to model 1 under label "fixed effects model" and model similar to model 2 dubbed "random effects model". They yielded very similar results. The models were obviously more complex, with the biggest difference in presence of standard errors along the `x` values for each `i`.

Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: What is the difference between these very simple normal and hierarchical normal BUGS models?

John K. Kruschke
Administrator

Your "Model 2" seems to be much like Figure 15.8, p. 406, of the book, with this correspondence:

y[ij] is your x[i], with the j in y[ij] playing the role of i in x[i].

mu[j] is your a[i],

tau[j] is your tau2, with your tau2 fixed at 100 instead of estimated,

muG is your mea


tauG is your tau

Right?



On Thu, Apr 25, 2013 at 4:24 PM, Adam Ryczkowski [via Doing Bayesian Data Analysis] <[hidden email]> wrote:
Thank you for your answer!

>When you say that the two models "yield the same results" what do you mean?
I ran several combinations of tau/tau2 with x being drawn from normal distribution (mean 100, sd 15) with length 10 and 40 samples. The post-hoc distribution of mea and tau was consistent with 100 and 15.

For the same 40 data points sampled from N(100; 15), from model 1 I get [<mean> (<HDI, as defined in your book>) ): mea = 101 (95.3; 106), sd = 16.7 (13.1; 20.6)
Model 2: mea = 101 (95.2; 106), sd = 16.7 (13.2; 20.7)

> the a[i] values should be shrunken versions of the x[i] values, and so tau in the complex model (the precision of the a[i]) should be larger than tau in the simpler model (the precision of the x[i]).

Well, they aren't. Really.

But in this, 3rd model, they are:
Model 3: mea = 99,1 (93,7; 105), sd = 13,1 (8,8; 17,6)

** Model 3 **

var
   ns,
   x[ns],
   a[ns],
   tau[ns],
   tau2,
   mea,
   sd2;

model{
for(i in 1:ns) {
   a[i] ~ dnorm(mea,tau2)
   x[i] ~ dnorm(a[i],tau[i])
   tau[i]~dgamma(0.001,0.001)
}
mea ~ dnorm(0, 0.001)
tau2~dgamma(0.001,0.001)
sd2<-1/tau2^0.5
}

There is also a second variation on the subject, let's call it Model 4:

** Model 4 **
Model 4: mea = 100 (95.1; 106), sd2 = 10 (0.014; 19.2), sd = 8,69 (0,0154; 18,9)

var
   ns,
   x[ns],
   a[ns],
   tau,
   tau2,
   mea,
   sd2,
   sd;

model{
for(i in 1:ns) {
   a[i] ~ dnorm(mea,tau2)
   x[i] ~ dnorm(a[i],tau)
}

mea ~ dnorm(0, 0.001)
tau~dgamma(0.001,0.001)
tau2~dgamma(0.001,0.001)
sd<-1/tau^0.5
sd2<-1/tau2^0.5
}

I just thought, that maybe the answer to this question is immediate simple and I just overlooked something. (My intuition tells me, that adding another layer of parameter hierarchy in normal model, when each meta-parameter is just a mean of an input variable, changes very little in the model).

I think, that the only difference between the models that may be relevant may surface when input data are heteroskedastic, but I'm still not certain what exactly will happen.

The problem arose, when I was reviewing various BUGS models for network meta-analysis. In one publication found BUGS model similar to model 1 under label "fixed effects model" and model similar to model 2 dubbed "random effects model". They yielded very similar results. The models were obviously more complex, with the biggest difference in presence of standard errors along the `x` values for each `i`.




To start a new topic under Doing Bayesian Data Analysis, email [hidden email]
To unsubscribe from Doing Bayesian Data Analysis, click here.
NAML

Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: What is the difference between these very simple normal and hierarchical normal BUGS models?

Adam Ryczkowski
I don't understand, why you've written "y[ij]". Neither in your book (ch. 15.2.2, p. 407 "SystemsBrugs.R") nor in the published repository of your code, y is indexed by two indices. Assuming that by "ij" you have meant "i", there is still a fundamental difference between models. My models 2, 3 and 4 assume separate random variable "a" _for each_ data point separately; it is as if you had only one measurement for each airplane.

In that case, i.e. if I would take your model SystemsJags.R and assumed subj[i] == i, I'd get my Model 3.

But it isn't the model, I'm asking about. Please, take another look at my Model 1 and Model 2. My intuition and limited numeric evidence tells me, that for all practical purposes they are equivalent, despite the fact that one is hierarchic, and the other one is not. On the other hand, it seems unlikely that adding extra layer of hierarchy into model changes so little.

**Model 1 ("normal")**

    var
       Ndata, #number of values
       y[Ndata], #individual normal measurments
       tauG, #precision
       muG, #fitted mean
       sG; #fitted sd
   
    model{
    for(i in 1:Ndata) {
       y[i] ~ dnorm(muG,tauG)
    }
    muG ~ dnorm(0, 0.001)
    tauG~dgamma(0.001,0.001)
    sG<-1/tauG^0.5
    }


**Model 2 ("meta-normal")**

    var
       Ndata, #number of values
       y[Ndata], #individual normal measurments
       tauG, #precision
       mu[Ndata], #meta-value for x
       tauG2, #superfluous(?) precision parameter
       muG, #fitted mean
       sG; #fitted sd
   
    model{
    for(i in 1:Ndata) {
       y[i] ~ dnorm(mu[i],tauG2)
       mu[i] ~ dnorm(muG,tauG)
    }
   
    muG ~ dnorm(0, 0.001)
    tauG~dgamma(0.001,0.001)
    tauG2<-100 #or whatever else, but fixed
    sG<-1/tauG^0.5
    }


** Model 3 **

var
   Ndata,
   y[Ndata],
   mu[Ndata],
   tau[nData],
   tauG,
   muG,
   sG;

model{
for(i in 1:Ndata) {
   y[i] ~ dnorm(mu[i],tau[i])
   mu[i] ~ dnorm(muG,tauG)
   tau[i] ~ dgamma(0.001,0.001)
}
muG ~ dnorm(0, 0.001)
tauG~dgamma(0.001,0.001)
sG<-1/tauG^0.5
}

** Model 4 **

var
   Ndata,
   y[Ndata],
   mu[Ndata],
   tauG,
   tauG2,
   muG,
   sG2,
   sG;

model{
for(i in 1:Ndata) {
   y[i] ~ dnorm(mu[i],tauG)
   mu[i] ~ dnorm(muG,tauG2)
}

muG ~ dnorm(0, 0.001)
tauG~dgamma(0.001,0.001)
tauG2~dgamma(0.001,0.001)
sG<-1/tauG^0.5
sG2<-1/tauG2^0.5
}
Loading...