ANÁLISIS ESTADÍSTICO DE VOLCANES ACTIVOS A NIVEL GLOBAL
MODELO EXPONENCIAL

CARRERA DE GEOLOGÍA

GRUPO N°2

ANÁLISIS ESTADÍSTICO INFERENCIAL

CARGA DE LIBRERIAS Y DATOS

Carga de librerias

## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.6'
## (as 'lib' is unspecified)
## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.6'
## (as 'lib' is unspecified)
## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.6'
## (as 'lib' is unspecified)
## Installing package into '/cloud/lib/x86_64-pc-linux-gnu-library/4.6'
## (as 'lib' is unspecified)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

lectura del dataset

#CARGA DE DATASET
Volcanes_Globales <- read.csv("global_volcano_eruption_intelligence.csv", header = T, sep = ";", dec = ".")

SELECCIÓN DE VARIABLE

retornoV<- Volcanes_Globales$avg_eruption_return_period_years

# Eliminar NA
retorno <- na.omit(retornoV)

TABLA DE DISTRIBUCIÓN DE FRECUENCIAS

## [1] 898
## [1] 8
## [1] 200
# Frecuencia absoluta
ni <- numeric(k)

for (i in 1:k) {
  
  if (i == 1) {
    
    # Primer intervalo: [0, 200]
    ni[i] <- sum(
      retorno >= Li[i] &
        retorno <= Ls[i],
      na.rm = TRUE
    )
    
  } else {
    
    # Resto de intervalos: (Li, Ls]
    ni[i] <- sum(
      retorno > Li[i] &
        retorno <= Ls[i],
      na.rm = TRUE
    )
    
  }
}

# Frecuencia relativa
hi <- ni / n

# ==============================================================================
# TABLA DE DISTRIBUCIÓN DE FRECUENCIAS
# ==============================================================================

TDFretorno <- data.frame(
  Intervalo = etiquetas,
  ni = ni,
  hi = round(hi * 100, 2)
)

# Agregar fila TOTAL
TDFretorno <- rbind(
  TDFretorno,
  data.frame(
    Intervalo = "TOTAL",
    ni = sum(ni),
    hi = 100
  )
)

# Mostrar tabla
TDFretorno
##       Intervalo  ni     hi
## 1     [0 - 200] 626  69.71
## 2   (200 - 400] 176  19.60
## 3   (400 - 600]  46   5.12
## 4   (600 - 800]  24   2.67
## 5  (800 - 1000]  13   1.45
## 6 (1000 - 1200]   6   0.67
## 7 (1200 - 1400]   6   0.67
## 8 (1400 - 1600]   1   0.11
## 9         TOTAL 898 100.00

Tabla de frecuencias

Distribución de Frecuencias del Período de Retorno de los Volcanes
Análisis de frecuencias mediante ocho intervalos de amplitud de 200 años
Intervalo del Período de Retorno (años) Frecuencia Absoluta (ni) Frecuencia Relativa (hi) (%)
[0 - 200] 626 69.71
(200 - 400] 176 19.60
(400 - 600] 46 5.12
(600 - 800] 24 2.67
(800 - 1000] 13 1.45
(1000 - 1200] 6 0.67
(1200 - 1400] 6 0.67
(1400 - 1600] 1 0.11
TOTAL 898 100.00

GRÁFICA DE DISTRIBUCIÓN DE FRECUENCIA

Diagrama de frecuencia absoluta

hist(
  retorno,
  breaks = breaks_retorno,
  freq = TRUE,
  xaxt = "n",
  col = "lightblue",
  border = "black",
  main = "Histograma de Frecuencias Absolutas",
  xlab = "Período de Retorno (años)",
  ylab = "Cantidad"
)

axis(
  side = 1,
  at = breaks_retorno,
  labels = round(breaks_retorno, 0),
  las = 2
)

CONJETURA DEL MODELO

DEBIDO A LA SIMILITUD VISUAL SE CONJETURA UN MODELO EXPONENCIAL

PARAMETROS

