MODELO DE PROBABILIDAD PARA LONGITUDE

Estadística inferencial

Análisis inferencial de la variable cuantitativa continua Longitude, estructurado según los pasos solicitados.

1. CARGA DE LIBRERÍAS

# ====================================================
# Variable-Cuantitativa-Continua-Longitude.R
# Autor: r3532706
# Fecha: 2025-12-18
# ====================================================

# ====================================================
# LIBRERÍAS
# ====================================================
library(readxl)
library(dplyr)
library(gt)

# ====================================================
# CARGA DE LA BASE DE DATOS
# ====================================================

if (file.exists("datos_nuevoartes.xlsx")) {
  
  datos_nuevoartes <- read_excel("datos_nuevoartes.xlsx")
  
} else {
  
  archivos_xlsx <- list.files(pattern = "\\.xlsx$", ignore.case = TRUE)
  archivos_csv <- list.files(pattern = "\\.csv$", ignore.case = TRUE)
  
  if (length(archivos_xlsx) > 0) {
    
    datos_nuevoartes <- read_excel(archivos_xlsx[1])
    
  } else if (length(archivos_csv) > 0) {
    
    datos_nuevoartes <- read.csv(
      archivos_csv[1],
      header = TRUE,
      sep = ",",
      dec = ".",
      stringsAsFactors = FALSE
    )
    
    if (ncol(datos_nuevoartes) == 1) {
      datos_nuevoartes <- read.csv(
        archivos_csv[1],
        header = TRUE,
        sep = ";",
        dec = ",",
        stringsAsFactors = FALSE
      )
    }
    
  } else {
    
    stop("No se encontró datos_nuevoartes.xlsx ni ningún archivo CSV/XLSX en la carpeta del RMD.")
  }
}

datos_nuevoartes <- as.data.frame(datos_nuevoartes)
names(datos_nuevoartes) <- trimws(names(datos_nuevoartes))

2. SELECCIÓN DE LA VARIABLE

# ====================================================
# VARIABLE: LONGITUDE
# ====================================================
longitude <- datos_nuevoartes$longitude
longitude <- longitude[!is.na(longitude)]
longitude <- as.numeric(longitude)
longitude <- longitude[is.finite(longitude)]

3. TABLA DE FRECUENCIAS

# ====================================================
# PARÁMETROS DE CLASIFICACIÓN
# ====================================================
k_long <- 12
n_long <- length(longitude)

min_long <- min(longitude)
max_long <- max(longitude)

R_long <- max_long - min_long
A_real <- R_long / k_long

# ====================================================
# AJUSTE DE AMPLITUD DE CLASE
# ====================================================
A_long <- ifelse(
  A_real <= 2, 2,
  ifelse(
    A_real <= 5, 5,
    ifelse(A_real <= 10, 10, ceiling(A_real / 10) * 10)
  )
)

# ====================================================
# LÍMITES DE CLASE
# ====================================================
Li0 <- floor(min_long / A_long) * A_long

Li_long <- seq(Li0, by = A_long, length.out = k_long)
Ls_long <- Li_long + A_long
MC_long <- round((Li_long + Ls_long) / 2, 2)

# ====================================================
# FRECUENCIAS
# ====================================================
ni_long <- numeric(k_long)

for (i in 1:k_long) {
  if (i < k_long) {
    ni_long[i] <- sum(longitude >= Li_long[i] & longitude < Ls_long[i])
  } else {
    ni_long[i] <- sum(longitude >= Li_long[i] & longitude <= max_long)
  }
}

hi_long <- round((ni_long / sum(ni_long)) * 100, 2)

Ni_asc_long <- cumsum(ni_long)
Ni_dsc_long <- rev(cumsum(rev(ni_long)))

Hi_asc_long <- round(cumsum(hi_long), 2)
Hi_dsc_long <- round(rev(cumsum(rev(hi_long))), 2)

# ====================================================
# TABLA DE FRECUENCIAS
# ====================================================
TDF_longitude <- data.frame(
  Li = Li_long,
  Ls = Ls_long,
  MC = MC_long,
  ni = ni_long,
  hi = hi_long,
  Ni_asc = Ni_asc_long,
  Ni_dsc = Ni_dsc_long,
  Hi_asc = Hi_asc_long,
  Hi_dsc = Hi_dsc_long
)

TDF_longitude <- rbind(
  TDF_longitude,
  data.frame(
    Li = "TOTAL",
    Ls = "",
    MC = "",
    ni = sum(ni_long),
    hi = 100,
    Ni_asc = "",
    Ni_dsc = "",
    Hi_asc = "",
    Hi_dsc = ""
  )
)

