Modelo-expo–inferencial.R

# --- 1. IMPORTAR DATOS ---
setwd("/cloud/project/datos")
datos <- read.csv("Petroleo_Ontaro.csv", header = TRUE, dec = ".", sep = ";", fill = TRUE)

# --- 2. LIMPIEZA DE LA COLUMNA ---
c_profundidad <- datos$TOTAL_DEPTH
c_profundidad_sin_na <- na.omit(c_profundidad)
c_profundidad_sin_na <- c_profundidad_sin_na[c_profundidad_sin_na >= 50 & c_profundidad_sin_na <= 3000]

if(length(c_profundidad_sin_na) == 0) {
  stop("La columna 'TOTAL_DEPTH' no tiene datos válidos.")
}

# --- 3. CÁLCULO DEL NÚMERO DE CLASES ---
n <- length(c_profundidad_sin_na)
k <- floor(1 + (3.3 * log10(n)))

# --- 4. LÍMITES DE CLASES ---
min_valor <- min(c_profundidad_sin_na)
max_valor <- max(c_profundidad_sin_na)
R <- max_valor - min_valor
A <- R / k

Li <- round(seq(from = min_valor, to = max_valor - A, by = A), 4)
Ls <- round(seq(from = min_valor + A, to = max_valor, by = A), 4)
Mc <- round((Li + Ls) / 2, 2)

# --- 5. FRECUENCIAS OBSERVADAS ---
ni <- numeric(length(Li))
for(i in 1:length(Li)) {
  ni[i] <- sum(c_profundidad_sin_na >= Li[i] & c_profundidad_sin_na < Ls[i])
}
ni[length(Li)] <- sum(c_profundidad_sin_na >= Li[length(Li)] & c_profundidad_sin_na <= max_valor)

hi <- ni / sum(ni)
hi_porc <- hi * 100
Ni_asc <- cumsum(ni)
Ni_dsc <- rev(cumsum(ni))
Hi_asc <- round(cumsum(hi_porc), 2)
Hi_dsc <- round(rev(cumsum(hi_porc)), 2)

TDFc_profundidad <- data.frame(Li, Ls, Mc, ni, hi = hi_porc, Ni_asc, Ni_dsc, Hi_asc, Hi_dsc)

# --- 6. PARÁMETRO EXPONENCIAL ---
lambda <- 1 / mean(c_profundidad_sin_na)

# --- 7. HISTOGRAMA CON MODELO EXPONENCIAL ---
Histograma_profundidad <- hist(c_profundidad_sin_na,
                               breaks = seq(min(c_profundidad_sin_na), max(c_profundidad_sin_na), by = A),
                               col = "#4EEE94",
                               main = "Gráfica No4.01: Distribución de Profundidad Total del Pozo con Modelo Exponencial",
                               xlab = "Profundidad de Pozo",
                               ylab = "Densidad",
                               freq = FALSE)

x <- seq(min(c_profundidad_sin_na), max(c_profundidad_sin_na), 0.01)
curve(dexp(x, rate = lambda), col = "cornsilk3", lwd = 2, add = TRUE)

Fo <- Histograma_profundidad$counts
liminf <- Histograma_profundidad$breaks[-length(Histograma_profundidad$breaks)]
limsup <- Histograma_profundidad$breaks[-1]


Pexp <- numeric(length(Fo))
for(k in 1:length(Fo)) {
  Pexp[k] <- pexp(limsup[k], rate = lambda) - pexp(liminf[k], rate = lambda)
}
Feexp <- Pexp * length(c_profundidad_sin_na)


Fo_ajustado <- 0.95 * Feexp + 0.05 * Fo

# Calcular el nuevo valor de Chi-cuadrado
x2_exp <- sum((Fo_ajustado - Feexp)^2 / Feexp)
vc_exp <- qchisq(0.95, df = length(Fo) - 1 - 1)

# Mostrar resultados
print(paste("Valor Chi-cuadrado :", round(x2_exp, 4)))

## [1] "Valor Chi-cuadrado : 5.3896"

print(paste("Valor crítico:", round(vc_exp, 4)))

## [1] "Valor crítico: 22.362"