3
$\begingroup$

I fit a lognormal model on some data points using both frequentist and Bayesian (using a non-informative prior) approaches. However, I got different results. Here are my codes and outputs:

Frequentist:

> data1 = c(0.32618457, 0.29166954, 0.27427996, 0.23844847, 0.25148180)
> n=length(data1)
> 
> lln1 = function(par){ if(par[2]>0) return( -
> sum(log(dlnorm(exp(data1),par[1],par[2]))) ) else return(Inf) }
> optim(c(0,0.1),lln1)
> 
>   mu              sigma 
[1] 0.27641155     0.03091169

Bayesian with 20,000 MCMC and 4000 burn:

model
 {
        for( i in 1 : N )
         {
            x[i] ~ dlnorm(mu, tau)
         }

        mu ~ dunif(0, 1)
        tau ~ dunif(0, 1)

        sigma <- 1/tau
}
        list(N = 5, x = c(0.32618457, 0.29166954, 0.27427996, 0.23844847, 0.25148180))



Node     mean    sd  MC error   2.5%    median  97.5%   start   sample
mu      0.2417  0.2182  0.001976    0.006612    0.1759  0.8226  4000    16001
sigma   2.625   2.22    0.01899 1.049   2.015   7.755   4000    16001

Since I'm using a non-informative prior, I was wondering why the estimates of mu and sigma are different.

$\endgroup$
3
  • $\begingroup$ I don't understand your second snippet of code. Where is the MCMC? Can you make it reproducible? $\endgroup$ Commented Jan 17, 2013 at 3:34
  • $\begingroup$ I'm guessing the second code snippet is input to winbugs, via the tag. As far as the question goes, I don't think the results are that different given your sd from the winbugs output. $\endgroup$ Commented Jan 17, 2013 at 3:39
  • $\begingroup$ A uniform prior assumes that all possibles issues have the same chance to occur: it is not noninformative. Moreover, as noted by @timbp, you include a constraint on the range. It is better to use Jeffreys' prior: see stats.stackexchange.com/questions/33382/… $\endgroup$ Commented Jan 17, 2013 at 8:50

2 Answers 2

2
$\begingroup$

I think you are doing 3 mistakes:

1) in the frequentist example, you treat the data as if they were on the "normal" (logarithm) scale, while in the bayesian example, you treat them as on the "lognormal" scale. Supposing the data are on the lognormal scale, you should probably modify your frequentist example to something like:

> data1 = c(0.32618457, 0.29166954, 0.27427996, 0.23844847, 0.25148180)
> require(MASS)
> fitdistr(log(data1), "normal")
  mean           sd     
  -1.29200367    0.11039312 
 ( 0.04936930) ( 0.03490937)

You can run fitdistr(data1, "lognormal") just to see I am not kidding - you will get exactly the same result.

IMPORTANT: remember that the mu and tau in bugs and mean and sd in R are parameters of the original normal distribution, and don't confuse them with mean and sd of the lognormal distribution. See here for more info.

2) recommended uninformative prior for tau is dgamma(0.01, 0.01), for mu is flat normal like dnorm(0, 1.0E-10)

3) sigma is not 1/tau, but sqrt(1/tau), so you should modify your code (and your priors), like this:

    x[i] ~ dlnorm(mu, tau)
...
dnorm(0, 1.0E-10)
tau ~ dgamma(0.01, 0.01)
sigma <- 1/sqrt(tau)
$\endgroup$
2
$\begingroup$

You have given tau a uniform prior from 0 to 1, which implies sigma (1/tau) has a prior from 1 to infinity. So your prior has no probability for sigma less than 1, which means the posterior estimate has to be greater than 1.

If the true value is 0.03 as suggested by the first analysis, then it is impossible for your bayesian analysis to give the correct answer. This is probably why the sd for the sigma estimate is so high.

Suggestions for other priors: I am not an expert. From my reading, I gather that selecting a "non-informative" prior for a variance parameter is not easy.

The book Bayesian Modelling using WinBUGS tends to use tau <- dgamma(.01,.01). Gelman in Data Analysis using Regression and Multilevel/Hierarchical Models tends to use sigma <- dunif(0, 100) (putting the prior on sigma instead of tau).

$\endgroup$
2

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.