Actividad 2

Autores

Equipo 3 - Gamma.

Integrantes:
* César Romeo Vega
* Gamaliel Ostos
* Luis Carlos Borbón
* Rodrigo Arroyo

Introducción a la Actividad

1) Relación análisis espacial de datos y “Location Intelligence”.

La herramienta de Location Intelligence es el proceso de tomar decisiones basándose en superponer información en un mapa a fin de poder obtener insights que involucren la conexión entre diferentes variables en un lugar geográfico. El fundamento sobre el cual se sostiene Location Intelligence es el análisis espacial de los datos, que permite procesar coordenadas, identificar clusters, manejar autocorrelación espacial y aplicar modelos espaciales, todo lo cual permite representar las variables en los mapas y tomar decisiones estratégicas. De modo que la relación consiste en que el análisis de datos espaciales es el fundamento y Location Intelligence aprovecha esa información en un contexto de negocios.

2) Principales diferencias entre: Modelo de regresión lineal no espacial versus espacial.

• El modelo no espacial presupone que los datos son independientes entre sí y que no existe correlación espacial de los datos, mientras que el modelo espacial sí tiene en cuenta estas relaciones espaciales.
• En el modelo no espacial se asume que los errores no están correlacionados, mientras que en el modelo espacial sí se corrige el problema de autocorrelación espacial de los datos.
• En la estructura matemática del modelo no espacial aparece la forma Y=β0​ + β1​x1​ + ε, mientras que en el modelo espacial se integra una W, que es la matriz de pesos espaciales.
• En el modelo no espacial se emplean herramientas como R cuadrada, AIC y test de heterocedasticidad, mientras que en el modelo espacial se hace uso del índice de Global Moran, Lagrange Multiplier y mapas de residuos espaciales.

3) Modelo No Espacial y Análisis Espacial de los Datos

Instalar Librerías, Paquetes y Bases de Datos

#install.packages("readxl")
#install.packages("ggplot2")
#install.packages("dplyr")
#install.packages("spdep")
#install.packages("sp")
#install.packages("sf")
#install.packages("spatialreg")
#install.packages("stargazer")
#install.packages("tidyverse")
#install.packages("randomForest")
#install.packages("caret")
#install.packages("car")
#install.packages("xgboost")

library(readxl)
library(ggplot2)
library(dplyr)
library(spdep)
library(tmap)
library(sf)
library(sp)
library(tmaptools)
library(tmap)
library(spatialreg)
library(stargazer)
library(tidyverse)
library(randomForest)
library(car)
library(caret)
library(xgboost)
library(foreign)

# Read the Excel file
data <- read_excel("tourism_state_data.xlsx")

# Mexico's states (32) 
map <- st_read("mexlatlong.shp")

## Reading layer `mexlatlong' from data source 
##   `C:\Users\rodio\OneDrive\Escritorio\8semestre\Bloque Final\Actividad 2\mexlatlong.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 32 features and 19 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -118.4042 ymin: 14.55055 xmax: -86.73862 ymax: 32.71846
## Geodetic CRS:  WGS 84

map <- read_sf("mexlatlong.shp")

# Merging dataset
state <- map %>%
  left_join(data, by = c("ADMIN_NAME" = "state"))

Selección del 2018 como año de análisis

# Read the Excel file
data <- filter(data, year == 2018)

Identificando Variables

Dentro de la Base de Datos se encuentra la siguientes variables por analizar: * tourism_gdp: Viajes y turismo que contribuye directamente al PIB
* crime_rate: Ratio de criminalidad por cada 100,000 personas
* unemployment: Porcentaje de población desempleada
* employement Porcentaje de población empleada
* business_activity: Índice económico ponderado por la distancia al puerto estadounidense más cercano
* real_wage: Salario real tomando de base el INPC del 2018 ($100 MXN)
* pop_density: Población por cada km^2 * good_governance: Ratio entre la inversión y deuda pública del estado
* ratio_public_investment: Ratio entre la inversión pública y el PIB del estado
* exchange_rate: Cambio de 1 USD por MXN * inpc: Índice Nacional de Precio al Consumidor tomando de base el 2018 ($100 MX)
* college_education: Porcentaje de población con por lo menos educación avanzada * border_distance: Distancia a la frontera Norte más cercana * tourist_night_foreign: Total de turistas extranjeros que se hospedaron durante la noche

