Assignment 3

Author

Aki Vehtari et al.

1 General information

The exercises here refer to the lecture 3/BDA chapters 2-3 content.

The exercises constitute 90% of the Quiz 3 grade.

We prepared a quarto notebook specific to this assignment to help you get started. You still need to fill in your answers on Mycourses! You can inspect this and future templates

  • as a rendered html file (to access the qmd file click the “</> Code” button at the top right hand corner of the template)
General Instructions for Answering the Assignment Questions
  • Questions below are exact copies of the text found in the MyCourses quiz and should serve as a notebook where you can store notes and code.
  • We recommend opening these notebooks in the Aalto JupyterHub, see how to use R and RStudio remotely.
  • For inspiration for code, have a look at the BDA R Demos and the specific Assignment code notebooks
  • Recommended additional self study exercises for each chapter in BDA3 are listed in the course web page. These will help to gain deeper understanding of the topic.
  • Common questions and answers regarding installation and technical problems can be found in Frequently Asked Questions (FAQ).
  • Deadlines for all assignments can be found on the course web page and in MyCourses.
  • You are allowed to discuss assignments with your friends, but it is not allowed to copy solutions directly from other students or from internet.
  • Do not share your answers publicly.
  • Do not copy answers from the internet or from previous years. We compare the answers to the answers from previous years and to the answers from other students this year.
  • Use of AI is allowed on the course, but the most of the work needs to by the student, and you need to report whether you used AI and in which way you used them (See points 5 and 6 in Aalto guidelines for use of AI in teaching).
  • All suspected plagiarism will be reported and investigated. See more about the Aalto University Code of Academic Integrity and Handling Violations Thereof.
  • If you have any suggestions or improvements to the course material, please post in the course chat feedback channel, create an issue, or submit a pull request to the public repository!

1.1 Assignment questions

For convenience the assignment questions are copied below. Answer the questions in MyCourses.

1. Inference for normal mean and deviation


A factory has a production line for manufacturing car windshields. A sample of windshields has been taken for testing hardness. The observed hardness values \( y_1 \) can be found in the dataset windshieldy1 in the aaltobda package.

We may assume that the observations follow a normal distribution with an unknown standard deviation \( \sigma \). We wish to obtain information about the unknown average hardness \( \mu \). For simplicity we assume standard uninformative prior discussed in the book, that is, \( p(\mu, \sigma^2) \propto (\sigma^2)^{-1} \). It is not necessary to derive the posterior distribution in this quiz, as it has already been done in the book (see section 3.2). As in the lecture, define \(n\) as the number of observations, \(\bar{y} \) as the arithmatic mean estimate and \(s^2 = \frac{1}{n-1}\sum^n_{i=1} (y_i-\bar{y})^2\).

1.1 The likelihood \( p(y|\mu,\sigma) \) can be expressed as:


1.2 The resulting joint posterior for \(\mu\) and \(\sigma^2\) may be expressed as:

1.3 The resulting marginal posterior for \(\mu \) may be expressed as:

1.4 The resulting marginal posterior for \(\sigma^2 \) may be expressed as:

What can you say about the unknown \( \mu \)?

1.5 Compute and report the point estimate \( E(\mu|y) \) (report your answer with 3 decimal digits):

1.6 Compute and report the central 95% posterior interval (report your answer with 3 decimal digits). R does not have a built in quantile function for t-distributions with non-zero mean and scale different from 1. The code template gives instructions how you can make your own function:
NB: Posterior intervals are also called credible intervals and are different from confidence intervals.

1.7 Using the code in the template, what can you say about the joint posterior of \(\mu,\sigma\)? 

What can you say about the hardness of the next windshield coming from the production line before actually measuring the hardness?

1.7 Compute and report the point estimate \(E(\tilde{y}|y)\)(report your answer with 3 decimal digits) :

1.8 Compute and report a posterior predictive 95%-interval (report your answer with 3 decimal digits):
NB: Posterior predictive intervals are different from posterior intervals of parameters of the model.

2. Inference for the difference between proportions


