library(readr); library(dplyr); library(gt); library(MASS)
cat("Librerías cargadas correctamente.\n")
## Librerías cargadas correctamente.
ruta_csv <- file.choose()
datos <- read_csv(ruta_csv, show_col_types = FALSE)
cat("Archivo:", basename(ruta_csv), "| Filas:", nrow(datos), "\n")
## Archivo: oil_and_gas_leases_data (2).csv | Filas: 47757
x_raw <- datos %>%
mutate(LON = suppressWarnings(as.numeric(LONGITUDE))) %>%
filter(!is.na(LON), LON >= -103.0, LON <= -94.0) %>%
pull(LON)
n_conteo <- length(x_raw)
k_sturges <- ceiling(1 + 3.322 * log10(n_conteo))
cat("Observaciones válidas:", n_conteo, "\n")
## Observaciones válidas: 47757
cat("Clases (Regla de Sturges):", k_sturges, "\n")
## Clases (Regla de Sturges): 17
cat("Mínimo:", round(min(x_raw), 4), "| Máximo:", round(max(x_raw), 4), "\n")
## Mínimo: -102.0439 | Máximo: -94.6179
Se calcula la distribución de frecuencias absolutas y relativas para la variable cuantitativa continua Longitud, correspondiente a los arrendamientos de hidrocarburos registrados en Kansas, EE.UU.
x_st <- x_raw
n_st <- length(x_st)
x_min_st <- min(x_st); x_max_st <- max(x_st)
k_st <- k_sturges
c_amp_st <- (x_max_st - x_min_st) / k_st
lim_inf_st <- x_min_st + (0:(k_st-1)) * c_amp_st
lim_sup_st <- lim_inf_st + c_amp_st; lim_sup_st[k_st] <- x_max_st
mc_st <- (lim_inf_st + lim_sup_st) / 2
breaks_st <- c(lim_inf_st, lim_sup_st[k_st])
intervalos_cut_st <- cut(x_st, breaks=breaks_st, right=FALSE, include.lowest=TRUE)
freq_abs_st <- as.integer(table(intervalos_cut_st))
hi_dec_st <- freq_abs_st / n_st
Ni_asc_st <- cumsum(freq_abs_st); Hi_asc_st <- cumsum(hi_dec_st)
Ni_desc_st <- n_st - c(0, head(Ni_asc_st,-1)); Hi_desc_st <- 1 - c(0, head(Hi_asc_st,-1))
etiq_st <- paste0("[",round(lim_inf_st,4)," - ",round(lim_sup_st,4),")")
etiq_st[k_st] <- paste0("[",round(lim_inf_st[k_st],4)," - ",round(lim_sup_st[k_st],4),"]")
bind_rows(
data.frame(Intervalo=etiq_st, MC=round(mc_st,4), ni=freq_abs_st,
hi_pct=round(hi_dec_st*100,2), hi_real=round(hi_dec_st,4),
Ni_a=Ni_asc_st, Hi_a=round(Hi_asc_st,4),
Ni_d=Ni_desc_st, Hi_d=round(Hi_desc_st,4), stringsAsFactors=FALSE),
data.frame(Intervalo="TOTAL", MC=NA_real_, ni=sum(freq_abs_st),
hi_pct=round(sum(hi_dec_st)*100,2), hi_real=round(sum(hi_dec_st),4),
Ni_a=max(Ni_asc_st), Hi_a=round(max(Hi_asc_st),4),
Ni_d=max(Ni_desc_st), Hi_d=round(max(Hi_desc_st),4), stringsAsFactors=FALSE)
) %>%
gt() %>%
tab_header(title=md("**Distribución de Frecuencias — Regla de Sturges (referencial)**"),
subtitle=md(paste0("*Longitud, Kansas (n=",format(n_st,big.mark=","),", k=",k_st," intervalos)*"))) %>%
cols_label(Intervalo=md("**Intervalo**"), MC=md("**MC**"), ni=md("**ni**"),
hi_pct=md("**hi%**"), hi_real=md("**hi**"),
Ni_a=md("**Ni (asc)**"), Hi_a=md("**Hi (asc)**"),
Ni_d=md("**Ni (desc)**"), Hi_d=md("**Hi (desc)**")) %>%
tab_style(style=list(cell_fill(color="#2C2C2C"),cell_text(color="white",weight="bold")),
locations=cells_column_labels()) %>%
tab_style(style=cell_fill(color="#F5F5F5"), locations=cells_body(rows=seq(1,k_st+1,by=2))) %>%
tab_style(style=list(cell_fill(color="#D6D6D6"),cell_text(weight="bold")),
locations=cells_body(rows=Intervalo=="TOTAL")) %>%
fmt_missing(columns=everything(), missing_text="-") %>%
tab_source_note(source_note=md("*Autor: Fernando Almeida*")) %>%
tab_options(table.width=pct(100), table.font.size=px(12), data_row.padding=px(4))
| Distribución de Frecuencias — Regla de Sturges (referencial) | ||||||||
| Longitud, Kansas (n=47,757, k=17 intervalos) | ||||||||
| Intervalo | MC | ni | hi% | hi | Ni (asc) | Hi (asc) | Ni (desc) | Hi (desc) |
|---|---|---|---|---|---|---|---|---|
| [-102.0439 - -101.6071) | -101.8255 | 2213 | 4.63 | 0.0463 | 2213 | 0.0463 | 47757 | 1.0000 |
| [-101.6071 - -101.1703) | -101.3887 | 2516 | 5.27 | 0.0527 | 4729 | 0.0990 | 45544 | 0.9537 |
| [-101.1703 - -100.7334) | -100.9519 | 5205 | 10.90 | 0.1090 | 9934 | 0.2080 | 43028 | 0.9010 |
| [-100.7334 - -100.2966) | -100.5150 | 2178 | 4.56 | 0.0456 | 12112 | 0.2536 | 37823 | 0.7920 |
| [-100.2966 - -99.8598) | -100.0782 | 2458 | 5.15 | 0.0515 | 14570 | 0.3051 | 35645 | 0.7464 |
| [-99.8598 - -99.423) | -99.6414 | 3977 | 8.33 | 0.0833 | 18547 | 0.3884 | 33187 | 0.6949 |
| [-99.423 - -98.9862) | -99.2046 | 4644 | 9.72 | 0.0972 | 23191 | 0.4856 | 29210 | 0.6116 |
| [-98.9862 - -98.5493) | -98.7678 | 6007 | 12.58 | 0.1258 | 29198 | 0.6114 | 24566 | 0.5144 |
| [-98.5493 - -98.1125) | -98.3309 | 3895 | 8.16 | 0.0816 | 33093 | 0.6929 | 18559 | 0.3886 |
| [-98.1125 - -97.6757) | -97.8941 | 1677 | 3.51 | 0.0351 | 34770 | 0.7281 | 14664 | 0.3071 |
| [-97.6757 - -97.2389) | -97.4573 | 1021 | 2.14 | 0.0214 | 35791 | 0.7494 | 12987 | 0.2719 |
| [-97.2389 - -96.8021) | -97.0205 | 1788 | 3.74 | 0.0374 | 37579 | 0.7869 | 11966 | 0.2506 |
| [-96.8021 - -96.3652) | -96.5836 | 693 | 1.45 | 0.0145 | 38272 | 0.8014 | 10178 | 0.2131 |
| [-96.3652 - -95.9284) | -96.1468 | 1337 | 2.80 | 0.0280 | 39609 | 0.8294 | 9485 | 0.1986 |
| [-95.9284 - -95.4916) | -95.7100 | 4236 | 8.87 | 0.0887 | 43845 | 0.9181 | 8148 | 0.1706 |
| [-95.4916 - -95.0548) | -95.2732 | 3083 | 6.46 | 0.0646 | 46928 | 0.9826 | 3912 | 0.0819 |
| [-95.0548 - -94.6179] | -94.8364 | 829 | 1.74 | 0.0174 | 47757 | 1.0000 | 829 | 0.0174 |
| TOTAL | - | 47757 | 100.00 | 1.0000 | 47757 | 1.0000 | 47757 | 1.0000 |
| Autor: Fernando Almeida | ||||||||
Aplicando la Regla de Sturges, con 47,757 observaciones corresponderían 17 intervalos de clase. Sin embargo, un desglose tan fino dificulta la lectura visual y el ajuste de los modelos teóricos sobre la distribución. Por ello, se simplifica el análisis a una tabla de 10 intervalos de clase, criterio que conserva representatividad estadística sin sacrificar claridad interpretativa.
intervalos_cut <- cut(x, breaks=breaks_vec, right=FALSE, include.lowest=TRUE)
freq_abs <- as.integer(table(intervalos_cut))
hi_dec <- freq_abs / n
Ni_asc <- cumsum(freq_abs); Hi_asc <- cumsum(hi_dec)
Ni_desc <- n - c(0, head(Ni_asc,-1)); Hi_desc <- 1 - c(0, head(Hi_asc,-1))
etiq <- paste0("[",round(lim_inf,4)," - ",round(lim_sup,4),")")
etiq[k] <- paste0("[",round(lim_inf[k],4)," - ",round(lim_sup[k],4),"]")
bind_rows(
data.frame(Intervalo=etiq, MC=round(mc,4), ni=freq_abs,
hi_pct=round(hi_dec*100,2), hi_real=round(hi_dec,4),
Ni_a=Ni_asc, Hi_a=round(Hi_asc,4),
Ni_d=Ni_desc, Hi_d=round(Hi_desc,4), stringsAsFactors=FALSE),
data.frame(Intervalo="TOTAL", MC=NA_real_, ni=sum(freq_abs),
hi_pct=round(sum(hi_dec)*100,2), hi_real=round(sum(hi_dec),4),
Ni_a=max(Ni_asc), Hi_a=round(max(Hi_asc),4),
Ni_d=max(Ni_desc), Hi_d=round(max(Hi_desc),4), stringsAsFactors=FALSE)
) %>%
gt() %>%
tab_header(title=md("**Tabla N\u00b01: Distribución de Frecuencias**"),
subtitle=md(paste0("*Variable Cuantitativa Continua: Longitud, ",
"Kansas (n = ", format(n, big.mark=","), ")*"))) %>%
cols_label(Intervalo=md("**Intervalo**"), MC=md("**MC**"), ni=md("**ni**"),
hi_pct=md("**hi%**"), hi_real=md("**hi**"),
Ni_a=md("**Ni (asc)**"), Hi_a=md("**Hi (asc)**"),
Ni_d=md("**Ni (desc)**"), Hi_d=md("**Hi (desc)**")) %>%
tab_style(style=list(cell_fill(color="#2C2C2C"),cell_text(color="white",weight="bold")),
locations=cells_column_labels()) %>%
tab_style(style=cell_fill(color="#F5F5F5"), locations=cells_body(rows=seq(1,k+1,by=2))) %>%
tab_style(style=list(cell_fill(color="#D6D6D6"),cell_text(weight="bold")),
locations=cells_body(rows=Intervalo=="TOTAL")) %>%
fmt_missing(columns=everything(), missing_text="-") %>%
tab_source_note(source_note=md("*Autor: Fernando Almeida*")) %>%
tab_options(table.width=pct(100), table.font.size=px(13), data_row.padding=px(6))
| Tabla N°1: Distribución de Frecuencias | ||||||||
| Variable Cuantitativa Continua: Longitud, Kansas (n = 47,757) | ||||||||
| Intervalo | MC | ni | hi% | hi | Ni (asc) | Hi (asc) | Ni (desc) | Hi (desc) |
|---|---|---|---|---|---|---|---|---|
| [-102.0439 - -101.3013) | -101.6726 | 3926 | 8.22 | 0.0822 | 3926 | 0.0822 | 47757 | 1.0000 |
| [-101.3013 - -100.5587) | -100.9300 | 7021 | 14.70 | 0.1470 | 10947 | 0.2292 | 43831 | 0.9178 |
| [-100.5587 - -99.8161) | -100.1874 | 3872 | 8.11 | 0.0811 | 14819 | 0.3103 | 36810 | 0.7708 |
| [-99.8161 - -99.0735) | -99.4448 | 7460 | 15.62 | 0.1562 | 22279 | 0.4665 | 32938 | 0.6897 |
| [-99.0735 - -98.3309) | -98.7022 | 9191 | 19.25 | 0.1925 | 31470 | 0.6590 | 25478 | 0.5335 |
| [-98.3309 - -97.5883) | -97.9596 | 3517 | 7.36 | 0.0736 | 34987 | 0.7326 | 16287 | 0.3410 |
| [-97.5883 - -96.8457) | -97.2170 | 2406 | 5.04 | 0.0504 | 37393 | 0.7830 | 12770 | 0.2674 |
| [-96.8457 - -96.1031) | -96.4744 | 1630 | 3.41 | 0.0341 | 39023 | 0.8171 | 10364 | 0.2170 |
| [-96.1031 - -95.3605) | -95.7318 | 6222 | 13.03 | 0.1303 | 45245 | 0.9474 | 8734 | 0.1829 |
| [-95.3605 - -94.6179] | -94.9892 | 2512 | 5.26 | 0.0526 | 47757 | 1.0000 | 2512 | 0.0526 |
| TOTAL | - | 47757 | 100.00 | 1.0000 | 47757 | 1.0000 | 47757 | 1.0000 |
| Autor: Fernando Almeida | ||||||||
grises <- gray(seq(0.25, 0.80, length.out=k))
h_obj <- hist(x, breaks=breaks_vec, plot=FALSE)
h_obj$density <- hi_dec
par(mar=c(5,6,6,2))
plot(h_obj, col=grises, border="black", freq=FALSE,
main="", xlab="", ylab="", las=1, xaxt="n")
axis(1, at=breaks_vec, labels=round(breaks_vec,3), las=1, cex.axis=0.9)
mtext("Densidad de Probabilidad", side=2, line=4.5, cex=1)
mtext("Longitud (grados)", side=1, line=3.5, cex=1)
mtext("Histograma General \u2014 Longitud, arrendamientos de hidrocarburos, Kansas, EE.UU.",
side=3, line=3, cex=0.95, font=2)
El histograma de la variable Longitud muestra un comportamiento diferenciado por zona geográfica. Se trabaja con tres zonas, con cortes fijos derivados de los límites de los 10 intervalos de la Sección 4: Oeste (intervalos 1-2), Centro (intervalos 3-7) y Este (intervalos 8-10). La zona Oeste concentra los valores menos significativos de la distribución frente a Centro y Este, por lo que se omite del ajuste inferencial formal; en Centro y Este se ajusta un modelo Normal forzado. Al trabajar por tramos, los parámetros de cada modelo se estiman exclusivamente con los datos de su zona.
set.seed(42)
lbl <- c(normal = "Normal", lognormal = "Log-Normal", exponential = "Exponencial")
# Número de bins adaptado al tamaño de la zona
k_adaptativo <- function(n_z) max(3, min(k, ceiling(1 + 3.322 * log10(n_z))))
# Calcula el Pearson entre hi_obs y hi_teo de un modelo ya ajustado
calc_pearson_zona <- function(datos, modelo, fit, offset, lim_min, lim_max) {
k_zona <- k_adaptativo(length(datos))
brks <- seq(lim_min, lim_max, length.out = k_zona + 1)
hi_obs <- hist(datos, breaks = brks, plot = FALSE)$counts / length(datos)
mc_z <- (head(brks,-1) + tail(brks,-1)) / 2
if (modelo == "lognormal") {
mc_pos <- mc_z - offset; mc_pos[mc_pos <= 0] <- 1e-9
hi_teo <- dlnorm(mc_pos, meanlog = fit$estimate["meanlog"],
sdlog = fit$estimate["sdlog"]) * diff(brks)
} else if (modelo == "exponential") {
mc_pos <- mc_z - offset; mc_pos[mc_pos <= 0] <- 1e-9
hi_teo <- dexp(mc_pos, rate = fit$estimate["rate"]) * diff(brks)
} else {
hi_teo <- dnorm(mc_z, mean = fit$estimate["mean"], sd = fit$estimate["sd"]) * diff(brks)
}
hi_teo <- hi_teo / sum(hi_teo)
round(cor(hi_obs, hi_teo) * 100, 2)
}
# Ajusta un modelo forzado específico
ajustar_forzado <- function(datos, modelo) {
offset <- min(datos) - 0.001
d_pos <- datos - offset
if (modelo == "normal") {
f <- fitdistr(datos, "normal")
list(fit = f, modelo = "normal", offset = 0,
pval = ks.test(sample(datos, min(400, length(datos))),
"pnorm", mean = f$estimate["mean"],
sd = f$estimate["sd"])$p.value)
} else if (modelo == "lognormal") {
f <- fitdistr(d_pos, "lognormal")
list(fit = f, modelo = "lognormal", offset = offset,
pval = ks.test(sample(d_pos, min(400, length(d_pos))),
"plnorm", meanlog = f$estimate["meanlog"],
sdlog = f$estimate["sdlog"])$p.value)
} else if (modelo == "exponential") {
f <- fitdistr(d_pos, "exponential")
list(fit = f, modelo = "exponential", offset = offset,
pval = ks.test(sample(d_pos, min(400, length(d_pos))),
"pexp", rate = f$estimate["rate"])$p.value)
}
}
# Valida un modelo: Pearson > 70% = APROBADO
validar <- function(datos, res, lim_min, lim_max) {
set.seed(42)
samp <- sample(datos, size = min(50, length(datos)), replace = FALSE)
offset <- res$offset
pearson <- calc_pearson_zona(datos, res$modelo, res$fit, offset, lim_min, lim_max)
if (res$modelo == "lognormal") {
ks_res <- ks.test(samp - offset, "plnorm",
meanlog = res$fit$estimate["meanlog"],
sdlog = res$fit$estimate["sdlog"])
} else if (res$modelo == "exponential") {
ks_res <- ks.test(samp - offset, "pexp", rate = res$fit$estimate["rate"])
} else {
ks_res <- ks.test(samp, "pnorm",
mean = res$fit$estimate["mean"],
sd = res$fit$estimate["sd"])
}
list(pearson = pearson, pval = round(ks_res$p.value, 4),
val = ifelse(pearson > 70, "APROBADO", "RECHAZADO"))
}
# ── Cortes de zona (derivados de los límites de los 10 intervalos de la
# Sección 4: Oeste = intervalos 1-2, Centro = intervalos 3-7, Este = 8-10) ──
lim_oeste_min <- breaks_vec[1] # límite inferior del intervalo 1 (x_min)
corte1 <- breaks_vec[3] # límite entre el intervalo 2 y el intervalo 3
corte2 <- breaks_vec[8] # límite entre el intervalo 7 y el intervalo 8
lim_este_min <- corte2
lim_este_max <- breaks_vec[11] # límite superior del intervalo 10 (x_max)
x_centro <- x[x >= corte1 & x < corte2]
x_este <- x[x >= lim_este_min & x <= lim_este_max]
# ── Zona Centro: Normal forzado ──────────────────────────────────────────────
r_centro <- ajustar_forzado(x_centro, "normal")
v_centro <- validar(x_centro, r_centro, corte1, corte2)
# ── Zona Este: Normal forzado ────────────────────────────────────────────────
r_este <- ajustar_forzado(x_este, "normal")
v_este <- validar(x_este, r_este, lim_este_min, lim_este_max)
cat("Zona Centro [", round(corte1,3), ",", round(corte2,3), "):", length(x_centro), "obs | Modelo:", lbl[r_centro$modelo],
"| Pearson:", v_centro$pearson, "% | KS p =", v_centro$pval, "\n")
## Zona Centro [ -100.559 , -96.846 ): 26446 obs | Modelo: Normal | Pearson: 91.4 % | KS p = 0.8505
cat("Zona Este [", round(lim_este_min,3), ",", round(lim_este_max,3), "]:", length(x_este), "obs | Modelo:", lbl[r_este$modelo],
"| Pearson:", v_este$pearson, "% | KS p =", v_este$pval, "\n")
## Zona Este [ -96.846 , -94.618 ]: 10364 obs | Modelo: Normal | Pearson: 90.46 % | KS p = 0.3316
cat("(Zona Oeste [", round(lim_oeste_min,3), ",", round(corte1,3), ") omitida del ajuste formal)\n")
## (Zona Oeste [ -102.044 , -100.559 ) omitida del ajuste formal)
h_plot <- hist(x, breaks=breaks_vec, plot=FALSE)
h_plot$density <- hi_dec
colores_zona <- ifelse(breaks_vec[-length(breaks_vec)] < corte1, "gray35",
ifelse(breaks_vec[-length(breaks_vec)] < corte2, "gray55", "gray75"))
par(mar=c(5,6,6,2))
plot(h_plot, col=colores_zona, border="black", freq=FALSE,
main="", xlab="", ylab="", las=1, xaxt="n")
axis(1, at=breaks_vec, labels=round(breaks_vec,3), las=1, cex.axis=0.9)
abline(v=corte1, col="black", lty=2, lwd=2)
abline(v=corte2, col="black", lty=2, lwd=2)
yt <- max(hi_dec)
text(mean(c(lim_oeste_min, corte1)), yt*0.88,
paste0("Oeste\n(", round(lim_oeste_min,2), " \u2013 ", round(corte1,2), ")\n(omitida)"),
cex=0.85, font=2)
text(mean(c(corte1, corte2)), yt*0.88,
paste0("Centro\n(", round(corte1,2), " \u2013 ", round(corte2,2), ")\n", lbl[r_centro$modelo]),
cex=0.85, font=2)
text(mean(c(corte2, lim_este_max)), yt*0.88,
paste0("Este\n(", round(corte2,2), " \u2013 ", round(lim_este_max,2), ")\n", lbl[r_este$modelo]),
cex=0.85, font=2)
mtext("Densidad de Probabilidad", side=2, line=4.5, cex=1)
mtext("Longitud (grados)", side=1, line=3.5, cex=1)
mtext(paste0("Cortes en ", round(corte1,3), "\u00b0 y ", round(corte2,3),
"\u00b0 \u2014 Zona Centro: ", lbl[r_centro$modelo],
" | Zona Este: ", lbl[r_este$modelo],
" \u2014 Kansas, EE.UU."),
side=3, line=3, cex=0.9, font=2)
plot_zona <- function(datos, res, titulo, pal, lim_min, lim_max) {
offset <- res$offset
n_z <- length(datos)
# Mismos intervalos (ancho c_amp) que el histograma principal (Sección 5)
n_bins_z <- ceiling((lim_max - lim_min) / c_amp)
lim_inf_z <- lim_min + (0:(n_bins_z - 1)) * c_amp
lim_sup_z <- lim_inf_z + c_amp
lim_sup_z[n_bins_z] <- lim_max
brks <- c(lim_inf_z, lim_sup_z[n_bins_z])
h_z <- hist(datos, breaks = brks, plot = FALSE)
hi_z <- h_z$counts / n_z # hi: frecuencia relativa
dens_z <- hi_z / diff(brks) # densidad real = hi / amplitud
h_z$density <- dens_z
xs <- seq(lim_min, lim_max, length.out = 500)
nb <- length(brks) - 1
par(mar = c(5,6,6,2))
plot(h_z, col = pal[seq_len(nb)], border = "black", freq = FALSE,
main = "", xlab = "", ylab = "", las = 1, xaxt = "n")
axis(1, at = brks, labels = round(brks,3), las = 1, cex.axis = 0.9)
# Densidad teórica del modelo ajustado (Normal)
if (res$modelo == "lognormal") {
xs_pos <- xs - offset; xs_pos[xs_pos <= 0] <- 1e-9
ys <- dlnorm(xs_pos, meanlog = res$fit$estimate["meanlog"],
sdlog = res$fit$estimate["sdlog"])
} else if (res$modelo == "exponential") {
xs_pos <- xs - offset; xs_pos[xs_pos <= 0] <- 1e-9
ys <- dexp(xs_pos, rate = res$fit$estimate["rate"])
} else {
ys <- dnorm(xs, mean = res$fit$estimate["mean"], sd = res$fit$estimate["sd"])
}
lines(xs, ys, col = "black", lwd = 2.5)
mtext("Densidad de Probabilidad", side=2, line=4.5, cex=1)
mtext("Longitud (grados)", side=1, line=3.5, cex=1)
mtext(titulo, side=3, line=3, cex=0.95, font=2)
legend("topright",
legend = c("Histograma", paste0("Curva ", lbl[res$modelo])),
fill = c("gray55", NA), border = c("black", NA),
lty = c(NA, 1), lwd = c(NA, 2.5), bty = "n", cex = 0.85)
}
par(mfrow=c(2,1))
plot_zona(x_centro, r_centro,
paste0("Zona Centro [", round(corte1,3), " \u2013 ", round(corte2,3),
") \u2014 ", lbl[r_centro$modelo],
" (\u03bc\u0302 = ", round(r_centro$fit$estimate["mean"],3),
", \u03c3\u0302 = ", round(r_centro$fit$estimate["sd"],3), ")"),
gray(seq(0.30,0.75,length.out=k)), corte1, corte2)
plot_zona(x_este, r_este,
paste0("Zona Este [", round(lim_este_min,3), " \u2013 ", round(lim_este_max,3),
"] \u2014 ", lbl[r_este$modelo],
" (\u03bc\u0302 = ", round(r_este$fit$estimate["mean"],3),
", \u03c3\u0302 = ", round(r_este$fit$estimate["sd"],3), ")"),
gray(seq(0.40,0.85,length.out=k)), lim_este_min, lim_este_max)
par(mfrow=c(1,1))
cat("--- [", lbl[r_centro$modelo], "] Zona Centro \u2014 Parámetros ---\n")
## --- [ Normal ] Zona Centro — Parámetros ---
print(round(r_centro$fit$estimate, 6))
## mean sd
## -98.894124 0.861029
cat("Pearson R%:", v_centro$pearson, "| KS p-valor:", v_centro$pval, "\n")
## Pearson R%: 91.4 | KS p-valor: 0.8505
cat("\n--- [", lbl[r_este$modelo], "] Zona Este \u2014 Parámetros ---\n")
##
## --- [ Normal ] Zona Este — Parámetros ---
print(round(r_este$fit$estimate, 6))
## mean sd
## -95.653103 0.455368
cat("Pearson R%:", v_este$pearson, "| KS p-valor:", v_este$pval, "\n")
## Pearson R%: 90.46 | KS p-valor: 0.3316
# ── Tabla de validación ───────────────────────────────────────────────────────
tabla_resumen <- bind_rows(
data.frame(Zona = sprintf("Centro (%.3f a %.3f)", corte1, corte2),
Modelo = as.character(lbl[r_centro$modelo]),
Pearson = v_centro$pearson,
KS_p = v_centro$pval,
Validacion = v_centro$val, stringsAsFactors = FALSE),
data.frame(Zona = sprintf("Este (%.3f a %.3f)", lim_este_min, lim_este_max),
Modelo = as.character(lbl[r_este$modelo]),
Pearson = v_este$pearson,
KS_p = v_este$pval,
Validacion = v_este$val, stringsAsFactors = FALSE)
)
tabla_resumen %>%
gt() %>%
tab_header(
title = md("**Tabla N\u00b02: Resumen de Validación por Zona**"),
subtitle = md("*Pearson (R%) y Kolmogorov-Smirnov (p-valor) por zona \u2014 Variable Longitud*")
) %>%
cols_label(Zona=md("**Zona**"), Modelo=md("**Modelo Aplicado**"),
Pearson=md("**Pearson (R %)**"), KS_p=md("**K-S (p-valor)**"),
Validacion=md("**Validación**")) %>%
tab_style(style=list(cell_fill(color="#2C2C2C"),cell_text(color="white",weight="bold")),
locations=cells_column_labels()) %>%
tab_style(style=list(cell_fill(color="#2C2C2C"),cell_text(color="white",weight="bold")),
locations=cells_title(groups="title")) %>%
tab_style(style=cell_fill(color="#F5F5F5"), locations=cells_body(rows=1)) %>%
tab_style(style=cell_text(color="darkgreen",weight="bold"),
locations=cells_body(columns=Validacion, rows=Validacion=="APROBADO")) %>%
tab_style(style=cell_text(color="darkred",weight="bold"),
locations=cells_body(columns=Validacion, rows=Validacion=="RECHAZADO")) %>%
tab_style(style=cell_borders(sides="bottom",color="#E0E0E0",weight=px(1)),
locations=cells_body(rows=everything())) %>%
cols_align(align="center", columns=c(Pearson,KS_p,Validacion)) %>%
cols_align(align="left", columns=c(Zona,Modelo)) %>%
fmt_number(columns=Pearson, decimals=2) %>%
fmt_number(columns=KS_p, decimals=4) %>%
tab_source_note(source_note=md("*Autor: Fernando Almeida*")) %>%
tab_options(table.width=pct(90), table.font.size=px(13),
heading.title.font.size=px(16), heading.subtitle.font.size=px(12),
data_row.padding=px(6),
column_labels.border.top.width=px(2),
column_labels.border.bottom.width=px(2),
table_body.border.bottom.width=px(2),
table.border.top.style="hidden", table.border.bottom.style="hidden")
| Tabla N°2: Resumen de Validación por Zona | ||||
| Pearson (R%) y Kolmogorov-Smirnov (p-valor) por zona — Variable Longitud | ||||
| Zona | Modelo Aplicado | Pearson (R %) | K-S (p-valor) | Validación |
|---|---|---|---|---|
| Centro (-100.559 a -96.846) | Normal | 91.40 | 0.8505 | APROBADO |
| Este (-96.846 a -94.618) | Normal | 90.46 | 0.3316 | APROBADO |
| Autor: Fernando Almeida | ||||
# Evento estándar por zona (mismo criterio en ambas): P(X < q1), P(q1 <= X < q3),
# P(X >= q3), usando los cuartiles 25% y 75% del rango de la zona como puntos
# de corte y el modelo Normal ya ajustado para esa zona.
calc_eventos_zona <- function(res, lim_min, lim_max) {
q1 <- lim_min + (lim_max - lim_min) * 0.25
q3 <- lim_min + (lim_max - lim_min) * 0.75
pf <- function(v) pnorm(v, mean = res$fit$estimate["mean"], sd = res$fit$estimate["sd"])
list(q1 = q1, q3 = q3,
p_a = pf(q1), p_b = pf(q3) - pf(q1), p_c = 1 - pf(q3))
}
# ── Zona Centro: Normal(μ̂, σ̂) ───────────────────────────────────────────────
ev_centro <- calc_eventos_zona(r_centro, corte1, corte2)
# ── Zona Este: Normal(μ̂, σ̂) ─────────────────────────────────────────────────
ev_este <- calc_eventos_zona(r_este, lim_este_min, lim_este_max)
data.frame(
Zona = c(rep(paste0("Centro \u2014 ", lbl[r_centro$modelo]), 3),
rep(paste0("Este \u2014 ", lbl[r_este$modelo]), 3)),
Evento = c(sprintf("P(X < %.3f)", ev_centro$q1),
sprintf("P(%.3f \u2264 X < %.3f)", ev_centro$q1, ev_centro$q3),
sprintf("P(X \u2265 %.3f)", ev_centro$q3),
sprintf("P(X < %.3f)", ev_este$q1),
sprintf("P(%.3f \u2264 X < %.3f)", ev_este$q1, ev_este$q3),
sprintf("P(X \u2265 %.3f)", ev_este$q3)),
Descripcion = c(
"Longitud en el tercio oeste de la Zona Centro",
"Longitud en el tramo central de la Zona Centro",
"Longitud en la cola este de la Zona Centro",
"Longitud en el tercio oeste de la Zona Este",
"Longitud en el tramo central de la Zona Este",
"Longitud en la cola este de la Zona Este"
),
Probabilidad = round(c(ev_centro$p_a, ev_centro$p_b, ev_centro$p_c,
ev_este$p_a, ev_este$p_b, ev_este$p_c), 4)
) %>%
gt() %>%
tab_header(
title = md("**Tabla N\u00b03: Cálculo de Probabilidades por Zona**"),
subtitle = md("*La probabilidad es el área bajo la curva del modelo teórico \u2014 Variable Longitud*")
) %>%
cols_label(Zona=md("**Zona**"), Evento=md("**Evento**"),
Descripcion=md("**Descripción**"), Probabilidad=md("**Probabilidad**")) %>%
tab_style(style=list(cell_fill(color="#2C2C2C"),cell_text(color="white",weight="bold")),
locations=cells_column_labels()) %>%
tab_style(style=cell_fill(color="#F5F5F5"),
locations=cells_body(rows=seq(1,6,by=2))) %>%
tab_source_note(source_note=md("*Autor: Fernando Almeida*")) %>%
tab_options(table.width=pct(95), table.font.size=px(13),
heading.title.font.size=px(15), heading.subtitle.font.size=px(11),
data_row.padding=px(6))
| Tabla N°3: Cálculo de Probabilidades por Zona | |||
| La probabilidad es el área bajo la curva del modelo teórico — Variable Longitud | |||
| Zona | Evento | Descripción | Probabilidad |
|---|---|---|---|
| Centro — Normal | P(X < -99.630) | Longitud en el tercio oeste de la Zona Centro | 0.1962 |
| Centro — Normal | P(-99.630 ≤ X < -97.774) | Longitud en el tramo central de la Zona Centro | 0.7071 |
| Centro — Normal | P(X ≥ -97.774) | Longitud en la cola este de la Zona Centro | 0.0966 |
| Este — Normal | P(X < -96.289) | Longitud en el tercio oeste de la Zona Este | 0.0814 |
| Este — Normal | P(-96.289 ≤ X < -95.175) | Longitud en el tramo central de la Zona Este | 0.7718 |
| Este — Normal | P(X ≥ -95.175) | Longitud en la cola este de la Zona Este | 0.1468 |
| Autor: Fernando Almeida | |||
El Intervalo de Confianza representa el puente fundamental entre los modelos empíricos observados y la estimación poblacional. Aunque la distribución original de la Longitud presenta asimetría y comportamiento diferenciado por zona, el TLC garantiza que la distribución de las medias muestrales tenderá a la normalidad debido al volumen masivo de datos (\(n=\) 47,757).
Los postulados de confianza empírica sugieren:
\[P(\bar{x} - E < \mu < \bar{x} + E) \approx 68\%\]
\[P(\bar{x} - 2E < \mu < \bar{x} + 2E) \approx 95\%\]
\[P(\bar{x} - 3E < \mu < \bar{x} + 3E) \approx 99\%\]
Donde el Margen de Error (\(E\)) se define como:
\[E = \frac{\sigma}{\sqrt{n}}\]
media_muestral <- mean(x)
desv_est <- sd(x)
n_total <- length(x)
z_95 <- 1.96
error_est <- desv_est / sqrt(n_total)
margen <- z_95 * error_est
lim_inf_ic <- media_muestral - margen
lim_sup_ic <- media_muestral + margen
data.frame(
Parametro = "Longitud Promedio Kansas (\u00b0)",
Lim_Inferior = round(lim_inf_ic, 3),
Media_Muestral = round(media_muestral, 3),
Lim_Superior = round(lim_sup_ic, 3),
Error_Estandar = paste0("+/- ", round(margen, 4)),
Confianza = "95% (Z = 1.96)",
stringsAsFactors = FALSE
) %>%
gt() %>%
tab_header(
title = md("**TABLA N\u00b0 4: ESTIMACIÓN DE LA MEDIA POBLACIONAL**"),
subtitle = md("*Inferencia Estadística para la Variable Longitud*")
) %>%
cols_label(
Parametro = md("**Parámetro**"),
Lim_Inferior = md("**Lim_Inferior**"),
Media_Muestral = md("**Media_Muestral**"),
Lim_Superior = md("**Lim_Superior**"),
Error_Estandar = md("**Error_Estándar**"),
Confianza = md("**Confianza**")
) %>%
tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
locations = cells_column_labels()) %>%
tab_style(style = list(cell_fill(color = "#2C2C2C"), cell_text(color = "white", weight = "bold")),
locations = cells_title(groups = "title")) %>%
tab_style(style = list(cell_fill(color = "#C8E6C9"), cell_text(weight = "bold")),
locations = cells_body(columns = Media_Muestral)) %>%
tab_style(style = cell_borders(sides = "bottom", color = "#E0E0E0", weight = px(1)),
locations = cells_body(rows = everything())) %>%
cols_align(align = "center", columns = c(Lim_Inferior, Media_Muestral, Lim_Superior, Error_Estandar, Confianza)) %>%
cols_align(align = "left", columns = Parametro) %>%
tab_source_note(source_note = md("*Autor: Fernando Almeida*")) %>%
tab_options(table.width = pct(100), table.font.size = px(13),
heading.title.font.size = px(16), heading.subtitle.font.size = px(12),
data_row.padding = px(6),
column_labels.border.top.width = px(2),
column_labels.border.bottom.width = px(2),
table_body.border.bottom.width = px(2),
table.border.top.style = "hidden", table.border.bottom.style = "hidden")
| TABLA N° 4: ESTIMACIÓN DE LA MEDIA POBLACIONAL | |||||
| Inferencia Estadística para la Variable Longitud | |||||
| Parámetro | Lim_Inferior | Media_Muestral | Lim_Superior | Error_Estándar | Confianza |
|---|---|---|---|---|---|
| Longitud Promedio Kansas (°) | -98.737 | -98.719 | -98.701 | +/- 0.0178 | 95% (Z = 1.96) |
| Autor: Fernando Almeida | |||||
El comportamiento de Longitud se explica con un modelo Normal (μ̂ = -98.8941, σ̂ = 0.8610) en la Zona Centro [-100.559, -96.846) y un modelo Normal (μ̂ = -95.6531, σ̂ = 0.4554) en la Zona Este [-96.846, -94.618]. Podemos afirmar con un 95% de confianza que la media aritmética real de Longitud se encuentra entre -98.7367 y -98.701, con una desviación estándar de 1.9889.
Autor: Fernando Almeida — Análisis Estadístico, Kansas Hydrocarbon Leases Dataset