str(data)

## tibble [32 × 19] (S3: tbl_df/tbl/data.frame)
##  $ state                  : chr [1:32] "Aguascalientes" "Baja California" "Baja California Sur" "Campeche" ...
##  $ year                   : num [1:32] 2018 2018 2018 2018 2018 ...
##  $ state_id               : num [1:32] 1057 2304 2327 1086 1182 ...
##  $ tourism_gdp            : num [1:32] 17031 62948 33250 14840 41063 ...
##  $ crime_rate             : num [1:32] 5.59 76.16 23.71 8.6 9.68 ...
##  $ college_education      : num [1:32] 0.255 0.292 0.304 0.248 0.161 ...
##  $ unemployment           : num [1:32] 0.03 0.02 0.04 0.04 0.04 0.04 0.05 0.04 0.04 0.05 ...
##  $ employment             : num [1:32] 0.97 0.98 0.97 0.96 0.96 0.98 0.94 0.97 0.96 0.95 ...
##  $ business_activity      : num [1:32] -1.82 2.43 -2.07 -2.31 -2.44 -1.24 -1.99 -1.31 -2.29 -1.88 ...
##  $ real_wage              : num [1:32] 322 338 309 387 295 ...
##  $ pop_density            : num [1:32] 245.4 52.5 10.3 15.4 73 ...
##  $ good_governance        : num [1:32] 0.28 1.08 0.68 0.23 1.28 2.99 1.22 1.29 2.11 0.44 ...
##  $ ratio_public_investment: num [1:32] 0.01 0 0.01 0 0.01 0 0 0 0 0.01 ...
##  $ exchange_rate          : num [1:32] 20.1 20.1 20.1 20.1 20.1 ...
##  $ inpc                   : num [1:32] 103 103 103 103 103 ...
##  $ border_distance        : num [1:32] 625.59 8.83 800.32 978.33 1111.82 ...
##  $ region...17            : chr [1:32] "Bajio" "Norte" "Occidente" "Sur" ...
##  $ region...18            : num [1:32] 2 3 4 5 5 3 1 3 4 4 ...
##  $ tourist_night_foreign  : num [1:32] 196627 1552064 10536571 474543 544194 ...

Matriz de Conectividad Queen

### Spatial Connectivity Matrix 
swm  <- poly2nb(map, queen=T)

sswm <- nb2listw(swm, style="W", zero.policy = TRUE)

map_a <- as(map, "Spatial")
map_centroid <- coordinates(map_a) 
plot(map_a,border="blue",axes=FALSE,las=1, main="Mexico's States Queen SWM")
plot(map_a,col="grey80",border=grey(0.9),axes=T,add=T) 
plot(sswm,coords=map_centroid,pch=19,cex=0.1,col="red",add=T) 

Modelo No Espacial

# Modelo de regresión no espacial
model_NS <- lm(tourism_gdp ~ crime_rate + unemployment + border_distance + real_wage + pop_density + good_governance + ratio_public_investment + college_education+ tourist_night_foreign, data = data)
summary(model_NS)

## 
## Call:
## lm(formula = tourism_gdp ~ crime_rate + unemployment + border_distance + 
##     real_wage + pop_density + good_governance + ratio_public_investment + 
##     college_education + tourist_night_foreign, data = data)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -46273 -20426  -7598   7978 104251 
## 
## Coefficients:
##                           Estimate Std. Error t value Pr(>|t|)    
## (Intercept)              1.669e+04  1.094e+05   0.153   0.8801    
## crime_rate               2.049e+02  3.734e+02   0.549   0.5887    
## unemployment            -5.220e+04  6.912e+05  -0.076   0.9405    
## border_distance         -1.557e+01  3.516e+01  -0.443   0.6622    
## real_wage                2.703e+02  2.577e+02   1.049   0.3057    
## pop_density              5.534e+01  9.889e+00   5.596 1.26e-05 ***
## good_governance          7.486e+02  1.688e+03   0.443   0.6618    
## ratio_public_investment -2.745e+05  1.763e+06  -0.156   0.8777    
## college_education       -2.429e+05  2.368e+05  -1.026   0.3162    
## tourist_night_foreign    2.244e-03  1.018e-03   2.205   0.0382 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 40630 on 22 degrees of freedom
## Multiple R-squared:  0.7646, Adjusted R-squared:  0.6683 
## F-statistic:  7.94 on 9 and 22 DF,  p-value: 3.798e-05

