\[-\triangledown^2l(\mu,\sigma)=\begin{bmatrix} \frac{n}{\sigma^2} & \frac{2}{\sigma^3}\sum_{i = 1}^{n}(x_i-\mu) \\ \frac{2}{\sigma^3}\sum_{i = 1}^{n}(x_i-\mu) & -\frac{n}{\sigma^2}+\frac{3}{\sigma^4}\sum_{i = 1}^{n}(x_i-\mu)^2 \\ \end{bmatrix} \]
\[I(\mu,\sigma)=E\{-\triangledown^2l(\mu,\sigma)\}=\begin{bmatrix} \frac{n}{\sigma^2} & 0 \\ 0 & \frac{2n}{\sigma^2} \\ \end{bmatrix}\]
\[I(\mu,\sigma^2)=\begin{bmatrix} \frac{n}{\sigma^2} & 0 \\ 0 & \frac{n}{2\sigma^4} \\ \end{bmatrix}\]
Since \[I^{-1}(\hat{\mu},\hat{\sigma})=\begin{bmatrix} \frac{\Sigma(x_i-\overline{x})^2}{n^2} & 0 \\ 0 & \frac{\Sigma(x_i-\overline{x})^2}{2n^2} \\ \end{bmatrix}\] and \[I^{-1}(\hat{\mu},\hat{\sigma}^2)=\begin{bmatrix} \frac{\Sigma(x_i-\overline{x})^2}{n^2} & 0 \\ 0 & \frac{2\{\Sigma(x_i-\overline{x})^2\}^2}{n^3} \\ \end{bmatrix}\] therefore, \[s.e(\hat{\mu})=\sqrt{\frac{\Sigma(x_i-\overline{x})^2}{n^2}}\] \[s.e(\hat{\sigma})=\sqrt{\frac{\Sigma(x_i-\overline{x})^2}{2n^2}}\] \[s.e(\hat{\sigma}^2)=\sqrt{\frac{2\{\Sigma(x_i-\overline{x})^2\}^2}{n^3}}\]
\[I^{-1}(\pi)=\begin{bmatrix} 0.0009375 & -0.0003125 & -0.000625 \\ -0.0003125 & 0.0009375 & -0.000625 \\ -0.000625 & -0.000625 & 0.00125 \\ \end{bmatrix} \]
General R function
set.seed(2022)
n_sim =
n =
p1 =
p2 =
T = matrix(NA, ncol=3, nrow=n_sim)
for (i in 1:n_sim){
sam = sample(c(1,2,3), size=n, replace= TRUE,
prob=list(p1,p2,1-p1-p2))
t = table(sam)
T[i,] = t
}
# Simulation(n_sim=1000)
set.seed(2022)
n_sim=1000
n=200
p1=p2=1/4
T=matrix(NA, ncol=3, nrow=n_sim)
for (i in 1:n_sim){
sam = sample(c(1,2,3), size=n, replace= TRUE, prob=list(p1,p2,1-p1-p2))
t = table(sam)
T[i,] = t
}
# MLE matrix
M <- matrix(T/n, nrow=1000, ncol=3)
head(M, 5)
## [,1] [,2] [,3]
## [1,] 0.300 0.275 0.425
## [2,] 0.225 0.255 0.520
## [3,] 0.235 0.315 0.450
## [4,] 0.225 0.230 0.545
## [5,] 0.205 0.270 0.525
# Covariance
cov(M)
## [,1] [,2] [,3]
## [1,] 0.0009994282 -0.0003154114 -0.0006840168
## [2,] -0.0003154114 0.0009641209 -0.0006487095
## [3,] -0.0006840168 -0.0006487095 0.0013327263
M_1 <- cov(M)
apprx_v <- c(M_1[1,1],M_1[2,2],M_1[3,3],M_1[1,2],M_1[1,3],M_1[2,3])
true_v <- c(0.0009375, 0.0009375, 0.00125, -0.0003125, -0.000625, -0.000625)
abs_v <- abs(true_v-apprx_v)
# Table
table <- cbind(true_v, apprx_v, abs_v)
rownames(table) <- c("p1&p1", "p2&p2", "p3&p3", "p1&p2", "p1&p3", "p2&p3")
table
## true_v apprx_v abs_v
## p1&p1 0.0009375 0.0009994282 6.192820e-05
## p2&p2 0.0009375 0.0009641209 2.662090e-05
## p3&p3 0.0012500 0.0013327263 8.272633e-05
## p1&p2 -0.0003125 -0.0003154114 2.911386e-06
## p1&p3 -0.0006250 -0.0006840168 5.901682e-05
## p2&p3 -0.0006250 -0.0006487095 2.370951e-05
# Simulation(n_sim=10,000)
set.seed(2022)
n_sim=10000
n=200
p1=p2=1/4
T=matrix(NA, ncol=3, nrow=n_sim)
for (i in 1:n_sim){
sam = sample(c(1,2,3), size=n, replace= TRUE, prob=list(p1,p2,1-p1-p2))
t = table(sam)
T[i,] = t
}
# MLE matrix
M <- matrix(T/n, nrow=10000, ncol=3)
head(M, 5)
## [,1] [,2] [,3]
## [1,] 0.300 0.275 0.425
## [2,] 0.225 0.255 0.520
## [3,] 0.235 0.315 0.450
## [4,] 0.225 0.230 0.545
## [5,] 0.205 0.270 0.525
# Covariance
cov(M)
## [,1] [,2] [,3]
## [1,] 0.0009365716 -0.0003031699 -0.0006334017
## [2,] -0.0003031699 0.0009277465 -0.0006245766
## [3,] -0.0006334017 -0.0006245766 0.0012579783
M_1 <- cov(M)
apprx_v <- c(M_1[1,1],M_1[2,2],M_1[3,3],M_1[1,2],M_1[1,3],M_1[2,3])
true_v <- c(0.0009375, 0.0009375, 0.00125, -0.0003125, -0.000625, -0.000625)
abs_v <- abs(true_v-apprx_v)
# Table
table <- cbind(true_v, apprx_v, abs_v)
rownames(table) <- c("p1&p1", "p2&p2", "p3&p3", "p1&p2", "p1&p3", "p2&p3")
table
## true_v apprx_v abs_v
## p1&p1 0.0009375 0.0009365716 9.283991e-07
## p2&p2 0.0009375 0.0009277465 9.753532e-06
## p3&p3 0.0012500 0.0012579783 7.978298e-06
## p1&p2 -0.0003125 -0.0003031699 9.330114e-06
## p1&p3 -0.0006250 -0.0006334017 8.401715e-06
## p2&p3 -0.0006250 -0.0006245766 4.234173e-07
=> As we increase the sample size, the approximate values of the covariance of MLE seem to get closer to their true values.
# Histogram
par(mfrow=c(1,3))
H1 <- hist(M[,1], col="blue", main="MLE of p1", xlab="p1")
H2 <- hist(M[,2], col="red", main="MLE of p2", xlab="p2")
H3 <- hist(M[,3], col="green", main="MLE of p3", xlab="p3")
=> As can be seen from the graphs, the shape of the distribution of
each column of M is similar to the shape of normal distribution. Based
on what we have learned in class about the asymptotic distribution of
MLE, we can confirm that MLE is asymptotically normally
distributed. Each graph shows the marginal function of each MLE.
For p1, p2, its mean is 0.25 and the s.d is around 0.1, and for p3, the
mean is 0.5 and the s.d is around 0.125. These are close to their true
values. Specifically, the s.d values are close to the values from the
diagonals of the inverse of the Fisher-Information matrix.