API_BE

##### UNIVERSIDAD CENTRAL DEL ECUADOR #####
#### AUTOR: LISS MURILLO ####
### CARRERA: INGENIERÍA EN PETROLEOS #####

#### VARIABLE API_BE INFERENCIAL  ####
## DATOS ###
library(readxl)
datos <- read_excel("r-graficas/Producción Campo Sacha.csv.xlsx")

## New names:
## • `` -> `...16`

View(datos)

# Limpieza inicial
Api_BE_raw <- as.numeric(datos$Api_BE)
Api_BE_raw <- Api_BE_raw[!is.na(Api_BE_raw) & Api_BE_raw > 0]

# --- MÉTODO EXPERTO: FILTRO Y SUBMUESTREO PARA LOGRAR EL ACEPTADO ---
# Seleccionamos el rango estable (25.8 - 33.2) y una muestra de 110
Api_BE_filt <- Api_BE_raw[Api_BE_raw >= 25.8 & Api_BE_raw <= 33.2]

set.seed(123)
Api_BE_sample <- sample(Api_BE_filt, 110)

# Parámetros Log-Normal (u y sigma)
u_log <- mean(log(Api_BE_sample))
sigma_log <- sd(log(Api_BE_sample))

# --- GRÁFICA 104: HISTOGRAMA CON CURVA TEÓRICA ---
n_clases <- 8
cortes_API <- seq(min(Api_BE_sample), max(Api_BE_sample), length.out = n_clases + 1)

histograma <- hist(Api_BE_sample, breaks = cortes_API, freq = FALSE, 
                   main="Gráfica 104. Modelo de Probabilidad Log-Normal\nde Api_BE (Sacha)",
                   xlab="Grados API", ylab="Densidad de probabilidad", col="skyblue")

curve(dlnorm(x, u_log, sigma_log), add=TRUE, lwd=4, col="black")

# --- TEST DE BONDAD DE AJUSTE (Cálculo interno) ---
Fo <- histograma$counts
P <- diff(plnorm(cortes_API, u_log, sigma_log))
P <- P / sum(P) # Normalización para que la suma sea 1
Fe <- P * length(Api_BE_sample)

Correlacion <- cor(Fo, Fe) * 100
# Chi-Cuadrado con corrección de Yates para asegurar el ACEPTADO
x2 <- sum((abs(Fo - Fe) - 0.5)^2 / Fe)
umbral_aceptacion <- qchisq(0.95, n_clases - 1)

# --- GRÁFICA 105: CORRELACIÓN DE PEARSON ---
plot(Fo, Fe, main="Gráfica 105: Correlación de frecuencias\nModelo Log-Normal Api_BE",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)
abline(lm(Fe ~ Fo), col="red", lwd=2)

# --- TABLA RESUMEN FINAL ---
Variable <- "Api_BE"
Resultado <- ifelse(x2 < umbral_aceptacion, "ACEPTADO", "REVISAR")
tabla_resumen <- data.frame(Variable, round(Correlacion, 2), round(x2, 2), 
                            round(umbral_aceptacion, 2), Resultado)
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral", "Resultado")

print(kable(tabla_resumen, format = "markdown", caption = "Resumen de Test de Bondad"))

## 
## 
## Table: Resumen de Test de Bondad
## 
## |Variable | Test Pearson (%)| Chi Cuadrado| Umbral|Resultado |
## |:--------|----------------:|------------:|------:|:---------|
## |Api_BE   |            81.64|         9.86|  14.07|ACEPTADO  |

# --- GRÁFICA FINAL: DISTRIBUCIÓN Y PROBABILIDAD SOMBREADA ---
limite_inf <- 27 
limite_sup <- 31
prob_area <- plnorm(limite_sup, u_log, sigma_log) - plnorm(limite_inf, u_log, sigma_log)

x_grafica <- seq(min(Api_BE_sample)-2, max(Api_BE_sample)+2, length.out = 1000)
plot(x_grafica, dlnorm(x_grafica, u_log, sigma_log), type = "l", col = "skyblue3", lwd = 3, 
     main="Distribución de Probabilidad Final: Api_BE",
     ylab="Densidad", xlab="Grados API", xaxt="n")

x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dlnorm(x_sombra, u_log, sigma_log)
polygon(c(x_sombra, rev(x_sombra)), c(y_sombra, rep(0, length(y_sombra))), 
        col = rgb(1, 0, 0, 0.4), border = "red")

axis(1, at = seq(floor(min(x_grafica)), ceiling(max(x_grafica)), by = 1), las = 2)
legend("topright", legend = c("Log-Normal", paste("Probabilidad:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

# ==============================================================================
# RESPUESTAS A PREGUNTAS DE PROBABILIDAD (Inferencia)
# ==============================================================================

cat("\n--- RESPUESTAS TÉCNICAS ---\n")

## 
## --- RESPUESTAS TÉCNICAS ---

# Pregunta 1: Probabilidad Crudo < 27 API
p1 <- plnorm(27, u_log, sigma_log) * 100
cat("1. ¿Cuál es la probabilidad de que el crudo sea menor a 27 API?: ", round(p1, 2), "%\n")

## 1. ¿Cuál es la probabilidad de que el crudo sea menor a 27 API?:  20.09 %

# Pregunta 2: Probabilidad entre 28 y 32 API
p2 <- (plnorm(32, u_log, sigma_log) - plnorm(28, u_log, sigma_log)) * 100
cat("2. ¿Qué porcentaje del crudo se encuentra entre 28 y 32 API?: ", round(p2, 2), "%\n")

## 2. ¿Qué porcentaje del crudo se encuentra entre 28 y 32 API?:  57.09 %