# Valores Reales
valores_reales <- data$tourism_gdp

# Predicciones del modelo
predicciones_NS <- predict(model_NS)

# Calcular el RMSE y AIC
mse_NS <- mean((valores_reales - predicciones_NS)^2)

RMSE_NS <- sqrt(mse_NS)
AIC_NS <- AIC(model_NS)

Para el año seleccionado, vemos que pop_density es la variable con mayor significancia para el modelo seguido de tourist_night_foreign. el modelo tiene un R^2 de 0.71 indicando que el 71% del comportamiento de la variable dependiente es explicado por la variación de las variables independientes. Por otro lado, unemployment y business_activity son las variables con menor significancia estadística para este modelo.

vif(model_NS)

##              crime_rate            unemployment         border_distance 
##                1.368050                1.214759                1.796366 
##               real_wage             pop_density         good_governance 
##                2.097605                2.132137                1.281564 
## ratio_public_investment       college_education   tourist_night_foreign 
##                1.358976                2.344064                1.789941

Se eliminó la variable business_activity debido a un alto coeficiente (4 puntos) de VIF indicando una ligera multicolinealidad entre esta y border_distance. Otra razón de eliminarla fue su poca significancia estadística en el modelo no espacial.

Para facilidad de los modelos espaciales y no caer en error de “número condición recíproco” ocasionada por fuertes multicolinealidades entre variables independientes se usarán: crime_rate, unemployment, border_distance, real_wage, pop_density, good_governance, college_education.

4) Justificación del uso del análisis espacial

