Problema 2: Propiedades de los estimadores

CONSIDERACIONES INICIALES:

En inferencia estadística, para estimar valores de una población (parámetro poblacional), se hace uso de una muestra cuyos estimadores cumplan ciertas propiedades esenciales que lo hagan representativo de la población (estadístico muestral). Estas propiedades son cruciales porque determinan la confiabilidad de los estimadores para hacer inferencias adecuadas. Las propiedades que deben cumplir los estimadores son las siguientes:

1. Insesgadez: un estimador es insesgado si el valor esperado de la muestra es igual que su valor real poblacional.

2. Eficiencia: un estimador es eficiente si su error estándar coincide con la varianza poblacional.

3. Consistencia: un estimador es consistente si a medida que se aumenta el tamaño de la muestra, el estimador converge a los valores reales de la población.


Desarrollo del problema:

Para empezar, se genera una función que contenga la fórmula de los 4 estimadores a los cuales se les va a comprobar si cumplen con las propiedades de los estimadores, esto es, si son insesgados, eficientes y consistentes. Esta función se ejecutará en diferentes escenarios de tamaño de muestra (20, 50, 100 y 1000), y se harán las conclusiones pertinentes. 

parametro_teta <- 4
tamano_muestra <- c(20, 50, 100, 1000)



teta <- function(x){
  x1 <- sample(length(x),1)
  x2 <- sample(length(x),1)
  x3 <- sample(length(x),1)
  x4 <- sample(length(x),1)
  resultado1 <- (x[x1] + x[x2])/6 + (x[x3] + x[x4])/3
  resultado2 <- (x[x1] + 2*x[x2] + 3*x[x3] + 4*x[x4])/5
  resultado3 <- (x[x1] + x[x2] + x[x3] + x[x4])/4
  resultado4 <- (min(x) + max(x))/2
  media <- mean(x)
  
  resultados <- c(resultado1, resultado2, resultado3, resultado4, media)
  
  return(resultados)
}


Una vez creada la función, se efectuarán 1000 repeticiones a la función pra 4 escenarios; una para una muestra de 20, 50, 100 y 1000 registros. 

resultados <- data.frame(Estimador1 = numeric(1000),
                         Estimador2 = numeric(1000),
                         Estimador3 = numeric(1000),
                         Estimador4 = numeric(1000))


resultados2 <- data.frame(Estimador1 = numeric(1000),
                          Estimador2 = numeric(1000),
                          Estimador3 = numeric(1000),
                          Estimador4 = numeric(1000))

resultados3 <- data.frame(Estimador1 = numeric(1000),
                          Estimador2 = numeric(1000),
                          Estimador3 = numeric(1000),
                          Estimador4 = numeric(1000))

resultados4 <- data.frame(Estimador1 = numeric(1000),
                          Estimador2 = numeric(1000),
                          Estimador3 = numeric(1000),
                          Estimador4 = numeric(1000))



# Escenario 1: muestra de 20 registros

set.seed(284)
for (i in 1:1000) {
  muestra <- rexp(n = 20, 
                  rate = 1/parametro_teta)
  
  resumen_resultados <- teta(muestra)
  
  resultados$Estimador1[i] <- resumen_resultados[1]
  resultados$Estimador2[i] <- resumen_resultados[2]
  resultados$Estimador3[i] <- resumen_resultados[3]
  resultados$Estimador4[i] <- resumen_resultados[4]
  
  resultados$MediaMuestra[i] <- mean(muestra)
}


# Escenario 2: muestra de 50 registros

set.seed(384)
for (i in 1:1000) {
  muestra <- rexp(n = 50, 
                  rate = 1/parametro_teta)
  
  resumen_resultados <- teta(muestra)
  
  resultados2$Estimador1[i] <- resumen_resultados[1]
  resultados2$Estimador2[i] <- resumen_resultados[2]
  resultados2$Estimador3[i] <- resumen_resultados[3]
  resultados2$Estimador4[i] <- resumen_resultados[4]
  
  resultados2$MediaMuestra[i] <- mean(muestra)
}