An experiment was performed to estimate the effect of beta-blockers on mortality of cardiac patients. A group of patients was randomly assigned to treatment and control groups: out of 674 patients receiving the control, 39 died, and out of 680 receiving the treatment, 22 died. Assume that the outcomes are independent and binomially distributed, with probabilities of death of \( p_0 \) and \( p_1 \) under the control and treatment, respectively. Set up a noninformative or weakly informative prior distribution on \( (p_0,p_1) \). In the below, n refers to the number of trials, y to the number of successes and \(\theta\) to the probability within a binomial model.

Formulate model below.

2.1 Take \(y_C\) as the number of deaths in the control group and \(y_T\) as the number of deaths in the treatment group. The data model can be written as:

2.2 In the context of the data and model, the prior:

2.3 The resulting posterior with independent \(Beta(1,1) \) priors for \(p_1\) and \(p_2\) can be expressed as:

Using the \(Beta(1,1)\) prior for \(p_0\) and \(p_1\) independently, summarize the posterior distribution for the odds ratio, \( \mathrm{OR} = (p_1/(1-p_1))/(p_0/(1-p_0)) \). With this (and any Beta prior), a the posteriors for \((p_0,p_1)\) are also Beta distributions, respectively (see equations in the book). First, use rbeta() to obtain posterior draws for the posterior distributions of \(p_0\) and \(p_1\). Then use these draws and odds ratio equation above to get the posterior of the odds ratio. This is called a push-forward distribution.  Obtain at least 1000 draws. If you are unsure how to get started, check out the code template.

2.4 Compute a point estimate for \(E(\mathrm{OR}|y_0,y_1)\).Report the result in decimals with two decimal digits:  

2.5 Compute an estimate for the posterior 95% central interval. Report the result in decimal with three decimal digits:

2.6 Now use what might be deemed a weakly informative \(Beta(2,2)\) prior on both probabilities \(p_0\) and \(p_1\) and discuss the sensitivity of your inference to your choice of prior density:

2.7 Use Frank Harrell's recommendations on how to state results in Bayesian two group comparison:


3. Inference for the difference between normal means


Consider a case where the same factory has two production lines for manufacturing car windshields. Independent samples from the two production lines were tested for hardness. The hardness measurements for the two samples \( \mathbf{y}_1 \) and \( \mathbf{y}_2 \) be found in  the datasets windshieldy1 and windshieldy2 in the aaltobda package and in this exercise we will model the measurements with two independent normal distributions respectively.

For the model, we assume that standard deviations \(\sigma_1\) and \(\sigma_2\) of the normal models are unknown. Let \(\bar{y_1}\) and \(\bar{y_2}\) denote averages for  \( \mathbf{y}_1 \) and \( \mathbf{y}_2 \), respectively and \(s_1^2\), \(s_2^2\) denote corresponding sample standard deviations for  \( \mathbf{y}_1 \) and \( \mathbf{y}_2 \), computed in a same manner as in exercise 1. Also, \(n\) denotes number of samples in each dataset(same number for both). Use the uninformative prior and answer the following questions:

3.1 The data model may be expressed as:

3.2 The prior can be expressed as:

3.3 The resulting marginal posterior for \(\mu_1\) can be expressed as:

3.4 The resulting joint posterior for \((\mu_1,{\sigma_1}^2)\) can be expressed as:

3.5 The resulting joint posterior for \((\mu_2,{\sigma_2}^2)\) can be expressed as:


What can you say about \( \mu_d = \mu_1 - \mu_2 \)?
You can use the rtnew() function to sample from the posterior distributions of \(\mu_1\) and \(\mu_2\), and use these samples to get a sample from the distribution of the difference \(\mu_d = \mu_1 - \mu_2\). Draw at least 10000 samples.

3.6 Compute the point estimate \( E(\mu_d|y_1, y_2) \). Report the result in decimals with three decimal digits:

3.7 Compute a posterior 95%-interval. Report the result in decimals with three decimal digits:

3.8 Given this specific model, what is the probability that the means are exactly the same (\(\mu_1 = \mu_2\)) (Use Frank Harrell's recommendations)? Explain your reasoning. 

3.9 Compute the probability that \(\mu_1 < \mu_2\) (same as the probability that \(\mu_1 - \mu_2 < 0\)). Report the result in decimals with two decimal digits: