2
$\begingroup$

This is a homework problem. I am very new to Bayesian conjugate analysis, so hang on with me.

I have a sample $(x_1...x_n)$ of $n=20$ observations from an experiment. These observations are breakdown times for some machinery. The average of these observations is 495 minutes.

As information (given by the problem, an expert's opinion) I am suggested that the breakdown times are distributed exponentially with an average breakdown time of 400 minutes and one would be surprised if the machinery lasted more than 600 minutes.

So, by the text, my sample distribution would be $f_{X|\Theta }(x|\theta)= \frac{1}{\theta} e^{\frac{-x}{\theta}}$,$x\geq0$, and since I'm assuing the observations to be iid, my likelihood function would turn out to be $L(\theta;\vec{x_n})= \frac{1}{\theta^{n}} e^{\frac{-\sum x_i}{\theta}}$. (I would really like to use this parametrization, since I can interpret the expected value of my random variable as a time average).

My goal is to find the parameters of a conjugate prior (to then find the probability that the average machinery lifetime exceeds some value, through the posterior distribution). I know that an inverse Gamma distribution is a conjugate prior for my sample distribution. For it to be so, I must use the following parametrization: $$f_\Theta(\theta)=\frac{\beta ^{\alpha }}{\Gamma (\alpha )}\theta^{-\alpha -1}e^{-\frac{\beta}{\theta}}, \theta \geq 0$$ Using Bayes rule, I know that the posterior distribution must have the form of $$\Theta|Xn\sim IG(\alpha + n, \beta + \sum_{i=1}^{n}x_i)$$.

Attempt:

a) I interpret my hyperparameters $\alpha$ and $\beta$ as "number of prior observations" and "sum of prior observations" as the posterior distribution parameters suggest. I also suppose that the extra information (would be surprised if machinery lasted more than 600 minutes) as $P(\Theta > 600)=1-P(\Theta<600)=\gamma$ where $\gamma$ is some arbitrarily small value. This would lead me to saying that the subjective information on the average life of the machinery (average breakdown time of 400 minutes) can be translated as $\frac{\beta}{\alpha}=400$, as this would be some sort of "prior sample mean".

b2) I say that $E(\Theta)=\frac{\beta}{\alpha-1}=400$, while the extra information stays the same.

Both ways, when I try to find the parameters $\alpha$ and $\beta$ that satisfy the two conditions, I get into real trouble.

###Assuming beta/alpha = 400###

library(invgamma)
values.alfa <- seq.(1,500, by = 1)

a <- pinvgamma(600, shape = values.alfa, scale = 400*values.alfa)

The results are the following

    > a <- pinvgamma(600, shape = values.alfa, scale = (400*(values.alfa)))
> a
  [1] 0.9999958 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
  [8] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
 [15] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
 [22] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000

suggesting that the combinations of $\alpha$ that satisfies the conditions $\frac{\beta}{\alpha}=400$ and $P(\Theta > 600)=1-P(\Theta<600)=\gamma$ (where $\gamma$ is, say, 0.01) are very very small. And of course, it gets worse if I make the assumptions at point 2).

So, questions:

1) Are the interpretations of the hyperparameters $\alpha$ and $\beta$ correct?

2) Which one (and why) would be the best rappresentation of my prior information? a) or b)?

3) Is the process of finding the hyperparameters correct? Why are the values so shifted?

$\endgroup$
  • $\begingroup$ The interpretation of the prior opinion is one of several. It could be on $\theta$ but it could also be on $X$ itself. Definitely not on the posterior probability in my opinion. Question: why do you use 400*values.alfa and not 400*(values.alfa-1). $\endgroup$ – Xi'an Oct 23 '16 at 17:23
  • $\begingroup$ @Xi'an I have omitted the output for the 400*(values.alfa-1) for brevity, since the results are equivalent. The reason why I used 400*values.alfa in this case is because I interpret the information (average breakdown times = 400) as a "prior sample mean", which would simply be $\frac{\beta}{\alpha}$, and therefore $\beta = 400\alpha$. I do (dare) this because of the interpretation I give to the hyperparameters, but don't know if this is a good move! $\endgroup$ – Easymode44 Oct 23 '16 at 18:40
2
$\begingroup$

If I was to find a prior such that $\mathbb{E}[\theta]=400$ and $\mathbb{P}(\theta>600)=.01$, I would also seek a conjugate prior, i.e., an inverse Gamma $\text{IG}(\alpha,\beta)$ distribution, leading to $\beta=400(\alpha-1)$ and solving $\mathbb{P}(\theta>600|\alpha,400(\alpha-1))=.01$. Using the R function

pral <-function(alf){
   pgamma(1/600,shap=alf,rat=400*(alf-1))}

that uses the inverse of the inverse Gamma probability, I obtain by sheer graphical identification that $\alpha$ is close to 40.

$\endgroup$
  • $\begingroup$ So it is after all a parametrization problem! The problem was that I was using a rate parameter as a scale parameter. So you would confirm that using $\frac{\beta}{\alpha} = 400$ would be ill-advised? $\endgroup$ – Easymode44 Oct 24 '16 at 14:53
  • $\begingroup$ Depending on how you define the Gamma scale/rate, there is no objection in picking $\frac{\beta}{\alpha} = 400$. $\endgroup$ – Xi'an Oct 24 '16 at 16:22

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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