# Escenario 3: muestra de 1000 registros

set.seed(484)
for (i in 1:1000) {
  muestra <- rexp(n = 100, 
                  rate = 1/parametro_teta)
  
  resumen_resultados <- teta(muestra)
  
  resultados3$Estimador1[i] <- resumen_resultados[1]
  resultados3$Estimador2[i] <- resumen_resultados[2]
  resultados3$Estimador3[i] <- resumen_resultados[3]
  resultados3$Estimador4[i] <- resumen_resultados[4]
  
  resultados3$MediaMuestra[i] <- mean(muestra)
}


# Escenario 4: muestra de 1000 registros

set.seed(584)
for (i in 1:1000) {
  muestra <- rexp(n = 1000, 
                  rate = 1/parametro_teta)
  
  resumen_resultados <- teta(muestra)
  
  resultados4$Estimador1[i] <- resumen_resultados[1]
  resultados4$Estimador2[i] <- resumen_resultados[2]
  resultados4$Estimador3[i] <- resumen_resultados[3]
  resultados4$Estimador4[i] <- resumen_resultados[4]
  
  resultados4$MediaMuestra[i] <- mean(muestra)
  
}

head(resultados)
##   Estimador1 Estimador2 Estimador3 Estimador4 MediaMuestra
## 1   3.436907   7.027203   3.813298   4.433281     3.681022
## 2   7.592894  15.930232   6.404311   5.821576     3.278980
## 3   3.705870   8.642350   2.936738   4.970180     3.356678
## 4   6.816214  13.622531   6.648395   6.455837     4.037034
## 5   2.777749   5.162941   3.023098   5.162125     2.835756
## 6   3.045326   5.981364   2.857658   6.099695     3.464277
head(resultados2)
##   Estimador1 Estimador2 Estimador3 Estimador4 MediaMuestra
## 1   3.844965   8.021153   3.676661  10.685868     4.717861
## 2   1.735686   3.520724   1.704841  10.102810     4.366889
## 3   3.047629   6.488501   2.946296   7.860598     4.253561
## 4   2.022803   4.130001   2.037948   5.428794     3.174995
## 5   5.761109  12.088454   4.607667  10.035741     4.344423
## 6   7.371722  17.433723   5.725053  10.075408     4.225709
head(resultados3)
##   Estimador1 Estimador2 Estimador3 Estimador4 MediaMuestra
## 1   5.941125  11.895954   6.319184   5.542772     3.197732
## 2   4.108332   8.409466   4.555854  11.249606     3.706859
## 3   4.499828  10.197512   3.895539  13.273481     3.953416
## 4   2.628792   5.372716   2.781306  12.715610     4.134261
## 5   2.038944   3.992368   2.034275   9.778718     3.697770
## 6   4.933288  10.816298   5.217409   9.116356     3.366105
head(resultados4)
##   Estimador1 Estimador2 Estimador3 Estimador4 MediaMuestra
## 1   4.898151   9.185406   4.865148   14.23200     4.230598
## 2   3.872990   7.254053   3.309225   27.14864     4.058981
## 3   9.537734  18.066335   7.396090   14.37802     3.970126
## 4   3.616557   7.088431   3.015725   17.94011     4.141459
## 5   3.238137   6.857952   3.714104   13.82166     4.086248
## 6   8.251613  15.328268   6.353147   11.95082     3.867627


A continuación, se cálcula la media y la dispersión de cada uno de los estiamadores para concluir acerca de sus propiedades. Pra esto se hace uso de la media, la varianza y graficamente mediante gráficos de caja y bigote. 

#Estimador 1
media1 <- data.frame(
  media_estimador1 = mean(resultados$Estimador1),
  media_estimador2 = mean(resultados$Estimador2),
  media_estimador3 = mean(resultados$Estimador3), 
  media_estimador4 = mean(resultados$Estimador4),
  media_muestra = mean(resultados$MediaMuestra))


#Estimador 2
media2 <- data.frame(
  media_estimador1 = mean(resultados2$Estimador1),
  media_estimador2 = mean(resultados2$Estimador2),
  media_estimador3 = mean(resultados2$Estimador3), 
  media_estimador4 = mean(resultados2$Estimador4),
  media_muestra = mean(resultados2$MediaMuestra))


#Estimador 3
media3 <- data.frame(
  media_estimador1 = mean(resultados3$Estimador1),
  media_estimador2 = mean(resultados3$Estimador2),
  media_estimador3 = mean(resultados3$Estimador3), 
  media_estimador4 = mean(resultados3$Estimador4),
  media_muestra = mean(resultados3$MediaMuestra))


#Estimador 4
media4 <- data.frame(
  media_estimador1 = mean(resultados4$Estimador1),
  media_estimador2 = mean(resultados4$Estimador2),
  media_estimador3 = mean(resultados4$Estimador3), 
  media_estimador4 = mean(resultados4$Estimador4),
  media_muestra = mean(resultados4$MediaMuestra))

media <- rbind(media1, media2, media3, media4)

media
##   media_estimador1 media_estimador2 media_estimador3 media_estimador4
## 1         3.966451         7.908956         4.001685         7.215284
## 2         4.029348         8.040542         4.029078         9.084109
## 3         3.814307         7.584535         3.833864        10.406901
## 4         4.010412         7.990124         4.010694        14.942150
##   media_muestra
## 1      3.956583
## 2      3.996410
## 3      3.993683
## 4      4.003377


error_estandar <- function(x) {
  return(sd(x) / sqrt(length(x)))
}

#Estimador 1
error_estandar_1 <- data.frame(
  error_estandar_estimador1 = error_estandar(resultados$Estimador1),
  error_estandar_estimador2 = error_estandar(resultados$Estimador2),
  error_estandar_estimador3 = error_estandar(resultados$Estimador3),
  error_estandar_estimador4 = error_estandar(resultados$Estimador4),
  varianza_muestra = var(resultados$MediaMuestra))

#Estimador 2
error_estandar_2 <- data.frame(
  error_estandar_estimador1 = error_estandar(resultados2$Estimador1),
  error_estandar_estimador2 = error_estandar(resultados2$Estimador2),
  error_estandar_estimador3 = error_estandar(resultados2$Estimador3),
  error_estandar_estimador4 = error_estandar(resultados2$Estimador4),
  varianza_muestra = var(resultados2$MediaMuestra))

#Estimador 3
error_estandar_3 <- data.frame(
  error_estandar_estimador1 = error_estandar(resultados3$Estimador1),
  error_estandar_estimador2 = error_estandar(resultados3$Estimador2),
  error_estandar_estimador3 = error_estandar(resultados3$Estimador3),
  error_estandar_estimador4 = error_estandar(resultados3$Estimador4),
  varianza_muestra = var(resultados3$MediaMuestra))

#Estimador 4
error_estandar_4 <- data.frame(
  error_estandar_estimador1 = error_estandar(resultados4$Estimador1),
  error_estandar_estimador2 = error_estandar(resultados4$Estimador2),
  error_estandar_estimador3 = error_estandar(resultados4$Estimador3),
  error_estandar_estimador4 = error_estandar(resultados4$Estimador4),
  varianza_muestra = var(resultados4$MediaMuestra))


error_estandar <- rbind(error_estandar_1, error_estandar_2, error_estandar_3, error_estandar_4)