# ====================================================
# TABLA DE PRESENTACIÓN (GT)
# ====================================================
tabla_longitude <- TDF_longitude %>%
  gt() %>%
  fmt_number(columns = MC, decimals = 2) %>%
  tab_header(
    title = md("**Tabla N° 1**"),
    subtitle = md("Distribución de frecuencias de Longitude (12 clases)")
  ) %>%
  tab_source_note(
    source_note = md("Autor: Grupo Geología")
  ) %>%
  tab_style(
    style = cell_text(weight = "bold"),
    locations = cells_body(rows = Li == "TOTAL")
  )

tabla_longitude
Tabla N° 1
Distribución de frecuencias de Longitude (12 clases)
Li Ls MC ni hi Ni_asc Ni_dsc Hi_asc Hi_dsc
-180 -150 -165 85 0.77 85 11033 0.77 100.02
-150 -120 -135 1893 17.16 1978 10948 17.93 99.25
-120 -90 -105 1221 11.07 3199 9055 29 82.09
-90 -60 -75 1526 13.83 4725 7834 42.83 71.02
-60 -30 -45 228 2.07 4953 6308 44.9 57.19
-30 0 -15 318 2.88 5271 6080 47.78 55.12
0 30 15 334 3.03 5605 5762 50.81 52.24
30 60 45 279 2.53 5884 5428 53.34 49.21
60 90 75 2034 18.44 7918 5149 71.78 46.68
90 120 105 1784 16.17 9702 3115 87.95 28.24
120 150 135 1083 9.82 10785 1331 97.77 12.07
150 180 165 248 2.25 11033 248 100.02 2.25
TOTAL 11033 100.00
Autor: Grupo Geología 
# ====================================================
# HISTOGRAMA – LONGITUDE
# ====================================================
hist(
  longitude,
  breaks = c(Li_long, max(Ls_long)),
  right = FALSE,
  freq = TRUE,
  col = "grey",
  border = "black",
  main = " Modelo de Probabilidad de la longitud ",
  xlab = "Longitude (°)",
  ylab = "Densidad de probabilidad"
)

4. CONJETURA DEL MODELO

Conjetura del modelo

Para la variable Longitude se trabajan tres agrupaciones de datos: un modelo normal para el intervalo de −108 a −60, un modelo uniforme para el intervalo de −60 a 30 y un modelo gamma para valores mayores a 30. 

#====================================== AGRUPACIÓN 1 ======================================
#======================= AGRUPACIÓN 1 (−108 a −60) =========================
# MODELO EXPONENCIAL CRECIENTE
#TABLA → HISTOGRAMA → PEARSON → CHI-CUADRADO

# -------------------------------------------------------------------------
# 1. FILTRADO
# -------------------------------------------------------------------------
lon_108_60 <- longitude[longitude >= -108 & longitude <= -60]

# -------------------------------------------------------------------------
# 2. INTERVALOS
# -------------------------------------------------------------------------
Li_108_60 <- seq(-108, -68, by = 10)
LS_108_60 <- Li_108_60 + 10
MC_108_60 <- (Li_108_60 + LS_108_60) / 2
k_108_60  <- length(Li_108_60)

# -------------------------------------------------------------------------
# 3. FRECUENCIAS
# -------------------------------------------------------------------------
ni_108_60 <- numeric(k_108_60)

for (i in 1:k_108_60) {
  if (i < k_108_60) {
    ni_108_60[i] <- sum(lon_108_60 >= Li_108_60[i] & lon_108_60 < LS_108_60[i])
  } else {
    ni_108_60[i] <- sum(lon_108_60 >= Li_108_60[i] & lon_108_60 <= -60)
  }
}

hi_108_60  <- (ni_108_60 / sum(ni_108_60)) * 100
NI_108_60  <- cumsum(ni_108_60)
NID_108_60 <- rev(cumsum(rev(ni_108_60)))
HI_108_60  <- cumsum(hi_108_60)
HID_108_60 <- rev(cumsum(rev(hi_108_60)))

# -------------------------------------------------------------------------
# 4. TABLA DE FRECUENCIAS
# -------------------------------------------------------------------------
Tabla_108_60 <- data.frame(
  Li  = Li_108_60,
  LS  = LS_108_60,
  MC  = MC_108_60,
  ni  = ni_108_60,
  hi  = round(hi_108_60, 2),
  NI  = NI_108_60,
  NID = NID_108_60,
  HI  = round(HI_108_60, 2),
  HID = round(HID_108_60, 2)
)

tabla_gt_108_60 <- Tabla_108_60 %>%
  gt() %>%
  tab_header(
    title = md("**Tabla N° 1**"),
    subtitle = md("Distribución de frecuencias Longitude (−108 a −60)")
  ) %>%
  fmt_number(columns = c(hi, HI, HID), decimals = 2) %>%
  tab_source_note(
    source_note = md("Autor: Grupo Geología")
  )

