Estimación de Intervalos de Confianza mediante Bootstrap

Una empresa dedicada a la generación de energía renovable desea evaluar el rendimiento de sus 20 generadores eléctricos portátiles, alimentados con biocombustible. Estos generadores se utilizan en proyectos remotos, donde la eficiencia energética es crucial para minimizar costos y reducir la huella ambiental. Para el estudio, se mide la eficiencia energética de cada generador, expresada en kilovatios-hora por litro de combustible (kWh/L). Los resultados son:

Muestra de eficiencia energética (kWh/L)

muestra <- c(6.12, 5.87, 5.45, 6.33, 5.71, 6.04, 5.92, 5.65, 6.18, 5.78, 5.95, 6.21, 5.63, 5.79, 6.11, 5.88, 6.02, 5.76, 5.85, 6.10)

El objetivo es estimar un intervalo de confianza al 95% para la eficiencia media.

Para lo cual se realizan las siguientes actividades:

Punto a: Normalidad:

# Cargar librerías necesarias
library(boot)
library(ggplot2)

# Cargar librerías necesarias
library(boot)
library(ggplot2)

set.seed(123)  # Para asegurar reproducibilidad

# Datos originales
X <- c(6.12, 5.87, 5.45, 6.33, 5.71, 6.04, 5.92, 5.65, 6.18,5.78, 5.95, 6.21, 5.63, 5.79, 6.11, 5.88, 6.02, 5.76, 5.85, 6.10)

# Función para calcular la estadística de interés (varianza)
boot_function <- function(data, indices) {
  return(var(data[indices]))
}

# Número de muestras bootstrap
B <- 1000  

# Generar muestras bootstrap usando la librería boot
boot_results <- boot(data = X, statistic = boot_function, R = B)

# Obtener la distribución de varianzas bootstrap
bootstrap_var <- boot_results$t

# Calcular la media y desviación estándar de las varianzas bootstrap
mean_bootstrap <- mean(bootstrap_var)
sd_bootstrap <- sd(bootstrap_var)

Intervalos de Confianza

alpha <- 0.05

# Intervalo de Confianza Percentil
IC.p <- boot.ci(boot_results, type = "perc")$percent[4:5]

# Intervalo de Confianza Normal
IC.N <- boot.ci(boot_results, type = "norm")$normal[2:3]

# Intervalo de Confianza BCa
IC.BCa <- boot.ci(boot_results, type = "bca")$bca[4:5]

# Crear un dataframe para ggplot
df <- data.frame(bootstrap_var)

# Crear datos para la curva normal ajustada
x_vals <- seq(min(bootstrap_var), max(bootstrap_var), length.out = 100)
normal_curve <- data.frame(
  x = x_vals,
  y = dnorm(x_vals, mean = mean_bootstrap, sd = sd_bootstrap) # Corrección: Basado en varianza
)

Histograma Densidad de la Muestra

plot.boot.c<-ggplot(df, aes(x = bootstrap_var)) +
  geom_histogram(aes(y = after_stat(density)), bins = 30, fill = "lightblue", color = "black", alpha = 0.6) + # Corrección en el eje x
  geom_line(data = normal_curve, aes(x = x, y = y), color = "darkorange", linewidth = 1.5) + # Corrección: Curva normal basada en varianza
  geom_vline(xintercept = IC.p, color = "red", linetype = "dashed", linewidth = 1.2) + # IC Percentil
  geom_vline(xintercept = IC.N, color = "blue", linetype = "dotted", linewidth = 1.2) + # IC Normal
  geom_vline(xintercept = IC.BCa, color = "black", linetype = "dotdash", linewidth = 1.2) + # IC BCa (cambiado a negro)
  labs(title = "Histograma de la Varianza Bootstrap con Ajuste Normal",
       x = "Varianza Bootstrap", y = "Densidad") +
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5, size = 14)) +
  annotate("text", x = mean(IC.p), y = 0.9, label = "IC Percentil", color = "red", angle = 90, vjust = -1) +
  annotate("text", x = mean(IC.N), y = 0.9, label = "IC Normal", color = "blue", angle = 90, vjust = -0.5) +
  annotate("text", x = mean(IC.BCa), y = 0.9, label = "IC BCa", color = "black", angle = 90, vjust = -1) # Cambiado a negro

print(plot.boot.c)

IC.p 
## [1] 0.02295695 0.07561393
IC.N
## [1] 0.02657162 0.07857666
IC.BCa
## [1] 0.03125763 0.09096194

A partir de la simulación bootstrap, se obtuvieron los siguientes intervalos de confianza al 95% para la varianza de la población:

IC PercentilICBCa=[0.02295695 , 0.07561393]= NormalIC=[0.02657162, 0.07857666] IC BCa= [0.03125763 , 0.09096194]

Aplica las pruebas Shapiro-Wilk y Kolmogorov-Smirnov para evaluar la normalidad

# Definir la muestra
muestra <- c(6.12, 5.87, 5.45, 6.33, 5.71, 6.04, 5.92, 5.65, 6.18, 
             5.78, 5.95, 6.21, 5.63, 5.79, 6.11, 5.88, 6.02, 5.76, 5.85, 6.10)

# Prueba de normalidad de Shapiro-Wilk
shapiro_test <- shapiro.test(muestra)
print(shapiro_test)
## 
##  Shapiro-Wilk normality test
## 
## data:  muestra
## W = 0.98845, p-value = 0.9955
# Prueba de normalidad de Kolmogorov-Smirnov contra una distribución normal
ks_test <- ks.test(muestra, "pnorm", mean = mean(muestra), sd = sd(muestra))
print(ks_test)
## 
##  Exact one-sample Kolmogorov-Smirnov test
## 
## data:  muestra
## D = 0.093219, p-value = 0.9882
## alternative hypothesis: two-sided

Explicación: Se determino con la aplicación del Shapiro-Wilk normality test la normalidad de la muestra y se determinan los siguientes resultados:

W = 0.98845, p-value = 0.9955 D = 0.093219, p-value = 0.9882

Por tanto como el resultado de p-value resultan mayores que 0.05, se aprueba la hipótesis de normalidad ya que los datos siguen una distribución normal.

Punto b: Procedimiento Bootstrap:

Calcula la media de la muestra original

# Definir la muestra
muestra <-c(6.12, 5.87, 5.45, 6.33, 5.71, 6.04, 5.92, 5.65, 6.18, 5.78, 5.95, 6.21, 5.63, 5.79, 6.11, 5.88, 6.02, 5.76, 5.85, 6.10)

# Calcular la media
media_muestra <- mean(muestra)

# Mostrar el resultado
print(media_muestra)
## [1] 5.9175

Resultado de la media: 5.9175

Genera 1000 muestras bootstrap de tamaño n=20 mediante muestreo con reemplazo.

# Definir la muestra original
muestra <- c(6.12, 5.87, 5.45, 6.33, 5.71, 6.04, 5.92, 5.65, 6.18, 5.78, 
             5.95, 6.21, 5.63, 5.79, 6.11, 5.88, 6.02, 5.76, 5.85, 6.10)

# Número de muestras bootstrap
num_bootstrap <- 1000
n <- 20  # Tamaño de cada muestra bootstrap

# Generar las muestras bootstrap
bootstrap_samples <- replicate(num_bootstrap, sample(muestra, size = n, replace = TRUE))

# Mostrar las primeras 5 muestras bootstrap
bootstrap_samples[, 1:5]
##       [,1] [,2] [,3] [,4] [,5]
##  [1,] 5.95 5.85 5.65 5.65 6.18
##  [2,] 5.78 6.04 5.85 5.88 5.71
##  [3,] 6.10 5.85 6.21 6.02 5.85
##  [4,] 6.21 5.85 5.92 5.87 5.63
##  [5,] 5.79 5.88 5.71 5.87 6.33
##  [6,] 5.65 5.45 5.76 5.95 5.71
##  [7,] 5.92 6.02 6.21 5.95 5.63
##  [8,] 5.85 5.95 6.04 6.04 5.65
##  [9,] 5.78 6.18 6.12 5.76 5.92
## [10,] 6.10 6.02 5.87 6.33 5.92
## [11,] 5.65 5.71 6.04 5.71 6.11
## [12,] 5.76 5.45 5.65 5.78 5.85
## [13,] 5.95 5.45 5.65 5.71 5.63
## [14,] 5.88 6.10 5.78 5.87 6.21
## [15,] 5.63 5.88 6.04 6.12 5.76
## [16,] 5.78 6.10 5.65 6.04 5.65
## [17,] 6.11 6.02 5.88 6.12 5.65
## [18,] 5.63 5.45 6.33 5.85 6.33
## [19,] 5.76 6.12 5.88 6.12 5.78
## [20,] 5.95 6.21 6.33 6.21 5.87

Calcula la media muestral para cada muestra bootstrap

# Calcular la media muestral para cada muestra bootstrap
medias_bootstrap <- apply(bootstrap_samples, 2, mean)

# Mostrar las primeras 10 medias bootstrap
print(medias_bootstrap[1:10])
##  [1] 5.8615 5.8790 5.9285 5.9425 5.8685 5.9140 5.8860 5.9770 5.9085 5.8935

Obtén los percentiles 2.5 y 97.5 de las medias muestrales bootstrap.

# Calcular los percentiles 2.5 y 97.5 de las medias bootstrap
percentiles <- quantile(medias_bootstrap, probs = c(0.025, 0.975))

# Mostrar los resultados
print(percentiles)
##     2.5%    97.5% 
## 5.830475 6.013525

Calcula 2 tipos de intervalos de confianza para la media mediante Bootstrap

# IC Percentil (Usando los percentiles 2.5 y 97.5 de las medias bootstrap)
ic_percentil <- quantile(medias_bootstrap, probs = c(0.025, 0.975))

# IC Normal (Usando la media y desviación estándar de las medias bootstrap)
media_boot <- mean(medias_bootstrap)
sd_boot <- sd(medias_bootstrap)
ic_normal <- c(media_boot - 1.96 * sd_boot, media_boot + 1.96 * sd_boot)  # IC del 95%

# Mostrar los resultados
print(paste("IC Percentil (95%):", ic_percentil[1], "a", ic_percentil[2]))
## [1] "IC Percentil (95%): 5.830475 a 6.013525"
print(paste("IC Normal (95%):", ic_normal[1], "a", ic_normal[2]))
## [1] "IC Normal (95%): 5.82515685045187 a 6.01407314954813"

Expliacion

IC Percentil: Se calculan los percentiles 2.5% y 97.5% de las medias bootstrap, asumiendo que la distribución muestral es adecuada.

Compara ambos intervalos: ¿Son similares o presentan diferencias significativas?

# Calcular las longitudes de los intervalos
longitud_percentil <- ic_percentil[2] - ic_percentil[1]
longitud_normal <- ic_normal[2] - ic_normal[1]

# Comparar los intervalos
print(paste("IC Percentil (95%):", ic_percentil[1], "a", ic_percentil[2], "Longitud:", longitud_percentil))
## [1] "IC Percentil (95%): 5.830475 a 6.013525 Longitud: 0.18305"
print(paste("IC Normal (95%):", ic_normal[1], "a", ic_normal[2], "Longitud:", longitud_normal))
## [1] "IC Normal (95%): 5.82515685045187 a 6.01407314954813 Longitud: 0.188916299096251"
# Visualización con un gráfico
library(ggplot2)
df <- data.frame(
  Metodo = c("Percentil", "Normal"),
  LI = c(ic_percentil[1], ic_normal[1]),  # Límite Inferior
  LS = c(ic_percentil[2], ic_normal[2])   # Límite Superior
)

ggplot(df, aes(x = Metodo, ymin = LI, ymax = LS)) +
  geom_errorbar(width = 0.2, color = "blue") +
  geom_point(aes(y = (LI + LS) / 2), color = "red", size = 3) +
  labs(title = "Comparacion de Intervalos de Confianza (Bootstrap)", y = "Valor", x = "Método") +
  theme_minimal()

EXPLIACIÓN

Ls intervalos son similares (es decir, tienen longitudes parecidas y valores cercanos), esto sugiere que la distribución de las medias bootstrap es aproximadamente normal y que ambos métodos son confiables.

Punto c:Método Paramétrico:

Si se cumplen las condiciones para su aplicación, calcula un intervalo de confianza paramétrico. Compara este intervalo con los obtenidos mediante métodos no paramétricos.

# Cargar librería necesaria
library(stats)

# Calcular media y desviación estándar de la muestra original
media_muestra <- mean(muestra)
sd_muestra <- sd(muestra)

# Tamaño de la muestra
n <- length(muestra)

# Nivel de confianza del 95%
alpha <- 0.05
t_critico <- qt(1 - alpha/2, df = n - 1)  # Valor crítico t para n-1 grados de libertad

# Calcular el margen de error
margen_error <- t_critico * (sd_muestra / sqrt(n))

# Calcular los límites del intervalo de confianza
ic_parametrico <- c(media_muestra - margen_error, media_muestra + margen_error)

# Mostrar el resultado
print(paste("IC Parametrico(95%):", ic_parametrico[1], "a", ic_parametrico[2]))
## [1] "IC Parametrico(95%): 5.81303761518061 a 6.02196238481939"

EXPLICACIÓN

Comparación con Bootstrap:

La muestra sigue una distribución normal, el IC paramétrico y el IC Bootstrap Normal son similares.

library(yaml)

config <- yaml::yaml.load("
param1: valor1
param2: 100
param3: texto con espacios
")
print(config)
## $param1
## [1] "valor1"
## 
## $param2
## [1] 100
## 
## $param3
## [1] "texto con espacios"