alt text

Introduction

This is the second part of the exercise Beyesian inference You can see the first part in https://rpubs.com/TaylorCastro/716095 . For more details of the theory you can see http://en.wikipedia.org/wiki/Bayesianinference.

To complete the study, the laboratory also wants to estimate the weight average of the patients who present this ailment. To do this, you decide to choose a sample randomization of patients and obtain the necessary conclusions from the data collected in the sample. The data is assumed to come from a population N (µ, sigma^2), and it is assumed a prior distribution N (µ0, Sigma^2 0). As a representative sample, the laboratory decides to collect data from 10 people, and the values obtained are: 82.3 65.3 73.2 89.1 85.5 80.7 69.4 95.2 91 68.3 The population standard deviation = 3, is assumed to be known. From this information, the laboratory wishes to answer the following questions:

1. Calculate probability distributions, posterior means, deviations a posteriori typical and 95% probability intervals for each of the cases

µ0 = 70, 0 = 8;

I execute the command results <-normgcp (x, sigma.x, density = “normal”, params = c (mu0, sigma0), n.mu), where x is the vector containing the data, sigma x = 3, mu0 = 70, sigma0 = 8, n.mu = 1e4 (The results are stored in the results variable, which has four fields: results $ likelihood, results $ back, results $ mu and results $ mu.prior).

library(Bolstad)
x<-c(82.3,65.3,73.2,89.1,85.5,80.7,69.4,95.2,91,68.3)
results<-normgcp(x,80,sigma.x=3,density="normal",params=c(mu0=70,sigma0=8),n.mu=1e4)
## Known standard deviation :3

head(results$posterior,10)
##  [1] 0 0 0 0 0 0 0 0 0 0
results$pi
## NULL
head(results$mu,8)
## [1] 42.0000 42.0056 42.0112 42.0168 42.0224 42.0280 42.0336 42.0392

To calculate the posterior mean and posterior standard deviation, I follow a Reasoning similar to that of the first task of the case, using the sintegral function. I call I calculate the posteior mean doing

p1<-results$mu
p2<-results$posterior #I do
dens<-p1*p2 # and obtain the posterior mean as 
post.mean<-sintegral(p1,dens)$value
post.mean
## [1] 79.86135
round(post.mean,6)
## [1] 79.86135

I calculate the posterior standard deviatio doing

dens2<-(p1-post.mean)^2*p2

post.var<-sintegral(p1,dens2)$value

round(post.var,4)
## [1] 0.8875
post.sd<-sqrt(post.var)

round(post.sd,4)
## [1] 0.9421
post.sd<-sqrt(post.var)
round(post.sd,4)
## [1] 0.9421

2. Contrast the hypothesis H0: µ = 77.8 against H0: µ 6 = 77.8 in each of the cases. To calculate the probability interval, I first obtain the empirical distribution function. To do this, I use the command

p1<-results$mu
p2<-results$posterior
#and I use the command
cdf<-sintegral(p1,p2,n.pts = length(p1))$cdf

#Lower limit
lcb<-cdf$x[with(cdf,which.max(x[y<=0.025]))]
lcb
## [1] 78.0108
#Upper limit
ucb<-cdf$x[with(cdf,which.max(x[y<=0.975]))]
ucb
## [1] 81.7012

Since the value U = 77.8 is not within the calculated probability interval, we reject the null hypothesis.