tabla_gt_108_60
Tabla N° 1
Distribución de frecuencias Longitude (−108 a −60)
Li LS MC ni hi NI NID HI HID
-108 -98 -103 211 10.69 211 1974 10.69 100.00
-98 -88 -93 288 14.59 499 1763 25.28 89.31
-88 -78 -83 881 44.63 1380 1475 69.91 74.72
-78 -68 -73 440 22.29 1820 594 92.20 30.09
-68 -58 -63 154 7.80 1974 154 100.00 7.80
Autor: Grupo Geología 
#====================================== AGRUPACIÓN 2 ======================================
#======================= AGRUPACIÓN 2 (−60 a 30) =========================
# MODELO UNIFORME
# TABLA → HISTOGRAMA → CHI-CUADRADO
# INTERVALOS FIJOS: −60 a −30 | −30 a 0 | 0 a 30

# -------------------------------------------------------------------------
# 1. FILTRADO
# -------------------------------------------------------------------------
lon_60_30 <- longitude[longitude >= -60 & longitude <= 30]
n_60_30 <- length(lon_60_30)

# -------------------------------------------------------------------------
# 2. INTERVALOS
# -------------------------------------------------------------------------
Li_60_30 <- c(-60, -30, 0)
LS_60_30 <- c(-30, 0, 30)
MC_60_30 <- (Li_60_30 + LS_60_30) / 2
k_60_30  <- length(Li_60_30)

# -------------------------------------------------------------------------
# 3. FRECUENCIAS OBSERVADAS
# -------------------------------------------------------------------------
ni_60_30 <- numeric(k_60_30)

for (i in 1:k_60_30) {
  if (i < k_60_30) {
    ni_60_30[i] <- sum(lon_60_30 >= Li_60_30[i] & lon_60_30 < LS_60_30[i])
  } else {
    ni_60_30[i] <- sum(lon_60_30 >= Li_60_30[i] & lon_60_30 <= 30)
  }
}

hi_60_30  <- (ni_60_30 / sum(ni_60_30)) * 100
NI_60_30  <- cumsum(ni_60_30)
NID_60_30 <- rev(cumsum(rev(ni_60_30)))
HI_60_30  <- cumsum(hi_60_30)
HID_60_30 <- rev(cumsum(rev(hi_60_30)))

# -------------------------------------------------------------------------
# 4. TABLA DE FRECUENCIAS
# -------------------------------------------------------------------------
Tabla_60_30 <- data.frame(
  Li  = Li_60_30,
  LS  = LS_60_30,
  MC  = MC_60_30,
  ni  = ni_60_30,
  hi  = round(hi_60_30, 2),
  NI  = NI_60_30,
  NID = NID_60_30,
  HI  = round(HI_60_30, 2),
  HID = round(HID_60_30, 2)
)

Tabla_60_30
##    Li  LS  MC  ni    hi  NI NID     HI    HID
## 1 -60 -30 -45 228 25.91 228 880  25.91 100.00
## 2 -30   0 -15 318 36.14 546 652  62.05  74.09
## 3   0  30  15 334 37.95 880 334 100.00  37.95
# ====================================================
# 4. TABLA DE FRECUENCIAS – AGRUPACIÓN 2
# ====================================================

Tabla_60_30 <- data.frame(
  Li  = Li_60_30,
  LS  = LS_60_30,
  MC  = MC_60_30,
  ni  = ni_60_30,
  hi  = round(hi_60_30, 2),
  NI  = NI_60_30,
  NID = NID_60_30,
  HI  = round(HI_60_30, 2),
  HID = round(HID_60_30, 2)
)

tabla_gt_60_30 <- Tabla_60_30 %>%
  gt() %>%
  tab_header(
    title = md("**Tabla N° 2**"),
    subtitle = md("Distribución de frecuencias Longitude (−60 a 30)")
  ) %>%
  fmt_number(columns = c(hi, HI, HID), decimals = 2) %>%
  tab_source_note(
    source_note = md("Autor: Grupo Geología")
  )

tabla_gt_60_30
Tabla N° 2
Distribución de frecuencias Longitude (−60 a 30)
Li LS MC ni hi NI NID HI HID
-60 -30 -45 228 25.91 228 880 25.91 100.00
-30 0 -15 318 36.14 546 652 62.05 74.09
0 30 15 334 37.95 880 334 100.00 37.95
Autor: Grupo Geología 
#====================================== AGRUPACIÓN 3 ======================================
#======================= AGRUPACIÓN 3 (Longitude > 30) =====================
# MODELO GAMMA
# TABLA → HISTOGRAMA → CURVA GAMMA
# INTERVALOS EQUIDISTANTES (30 en 30)

# -------------------------------------------------------------------------
# 1. FILTRADO
# -------------------------------------------------------------------------
lon_30_max <- longitude[longitude > 30]
n_30_max <- length(lon_30_max)

