---
title: "MLE Normal"
author: "Prof. Eric A. Suess"
format: 
  html:
    self-contained: true
---

## Here we fit the MLE for the normal

$X_1, X_2, ...,X_n$ iid $Normal(\mu, \sigma^2)$, find the MMEs.

$f(x|\mu,\sigma)= \frac{1}{\sigma \sqrt{2\pi}}e^{-\frac{(x-\mu)^2}{2\sigma^2}}$

Consider a random sample from the normal distribution


```{r}
n <- 30
mu <- 70
sigma <- 3
sigma2 <- 3^2

x <- rnorm(n, mu, sigma)

x.orig <- x

x.mean <- mean(x)
x.sd <- sd(x)
x.sigma <- ((n-1)/n)*x.sd

x.mean
x.sd
x.sigma
```

# make a histogram

use the density estimator since this is a continuous random variable

```{r}
hist(x.orig, probability=T)
lines(density(x.orig))
```

MLE fitted normal 

```{r}
mu.hat.mle <- mean(x)
mu.hat.mle
sigma.hat.mle <- sqrt((1/length(x))*sum((x-mean(x))^2))
sigma.hat.mle
```


```{r}
library(fitdistrplus)

fit.norm <- fitdist(x.orig, "norm", method = "mle")
summary(fit.norm)
```

plot the fitted model to the simulated data

```{r}
denscomp(list(fit.norm), legendtext = c("Normal"))
```


```{r}
cdfcomp(list(fit.norm), legendtext = c("Normal"))
```

```{r}
qqcomp(list(fit.norm), legendtext = c("Normal"))
```

```{r}
ppcomp(list(fit.norm), legendtext = c("Normal"))
```


# exact confidence interval for mu

```{r}
mu.ci.exact <- c(x.mean - 1.96*x.sd/sqrt(n), x.mean + 1.96*x.sd/sqrt(n))
mu.ci.exact
```

# large sample confidence interval for mu

```{r}
mu.ci.large <- c(mu.hat.mle - 1.96*sigma.hat.mle/sqrt(n), mu.hat.mle +
	1.96*sigma.hat.mle/sqrt(n))
mu.ci.large
```

# parametric bootstrap the MLE

The AIC and BIC are defined as follows, with $n$ the sample size and $k$ the number of parameters in the model.

$$AIC = -2\log(L(\hat{\theta})) + 2k$$
$$BIC = -2\log(L(\hat{\theta})) + k\log(n)$$




```{r}
fit.norm <- fitdist(x.orig, "norm", method = "mle")
summary(fit.norm)

b1 <- bootdist(fit.norm, niter=101)
print(b1)
plot(b1)
plot(b1, enhance=TRUE)
summary(b1)
quantile(b1)
CIcdfplot(b1, CI.output = "quantile")
```

Try fitting the Gamma model and compare.

```{r}
fit.gamma <- fitdist(x.orig, "gamma", method = "mle")
summary(fit.gamma)

b1 <- bootdist(fit.gamma, niter=101)
print(b1)
plot(b1)
plot(b1, enhance=TRUE)
summary(b1)
quantile(b1)
CIcdfplot(b1, CI.output = "quantile")
```