Calculo de parametros del modelo

#PARAMETROS
# Número total de volcanes
n <- sum(TDFretorno$ni[TDFretorno$Intervalo != "TOTAL"])

# Datos agrupados (frecuencias absolutas)
x <- TDFretorno$ni[TDFretorno$Intervalo != "TOTAL"]

# Marcas de clase
marcas_clase <- (Li + Ls) / 2

# Media observada mediante datos agrupados
media_observada <- sum(marcas_clase * x) / n

# Parámetro de la distribución Exponencial LAMBDA
lambda <- 1 / media_observada

# Mostrar resultados
media_observada
## [1] 203.5635
lambda
## [1] 0.004912473

MODELO Y REALIDAD

Frecuencias observadas y esperadas

#==============================================================================
# FRECUENCIAS OBSERVADAS
# ==============================================================================

Fo <- TDFretorno$hi[TDFretorno$Intervalo != "TOTAL"]


# ==============================================================================
# FRECUENCIAS ESPERADAS - MODELO EXPONENCIAL
# ==============================================================================

P_exp <- numeric(length(Li))

for(i in 1:length(Li)){
  
  P_exp[i] <- pexp(
    Ls[i],
    rate = lambda
  ) -
    pexp(
      Li[i],
      rate = lambda
    )
}

# Frecuencia esperada porcentual
Fe <- P_exp * 100


# Mostrar resultados
Fo
## [1] 69.71 19.60  5.12  2.67  1.45  0.67  0.67  0.11
Fe
## [1] 62.56239590 23.42186209  8.76858400  3.28274776  1.22898211  0.46010146
## [7]  0.17225096  0.06448663

Grafica modelo y realidad

# GRÁFICA
par(oma = c(1, 1, 1, 1))

hist(
  retorno,
  freq = FALSE,
  breaks = breaks_retorno,
  main = "Gráfica Nº3. Comparación de la realidad con el modelo Exponencial del período de retorno volcánico",
  xlab = "Período de retorno (años)",
  ylab = "Densidad de probabilidad",
  col = "cyan",
  border = "gray30"
)

box(which = "outer", col = "black")


# Histograma sin graficar para obtener frecuencias
histograma <- hist(
  retorno,
  plot = FALSE,
  breaks = breaks_retorno
)

h <- length(histograma$counts)
h
## [1] 8
# Curva del modelo exponencial
curve(
  dexp(
    x,
    rate = lambda
  ),
  from = min(retorno),
  to = max(retorno),
  add = TRUE,
  col = "red",
  lwd = 3
)

TESTS DE APROBACIÓN

Gráfica de correlación del modelo Exponencial y la realidad

# Frecuencias observadas (%)
Fo <- TDFretorno$hi[TDFretorno$Intervalo != "TOTAL"]

# Frecuencias esperadas (%)
Fe <- P_exp * 100


# Gráfica
plot(
  Fo,
  Fe,
  main = "Gráfica N°4: Correlación de frecuencias\nentre la realidad y el modelo exponencial del período de retorno",
  xlab = "Frecuencia observada (%)",
  ylab = "Frecuencia esperada (%)",
  pch = 19,
  col = "darkblue",
  xlim = c(0, max(c(Fo, Fe)) * 1.1),
  ylim = c(0, max(c(Fo, Fe)) * 1.1)
)

# Línea de ajuste perfecto
abline(
  a = 0,
  b = 1,
  col = "red",
  lwd = 2
)

Test de Pearson

# Correlación entre realidad y modelo Exponencial
Correlacion <- cor(Fo, Fe) * 100
Correlacion
## [1] 99.40248
# Evaluación del Test de Pearson (>80%)
if(Correlacion > 80){
  
  print("APRUEBA EL TEST DE PEARSON")
  
}else{
  
  print("NO APRUEBA EL TEST DE PEARSON")
  
}
## [1] "APRUEBA EL TEST DE PEARSON"

Test Chi cuadrado

# Frecuencias absolutas observadas
Fo_abs <- TDFretorno$ni[TDFretorno$Intervalo != "TOTAL"]


# Tamaño de la muestra
n <- sum(Fo_abs)
n
## [1] 898
# Frecuencias absolutas esperadas
Fe_abs <- P_exp * n


# Grados de libertad
# gl = número de intervalos - parámetros estimados - 1
# La distribución exponencial tiene 1 parámetro: lambda

grados_libertad <- length(Fo_abs) - 2
grados_libertad
## [1] 6
# Nivel de confianza
nivel_confianza <- 0.9999999


# Estadístico Chi-cuadrado

x2 <- sum(
  (Fo_abs - Fe_abs)^2 / Fe_abs
)

x2
## [1] 41.89031
# Valor crítico Chi-cuadrado

umbral_aceptacion <- qchisq(
  nivel_confianza,
  df = grados_libertad
)

umbral_aceptacion
## [1] 43.33776
# Evaluación del modelo

if (x2 < umbral_aceptacion) {
  
  print("APRUEBA EL TEST DE CHI-CUADRADO")
  
} else {
  
  print("NO APRUEBA EL TEST DE CHI-CUADRADO")
  
}
## [1] "APRUEBA EL TEST DE CHI-CUADRADO"

Tabla de resumen

Evaluación estadística del modelo exponencial del período de retorno volcánico
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Período de retorno (años) 99.4 41.89 43.34

CALCULO DE PROBABILIDAD

Pregunta

Calculo de valores

## Trabajamos con todo nuestro dataset para el cálculo de parámetros

retorno_completo <- as.numeric(retorno)

retorno_completo <- na.omit(retorno_completo)

# Tamaño de la población

n_completo <- length(retorno_completo)

n_completo
## [1] 898

Media aritmetica

media_completa <- mean(retorno_completo)

media_completa
## [1] 189.8882
# ==============================================================================
# Parámetro lambda del modelo exponencial
# ==============================================================================

lambda_completo <- 1 / media_completa

lambda_completo
## [1] 0.005266256

Calculo de probabilidad

limite_inferior <- 400
limite_superior <- 800


probabilidad_retorno <- pexp(
  limite_superior,
  rate = lambda_completo
) -
  pexp(
    limite_inferior,
    rate = lambda_completo
  )


# Probabilidad en porcentaje

probabilidad_retorno * 100
## [1] 10.68609

Area de probabilidad

Probabilidad 2

Calculo de probabilidad 2

total_volcanes <- 85
volcanes_esperados <- probabilidad_retorno * total_volcanes
volcanes_esperados
## [1] 9.08318

INTERVALOS DE CONFIANZA

Cálculo de parámetros e intervalos

# Datos completos
retorno_completo <- as.numeric(retorno_completo)
retorno_completo <- na.omit(retorno_completo)

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

# Media muestral (solo usada para calcular el intervalo)
media_muestral <- mean(retorno_completo)
media_muestral
## [1] 189.8882
# Desviación estándar muestral
desviacion <- sd(retorno_completo)
desviacion
## [1] 209.0615
# Error estándar
e <- desviacion / sqrt(n)
e
## [1] 6.976474
# Límites del intervalo de confianza aproximado 95%
li <- media_muestral - 2 * e
ls <- media_muestral + 2 * e

Tabla de resumen

Tabla Nro.3: Intervalo de confianza de la media poblacional del período de retorno volcánico
Límite inferior Media poblacional Límite superior Error estándar
175.94 μ 203.84 6.98

CONCLUSIÓN

La variable período de retorno de los volcanes activos se explica mediante un modelo probabilístico Exponencial, donde la media aritmética es de 189.89, con una desviación estándar de 209.06, evidenciando la variabilidad temporal existente entre los registros analizados. De esta manera, es posible realizar estimaciones del comportamiento de la variable, como el rango esperado del período de retorno promedio.

Mediante el teorema del límite central, se establece que la media aritmética poblacional del período de retorno volcánico se encuentra entre 175.94 y 203.84 con un 95% de confianza.

.