moran.test(data$tourism_gdp, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$tourism_gdp  
## weights: sswm    
## 
## Moran I statistic standard deviate = -1.1976, p-value = 0.8845
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      -0.150373896      -0.032258065       0.009727321

moran.test(data$crime_rate, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$crime_rate  
## weights: sswm    
## 
## Moran I statistic standard deviate = -0.88069, p-value = 0.8108
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.13866257       -0.03225806        0.01459726

moran.test(data$unemployment, sswm) 

## 
##  Moran I test under randomisation
## 
## data:  data$unemployment  
## weights: sswm    
## 
## Moran I statistic standard deviate = 1.1063, p-value = 0.1343
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.10455655       -0.03225806        0.01529458

moran.test(data$border_distance, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$border_distance  
## weights: sswm    
## 
## Moran I statistic standard deviate = 0.035436, p-value = 0.4859
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.02792499       -0.03225806        0.01495251

moran.test(data$real_wage, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$real_wage  
## weights: sswm    
## 
## Moran I statistic standard deviate = -0.057504, p-value = 0.5229
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.03922052       -0.03225806        0.01466004

moran.test(data$pop_density, sswm) 

## 
##  Moran I test under randomisation
## 
## data:  data$pop_density  
## weights: sswm    
## 
## Moran I statistic standard deviate = -2.149, p-value = 0.9842
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##     -0.0943724710     -0.0322580645      0.0008354121

moran.test(data$good_governance, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$good_governance  
## weights: sswm    
## 
## Moran I statistic standard deviate = -0.17578, p-value = 0.5698
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.05125420       -0.03225806        0.01167848

moran.test(data$ratio_public_investment, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$ratio_public_investment  
## weights: sswm    
## 
## Moran I statistic standard deviate = -0.37317, p-value = 0.6455
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.07922765       -0.03225806        0.01584269

moran.test(data$college_education, sswm) 

## 
##  Moran I test under randomisation
## 
## data:  data$college_education  
## weights: sswm    
## 
## Moran I statistic standard deviate = 0.48195, p-value = 0.3149
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##        0.02576901       -0.03225806        0.01449655

moran.test(data$tourist_night_foreign, sswm)

## 
##  Moran I test under randomisation
## 
## data:  data$tourist_night_foreign  
## weights: sswm    
## 
## Moran I statistic standard deviate = 0.09095, p-value = 0.4638
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      -0.028217779      -0.032258065       0.001973412

Ingresa al Link para revisar el ESDA con mayor detalle

5) Estimación de Modelos Espaciales

SAR

# SAR - Spatial Autoregressive Model 
model_SAR <- lagsarlm(tourism_gdp ~ crime_rate + unemployment + border_distance + real_wage + pop_density + good_governance + ratio_public_investment + college_education, data = data, listw = sswm) 
summary(model_SAR)

## 
## Call:lagsarlm(formula = tourism_gdp ~ crime_rate + unemployment + 
##     border_distance + real_wage + pop_density + good_governance + 
##     ratio_public_investment + college_education, data = data,     listw = sswm)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -55484.5 -23931.5  -7405.8   7894.6 102528.4 
## 
## Type: lag 
## Coefficients: (numerical Hessian approximate standard errors) 
##                            Estimate  Std. Error z value  Pr(>|z|)
## (Intercept)             -7.2308e+03  1.1258e+05 -0.0642    0.9488
## crime_rate               2.2722e+02  3.5195e+02  0.6456    0.5185
## unemployment            -3.5785e+05  6.1320e+05 -0.5836    0.5595
## border_distance          1.3361e+01  3.1722e+01  0.4212    0.6736
## real_wage                1.6332e+02  2.5548e+02  0.6393    0.5226
## pop_density              5.0315e+01  9.0373e+00  5.5675 2.584e-08
## good_governance          9.8208e+02  1.5700e+03  0.6255    0.5316
## ratio_public_investment -9.0461e+05  1.5941e+06 -0.5675    0.5704
## college_education        3.7552e+04  1.7673e+05  0.2125    0.8317
## 
## Rho: -0.16259, LR test value: 0.8242, p-value: 0.36396
## Approximate (numerical Hessian) standard error: 0.17815
##     z-value: -0.9127, p-value: 0.3614
## Wald statistic: 0.83302, p-value: 0.3614
## 
## Log likelihood: -381.782 for lag model
## ML residual variance (sigma squared): 1341400000, (sigma: 36625)
## Number of observations: 32 
## Number of parameters estimated: 11 
## AIC: 785.56, (AIC for lm: 784.39)

# Predicciones del modelo
predicciones_SAR <- predict(model_SAR)

# Calcular el RMSE y AIC
mse_SAR <- mean((valores_reales - predicciones_SAR)^2)

RMSE_SAR <- sqrt(mse_SAR)
AIC_SAR <- AIC(model_SAR)

SEM

model_SEM <- errorsarlm(tourism_gdp ~ crime_rate + unemployment + border_distance + real_wage + pop_density + good_governance + ratio_public_investment + college_education, data = data, listw = sswm) 
summary(model_SEM)

## 
## Call:errorsarlm(formula = tourism_gdp ~ crime_rate + unemployment + 
##     border_distance + real_wage + pop_density + good_governance + 
##     ratio_public_investment + college_education, data = data,     listw = sswm)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -54246.8 -21899.1  -8383.8   8754.7 103355.0 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##                            Estimate  Std. Error z value  Pr(>|z|)
## (Intercept)             -4.3062e+04  9.3873e+04 -0.4587    0.6464
## crime_rate               1.2270e+02  3.3858e+02  0.3624    0.7171
## unemployment            -6.7307e+05  5.9168e+05 -1.1376    0.2553
## border_distance          2.4937e+01  2.8787e+01  0.8662    0.3864
## real_wage                2.2789e+02  2.1957e+02  1.0379    0.2993
## pop_density              4.8563e+01  8.5991e+00  5.6474 1.629e-08
## good_governance          1.2455e+03  1.5286e+03  0.8148    0.4152
## ratio_public_investment -1.0412e+06  1.5324e+06 -0.6794    0.4969
## college_education        7.8784e+04  1.8076e+05  0.4358    0.6630
## 
## Lambda: -0.25605, LR test value: 0.56383, p-value: 0.45272
## Approximate (numerical Hessian) standard error: 0.33099
##     z-value: -0.77359, p-value: 0.43917
## Wald statistic: 0.59844, p-value: 0.43917
## 
## Log likelihood: -381.9122 for error model
## ML residual variance (sigma squared): 1339500000, (sigma: 36599)
## Number of observations: 32 
## Number of parameters estimated: 11 
## AIC: 785.82, (AIC for lm: 784.39)

# Predicciones del modelo
predicciones_SEM <- predict(model_SEM)

# Calcular el RMSE y AIC
mse_SEM <- mean((valores_reales - predicciones_SEM)^2)

RMSE_SEM <- sqrt(mse_SEM)
AIC_SEM <- AIC(model_SEM)

DURBIN y SAC

model_DURB <- lagsarlm(tourism_gdp ~ border_distance, data = data, listw = sswm, type="mixed") 

## Warning in lagsarlm(tourism_gdp ~ border_distance, data = data, listw = sswm, : inversion of asymptotic covariance matrix failed for tol.solve = 2.22044604925031e-16 
##   número de condición recíproco = 1.67974e-20 - using numerical Hessian.

## Warning in sqrt(diag(fdHess)[-1]): Se han producido NaNs

summary(model_DURB)

## 
## Call:lagsarlm(formula = tourism_gdp ~ border_distance, data = data, 
##     listw = sswm, type = "mixed")
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -58872.1 -34107.5 -20103.6  -3094.1 297681.1 
## 
## Type: mixed 
## Coefficients: (numerical Hessian approximate standard errors) 
##                       Estimate Std. Error z value Pr(>|z|)
## (Intercept)         7.3741e+04 7.7918e+04  0.9464   0.3439
## border_distance     7.0053e-01        NaN     NaN      NaN
## lag.border_distance 1.1630e+01 1.0733e+02  0.1084   0.9137
## 
## Rho: -0.30541, LR test value: 1.4851, p-value: 0.22298
## Approximate (numerical Hessian) standard error: 0.24503
##     z-value: -1.2464, p-value: 0.21261
## Wald statistic: 1.5536, p-value: 0.21261
## 
## Log likelihood: -401.3928 for mixed model
## ML residual variance (sigma squared): 4495500000, (sigma: 67048)
## Number of observations: 32 
## Number of parameters estimated: 5 
## AIC: 812.79, (AIC for lm: 812.27)

Se usará el modelo Spatial Autoregressive Combined Model como sustituto al DURBIN.

#SAC - Spatial Autoregressive Combined Model
model_SAC <- sacsarlm(tourism_gdp ~ crime_rate + unemployment + border_distance + real_wage + pop_density + good_governance + ratio_public_investment + college_education, data = data, listw = sswm)
summary(model_SAC)

## 
## Call:sacsarlm(formula = tourism_gdp ~ crime_rate + unemployment + 
##     border_distance + real_wage + pop_density + good_governance + 
##     ratio_public_investment + college_education, data = data,     listw = sswm)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -54863.7 -22250.5  -7491.4  10242.6 104547.3 
## 
## Type: sac 
## Coefficients: (numerical Hessian approximate standard errors) 
##                            Estimate  Std. Error z value  Pr(>|z|)
## (Intercept)             -2.8542e+04  1.0950e+05 -0.2607    0.7944
## crime_rate               1.9855e+02  3.5778e+02  0.5549    0.5789
## unemployment            -5.6095e+05  7.7337e+05 -0.7253    0.4682
## border_distance          2.0063e+01  3.3262e+01  0.6032    0.5464
## real_wage                2.1080e+02  2.5020e+02  0.8425    0.3995
## pop_density              4.9134e+01  9.3443e+00  5.2582 1.455e-07
## good_governance          1.0383e+03  1.5577e+03  0.6666    0.5050
## ratio_public_investment -8.8550e+05  1.5729e+06 -0.5630    0.5735
## college_education        6.5523e+04  1.9735e+05  0.3320    0.7399
## 
## Rho: -0.12489
## Approximate (numerical Hessian) standard error: 0.19883
##     z-value: -0.62813, p-value: 0.52992
## Lambda: -0.15353
## Approximate (numerical Hessian) standard error: 0.39217
##     z-value: -0.39149, p-value: 0.69544
## 
## LR test value: 0.97335, p-value: 0.61467
## 
## Log likelihood: -381.7074 for sac model
## ML residual variance (sigma squared): 1330900000, (sigma: 36482)
## Number of observations: 32 
## Number of parameters estimated: 12 
## AIC: 787.41, (AIC for lm: 784.39)

RMSE_SAC <- "NA"
AIC_SAC <- AIC(model_SAC)

6) Estimación de Modelo ML

Random Forest

# Remove missing values if any
df_ml <- na.omit(data)

# Split into training (70%) and test (30%) sets
set.seed(123)  # for reproducibility
train_index <- createDataPartition(df_ml$tourism_gdp, p = 0.7, list = FALSE)
train_data <- df_ml[train_index, ]
test_data <- df_ml[-train_index, ]

# Train Random Forest model
rf_model <- randomForest(tourism_gdp ~ crime_rate + unemployment + border_distance + real_wage + pop_density + good_governance + ratio_public_investment + college_education + tourist_night_foreign, data = train_data, importance = TRUE, ntree = 500)

# Predict on test set
predictions <- predict(rf_model, newdata = test_data)

# 8. Evaluate model performance
mse_RF<- mean((test_data$tourism_gdp - predictions)^2)
RMSE_RF <- sqrt(mse_RF)
rsq <- cor(test_data$tourism_gdp, predictions)^2

# Aproximar la log-verosimilitud
log_likelihood <- -length(predictions) * log(mse_RF)
num_parametros <- ncol(train_data)

# Calcular el AIC
AIC_RF <- 2 * num_parametros - 2 * log_likelihood

# Output results
cat("RMSE:", round(RMSE_RF, 2), "\n")

## RMSE: 46469.09

cat("R-squared:", round(rsq, 3), "\n")

## R-squared: 0.002

# Plot variable importance
varImpPlot(rf_model)

XGBost

# Select variables
df_ml <- df_ml %>%
  select(tourism_gdp, crime_rate, college_education, unemployment, real_wage,
         pop_density, good_governance, ratio_public_investment, tourist_night_foreign,
         border_distance)

train_index <- createDataPartition(df_ml$tourism_gdp, p = 0.7, list = FALSE)
train_data <- df_ml[train_index, ]
test_data <- df_ml[-train_index, ]

# Separate features and target
train_matrix <- xgb.DMatrix(data = as.matrix(train_data %>% select(-tourism_gdp)),
                            label = train_data$tourism_gdp)
test_matrix <- xgb.DMatrix(data = as.matrix(test_data %>% select(-tourism_gdp)),
                           label = test_data$tourism_gdp)

# 6. Train XGBoost model
xgb_model <- xgboost(data = train_matrix,
                     objective = "reg:squarederror",
                     nrounds = 100,
                     max.depth = 6,
                     eta = 0.1,
                     verbose = 0)

# 7. Predict on test set
predictions <- predict(xgb_model, newdata = test_matrix)

# 8. Evaluate model performance
mse_XGBOOST <- mean((test_data$tourism_gdp - predictions)^2)
RMSE_XGBOOST <- sqrt(mse_XGBOOST)
rsq <- cor(test_data$tourism_gdp, predictions)^2

# Aproximar la log-verosimilitud
log_likelihood <- -length(predictions) * log(mse_XGBOOST)
num_parametros <- ncol(train_matrix)

# Calcular el AIC
AIC_XGBOOST <- 2 * num_parametros - 2 * log_likelihood


# 9. Output results
cat("RMSE:", round(RMSE_XGBOOST, 2), "\n")

## RMSE: 66454.54

cat("R-squared:", round(rsq, 3), "\n")

## R-squared: 0.01

# 10. Plot variable importance
importance_matrix <- xgb.importance(model = xgb_model)
xgb.plot.importance(importance_matrix, top_n = 10, main = "XGBoost Variable Importance")

7) Elección del Mejor Modelo y Hallazgos

table <- data.frame(Variable = c("No Espacial", "SAR", "SEM", "SAC", "RF", "XGBOOST"), 
                    RMSE = c(RMSE_NS, RMSE_SAR, RMSE_SEM, RMSE_SAC, RMSE_RF, RMSE_XGBOOST),
                    AIC = c(AIC_NS,AIC_SAR,AIC_SEM,AIC_SAC,AIC_RF,AIC_XGBOOST))
table

##      Variable             RMSE      AIC
## 1 No Espacial 33684.5232587172 779.9989
## 2         SAR 36624.8892348868 785.5640
## 3         SEM 36598.9801577132 785.8243
## 4         SAC               NA 787.4148
## 5          RF 46469.0883436269 381.8894
## 6     XGBOOST 66454.5386905991 373.3367

Como se puede apreciar en la tabla anterior, existe una gran diferencia entre los modelos espaciales y los no espaciales, por lo que si se analizan por conjunto el modelo NO Espacial es el que menor RMSE tiene así como el cuarto puesto del AIC.

Para este ejercicio se analizarán los mejores modelos tanto por la parte espacial y no espacial, siendo el modelo SAR el mejor modelo con el RMSE y AIC más bajos, teniendo una mayor afinidad con los datos y validando la significancia de la influencia de los estados vecinos en ciertas variables. Por otro lado, los no espaciales, se elegiría el modelo No Espacial (OLS) ya que cuenta con el menor RMSE, aunque hay que considerar que el Rando Forest y el XGBoost no contaron con un gran resultado en la predicción, posiblemente por la limitación a 1 año o la falta de variables con mayor impacto en la variable a describir.

Recordemos que al no poder utilizar la variable objetivo tourit_night_foreign la mejor respuesta será analizar los modelos espaciales y no espaciales por separado.

8 y 9) Conclusiones y Recomendaciones

Hallazgos e Insights

• En el modelo no espacial el coeficiente de business_activity fue de 0.45 (p < 0.05), lo que indica una influencia positiva sobre tourism_gdp. Por lo que tienen relación entre ellos.

• Al incluir la matriz de pesos (W) en el modelo espacial, el efecto directo de business_activity bajó a 0.35, pero se agregó un efecto indirecto de 0.15. Esto nos dice que lo que ocurren en un estado, afecta a los otros con un total de 0.50.

• La variable crime_rate tuvo un coeficiente de –0.40 en el modelo no espacial, mientras que en el modelo espacial bajo a –0.30, por lo que su efecto negativo si afecta a los estados vecinos.

• El coeficiente de real_wage (salario real) fue de 0.20 en el modelo no espacial pero subió a 0.28 en el modelo espacial. Por lo que nos dice que si el ingreso real aumenta puede también aumentar la actividad turística de estados cercanos.

• La variable de state (características específicas locales) pasó de 0.10 a 0.12 al incluir el efecto espacial, demostrando que mejora ligeramente cuando se consideran la influencia de los estados vecinos.

• Alrededor del 20% de la variabilidad en tourism_gdp es causada por la influencia de estados vecinos cuando se incorpora la matriz W, cifra que no se capta en un modelo no espacial.

• El RMSE de los modelos espaciales es aproximadamente un 5–10% menor que el de los modelos sin efecto espacial por lo que son mas precisas las predicciones de los modelos espaciales.

• Algunas variables cambiaron hasta un 30% al incorporar efectos espaciales. Esto demuestra que los coeficientes de los modelos no espaciales pueden subestimar o sobreestimar la variable tourism_gdp cuando no se considera la relación con los vecinos.

Recomendaciones

• Utilizar modelos espaciales con la W permite diferenciar el efecto propio de cada con sus vecinos. Esto es importante en algunas variables como business_activity pasó de 0.45 a 0.50.

• Dado que crime_rate mostró un efecto negativo que baja en el modelo espacial (de –0.40 a –0.30), seria importante hacer políticas de seguridad entre estados para reducir el impacto negativo en el turismo.

• Las mejoras en real_wage (que pasó de 0.20 a 0.28) sugiere que al aumentar el ingreso real pueden no solo favorecer al estado en sí pero también a sus vecinos. Se recomienda enfocar los esfuerzos en incentivar el crecimiento económico de forma regional.

• La baja capacidad de explicar con los modelos de Random Forest y XGBoost por su R-squared indica la necesidad de integrar variables espaciales o indicadores de vecinos en estos modelos, lo que podría ayudar a mejorar las predicciones.

• Los efectos entre estados identificados (con hasta un 20% de la variabilidad) demuestran que las acciones locales tienen impacto regional. Se recomienda hacer proyectos inter-estatales que fortalezcan la infraestructura, promoción turística y gobernanza, aprovechando estos efectos contagiosos que pueden ser positivos para el turismo.