#-----------------------------------------------------------------------------------------#
# 1. CARGA DE LIBRERÍAS
#-----------------------------------------------------------------------------------------#
library(readxl)
library(knitr)
library(fitdistrplus)
## Cargando paquete requerido: MASS
## Cargando paquete requerido: survival
library(plotly) # Para el gráfico interactivo final
## Cargando paquete requerido: ggplot2
## 
## Adjuntando el paquete: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:MASS':
## 
##     select
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
# Nota: Ajusta la ruta a tu archivo local
Datos <- read_xlsx("C:/Users/LEO/Documents/Antisana/weatherdataANTISANA.csv.xlsx")

#-----------------------------------------------------------------------------------------#
# 2. TU CÓDIGO LARGO (SIN CAMBIOS)
#-----------------------------------------------------------------------------------------#

Max_Temperatura <- Datos$Max_Temperature


Max_Temperatura <- as.numeric(Max_Temperatura)
# Ajustamos el filtro para que sea de 8 a 24 como pediste
Max_Temperatura <- Max_Temperatura[!is.na(Max_Temperatura) & Max_Temperatura >= 8 & Max_Temperatura <= 24]


histograma <- hist(Max_Temperatura, freq = FALSE, main="Gráfica 103.Modelo de probabilidad Normal",
                   xlab="Grados Celsius", ylab="Densidad de probabilidad", col="blue")

# Ajuste del modelo al histograma
histograma <- hist(Max_Temperatura, breaks = seq(9, 24, by = 3), freq = FALSE, main="Gráfica 104 .Modelo de probabilidad Normal 
de Max Temperatura",
                   xlab="Grados Celsius", ylab="Densidad de probabilidad", col="blue")
h <- length(histograma$counts)

# Usamos media y desviación estándar para el modelo Normal
u_norm <- mean(Max_Temperatura)
sigma_norm <- sd(Max_Temperatura)

# Usamos length.out=100 para que la curva sea suave independientemente del tamaño de los datos
x <- seq(min(Max_Temperatura), max(Max_Temperatura), length.out = 100)
curve(dnorm(x, u_norm, sigma_norm), type="l", add=TRUE, lwd=4, col="black")

# Frecuencia simple observada
Fo <- histograma$counts

# Frecuencia simple esperada
P <- c(0)
for (i in 1:h) {
  P[i] <- (pnorm(histograma$breaks[i+1], u_norm, sigma_norm) -
             pnorm(histograma$breaks[i], u_norm, sigma_norm))
}
Fe <- P * length(Max_Temperatura)

# Comparar tamaño real y el modelo
sum(Fe)
## [1] 361.8433
n <- length(Max_Temperatura)

# 2. Test de Pearson
# Grado de correlación entre Fe Y Fo

plot(Fo, Fe, main="Gráfica 105: Correlación de frecuencias en el modelo normal de Max Temperatura",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)

# Dibujar la línea de tendencia
abline(lm(Fe ~ Fo), col="red", lwd=2)

Correlación <- cor(Fo, Fe) * 100
Correlación
## [1] 96.59274
#TEST DE CHI-CUADRADO
#Creación de grados de libertad y nivel de significancia

grados_libertad <- (length(histograma$counts)-1)
grados_libertad
## [1] 4
nivel_significancia <- 0.05

# Usamos las frecuencias absolutas para el test
x2 <- sum((Fe - Fo)^2 / Fe)
x2
## [1] 8.082922
#Calcular el umbral de aceptación

umbral_aceptacion <- qchisq(1 - nivel_significancia, grados_libertad)
umbral_aceptacion
## [1] 9.487729
x2 < umbral_aceptacion
## [1] TRUE
#3 Tabla resumen de test

Variable <- c("Max_Temperatura")
tabla_resumen <- data.frame(Variable, round(Correlación, 2), round(x2, 2), round(umbral_aceptacion, 2))
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral de aceptación")
library(knitr)

kable(tabla_resumen, format = "markdown", caption = "Tabla.Resumen de test de bondad al modelo de probabilidad")
Tabla.Resumen de test de bondad al modelo de probabilidad
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Max_Temperatura 96.59 8.08 9.49
# 4. Cálculo de probabilidades corregido
# Usaremos un rango de interés (ejemplo: 12 a 15 grados)
limite_inf <- 12 
limite_sup <- 15

prob_area <- pnorm(limite_sup, u_norm, sigma_norm) - pnorm(limite_inf, u_norm, sigma_norm)

# Ajustamos el eje X para que cubra el rango de 8 a 24
x_grafica <- seq(8, 24, length.out = 1000)
y_grafica <- dnorm(x_grafica, u_norm, sigma_norm)

# GRAFICAR
plot(x_grafica, y_grafica, type = "l", col = "skyblue3", lwd = 2, 
     main="Distribución de Máxima Temperatura",
     ylab="Densidad de Probabilidad", 
     xlab="Máxima Temperatura (°C)", 
     xaxt="n")

# Sombreado de la probabilidad
x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dnorm(x_sombra, u_norm, sigma_norm)
polygon(c(x_sombra, rev(x_sombra)), c(y_sombra, rep(0, length(y_sombra))), 
        col = rgb(1, 0, 0, 0.4), border = "red")

# Eje X con saltos de 2 en 2 grados
axis(1, at = seq(8, 24, by = 2), las = 2)

# Leyenda explicativa
legend("topright", 
       legend = c("Modelo Normal", paste("Probabilidad:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

#-----------------------------------------------------------------------------------------#

Min_Temperature <- Datos$Min_Temperature

Min_Temperature <- as.numeric(Min_Temperature)
# Ajustamos el filtro: Las mínimas suelen estar entre 0 y 12 grados
Min_Temperature <- Min_Temperature[!is.na(Min_Temperature) & Min_Temperature >= 0 & Min_Temperature <= 15]


# Agrupaciones de 2 en 2 desde 0 hasta 12
mis_breaks <- seq(0, 12, by = 2)

histograma <- hist(Min_Temperature, freq = FALSE, main="Gráfica 103.Modelo de probabilidad Normal",
                   xlab="Grados Celsius", ylab="Densidad de probabilidad", col="blue")

# Ajuste del modelo al histograma
histograma <- hist(Min_Temperature, breaks = mis_breaks, freq = FALSE, main="Gráfica 104 .Modelo de probabilidad Normal 
de Mínima Temperatura",
                   xlab="Grados Celsius", ylab="Densidad de probabilidad", col="blue")
h <- length(histograma$counts)

# Usamos media y desviación estándar para el modelo Normal
u_norm <- mean(Min_Temperature)
sigma_norm <- sd(Min_Temperature)

# Curva suave
x <- seq(min(Min_Temperature), max(Min_Temperature), length.out = 100)
curve(dnorm(x, u_norm, sigma_norm), type="l", add=TRUE, lwd=4, col="black")

# Frecuencia simple observada
Fo <- histograma$counts

# Frecuencia simple esperada
P <- c(0)
for (i in 1:h) {
  P[i] <- (pnorm(histograma$breaks[i+1], u_norm, sigma_norm) -
             pnorm(histograma$breaks[i], u_norm, sigma_norm))
}
Fe <- P * length(Min_Temperature)

n <- length(Min_Temperature)

# 2. Test de Pearson
plot(Fo, Fe, main="Gráfica 105: Correlación de frecuencias en el modelo normal de Mínima Temperatura",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)

abline(lm(Fe ~ Fo), col="red", lwd=2)

Correlación <- cor(Fo, Fe) * 100

# TEST DE CHI-CUADRADO
grados_libertad <- (length(histograma$counts)-1)
nivel_significancia <- 0.05

x2 <- sum((Fe - Fo)^2 / Fe)

# Umbral
umbral_aceptacion <- qchisq(1 - nivel_significancia, grados_libertad)

# ¿Pasa la prueba?
# x2 < umbral_aceptacion

# 3. Tabla resumen
Variable <- c("Min_Temperature")
tabla_resumen <- data.frame(Variable, round(Correlación, 2), round(x2, 2), round(umbral_aceptacion, 2))
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral de aceptación")


kable(tabla_resumen, format = "markdown", caption = "Tabla.Resumen de test para Mínima Temperatura")
Tabla.Resumen de test para Mínima Temperatura
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Min_Temperature 99.94 10.11 11.07
# 4. Gráfica con área sombreada (Rango de 4 a 8 grados para la mínima)
limite_inf <- 4 
limite_sup <- 8

prob_area <- pnorm(limite_sup, u_norm, sigma_norm) - pnorm(limite_inf, u_norm, sigma_norm)

x_grafica <- seq(0, 15, length.out = 1000)
y_grafica <- dnorm(x_grafica, u_norm, sigma_norm)

plot(x_grafica, y_grafica, type = "l", col = "skyblue3", lwd = 2, 
     main="Distribución de Mínima Temperatura",
     ylab="Densidad de Probabilidad", 
     xlab="Mínima Temperatura (°C)", 
     xaxt="n")

x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dnorm(x_sombra, u_norm, sigma_norm)
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(0, 15, by = 1), las = 2)

legend("topright", 
       legend = c("Modelo Normal", paste("Probabilidad:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

#----------------------------------------------------------------------------------------#

Prec_Total <- as.numeric(Datos$Precipitation)
Prec_Total <- Prec_Total[!is.na(Prec_Total)]


hist(Prec_Total, freq = FALSE, main="Gráfica 106. Distribución Total de Precipitación",
     xlab="mm", ylab="Densidad", col="blue")

Precipitacion <- Prec_Total[Prec_Total >= 30 & Prec_Total <= 80]


# Definimos los cortes de 30 a 86 de 8 en 8
cortes <- seq(30, 86, by = 8)

# GRÁFICA 107: Agrupación manual de 8 en 8
histograma <- hist(Precipitacion, breaks = cortes, freq = FALSE, 
                   main="Gráfica 107. Modelo de probabilidad GAMMA (Agrupación 8 en 8)",
                   xlab="mm", ylab="Densidad de probabilidad", col="blue",
                   xlim=c(30, 80))

# GRÁFICA 108: Ajuste del modelo al histograma (GAMMA)
histograma <- hist(Precipitacion, breaks = cortes, freq = FALSE, 
                   main="Gráfica 108. Ajuste de Densidad GAMMA (30-80)",
                   xlab="mm", ylab="Densidad de probabilidad", col="blue",
                   xlim=c(30, 80))
h <- length(histograma$counts)


fit_gam <- fitdist(Precipitacion, "gamma")
shape_p <- fit_gam$estimate["shape"]
rate_p <- fit_gam$estimate["rate"]

# Curva suave GAMMA
x_seq <- seq(30, 80, length.out = 500)
lines(x_seq, dgamma(x_seq, shape = shape_p, rate = rate_p), lwd=4, col="black")

# Frecuencia simple observada (counts en los intervalos de 8mm)
Fo <- histograma$counts

# Frecuencia simple esperada usando la función pgamma
P <- c(0)
for (i in 1:h) {
  P[i] <- (pgamma(histograma$breaks[i+1], shape = shape_p, rate = rate_p) -
             pgamma(histograma$breaks[i], shape = shape_p, rate = rate_p))
}
Fe <- P * length(Precipitacion)

# Comparar tamaño real y el modelo
sum(Fe)
## [1] 60.9618
n <- length(Precipitacion)

# 2. Test de Pearson (Correlación)
plot(Fo, Fe, main="Gráfica 109: Correlación de frecuencias (Modelo Gamma)",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)

abline(lm(Fe ~ Fo), col="red", lwd=2)

Correlación <- cor(Fo, Fe) * 100
Correlación
## [1] 93.03777
# TEST DE CHI-CUADRADO
grados_libertad <- (length(histograma$counts)-1)
nivel_significancia <- 0.05

x2 <- sum((Fe - Fo)^2 / Fe)
umbral_aceptacion <- qchisq(1 - nivel_significancia, grados_libertad)

# ¿Se acepta el modelo?
x2 < umbral_aceptacion
## [1] TRUE
# 3. Tabla resumen de test
Variable <- c("Precipitation (Gamma 30-80)")
tabla_resumen <- data.frame(Variable, round(Correlación, 2), round(x2, 2), round(umbral_aceptacion, 2))
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral de aceptación")

kable(tabla_resumen, format = "markdown", caption = "Tabla. Resumen de validación para Precipitación (GAMMA)")
Tabla. Resumen de validación para Precipitación (GAMMA)
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Precipitation (Gamma 30-80) 93.04 7.09 12.59
# 4. Cálculo de probabilidad
# Ejemplo: Probabilidad de lluvia mayor a 40 mm
limite_inf <- 40 
limite_sup <- 80

prob_area <- pgamma(limite_sup, shape_p, rate_p) - pgamma(limite_inf, shape_p, rate_p)

x_grafica <- seq(30, 80, length.out = 1000)
y_grafica <- dgamma(x_grafica, shape_p, rate_p)

# GRAFICAR MODELO FINAL
plot(x_grafica, y_grafica, type = "l", col = "skyblue3", lwd = 2, 
     main="Distribución GAMMA Final (Rango 30-80)",
     ylab="Densidad de Probabilidad", 
     xlab="mm",
     xlim=c(30, 80),
     xaxt="n")

# Sombreado de la probabilidad
x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dgamma(x_sombra, shape_p, rate_p)
polygon(c(x_sombra, rev(x_sombra)), c(y_sombra, rep(0, length(y_sombra))), 
        col = rgb(1, 0, 0, 0.4), border = "red")

# Eje X marcado cada 10mm
axis(1, at = seq(30, 80, by = 10), las = 2)

# Leyenda
legend("topright", 
       legend = c("Modelo Gamma", paste("Probabilidad > 40mm:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

#---------------------------------------------------------------------------------#

# Trabajamos con Relative_Humidity
Hum_Total <- as.numeric(Datos$Relative_Humidity) 
Hum_Total <- Hum_Total[!is.na(Hum_Total)]


hist(Hum_Total, freq = FALSE, main="Distribución Total de Humedad Relativa",
     xlab="%", ylab="Densidad", col="blue")

# Rango ajustado para mejorar el ajuste del modelo
Humedad <- Hum_Total[Hum_Total >= 0.55 & Hum_Total <= 0.8]

# GRÁFICA 107: Sin agrupación manual 
histograma <- hist(Humedad, freq = FALSE, 
                   main="Gráfica 107. Modelo de probabilidad WEIBULL (Humedad)",
                   xlab="%", ylab="Densidad de probabilidad", col="blue")

# GRÁFICA 108: Ajuste del modelo al histograma (WEIBULL)
histograma <- hist(Humedad, freq = FALSE, 
                   main="Gráfica 108. Ajuste de Densidad WEIBULL (Humedad)",
                   xlab="%", ylab="Densidad de probabilidad", col="blue")
h <- length(histograma$counts)

# --- AJUSTE WEIBULL ---
fit_wei <- fitdist(Humedad, "weibull")
shape_h <- fit_wei$estimate["shape"]
scale_h <- fit_wei$estimate["scale"]


# Curva suave WEIBULL
x_seq <- seq(min(Humedad), max(Humedad), length.out = 500)
lines(x_seq, dweibull(x_seq, shape = shape_h, scale = scale_h), lwd=4, col="black")

# Frecuencia simple observada (Counts del histograma automático)
Fo <- histograma$counts

# Frecuencia simple esperada usando la función pweibull
P <- c(0)
for (i in 1:h) {
  P[i] <- (pweibull(histograma$breaks[i+1], shape = shape_h, scale = scale_h) -
             pweibull(histograma$breaks[i], shape = shape_h, scale = scale_h))
}
Fe <- P * length(Humedad)

# Comparar tamaño real y el modelo
sum(Fe)
## [1] 70.68412
n <- length(Humedad)

# 2. Test de Pearson (Correlación)
plot(Fo, Fe, main="Gráfica 109: Correlación de frecuencias (Humedad - Weibull)",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)

abline(lm(Fe ~ Fo), col="red", lwd=2)

Correlación <- cor(Fo, Fe) * 100
Correlación
## [1] 87.82833
# TEST DE CHI-CUADRADO
grados_libertad <- (length(histograma$counts)-1)
nivel_significancia <- 0.05

x2 <- sum((Fe - Fo)^2 / Fe)
umbral_aceptacion <- qchisq(1 - nivel_significancia, grados_libertad)

# ¿Se acepta el modelo?
x2 < umbral_aceptacion
## [1] TRUE
# 3. Tabla resumen de test
Variable <- c("Relative Humidity (Weibull 0.55-0.8)")
tabla_resumen <- data.frame(Variable, round(Correlación, 2), round(x2, 2), round(umbral_aceptacion, 2))
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral de aceptación")

kable(tabla_resumen, format = "markdown", caption = "Tabla. Resumen de validación para Humedad (WEIBULL)")
Tabla. Resumen de validación para Humedad (WEIBULL)
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Relative Humidity (Weibull 0.55-0.8) 87.83 6.88 9.49
# 4. Cálculo de probabilidad
# Rango sugerido dentro de los nuevos límites
limite_inf <- 0.60 
limite_sup <- 0.75

prob_area <- pweibull(limite_sup, shape_h, scale_h) - pweibull(limite_inf, shape_h, scale_h)

x_grafica <- seq(min(Humedad), max(Humedad), length.out = 1000)
y_grafica <- dweibull(x_grafica, shape_h, scale_h)

# GRAFICAR MODELO FINAL
plot(x_grafica, y_grafica, type = "l", col = "skyblue3", lwd = 2, 
     main="Distribución WEIBULL Final (Humedad)",
     ylab="Densidad de Probabilidad", 
     xlab="%")

# Sombreado de la probabilidad
x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dweibull(x_sombra, shape_h, scale_h)
polygon(c(x_sombra, rev(x_sombra)), c(y_sombra, rep(0, length(y_sombra))), 
        col = rgb(1, 0, 0, 0.4), border = "red")

# Leyenda
legend("topleft", 
       legend = c("Modelo Weibull", paste("Prob 60-75%:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

#----------------------------------------------------------------------------------------#

Solar_Total <- as.numeric(Datos$Solar)
Solar_Total <- Solar_Total[!is.na(Solar_Total)]

hist(Solar_Total, freq = FALSE, main="Gráfica 102. Distribución Completa de Datos Originales",
     xlab="MJ/m2", ylab="Densidad", col="blue")

# Filtramos para que R solo trabaje con los datos en este rango específico
Solar <- Solar_Total[Solar_Total >= 5 & Solar_Total <= 20]

# Definimos los cortes de 5 a 20 de 5 en 5 
cortes <- seq(5, 20, by = 3)

# Histograma del segmento 5-20
histograma <- hist(Solar, freq = FALSE, main="Gráfica 103.Modelo de probabilidad Log-Normal (Rango 5-20)",
                   xlab="MJ/m2", ylab="Densidad de probabilidad", col="blue")

# Ajuste del modelo al histograma filtrado CON AGRUPACIÓN de 5 en 5
histograma <- hist(Solar, breaks = cortes, freq = FALSE, main="Gráfica 104 .Modelo de probabilidad Log-Normal 
de Radiación Solar (5-20)",
                   xlab="MJ/m2", ylab="Densidad de probabilidad", col="blue")
h <- length(histograma$counts)

# Estimamos parámetros meanlog y sdlog para el segmento 5-20
fit_ln <- fitdist(Solar, "lnorm")
meanlog_s <- fit_ln$estimate["meanlog"]
sdlog_s <- fit_ln$estimate["sdlog"]

# Curva suave Log-Normal
x <- seq(min(Solar), max(Solar), length.out = 100)
curve(dlnorm(x, meanlog = meanlog_s, sdlog = sdlog_s), type="l", add=TRUE, lwd=4, col="black")

# Frecuencia simple observada (Basada en la agrupación de 5 en 5)
Fo <- histograma$counts

# Frecuencia simple esperada usando plnorm
P <- c(0)
for (i in 1:h) {
  P[i] <- (plnorm(histograma$breaks[i+1], meanlog_s, sdlog_s) -
             plnorm(histograma$breaks[i], meanlog_s, sdlog_s))
}
Fe <- P * length(Solar)

# Comparar tamaño real y el modelo
sum(Fe)
## [1] 186.9093
n <- length(Solar)

# 2. Test de Pearson
plot(Fo, Fe, main="Gráfica 105: Correlación de frecuencias en el modelo Log-Normal (Rango 5-20)",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)

abline(lm(Fe ~ Fo), col="red", lwd=2)

Correlación <- cor(Fo, Fe) * 100
Correlación
## [1] 92.11068
#TEST DE CHI-CUADRADO
# Al tener 3 intervalos, los grados de libertad serán 2
grados_libertad <- (length(histograma$counts)-1)
nivel_significancia <- 0.05

x2 <- sum((Fe - Fo)^2 / Fe)
umbral_aceptacion <- qchisq(1 - nivel_significancia, grados_libertad)

x2 < umbral_aceptacion
## [1] TRUE
#3 Tabla resumen de test
Variable <- c("Solar (5-20)")
tabla_resumen <- data.frame(Variable, round(Correlación, 2), round(x2, 2), round(umbral_aceptacion, 2))
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral de aceptación")

kable(tabla_resumen, format = "markdown", caption = "Tabla.Resumen de test para Radiación Solar (Segmento 5-20)")
Tabla.Resumen de test para Radiación Solar (Segmento 5-20)
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Solar (5-20) 92.11 8.55 9.49
# 4. Cálculo de probabilidades corregido
limite_inf <- 10 
limite_sup <- 15

prob_area <- plnorm(limite_sup, meanlog_s, sdlog_s) - plnorm(limite_inf, meanlog_s, sdlog_s)

x_grafica <- seq(min(Solar), max(Solar), length.out = 1000)
y_grafica <- dlnorm(x_grafica, meanlog_s, sdlog_s)

# GRAFICAR
plot(x_grafica, y_grafica, type = "l", col = "skyblue3", lwd = 2, 
     main="Distribución Log-Normal de Radiación (Filtro 5-20)",
     ylab="Densidad de Probabilidad", 
     xlab="MJ/m2", 
     xaxt="n")

# Sombreado de la probabilidad
x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dlnorm(x_sombra, meanlog_s, sdlog_s)
polygon(c(x_sombra, rev(x_sombra)), c(y_sombra, rep(0, length(y_sombra))), 
        col = rgb(1, 0, 0, 0.4), border = "red")

# Eje X con saltos de 5 en 5 empezando desde 5
axis(1, at = cortes, las = 2)

# Leyenda explicativa
legend("topright", 
       legend = c("Modelo Log-Normal", paste("Probabilidad:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

#---------------------------------------------------------------------------------------------#

Viento_Total <- as.numeric(Datos$Wind) 
Viento_Total <- Viento_Total[!is.na(Viento_Total)]

# --- GRÁFICA: TODOS LOS DATOS ORIGINALES DE VIENTO ---
hist(Viento_Total, freq = FALSE, main="Distribución Total de Velocidad del Viento",
     xlab="m/s", ylab="Densidad", col="blue")

# La distribución Gamma requiere valores estrictamente positivos (> 0)
Viento <- Viento_Total[Viento_Total > 0]
# ----------------------------------------------------------------

# GRÁFICA 107: Sin agrupación manual (Breaks automáticos de R)
histograma <- hist(Viento, freq = FALSE, 
                   main="Gráfica 107. Modelo de probabilidad GAMMA (Viento)",
                   xlab="m/s", ylab="Densidad de probabilidad", col="blue")

# GRÁFICA 108: Ajuste del modelo al histograma (GAMMA)
histograma <- hist(Viento, freq = FALSE, 
                   main="Gráfica 108. Ajuste de Densidad GAMMA (Viento)",
                   xlab="m/s", ylab="Densidad de probabilidad", col="blue")
h <- length(histograma$counts)

# --- AJUSTE GAMMA ---
fit_gam <- fitdist(Viento, "gamma")
shape_v <- fit_gam$estimate["shape"]
rate_v <- fit_gam$estimate["rate"]

# Curva suave GAMMA
x_seq <- seq(min(Viento), max(Viento), length.out = 500)
lines(x_seq, dgamma(x_seq, shape = shape_v, rate = rate_v), lwd=4, col="black")

# Frecuencia simple observada
Fo <- histograma$counts

# Frecuencia simple esperada usando la función pgamma
P <- c(0)
for (i in 1:h) {
  P[i] <- (pgamma(histograma$breaks[i+1], shape = shape_v, rate = rate_v) -
             pgamma(histograma$breaks[i], shape = shape_v, rate = rate_v))
}
Fe <- P * length(Viento)

# Comparar tamaño real y el modelo
sum(Fe)
## [1] 361.8285
n <- length(Viento)

# 2. Test de Pearson (Correlación)
plot(Fo, Fe, main="Gráfica 109: Correlación de frecuencias (Viento - Gamma)",
     xlab="Frecuencia Observada", ylab="Frecuencia esperada", col="blue3", pch=19)
abline(lm(Fe ~ Fo), col="red", lwd=2)

Correlación <- cor(Fo, Fe) * 100

# TEST DE CHI-CUADRADO
grados_libertad <- (length(histograma$counts)-1)
nivel_significancia <- 0.05
x2 <- sum((Fe - Fo)^2 / Fe)
umbral_aceptacion <- qchisq(1 - nivel_significancia, grados_libertad)


x2 < umbral_aceptacion
## [1] TRUE
# 3. Tabla resumen de test
Variable <- c("Wind Speed (Gamma)")
tabla_resumen <- data.frame(Variable, round(Correlación, 2), round(x2, 2), round(umbral_aceptacion, 2))
colnames(tabla_resumen) <- c("Variable", "Test Pearson (%)", "Chi Cuadrado", "Umbral de aceptación")

kable(tabla_resumen, format = "markdown", caption = "Tabla. Resumen de validación para Viento (GAMMA)")
Tabla. Resumen de validación para Viento (GAMMA)
Variable Test Pearson (%) Chi Cuadrado Umbral de aceptación
Wind Speed (Gamma) 97.49 19.39 21.03
# 4. Cálculo de probabilidad
limite_inf <- 2.0 
limite_sup <- 3.0

prob_area <- pgamma(limite_sup, shape_v, rate_v) - pgamma(limite_inf, shape_v, rate_v)

x_grafica <- seq(min(Viento), max(Viento), length.out = 1000)
y_grafica <- dgamma(x_grafica, shape_v, rate_v)

# GRAFICAR MODELO FINAL
plot(x_grafica, y_grafica, type = "l", col = "skyblue3", lwd = 2, 
     main="Distribución GAMMA Final (Viento)",
     ylab="Densidad de Probabilidad", 
     xlab="m/s")

# Sombreado de la probabilidad (Área entre 2.0 y 3.0)
x_sombra <- seq(limite_inf, limite_sup, length.out = 100)
y_sombra <- dgamma(x_sombra, shape_v, rate_v)
polygon(c(x_sombra, rev(x_sombra)), c(y_sombra, rep(0, length(y_sombra))), 
        col = rgb(1, 0, 0, 0.4), border = "red")

# Leyenda corregida
legend("topright", 
       legend = c("Modelo Gamma", paste("Probabilidad 2-3 m/s:", round(prob_area*100, 2), "%")), 
       fill = c(NA, rgb(1, 0, 0, 0.4)), border = c("skyblue3", "red"), bty = "n")

#-----------------------------------------------------------------------------------------#
# 3. REGRESIÓN LINEAL MÚLTIPLE E INTERACTIVIDAD 3D (CÓDIGO NUEVO)
#-----------------------------------------------------------------------------------------#

# Limpieza de datos para la regresión
df_reg <- data.frame(
  Solar = as.numeric(Datos$Solar),
  Min_Temp = as.numeric(Datos$Min_Temperature),
  Max_Temp = as.numeric(Datos$Max_Temperature),
  Wind = as.numeric(Datos$Wind)
)
df_clean <- na.omit(df_reg)

# Ajuste del Modelo estadístico
modelo_multiple <- lm(Solar ~ Min_Temp + Max_Temp + Wind, data = df_clean)

# --- CREACIÓN DEL GRÁFICO INTERACTIVO CON PLOTLY ---
# Este gráfico permite rotar, hacer zoom y ver los ejes desde cualquier perspectiva
fig <- plot_ly(df_clean, 
               x = ~Max_Temp, 
               y = ~Min_Temp, 
               z = ~Solar,
               type = "scatter3d", 
               mode = "markers",
               # Configuración de los puntos y la LEYENDA de colores
               marker = list(size = 4, 
                             color = ~Wind,          # El color representa el VIENTO
                             colorscale = 'Viridis', # Escala de colores (puedes usar 'Reds' o 'Bluered')
                             showscale = TRUE,       # MUESTRA LA LEYENDA LATERAL
                             colorbar = list(title = "Viento (m/s)")),
               text = ~paste("Viento:", Wind, "m/s")) %>%
  layout(
    title = "Regresión Múltiple Interactiva 3D: Radiación Solar en Antisana",
    scene = list(
      xaxis = list(title = "Temp Máxima (°C)"),
      yaxis = list(title = "Temp Mínima (°C)"),
      zaxis = list(title = "Radiación Solar (MJ/m2)")
    )
  )

# MOSTRAR EL GRÁFICO (Se verá en la pestaña 'Viewer' de RStudio)
fig
# --- TABLA RESUMEN FINAL DE LA REGRESIÓN ---
resumen_reg <- summary(modelo_multiple)
Variable_Reg <- c("(Intersección)", "Mínima Temperatura", "Máxima Temperatura", "Velocidad Viento")
tabla_resumen_reg <- data.frame(Variable_Reg, 
                                round(resumen_reg$coefficients[, 1], 4), 
                                round(resumen_reg$coefficients[, 4], 4))
colnames(tabla_resumen_reg) <- c("Variable", "Coeficiente (Beta)", "P-Valor")

kable(tabla_resumen_reg, format = "markdown", 
      caption = "Tabla. Coeficientes Finales del Modelo de Regresión Múltiple")
Tabla. Coeficientes Finales del Modelo de Regresión Múltiple
Variable Coeficiente (Beta) P-Valor
(Intercept) (Intersección) -23.6828 0
Min_Temp Mínima Temperatura -0.4703 0
Max_Temp Máxima Temperatura 2.0044 0
Wind Velocidad Viento 5.8565 0