# -------------------------------------------------------------------------
# 2. INTERVALOS
# -------------------------------------------------------------------------
ancho_clase <- 30
max_lon <- ceiling(max(lon_30_max) / ancho_clase) * ancho_clase

Li_30_max <- seq(30, max_lon - ancho_clase, by = ancho_clase)
LS_30_max <- Li_30_max + ancho_clase
MC_30_max <- (Li_30_max + LS_30_max) / 2
k_30_max  <- length(Li_30_max)

breaks_lon_30_max <- c(Li_30_max, max_lon)
breaks_lon_30_max[1] <- min(breaks_lon_30_max[1], min(lon_30_max))
breaks_lon_30_max[length(breaks_lon_30_max)] <- max(
  breaks_lon_30_max[length(breaks_lon_30_max)],
  max(lon_30_max)
)

# -------------------------------------------------------------------------
# 3. FRECUENCIAS OBSERVADAS
# -------------------------------------------------------------------------
ni_30_max <- numeric(k_30_max)

for (i in 1:k_30_max) {
  if (i < k_30_max) {
    ni_30_max[i] <- sum(lon_30_max >= Li_30_max[i] & lon_30_max < LS_30_max[i])
  } else {
    ni_30_max[i] <- sum(lon_30_max >= Li_30_max[i] & lon_30_max <= max_lon)
  }
}

Fo_abs_30_max <- ni_30_max

# -------------------------------------------------------------------------
# 4. TABLA DE FRECUENCIAS (GT)
# -------------------------------------------------------------------------
hi_30_max  <- (Fo_abs_30_max / sum(Fo_abs_30_max)) * 100
NI_30_max  <- cumsum(Fo_abs_30_max)
NID_30_max <- rev(cumsum(rev(Fo_abs_30_max)))
HI_30_max  <- cumsum(hi_30_max)
HID_30_max <- rev(cumsum(rev(hi_30_max)))

Tabla_30_max <- data.frame(
  Li  = Li_30_max,
  LS  = LS_30_max,
  MC  = MC_30_max,
  ni  = Fo_abs_30_max,
  hi  = round(hi_30_max, 2),
  NI  = NI_30_max,
  NID = NID_30_max,
  HI  = round(HI_30_max, 2),
  HID = round(HID_30_max, 2)
)

tabla_gt_30_max <- Tabla_30_max %>%
  gt() %>%
  tab_header(
    title = md("**Tabla N° 3**"),
    subtitle = md("Distribución de frecuencias Longitude (> 30)")
  ) %>%
  fmt_number(columns = c(hi, HI, HID), decimals = 2) %>%
  tab_source_note(
    source_note = md("Autor: Grupo Geología")
  )

tabla_gt_30_max
Tabla N° 3
Distribución de frecuencias Longitude (> 30)
Li LS MC ni hi NI NID HI HID
30 60 45 279 5.14 279 5428 5.14 100.00
60 90 75 2034 37.47 2313 5149 42.61 94.86
90 120 105 1784 32.87 4097 3115 75.48 57.39
120 150 135 1083 19.95 5180 1331 95.43 24.52
150 180 165 248 4.57 5428 248 100.00 4.57
Autor: Grupo Geología

5. CÁLCULO DE PARÁMETROS

#========================================================================================
# AGRUPACIÓN 1: PARÁMETROS DEL MODELO NORMAL
#========================================================================================

breaks_lon_108_60 <- c(Li_108_60, -60)

HistoLon1 <- hist(
  lon_108_60,
  breaks = breaks_lon_108_60,
  right = FALSE,
  include.lowest = TRUE,
  plot = FALSE
)

Fo_abs_108_60 <- HistoLon1$counts

# Parámetros de la normal
mu_hat    <- mean(lon_108_60)
sigma_hat <- sd(lon_108_60)

# Dominio de la curva
x_norm <- seq(min(breaks_lon_108_60), max(breaks_lon_108_60), length.out = 400)

# Escalado correcto de la normal
y_norm <- dnorm(x_norm, mean = mu_hat, sd = sigma_hat)
y_norm <- y_norm * diff(breaks_lon_108_60)[1] * sum(Fo_abs_108_60)

#========================================================================================
# AGRUPACIÓN 2: PARÁMETROS DEL MODELO UNIFORME
#========================================================================================

a_60_30 <- min(lon_60_30)
b_60_30 <- max(lon_60_30)

breaks_lon_60_30 <- c(Li_60_30, 30)

#========================================================================================
# AGRUPACIÓN 3: PARÁMETROS DEL MODELO GAMMA
#========================================================================================

x_gamma_data <- lon_30_max - 30

media_g <- mean(x_gamma_data)
var_g   <- var(x_gamma_data)

shape_hat <- media_g^2 / var_g
rate_hat  <- media_g / var_g

x_gamma <- seq(min(x_gamma_data), max(x_gamma_data), length.out = 500)
y_gamma <- dgamma(x_gamma, shape = shape_hat, rate = rate_hat)

# Escalado correcto a frecuencias absolutas
y_gamma <- y_gamma * ancho_clase * sum(Fo_abs_30_max)

6. GRÁFICA: CURVA SOBRE EL HISTOGRAMA

hist(
  lon_108_60,
  breaks = breaks_lon_108_60,
  right = FALSE,
  include.lowest = TRUE,
  freq = TRUE,
  col = "grey",
  border = "black",
  main = "Histograma con modelo exponencial (−108 a −60)",
  xlab = "Longitude (°)",
  ylab = "Frecuencia absoluta"
)

# Dibujar curva
lines(x_norm, y_norm, col = "red", lwd = 2)

legend(
  "topright",
  legend = c("Histograma", "Modelo normal"),
  col = c("grey", "red"),
  lwd = c(10, 2),
  bty = "n"
)

hist(
  lon_60_30,
  breaks = breaks_lon_60_30,
  freq = TRUE,
  col = "grey",
  border = "black",
  main = "Histograma con modelo uniforme (−60 a 30)",
  xlab = "Longitude (°)",
  ylab = "Frecuencia absoluta"
)

curve(
  dunif(x, min = a_60_30, max = b_60_30) * n_60_30 * diff(breaks_lon_60_30)[1],
  add = TRUE,
  lwd = 2,
  col = "red"
)

legend(
  "topright",
  legend = c("Histograma", "Modelo uniforme"),
  col = c("grey", "red"),
  lwd = c(10, 2),
  bty = "n"
)

y_max <- max(c(Fo_abs_30_max, y_gamma), na.rm = TRUE) * 1.30

hist(
  lon_30_max,
  breaks = breaks_lon_30_max,
  right = FALSE,
  include.lowest = TRUE,
  freq = TRUE,
  col = "grey",
  border = "black",
  ylim = c(0, y_max),
  main = "Histograma con modelo Gamma (Longitude > 30)",
  xlab = "Longitude (°)",
  ylab = "Frecuencia absoluta"
)

# Curva gamma
lines(x_gamma + 30, y_gamma, col = "red", lwd = 2)

legend(
  "topright",
  legend = c("Histograma", "Modelo Gamma"),
  col = c("grey", "red"),
  lwd = c(10, 2),
  bty = "n"
)

7. TEST

#-------------------------------------------------------------------------
# FRECUENCIAS ESPERADAS (MODELO NORMAL)
#-------------------------------------------------------------------------

k <- length(Fo_abs_108_60)

P_108_60 <- numeric(k)

for (i in 1:k) {
  P_108_60[i] <- pnorm(
    breaks_lon_108_60[i + 1],
    mean = mu_hat,
    sd = sigma_hat
  ) - pnorm(
    breaks_lon_108_60[i],
    mean = mu_hat,
    sd = sigma_hat
  )
}

P_108_60 <- P_108_60 / sum(P_108_60)

Fe_abs_108_60 <- P_108_60 * sum(Fo_abs_108_60)

Fe_abs_108_60
## [1] 147.3760 483.8968 720.1396 486.6754 135.9122
pearson_r_108_60 <- cor(Fo_abs_108_60, Fe_abs_108_60)
pearson_r_108_60
## [1] 0.8910311
pearson_pct_108_60 <- pearson_r_108_60 * 100
pearson_pct_108_60
## [1] 89.10311
limite_grafico_108_60 <- max(c(Fo_abs_108_60, Fe_abs_108_60), na.rm = TRUE)

plot(
  Fo_abs_108_60,
  Fe_abs_108_60,
  xlim = c(0, limite_grafico_108_60),
  ylim = c(0, limite_grafico_108_60),
  main = "Correlación de Pearson Fo vs Fe\nModelo normal – Longitude (−108 a −60)",
  xlab = "Frecuencia observada (Fo)",
  ylab = "Frecuencia esperada (Fe)",
  pch = 19,
  col = "blue3"
)

lines(
  x = c(0, limite_grafico_108_60),
  y = c(0, limite_grafico_108_60),
  col = "red",
  lwd = 2
)

Correlacion_Pearson_108_60 <- pearson_pct_108_60
Correlacion_Pearson_108_60
## [1] 89.10311
# ====================================================
# TEST DE CHI-CUADRADO
# ====================================================

nivel_significancia <- 0.05

Fo_108_60 <- (Fo_abs_108_60 / sum(Fo_abs_108_60)) * 100
Fo_108_60
## [1] 10.688956 14.589666 44.630193 22.289767  7.801418
Fe_108_60 <- P_108_60 * 100
Fe_108_60
## [1]  7.465856 24.513517 36.481234 24.654278  6.885115
x2_108_60 <- sum((Fe_108_60 - Fo_108_60)^2 / Fe_108_60)
x2_108_60
## [1] 7.577925
grados_libertad_108_60 <- length(Fo_108_60) - 1
grados_libertad_108_60
## [1] 4
umbral_aceptacion_108_60 <- qchisq(
  1 - nivel_significancia,
  grados_libertad_108_60
)

umbral_aceptacion_108_60
## [1] 9.487729
Decision_Chi_108_60 <- ifelse(
  x2_108_60 < umbral_aceptacion_108_60,
  "SE ACEPTA el modelo de distribución",
  "SE RECHAZA el modelo de distribución"
)

Decision_Chi_108_60
## [1] "SE ACEPTA el modelo de distribución"
#-------------------------------------------------------------------------
# FRECUENCIAS ESPERADAS (MODELO UNIFORME)
#-------------------------------------------------------------------------

Fo_abs_60_30 <- ni_60_30

P_60_30 <- punif(LS_60_30, min = a_60_30, max = b_60_30) -
  punif(Li_60_30, min = a_60_30, max = b_60_30)

P_60_30 <- P_60_30 / sum(P_60_30)

Fe_abs_60_30 <- P_60_30 * sum(Fo_abs_60_30)
Fe_abs_60_30
## [1] 290.3694 294.8910 294.7396
pearson_r_60_30 <- cor(Fo_abs_60_30, Fe_abs_60_30)
pearson_r_60_30
## [1] 0.9855951
pearson_pct_60_30 <- pearson_r_60_30 * 100
pearson_pct_60_30
## [1] 98.55951
limite_grafico_60_30 <- max(c(Fo_abs_60_30, Fe_abs_60_30), na.rm = TRUE)

plot(
  Fo_abs_60_30,
  Fe_abs_60_30,
  xlim = c(0, limite_grafico_60_30),
  ylim = c(0, limite_grafico_60_30),
  main = "Correlación de Pearson Fo vs Fe\nModelo uniforme – Longitude (−60 a 30)",
  xlab = "Frecuencia observada (Fo)",
  ylab = "Frecuencia esperada (Fe)",
  pch = 19,
  col = "blue3"
)

lines(
  x = c(0, limite_grafico_60_30),
  y = c(0, limite_grafico_60_30),
  col = "red",
  lwd = 2
)

Correlacion_Pearson_60_30 <- pearson_pct_60_30
Correlacion_Pearson_60_30
## [1] 98.55951
# ====================================================
# TEST DE CHI-CUADRADO
# ====================================================

Fo_60_30 <- (Fo_abs_60_30 / sum(Fo_abs_60_30)) * 100
Fe_60_30 <- P_60_30 * 100

x2_60_30 <- sum((Fo_60_30 - Fe_60_30)^2 / Fe_60_30)
x2_60_30
## [1] 2.322392
grados_libertad_60_30 <- length(Fo_60_30) - 1
grados_libertad_60_30
## [1] 2
umbral_aceptacion_60_30 <- qchisq(
  1 - nivel_significancia,
  grados_libertad_60_30
)

umbral_aceptacion_60_30
## [1] 5.991465
Decision_Chi_60_30 <- ifelse(
  x2_60_30 < umbral_aceptacion_60_30,
  "SE ACEPTA el modelo UNIFORME",
  "SE RECHAZA el modelo UNIFORME"
)

Decision_Chi_60_30
## [1] "SE ACEPTA el modelo UNIFORME"
# -------------------------------------------------------------------------
# FRECUENCIAS ESPERADAS (MODELO GAMMA)
# -------------------------------------------------------------------------

P_30_max <- numeric(k_30_max)

for (i in 1:k_30_max) {
  P_30_max[i] <- pgamma(LS_30_max[i] - 30, shape_hat, rate_hat) -
    pgamma(Li_30_max[i] - 30, shape_hat, rate_hat)
}

P_30_max <- P_30_max / sum(P_30_max)
Fe_abs_30_max <- P_30_max * sum(Fo_abs_30_max)

# -------------------------------------------------------------------------
# AJUSTE DE FRECUENCIAS
# -------------------------------------------------------------------------

validos <- which(
  Fe_abs_30_max > 0 &
    Fo_abs_30_max > 0 &
    is.finite(Fe_abs_30_max) &
    is.finite(Fo_abs_30_max)
)

Fo_val <- Fo_abs_30_max[validos]
Fe_val <- Fe_abs_30_max[validos]

# -------------------------------------------------------------------------
# TEST DE PEARSON
# -------------------------------------------------------------------------

pearson_r_30_max <- cor(Fo_val, Fe_val, method = "pearson")
pearson_pct_30_max <- pearson_r_30_max * 100

limite_grafico_30_max <- max(c(Fo_val, Fe_val), na.rm = TRUE)

plot(
  Fo_val,
  Fe_val,
  xlim = c(0, limite_grafico_30_max),
  ylim = c(0, limite_grafico_30_max),
  main = "Correlación de Pearson Fo vs Fe\nModelo Gamma – Longitude (> 30)",
  xlab = "Frecuencia observada (Fo)",
  ylab = "Frecuencia esperada (Fe)",
  pch = 19,
  col = "blue3"
)

lines(
  x = c(0, limite_grafico_30_max),
  y = c(0, limite_grafico_30_max),
  col = "red",
  lwd = 2
)

Correlacion_Pearson_30_max <- pearson_pct_30_max
Correlacion_Pearson_30_max
## [1] 98.35356
# -------------------------------------------------------------------------
# TEST DE CHI-CUADRADO
# -------------------------------------------------------------------------

Fo_30_max <- (Fo_val / sum(Fo_val)) * 100
Fe_30_max <- (Fe_val / sum(Fe_val)) * 100

x2_30_max <- sum((Fo_30_max - Fe_30_max)^2 / Fe_30_max)

k <- length(Fo_30_max)
p <- 2
grados_libertad_30_max <- k - p - 1

umbral_aceptacion_30_max <- qchisq(
  1 - nivel_significancia,
  grados_libertad_30_max
)

Decision_Chi_30_max <- ifelse(
  x2_30_max < umbral_aceptacion_30_max,
  "SE ACEPTA el modelo GAMMA",
  "SE RECHAZA el modelo GAMMA"
)

x2_30_max
## [1] 2.155173
grados_libertad_30_max
## [1] 2
umbral_aceptacion_30_max
## [1] 5.991465
Decision_Chi_30_max
## [1] "SE ACEPTA el modelo GAMMA"

8. CÁLCULO DE PROBABILIDADES

# =========================================================
# GRÁFICA DE PROBABILIDAD – MODELO NORMAL (−108 a −60)
# =========================================================

x_norm_prob <- seq(-108, -60, length.out = 500)

y_norm_prob <- dnorm(x_norm_prob, mean = mu_hat, sd = sigma_hat)

plot(
  x_norm_prob,
  y_norm_prob,
  type = "l",
  lwd = 2,
  col = "blue3",
  main = "Gráfica de probabilidad – Modelo normal\nLongitude (−108 a −60)",
  xlab = "Longitude (°)",
  ylab = "Densidad de probabilidad"
)

x_sec <- seq(-95, -75, by = 0.1)
y_sec <- dnorm(x_sec, mean = mu_hat, sd = sigma_hat)

polygon(
  c(x_sec, rev(x_sec)),
  c(y_sec, rep(0, length(y_sec))),
  col = rgb(1, 0, 0, 0.5),
  border = NA
)

lines(x_sec, y_sec, col = "red", lwd = 2)

legend(
  "topright",
  legend = c("Modelo normal", "Área de probabilidad"),
  col = c("blue3", "red"),
  lwd = 2,
  bty = "n"
)

P_norm <- pnorm(-75, mean = mu_hat, sd = sigma_hat) -
  pnorm(-95, mean = mu_hat, sd = sigma_hat)

P_norm_pct <- P_norm * 100

P_norm
## [1] 0.6342737
P_norm_pct
## [1] 63.42737
# =========================================================
# GRÁFICA DE PROBABILIDAD – MODELO UNIFORME (−60 a 30)
# =========================================================

a_prob_60_30 <- -60
b_prob_60_30 <- 30

x_uni <- seq(a_prob_60_30, b_prob_60_30, length.out = 500)
y_uni <- dunif(x_uni, min = a_prob_60_30, max = b_prob_60_30)

plot(
  x_uni,
  y_uni,
  type = "l",
  lwd = 2,
  col = "blue3",
  main = "Gráfica de probabilidad – Modelo uniforme\nLongitude (−60 a 30)",
  xlab = "Longitude (°)",
  ylab = "Densidad de probabilidad"
)

x_sec <- seq(-30, 0, by = 0.1)
y_sec <- dunif(x_sec, min = a_prob_60_30, max = b_prob_60_30)

polygon(
  c(x_sec, rev(x_sec)),
  c(y_sec, rep(0, length(y_sec))),
  col = rgb(1, 0, 0, 0.5),
  border = NA
)

lines(x_sec, y_sec, col = "red", lwd = 2)

legend(
  "topright",
  legend = c("Modelo uniforme", "Área de probabilidad"),
  col = c("blue3", "red"),
  lwd = 2,
  bty = "n"
)

P_uni <- punif(0, min = a_prob_60_30, max = b_prob_60_30) -
  punif(-30, min = a_prob_60_30, max = b_prob_60_30)

P_uni_pct <- P_uni * 100

