La simulación ayuda a entender y validad las propiedades de los estimadores estadísticos como son. insesgadez, eficiencia y la consistencia principalmente. El siguiente problema permite evidenciar las principales características de un grupo de estimadores propuestos para la estimación de un parámetro asociado a un modelo de probabilidad.
Se denomina sesgo de un estimador a la diferencia entre la esperanza (valor esperado) del estimador y el verdadero valor del parámetro a estimar. Es deseable que un estimador sea insesgado o centrado, esto es, que el sesgo sea nulo para que la esperanza del estimador sea igual al valor del parámetro que se desea estimar.
Un estimador se considera eficiente, cuando dentro de todos los estimadores insesgados presenta la menor varianza. La eficiencia relativa se mide por el ratio: Var(estimador1)/Var(estimador2). La eficiencia de los estimadores está limitada por las características de la distribución de probabilidad de la muestra de la que proceden. Un estimador es eficiente cuando verifica: Es insesgado y posee varianza mínima.
La consistencia de un estimador se refiere a la propiedad de que el estimador tiende a converger al verdadero valor del parámetro a medida que el tamaño de la muestra se acerca a infinito. En otras palabras, un estimador es consistente si su error tiende a cero a medida que la muestra crece. La consistencia es una propiedad fundamental de los buenos estimadores.
lambda=2
num_simulations <- 10000
estimates <- numeric(num_simulations)
set.seed(123)
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(20, rate = lambda)
theta_hat <- (theta[1] + theta[2] / 6) + (theta[3] + theta[4] / 3)
estimates[i] <- theta_hat
}
bias <- mean(estimates) - lambda
variance <- var(estimates)
estimates1<-estimates
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(20, rate = lambda)
theta_hat1 <- theta[1] + 2*theta[2] + 3*theta[3] + 4*theta[4] / 5
estimates[i] <- theta_hat1
}
bias1 <- mean(estimates) - lambda
variance1 <- var(estimates)
estimates2<-estimates
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(20, rate = lambda)
theta_hat2 <- theta[1] + theta[2] + theta[3] + theta[4] / 4
estimates[i] <- theta_hat2
}
bias2 <- mean(estimates) - lambda
variance2 <- var(estimates)
estimates3<-estimates
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(20, rate = lambda)
theta_hat3 <- min(theta) + max(theta) / 2
estimates[i] <- theta_hat3
}
bias3 <- mean(estimates) - lambda
variance3 <- var(estimates)
estimates4<-estimates
#################################################################################################################
sesgo20<-c(bias,bias1,bias2,bias3)
sesgo20
## [1] -0.7553894 1.4041177 -0.3865773 -1.0752675
varianza20<-c(variance,variance1,variance2,variance3)
varianza20
## [1] 0.5265128 3.6975056 0.7483613 0.1048577
estimadores20<-data.frame(estimates1,estimates2,estimates3,estimates4)
boxplot(estimadores20, main='Comparacion de estimadores con n=20')
abline(h=1/2,col='red')
#################################################################################################################
lambda=2
num_simulations <- 10000
estimates50 <- numeric(num_simulations)
set.seed(123)
for (i in 1:num_simulations) {
theta<- rexp(50, rate = lambda)
theta_hat4 <- (theta[1] + theta[2] / 6) + (theta[3] + theta[4] / 3)
estimates50[i] <- theta_hat4
}
bias4 <- mean(estimates50) - lambda
variance4 <- var(estimates50)
estimates5<-estimates50
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(50, rate = lambda)
theta_hat5 <- theta[1] + 2*theta[2] + 3*theta[3] + 4*theta[4] / 5
estimates50[i] <- theta_hat5
}
bias5 <- mean(estimates50) - lambda
variance5 <- var(estimates50)
estimates6<-estimates50
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(50, rate = lambda)
theta_hat6 <- theta[1] + theta[2] + theta[3] + theta[4] / 4
estimates50[i] <- theta_hat6
}
bias6 <- mean(estimates50) - lambda
variance6 <- var(estimates50)
estimates7<-estimates50
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(50, rate = lambda)
theta_hat7 <- min(theta) + max(theta) / 2
estimates50[i] <- theta_hat7
}
bias7 <- mean(estimates50) - lambda
variance7 <- var(estimates50)
estimates8<-estimates50
#################################################################################################################
sesgo50<-c(bias4,bias5,bias6,bias7)
sesgo50
## [1] -0.7617518 1.4156660 -0.3762028 -0.8711554
varianza50<-c(variance4,variance5,variance6,variance7)
varianza50
## [1] 0.5328618 3.7467212 0.7740415 0.1000852
estimadores50<-data.frame(estimates5,estimates6,estimates7,estimates8)
boxplot(estimadores50, main='Comparacion de estimadores con n=50')
abline(h=1/2,col='red')
#################################################################################################################
lambda=2
num_simulations <- 10000
estimates100 <- numeric(num_simulations)
for (i in 1:num_simulations) {
theta<- rexp(100, rate = lambda)
theta_hat8 <- (theta[1] + theta[2] / 6) + (theta[3] + theta[4] / 3)
estimates100[i] <- theta_hat8
}
bias8 <- mean(estimates100) - lambda
variance8 <- var(estimates100)
estimates9<-estimates100
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(100, rate = lambda)
theta_hat9 <- theta[1] + 2*theta[2] + 3*theta[3] + 4*theta[4] / 5
estimates100[i] <- theta_hat9
}
bias9 <- mean(estimates100) - lambda
variance9 <- var(estimates100)
estimates10<-estimates100
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(100, rate = lambda)
theta_hat10 <- theta[1] + theta[2] + theta[3] + theta[4] / 4
estimates100[i] <- theta_hat10
}
bias10 <- mean(estimates100) - lambda
variance10 <- var(estimates100)
estimates11<-estimates100
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(100, rate = lambda)
theta_hat11 <- min(theta) + max(theta) / 2
estimates100[i] <- theta_hat11
}
bias11 <- mean(estimates100) - lambda
variance11 <- var(estimates100)
estimates12<-estimates100
#################################################################################################################
sesgo100<-c(bias8,bias9,bias10,bias11)
sesgo100
## [1] -0.7618762 1.4063307 -0.3673588 -0.7003478
varianza100<-c(variance8,variance9,variance10,variance11)
varianza100
## [1] 0.52940709 3.65560105 0.79177777 0.09923152
estimadores100<-data.frame(estimates9,estimates10,estimates11,estimates12)
boxplot(estimadores100, main='Comparacion de estimadores con n=100')
abline(h=1/2,col='red')
#################################################################################################################
lambda=2
num_simulations <- 10000
estimates1000 <- numeric(num_simulations)
for (i in 1:num_simulations) {
theta<- rexp(1000, rate = lambda)
theta_hat12 <- (theta[1] + theta[2] / 6) + (theta[3] + theta[4] / 3)
estimates1000[i] <- theta_hat12
}
bias12 <- mean(estimates1000) - lambda
variance12 <- var(estimates1000)
estimates13<-estimates1000
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(1000, rate = lambda)
theta_hat13 <- theta[1] + 2*theta[2] + 3*theta[3] + 4*theta[4] / 5
estimates1000[i] <- theta_hat13
}
bias13 <- mean(estimates1000) - lambda
variance13 <- var(estimates1000)
estimates14<-estimates1000
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(1000, rate = lambda)
theta_hat14 <- theta[1] + theta[2] + theta[3] + theta[4] / 4
estimates1000[i] <- theta_hat14
}
bias14 <- mean(estimates1000) - lambda
variance14 <- var(estimates1000)
estimates15<-estimates1000
#################################################################################################################
for (i in 1:num_simulations) {
theta<- rexp(1000, rate = lambda)
theta_hat15 <- min(theta) + max(theta) / 2
estimates1000[i] <- theta_hat15
}
bias15 <- mean(estimates1000) - lambda
variance15 <- var(estimates1000)
estimates16<-estimates1000
#################################################################################################################
sesgo1000<-c(bias12,bias13,bias14,bias15)
sesgo1000
## [1] -0.7488828 1.4024907 -0.3793535 -0.1305345
varianza1000<-c(variance12,variance13,variance14,variance15)
varianza1000
## [1] 0.5330126 3.6899271 0.7354250 0.1035813
estimadores1000<-data.frame(estimates13,estimates14,estimates15,estimates16)
boxplot(estimadores1000, main='Comparacion de estimadores con n=1000')
abline(h=1/2,col='red')
#################################################################################################################