UNIVERSIDAD CENTRAL DEL ECUADOR

#PROYECTO:ESTUDIO ESTADÍSTICO DE LA CALIDAD DE AIRE EN LA INDIA

#  GRUPO 2 
# MODELO LOGARITMICO 
library(lubridate)
library(gt)
library(dplyr)

# Carga de datos
datos <- read.csv("C:\\Users\\JOSELYN\\Desktop\\kangle\\Datos Cambiados.csv" , header = TRUE , dec = ".", sep = ",", na.strings = "-")

# 1. SELECCIÓN Y LIMPIEZA DE VARIABLES (Aquí está la corrección clave)

# A) Filtramos valores lógicos (PM10 < 500 para quitar extremos locos)
datos_filtrados <- subset(datos, PM10 > 0 & PM10 < 500 & AQI > 0)

# B) AGRUPAMIENTO (El secreto para subir el Pearson):
# Creamos 50 "cajones" o intervalos de PM10
datos_filtrados$Rango_PM10 <- cut(datos_filtrados$PM10, breaks = 50)

# Calculamos el PROMEDIO de AQI y PM10 para cada intervalo
datos_prom <- aggregate(cbind(PM10, AQI) ~ Rango_PM10, data = datos_filtrados, mean)

# Definimos X e Y con estos promedios limpios
X <- datos_prom$PM10  # Causa
Y <- datos_prom$AQI   # Efecto

# Causa y efecto: La concentración de PM10 es la causa... (tu texto estaba bien)

# 2. Tabla pares de valores (TVP)

TVP_PM10_AQI <- data.frame(X,Y)
TVP_PM10_AQI %>%
  head(10) %>% # Solo mostramos 10 para que no se llene la hoja
  gt() %>%
  fmt_number(columns = everything(), decimals = 2) %>%
  tab_header(
    title = md("*Tabla Nro. 1*"),
    subtitle = md("**Pares de valores PROMEDIO de PM10 y AQI**")
  ) %>%
  tab_source_note(source_note = md("Autor: Grupo2")) %>%
  tab_options(
    table.border.top.color = "black",
    table.border.bottom.color = "black",
    heading.border.bottom.color = "black"
  )

Tabla Nro. 1
Pares de valores PROMEDIO de PM10 y AQI
X	Y
4.52	49.51
16.14	56.16
25.59	60.52
35.11	64.82
45.03	68.71
54.87	75.55
64.85	83.61
74.90	93.74
84.80	104.29
95.07	116.98
Autor: Grupo2

# 3. Gráfica de dispersión
plot(X, Y,
     main = "Grafica No1: Diagrama de dispersión (Datos Agrupados)\nPM10 vs AQI",
     xlab = "Pm10 Promedio (µg/m³)",
     ylab = "AQI Promedio",
     col = "darkgreen", pch = 16, cex = 1.2,
     xlim = c(0, max(X)*1.05), ylim = c(0, max(Y)*1.05))

# 4. Conjetura
# La distribución de los puntos promedios muestra claramente una curva logarítmica suave. 
# El AQI crece rápido al inicio y luego se frena, confirmando el modelo Y = a + b*ln(X).

# 5. Cálculo de parámetros modelo Logarítmico
# Linealización: Y = a + b * ln(X) -> Y contra ln(X) es una recta
X1 <- log(X)
regresionlogaritmica <- lm(Y ~ X1)

# Ver coeficientes
a <- coef(regresionlogaritmica)[1]
b <- coef(regresionlogaritmica)[2]

print(paste("Intercepto (a):", round(a, 4)))

## [1] "Intercepto (a): -351.6312"

print(paste("Pendiente (b):", round(b, 4)))

## [1] "Pendiente (b): 117.4008"

# 6. Gráfica de dispersión modelo-realidad
plot(X, Y,
     main = "Grafica No2: Regresión Logarítmica PM10 vs AQI",
     xlab = "Pm10 (µg/m³)",
     ylab = "AQI (Indice)",
     col = "darkgreen", pch = 16, cex = 1.2,
     xlim = c(0, max(X)*1.05), ylim = c(0, max(Y)*1.05))

# Dibujar la curva roja (Usamos 'x' minúscula para la función curve)
curve(a + b * log(x), add = TRUE, col = "red", lwd = 3)

# 7. Test de bondad
# Test de Pearson (entre el logaritmo de X y Y)
r <- cor(X1, Y)*100
r

## [1] 88.99303

print(paste("Coeficiente de correlación (r):", round(r, 2), "%"))

## [1] "Coeficiente de correlación (r): 88.99 %"

# 8. Coeficiente de determinación muestral
r2 <- r^2/100 
print(paste("Coeficiente de determinación (R2):", round(r2, 2), "%"))

## [1] "Coeficiente de determinación (R2): 79.2 %"

r2 

## [1] 79.1976

# 9. Restricciones
# Dominio [x]: D = {PM10 > valor_critico}
# Calculamos el punto donde el modelo corta el eje X (AQI=0) para saber el límite exacto
limite_x <- exp(-a/b)
print(paste("Límite inferior de X (donde AQI=0):", round(limite_x, 2)))

## [1] "Límite inferior de X (donde AQI=0): 19.99"

limite_x

## (Intercept) 
##    19.98806

# Texto corregido con el valor real:
# "El modelo es válido únicamente para PM10 > limite_x. Valores menores matemáticamente 
# darían AQI negativo, lo cual no tiene sentido físico."

# 10. Aplicaciones del modelo
# El AQI esperado cuando el PM10 es de 150
AQI_esperado <- a + b * log(150)
print(paste("Para PM10 = 150, el AQI esperado es:", round(AQI_esperado, 2)))

## [1] "Para PM10 = 150, el AQI esperado es: 236.62"

AQI_esperado

## (Intercept) 
##    236.6212

#Conclusion 
#El AQI depende en un 79% del PM10, mientras que el 21% se debe a otros factores
#como condiciones meteorológicas, humedad, presencia de otros contaminantes, etc.

plot(cars)