error_estandar
##   error_estandar_estimador1 error_estandar_estimador2 error_estandar_estimador3
## 1                0.07208275                 0.1478017                0.06989542
## 2                0.06751702                 0.1387353                0.06480907
## 3                0.06739853                 0.1372307                0.06318809
## 4                0.06500511                 0.1342204                0.06242749
##   error_estandar_estimador4 varianza_muestra
## 1                0.07856423       0.77812771
## 2                0.08640728       0.30440884
## 3                0.07941324       0.15003561
## 4                0.08459674       0.01630068
y_min <- min(resultados$Estimador1, resultados$Estimador2, resultados$Estimador3, resultados$Estimador4)
y_max <- max(resultados$Estimador1, resultados$Estimador2, resultados$Estimador3, resultados$Estimador4)

y_min2 <- min(resultados2$Estimador1, resultados2$Estimador2, resultados2$Estimador3, resultados2$Estimador4)
y_max2 <- max(resultados2$Estimador1, resultados2$Estimador2, resultados2$Estimador3, resultados2$Estimador4)

y_min3 <- min(resultados3$Estimador1, resultados3$Estimador2, resultados3$Estimador3, resultados3$Estimador4)
y_max3 <- max(resultados3$Estimador1, resultados3$Estimador2, resultados3$Estimador3, resultados3$Estimador4)

y_min4 <- min(resultados4$Estimador1, resultados4$Estimador2, resultados4$Estimador3, resultados4$Estimador4)
y_max4 <- max(resultados4$Estimador1, resultados4$Estimador2, resultados4$Estimador3, resultados4$Estimador4)


par(mfrow = c(1, 5)) 

boxplot(resultados$Estimador1, ylim = c(y_min, y_max))
title("Estimador 1")
boxplot(resultados$Estimador2, ylim = c(y_min, y_max))
title("Estimador 2")
boxplot(resultados$Estimador3, ylim = c(y_min, y_max)) 
title("Estimador 3")
boxplot(resultados$Estimador4, ylim = c(y_min, y_max))
title("Estimador 4")
boxplot(resultados$MediaMuestra, ylim = c(y_min, y_max))
title("Media muestral")
mtext("BOXPLOT PARA N = 20", outer = TRUE, cex = 0.8, line = -1)

par(mfrow = c(1, 5))

boxplot(resultados2$Estimador1, ylim = c(y_min, y_max))
title("Estimador 1")
boxplot(resultados2$Estimador2, ylim = c(y_min, y_max))
title("Estimador 2")
boxplot(resultados2$Estimador3, ylim = c(y_min, y_max)) 
title("Estimador 3")
boxplot(resultados2$Estimador4, ylim = c(y_min, y_max))
title("Estimador 4")
boxplot(resultados2$MediaMuestra, ylim = c(y_min, y_max))
title("Media muestral")
mtext("BOXPLOT PARA N = 50", outer = TRUE, cex = 0.8, line = -1)

par(mfrow = c(1, 5))

boxplot(resultados3$Estimador1, ylim = c(y_min, y_max))
title("Estimador 1")
boxplot(resultados3$Estimador2, ylim = c(y_min, y_max))
title("Estimador 2")
boxplot(resultados3$Estimador3, ylim = c(y_min, y_max)) 
title("Estimador 3")
boxplot(resultados3$Estimador4, ylim = c(y_min, y_max))
title("Estimador 4")
boxplot(resultados3$MediaMuestra, ylim = c(y_min, y_max))
title("Media muestral")
mtext("BOXPLOT PARA N = 100", outer = TRUE, cex = 0.8, line = -1)

par(mfrow = c(1, 5))

boxplot(resultados4$Estimador1, ylim = c(y_min, y_max))
title("Estimador 1")
boxplot(resultados4$Estimador2, ylim = c(y_min, y_max))
title("Estimador 2")
boxplot(resultados4$Estimador3, ylim = c(y_min, y_max)) 
title("Estimador 3")
boxplot(resultados4$Estimador4, ylim = c(y_min, y_max))
title("Estimador 4")
boxplot(resultados4$MediaMuestra, ylim = c(y_min, y_max))
title("Media muestral")
mtext("BOXPLOT PARA N = 1000", outer = TRUE, cex = 0.8, line = -1)

CONCLUSIONES FINALES: