How should I compute the Monte Carlo mean square error in R?

I try to get the Monte Carlo mean square error of the maximal likelihood estimators in R. I can write the calculation for the MLE that is repeated once, but I need to repeat the Monte Carlo calculation many times. How should I write this in R?

Firstly, we have sample size n=20 and the repetition N=30. The idea is that for the first repetition, we sample n=20 data from the Gumbel distribution with location 0 and scale 1 (true value mu=0 and sigma=1).

library(VGAM)
set.seed(101)
x <- rgumbel(20, loc = 0, scale = 1)

Then we compute the MLE of two parameters mu and sigma by the likelihood function

# Log-likelihood function of Gumbel distribution
ll_fn <- function(par, z, m){
  mu <- par[1]
  sigma <- par[2]
  -m * log(sigma) - sum((z - mu)/sigma) - sum(exp(-((z - mu)/sigma)))
}

m <- length(x)
fit<- optim(par = c(1, 1), 
              fn = ll_fn, z = y, m = m,
              hessian = TRUE, control = list(fnscale = -1))
fit3$par
#[1] 0.1718418 0.9813649

The Monte Carlo mean square error is that for the MLE of each i-th repetition

enter image description here



Comments

Popular posts from this blog

Spring Elasticsearch Operations

Network Error and Timeout on Authorize.net JS

Object oriented programming concepts (OOPs)