R_logaritmica

Regresión Logaritmica

UNIVERSIDAD CENTRAL DEL ECUADOR

PROYECTO: FOCOS DE CALOR EN EL ECUADOR

AUTORES: GUERRERO MARIA GABRIELA,PUCHAICELA MONICA, ZURITA JOHANNA

FECHA: 14/05/2025

library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.1.0     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(fitdistrplus)

## Loading required package: MASS
## 
## Attaching package: 'MASS'
## 
## The following object is masked from 'package:dplyr':
## 
##     select
## 
## Loading required package: survival

# Leer el archivo CSV
datos <- read.csv("focosdecalor_bdd_2021diciembre (1).csv", header = TRUE, sep = ",", dec = ".")
str(datos)

## 'data.frame':    22476 obs. of  17 variables:
##  $ MES_REPORT: int  11 11 8 6 5 6 11 9 3 3 ...
##  $ DIA_REPORT: int  20 20 6 10 28 10 20 29 22 22 ...
##  $ DPA_DESPRO: chr  "ZAMORA CHINCHIPE" "ZAMORA CHINCHIPE" "ZAMORA CHINCHIPE" "ZAMORA CHINCHIPE" ...
##  $ DPA_DESCAN: chr  "CHINCHIPE" "CHINCHIPE" "CHINCHIPE" "CHINCHIPE" ...
##  $ DPA_DESPAR: chr  "CHITO" "CHITO" "PUCAPAMBA" "PUCAPAMBA" ...
##  $ TXT_1     : chr  "PARROQUIA RURAL" "PARROQUIA RURAL" "PARROQUIA RURAL" "PARROQUIA RURAL" ...
##  $ LATITUDE  : chr  "-4,981720000000000" "-4,969160000000000" "-4,958520000000000" "-4,957820000000000" ...
##  $ LONGITUDE : chr  "-79,041280000000000" "-79,049490000000006" "-79,118430000000004" "-79,111859999999993" ...
##  $ BRIGHTNESS: chr  "354,759999999999990" "342,009999999999990" "331,860000000000010" "331,399999999999980" ...
##  $ SCAN      : chr  "0,510000000000000" "0,510000000000000" "0,150000000000000" "0,540000000000000" ...
##  $ TRACK     : chr  "0,490000000000000" "0,490000000000000" "0,380000000000000" "0,420000000000000" ...
##  $ SATELLITE : chr  "1" "1" "1" "1" ...
##  $ CONFIDENCE: chr  "n" "n" "n" "n" ...
##  $ VERSION   : chr  "2.0NRT" "2.0NRT" "2.0NRT" "2.0NRT" ...
##  $ BRIGHT_T31: chr  "299,420000000000020" "298,149999999999980" "299,160000000000030" "296,800000000000010" ...
##  $ FRP       : chr  "12,100000000000000" "6,870000000000000" "3,770000000000000" "5,500000000000000" ...
##  $ DAYNIGHT  : chr  "D" "D" "D" "D" ...

# Reemplazar comas por puntos y convertir columnas numéricas
datos$BRIGHTNESS <- as.numeric(gsub(",", ".", datos$BRIGHTNESS))

# Filtrar: meses 1–5 y BRIGHTNESS entre 250 y 300
datos <- datos %>%
  filter(MES_REPORT <= 5, BRIGHTNESS >= 250, BRIGHTNESS <= 300)

# Calcular promedios mensuales
monthly_avg <- datos %>%
  group_by(MES_REPORT) %>%
  summarise(BRIGHTNESS = mean(BRIGHTNESS, na.rm = TRUE)) %>%
  arrange(MES_REPORT)

# Definir función para calcular R²
r_squared_nls <- function(model, data, y_col) {
  y <- data[[y_col]]
  y_pred <- predict(model)
  ss_res <- sum((y - y_pred)^2)
  ss_tot <- sum((y - mean(y))^2)
  r2 <- 1 - (ss_res / ss_tot)
  return(r2)
}

# Regresión logarítmica con nls
model_brightness <- nls(
  BRIGHTNESS ~ a * log(MES_REPORT) + b, 
  data = monthly_avg, 
  start = list(a = 1, b = 290)
)

# Agregar predicciones al DataFrame
monthly_avg <- monthly_avg %>%
  mutate(pred_brightness = predict(model_brightness))

# Calcular R² antes de usarlo
r2_brightness <- r_squared_nls(model_brightness, monthly_avg, "BRIGHTNESS")

# Extraer coeficientes y redondear
coefs <- coef(model_brightness)
a <- round(coefs["a"], 3)
b <- round(coefs["b"], 3)
r2 <- round(r2_brightness, 4)

# Crear etiqueta para mostrar en el gráfico
equation_label <- paste0("BRIGHTNESS = ", a, " * log(MES) + ", b, 
                         "\nR² = ", r2)

# Graficar resultados
ggplot(monthly_avg, aes(x = MES_REPORT)) +
  geom_point(aes(y = BRIGHTNESS), color = "black") +
  geom_line(aes(y = pred_brightness), color = "blue") +
  labs(
    title = "Regresión logarítmica del BRIGHTNESS (Ene–May 2021)",
    subtitle = equation_label,
    x = "Mes",
    y = "BRIGHTNESS"
  )

# Imprimir R² en consola
cat("R² para BRIGHTNESS: ", r2, "\n")

## R² para BRIGHTNESS:  0.8975

#Conclusión
#Se aplicó un modelo de regresión logarítmica
#para analizar la relación entre BRIGHTNESS y el mes. 
#Los datos fueron depurados considerando solo los meses de enero a mayo de temporada seca 
#y valores de BRIGHTNESS entre 250 y 300 donde se encuentran focos moderados a altos, 
#eliminando outliers y ruido. Las ecuaciones obtenidas mostraron buen ajuste, 
#con un R² de 0.9 aproximadamente, reflejando tendencias térmicas consistentes.