P_uni
## [1] 0.3333333
P_uni_pct
## [1] 33.33333
# =========================================================
# GRÁFICA DE PROBABILIDAD – MODELO GAMMA (Longitude > 30)
# =========================================================

x_gamma_prob <- seq(0, max(lon_30_max) - 30, length.out = 500)

y_gamma_prob <- dgamma(
  x_gamma_prob,
  shape = shape_hat,
  rate = rate_hat
)

x_plot <- x_gamma_prob + 30

plot(
  x_plot,
  y_gamma_prob,
  type = "l",
  lwd = 2,
  col = "blue3",
  main = "Gráfica de probabilidad – Modelo Gamma\nLongitude (> 30)",
  xlab = "Longitude (°)",
  ylab = "Densidad de probabilidad"
)

x_sec <- seq(60, 120, by = 0.5)
x_sec_gamma <- x_sec - 30

y_sec <- dgamma(
  x_sec_gamma,
  shape = shape_hat,
  rate = rate_hat
)

polygon(
  c(x_sec, rev(x_sec)),
  c(y_sec, rep(0, length(y_sec))),
  col = rgb(1, 0, 0, 0.5),
  border = NA
)

lines(x_sec, y_sec, col = "red", lwd = 2)

legend(
  "topright",
  legend = c("Modelo Gamma", "Área de probabilidad"),
  col = c("blue3", "red"),
  lwd = 2,
  bty = "n"
)

x1 <- 60 - 30
x2 <- 120 - 30

P_gamma <- pgamma(x2, shape = shape_hat, rate = rate_hat) -
  pgamma(x1, shape = shape_hat, rate = rate_hat)

P_gamma_pct <- P_gamma * 100

P_gamma
## [1] 0.7518324
P_gamma_pct
## [1] 75.18324
Resumen_Probabilidades <- data.frame(
  Agrupacion = c("Normal (-108 a -60)", "Uniforme (-60 a 30)", "Gamma (> 30)"),
  Intervalo  = c("[-95, -75]", "[-30, 0]", "[60, 120]"),
  Probabilidad = c(P_norm, P_uni, P_gamma),
  Porcentaje   = c(P_norm_pct, P_uni_pct, P_gamma_pct)
)

Resumen_Probabilidades
##            Agrupacion  Intervalo Probabilidad Porcentaje
## 1 Normal (-108 a -60) [-95, -75]    0.6342737   63.42737
## 2 Uniforme (-60 a 30)   [-30, 0]    0.3333333   33.33333
## 3        Gamma (> 30)  [60, 120]    0.7518324   75.18324
# ====================================================
# TABLA RESUMEN DE PROBABILIDADES POR AGRUPACIÓN
# ====================================================

Resumen_Probabilidades <- data.frame(
  Agrupación = c(
    "Normal (−108 a −60)",
    "Uniforme (−60 a 30)",
    "Gamma (> 30)"
  ),
  Intervalo = c(
    "[−95, −75]",
    "[−30, 0]",
    "[60, 120]"
  ),
  Probabilidad = c(
    P_norm,
    P_uni,
    P_gamma
  ),
  Porcentaje = c(
    P_norm_pct,
    P_uni_pct,
    P_gamma_pct
  )
)

tabla_gt_prob <- Resumen_Probabilidades %>%
  gt() %>%
  tab_header(
    title = md("**Tabla N° X**"),
    subtitle = md("Resumen de probabilidades por agrupación y modelo de distribución")
  ) %>%
  fmt_number(
    columns = c(Probabilidad),
    decimals = 4
  ) %>%
  fmt_number(
    columns = c(Porcentaje),
    decimals = 2
  ) %>%
  cols_label(
    Agrupación  = "Agrupación / Modelo",
    Intervalo   = "Intervalo analizado (°)",
    Probabilidad = "Probabilidad",
    Porcentaje   = "Probabilidad (%)"
  ) %>%
  tab_source_note(
    source_note = md("Autor: Grupo Geología")
  )

tabla_gt_prob
Tabla N° X
Resumen de probabilidades por agrupación y modelo de distribución
Agrupación / Modelo Intervalo analizado (°) Probabilidad Probabilidad (%)
Normal (−108 a −60) [−95, −75] 0.6343 63.43
Uniforme (−60 a 30) [−30, 0] 0.3333 33.33
Gamma (> 30) [60, 120] 0.7518 75.18
Autor: Grupo Geología

9. CONCLUSIÓN

Conclusión

El análisis inferencial de la variable Longitude permitió organizar la información en tres agrupaciones principales. Para cada una se generó una tabla de frecuencias, se propuso un modelo de probabilidad, se representó la curva sobre el histograma, se aplicaron pruebas de ajuste mediante Pearson y chi-cuadrado, y finalmente se calcularon probabilidades para intervalos específicos.

FIN DEL MODELO INFERENCIAL