Análisis espacial y preparación de datos (California)
Para hacer la tabla de variables que utilizaremos ya en el modelo logístico, utilizaremos datos espaciales, mapas, etc. que no sólo nos ayudaran a determinar de forma numérica ciertos eventos que creemos cruciales para analizar si es viable poner un restaurante o no.
Sino también nos ayudarán a visualizar la información a traves de mapas etc.
Todo se hará en tramos censales, la unidad de medida utilizada por el US Census. Utilizando los TIGER files de la base de datos estadounidense la cual contiene datos en formato sf (spatial data) los cuales son datos geográficos donde están los edificios al igual que vias de acceso, rutas, etc.
Ocuparemos los tramos censales en lugar de los bloques censales, la diferencia siendo que un bloque censal es muy pequeño, un tramo representa lo que sería una colonia mientras que un bloque censal puede representar una manzana o cuadra, por lo que nuestras posibilidades de clasificar bien nuestras variables se disminuye.
También dentro de los TIGER files utilizaremos la base de datos de “landmarks”, que representa edificios importantes dentro de cada estado, escuelas, hospitales, iglesias ,edificios administrativos, Malls, aeropuertos entre otros.
#Reduciremos la base de datos al estado que necesitamosraw_finalrestaurant_df<-raw_finalrestaurant_df %>%filter(state=="California")cat("Las variables útiles de la base de datos exploratoria son: \n")
Las variables útiles de la base de datos exploratoria son:
Reading layer `tl_2018_06_pointlm' from data source
`C:\Users\wars_\Documents\Universidad\Actuaría\2025-A\Computo_cientifico\Proyecto_final\sf_california_files\tl_2018_06_pointlm\tl_2018_06_pointlm.shp'
using driver `ESRI Shapefile'
Simple feature collection with 82815 features and 5 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: -124.407 ymin: 32.53756 xmax: -114.1762 ymax: 42.0025
Geodetic CRS: NAD83
Creación de la tabla final para modelo LOGIT
#Conversion de coordenadas de restaurantes a sfrestaurants_sf<-st_as_sf(raw_finalrestaurant_df, coords =c("longitude","latitude"), crs=4326)#Checar sistema de coordenada de datos TRIGER/Line filesst_crs(restaurants_sf)
Coordinate Reference System:
User input: EPSG:4326
wkt:
GEOGCRS["WGS 84",
ENSEMBLE["World Geodetic System 1984 ensemble",
MEMBER["World Geodetic System 1984 (Transit)"],
MEMBER["World Geodetic System 1984 (G730)"],
MEMBER["World Geodetic System 1984 (G873)"],
MEMBER["World Geodetic System 1984 (G1150)"],
MEMBER["World Geodetic System 1984 (G1674)"],
MEMBER["World Geodetic System 1984 (G1762)"],
MEMBER["World Geodetic System 1984 (G2139)"],
MEMBER["World Geodetic System 1984 (G2296)"],
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ENSEMBLEACCURACY[2.0]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
AXIS["geodetic latitude (Lat)",north,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["geodetic longitude (Lon)",east,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]],
USAGE[
SCOPE["Horizontal component of 3D system."],
AREA["World."],
BBOX[-90,-180,90,180]],
ID["EPSG",4326]]
st_crs(california_censustracts)
Coordinate Reference System:
User input: NAD83
wkt:
GEOGCRS["NAD83",
DATUM["North American Datum 1983",
ELLIPSOID["GRS 1980",6378137,298.257222101,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
AXIS["latitude",north,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["longitude",east,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4269]]
#Coincidencia de sistema de coordenadascalifornia_censustracts<-st_transform(california_censustracts,4326)california_landmarks<-st_transform(california_landmarks,4326)#Contar numero de restaurantes por tramo censal:no_restaurantes<-st_join(restaurants_sf,california_censustracts["GEOID"]) %>%count(GEOID) %>%rename("restaurants_count"="n")no_restaurantes<-st_drop_geometry(no_restaurantes)#Contar numero de edificios importantes por tramo censal:no_edificios<-st_join(california_landmarks,california_censustracts["GEOID"]) %>%count(GEOID) %>%rename("landmarks_count"="n")no_edificios<-st_drop_geometry(no_edificios)#Juntar ya en la base de datos para el modelofinalrestaurant_california_df<-left_join(california_censustracts,no_restaurantes, by="GEOID")finalrestaurant_california_df<-left_join(finalrestaurant_california_df,no_edificios, by="GEOID")finalrestaurant_california_df<-finalrestaurant_california_df %>%mutate(restaurants_count=ifelse(is.na(restaurants_count),0,restaurants_count),landmarks_count=ifelse(is.na(landmarks_count),0,landmarks_count))
Ahora juntaremos esta informacion con datos censales del tidycensus package, como hemos mencionado, utilizaremos diferentes variables sociales y económicas y algunas serán calculadas como una tasa de relación. Por ejemplo, la tasa de pobreza es el número de personas en pobreza de ese tramo censal entre las personas totales del tramo. Para la tasa de empleo, es el número de personas en el momento con un empleo entre la fuerza laboral total.
#Juntar con la base de datos normal:finalrestaurant_california_df<-left_join(finalrestaurant_california_df,acs_data_2, by="GEOID") %>%mutate(poverty_rate=poverty/poverty_status_determination,employment_rate=employment/labor_force)#Generar base de datos sólo con información numéricafinalrestaurant_california_df_nogeometry<-st_drop_geometry(finalrestaurant_california_df)finalrestaurant_california_df_nogeometry<-finalrestaurant_california_df_nogeometry %>%select(STATEFP,GEOID,13:22) %>%mutate(viable=ifelse(restaurants_count>=1,1,0))
La variable “viable” es nuestra variable dependiente, 0 es no viable, 1 sí lo es.
Análisis Gráfico California
Veremos cuál es la relacion entre las condiciones socio-económicas que hay en los tramos censales y la ubicación de los restaurantes posicionados en ese estado. Esto como un análisis previo visual antes de hacer el modelo logit y ver si hay algún patrón a identificar fácilmente.
california_restaurantes_censustract_graph<-finalrestaurant_california_df %>%ggplot() +geom_sf(color="black") +geom_sf(data=restaurants_sf, color="blue", size=0.5) +labs(title ="Restaurantes de Fast Food en California, EU")california_restaurantes_censustract_graph
california_restaurantes_median_income<-finalrestaurant_california_df %>%mutate(median_income=ifelse(is.na(median_income),0,median_income)) %>%ggplot() +geom_sf(aes(fill=median_income),linewidth =0.05) +geom_sf(data=restaurants_sf, color="green", size=0.35) +scale_fill_scico(palette ="lajolla", direction=1) +labs(title ="Ingreso medio y restaurantes en California, EU")california_restaurantes_median_income
california_sanfrancisco_restaurantes_median_income<-finalrestaurant_california_df %>%mutate(median_income=ifelse(is.na(median_income),0,median_income)) %>%ggplot() +geom_sf(aes(fill=median_income),linewidth =0.05) +geom_sf(data=restaurants_sf, color="green", size=0.75) +coord_sf(xlim =c(-122.5, -121.5), ylim =c(37.2, 38)) +scale_fill_scico(palette="lajolla", direction=1) +labs(title ="Ingreso medio y restaurantes en San Francisco, California, EU")california_sanfrancisco_restaurantes_median_income
california_restaurantes_poverty_rate<-finalrestaurant_california_df %>%drop_na(poverty_rate) %>%ggplot() +geom_sf(aes(fill=poverty_rate),linewidth =0.05) +geom_sf(data=restaurants_sf, color="green", size=0.5) +scale_fill_scico(palette ="oslo", direction=1) +labs(title ="Tasa desempleo y restaurantes en California, EU")california_restaurantes_poverty_rate
california_sanfrancisco_restaurantes_poverty_rate<-finalrestaurant_california_df %>%drop_na(poverty_rate) %>%ggplot() +geom_sf(aes(fill=poverty_rate),linewidth =0.05) +geom_sf(data=restaurants_sf, color="green", size=0.35) +coord_sf(xlim =c(-122.5, -121.5), ylim =c(37.2, 38)) +scale_fill_scico(palette="oslo", direction=1) +labs(title ="Tasa desempleo y restaurantes en San Francisco, California, EU")california_sanfrancisco_restaurantes_poverty_rate
california_sanfrancisco_restaurantes_poverty<-finalrestaurant_california_df %>%drop_na(poverty) %>%ggplot() +geom_sf(aes(fill=poverty),linewidth =0.05) +geom_sf(data=restaurants_sf, color="green", size=0.35) +coord_sf(xlim =c(-122.5, -121.5), ylim =c(37.2, 38)) +scale_fill_scico(palette="lipari", direction=1) +labs(title ="Población desempleada y restaurantes en San Francisco, California, EU")california_sanfrancisco_restaurantes_poverty
california_sanfrancisco_restaurantes_employment<-finalrestaurant_california_df %>%drop_na(employment) %>%ggplot() +geom_sf(aes(fill=employment),linewidth =0.05) +geom_sf(data=restaurants_sf, color="red", size=0.75) +coord_sf(xlim =c(-122.5, -121.5), ylim =c(37.2, 38)) +scale_fill_scico(palette="devon", direction=1,limits =c(0, 7500)) +labs(title ="Población con empleo y restaurantes en San Francisco, California, EU")california_sanfrancisco_restaurantes_employment
california_sanfrancisco_restaurantes_population<-finalrestaurant_california_df %>%drop_na(population) %>%ggplot() +geom_sf(aes(fill=population),linewidth =0.05) +geom_sf(data=restaurants_sf, color="red", size=0.5) +coord_sf(xlim =c(-122.5, -121.5), ylim =c(37.2, 38)) +scale_fill_scico(palette="bamako", direction=1,limits =c(0, 20000)) +labs(title ="Población total y restaurantes en San Francisco, California, EU")california_sanfrancisco_restaurantes_population
california_restaurantes_population<-finalrestaurant_california_df %>%#drop_na(population) %>%ggplot() +geom_sf(aes(fill=population),linewidth =0.05) +geom_sf(data=restaurants_sf, color="red", size=0.35) +# coord_sf(xlim = c(-122.5, -121.5), ylim = c(37.2, 38)) +scale_fill_scico(palette="bamako", direction=1) +labs(title ="Población total y restaurantes en San Francisco, California, EU")california_restaurantes_population
Por otro lado generaremos la base de datos del Estado de Nueva York, California nos servirá como la base de datos para alimentar al modelo, y el estado de Nueva York para hacer comparaciones predictivas y análisis de conclusión.
Coordinate Reference System:
User input: EPSG:4326
wkt:
GEOGCRS["WGS 84",
ENSEMBLE["World Geodetic System 1984 ensemble",
MEMBER["World Geodetic System 1984 (Transit)"],
MEMBER["World Geodetic System 1984 (G730)"],
MEMBER["World Geodetic System 1984 (G873)"],
MEMBER["World Geodetic System 1984 (G1150)"],
MEMBER["World Geodetic System 1984 (G1674)"],
MEMBER["World Geodetic System 1984 (G1762)"],
MEMBER["World Geodetic System 1984 (G2139)"],
MEMBER["World Geodetic System 1984 (G2296)"],
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ENSEMBLEACCURACY[2.0]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
AXIS["geodetic latitude (Lat)",north,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["geodetic longitude (Lon)",east,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]],
USAGE[
SCOPE["Horizontal component of 3D system."],
AREA["World."],
BBOX[-90,-180,90,180]],
ID["EPSG",4326]]
Coordinate Reference System:
User input: NAD83
wkt:
GEOGCRS["NAD83",
DATUM["North American Datum 1983",
ELLIPSOID["GRS 1980",6378137,298.257222101,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
AXIS["latitude",north,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["longitude",east,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4269]]
Seccion de modelo logit
Desarrollaremos el modelo al igual que haremos sus pruebas estadísticas para ver la veracidad del modelo conforma a la información estadística, cabe mencionar que los datos que usamos para California son del año 2018 para que coincida con la última actualización de nuestra base de datos de restaurantes pública de Datafini.
Haremos nuevas predicciones y haremos comparaciones a través de mapas de calor sobre la viabilidad vs la ubicación real de los restaurantes en el estado y también en le ciudad.
# Predicciones para NY Citynuevos_datos <- finalrestaurant_new_york_df_nogeometry#Imputar faltantes si es necesarionuevos_datos$median_income <-impute_median(nuevos_datos$median_income)nuevos_datos$poverty_rate <-impute_median(nuevos_datos$poverty_rate)nuevos_datos$employment_rate <-impute_median(nuevos_datos$employment_rate)colnames(nuevos_datos)
# Seleccionar sólo variables para el modelonuevos_datos <- nuevos_datos[, c("GEOID","restaurants_count", "landmarks_count", "poverty", "median_income", "poverty_rate", "employment_rate")]# Predicciónprob_pred <-predict(modelo_cv, newdata = nuevos_datos, type ="prob")class_pred <-predict(modelo_cv, newdata = nuevos_datos)# Conversión a factor de predicciónpred_binarias <-ifelse(class_pred =="Yes", 1, 0)# Crear dataframe con resultadosresultados <-data.frame(nuevos_datos,prob_viable = prob_pred$Yes,pred_class = class_pred,pred_binaria = pred_binarias)# Guardar resultadoswrite.csv(resultados, "C:/Users/wars_/Documents/Universidad/Actuaría/2025-A/Computo_cientifico/Proyecto_final/Tablas_modelo_logit_newyork/resultados_prediccion_ny.csv", row.names =FALSE)resultados$GEOID <-as.character(resultados$GEOID)map_data <- new_york_censustracts %>%left_join(resultados, by ="GEOID")
Gráfico de viabilidad, en términos probabilísticos (entre 0 y 1) del Estado de Nueva York.
La viabilidad determinada como una matriz de probabilidades dependiendo del tramo censal, identificado por un “GEOID” que es el identificador dentro de las bases de datos de US CENSUS, la mostraremos en seguida a través de un mapa de calor, probabilidades más cercanas a 1 indican una viabilidad mayor siendo entonces más conveniente intentar colocar una franquicia en esa zona, una investigación más a fondo podría ingresar un umbral de “alta competencia”, donde si es un lugar con ya muchas franquicias establecidas será más conveniente no colocar el restaurante en esa zona.
library(ggplot2)ggplot(map_data) +geom_sf(aes(fill = prob_viable), color =NA) +scale_fill_viridis_c(option ="plasma") +labs(title ="Viabilidad por tracto censal en NY",fill ="Probabilidad") +theme_minimal()
También es conveniente hacer un análisis de otras variables como el ingreso medio, que es una de nuestras variables con mayor significancia para determinar la viabilidad y el mapa de calor de viabilidad que acabamos de mostrar.
NewYork_income_vs_viability<-plot_grid(NewYork_restaurantes_income,NewYork_Viability,ggdraw() +draw_label("Viabilidad vs Ingreso Promedio - Estado de Nueva York", fontface ='bold', size =14,x =0,hjust =0),labels =c("Ingreso promedio y restaurantes","Viabilidad de ubicación de restaurantes"),ncol=2,label_x =0, # 0 = far left, 1 = far rightlabel_y =1, # Optional: top of the plothjust =0, # Align text to the left edge of the label boxalign ="hv",axis ="tblr")ggsave("C:/Users/wars_/Documents/Universidad/Actuaría/2025-A/Computo_cientifico/Proyecto_final/Tablas_modelo_logit_newyork/Plots_finales/Plots_tramocensal/NewYork_income_vs_viability_plot.png", plot = NewYork_income_vs_viability, width =10, height =10, dpi =300)NewYork_income_vs_viability
Gráfico de viabilidad, en términos probabilísticos (entre 0 y 1) de la ciudad de Nueva York
Análogo al análisis previo sólo que enfocado en la ciudad de Nueva York, para ello, investigamos los condados que delimitan a la ciudad y en función a ello, filtrar nuestra base de datos.
#Filtracion para NY Citynyc_counties <-c("005", "047", "061", "081", "085")new_york_censustracts <- new_york_censustracts %>%filter(COUNTYFP %in% nyc_counties)map_data_nyc <- new_york_censustracts %>%left_join(resultados, by ="GEOID")NewYorkCity_viability<-ggplot(map_data_nyc) +geom_sf(aes(fill = prob_viable), color =NA) +coord_sf(xlim =c(-74.25884, -73.70017), ylim =c(40.47658, 40.91770)) +scale_fill_viridis_c(option ="plasma") +labs(fill ="Probabilidad") +theme_minimal()NewYorkCity_restaurantes_median_income<-finalrestaurant_new_york_df %>%mutate(median_income=ifelse(is.na(median_income),0,median_income)) %>%filter(COUNTYFP %in% nyc_counties) %>%ggplot() +geom_sf(aes(fill=median_income),linewidth =0.05) +geom_sf(data=restaurants_sf, color="green", size=0.75) +coord_sf(xlim =c(-74.25884, -73.70017), ylim =c(40.47658, 40.91770)) +scale_fill_scico(palette="lajolla", direction=1)+theme_minimal()NewYorkCity_income_vs_viability<-plot_grid(NewYork_restaurantes_income,NewYork_Viability,ggdraw(),labels =c("Ingreso promedio y restaurantes","Viabilidad de ubicación de restaurantes"),ncol=2,rel_widths =c(2, 2),label_x =0, # 0 = far left, 1 = far rightlabel_y =1, # Optional: top of the plothjust =0, # Align text to the left edge of the label boxalign ="hv",axis ="tblr")ggsave("C:/Users/wars_/Documents/Universidad/Actuaría/2025-A/Computo_cientifico/Proyecto_final/Tablas_modelo_logit_newyork/Plots_finales/Plots_tramocensal/NewYorkCity_income_vs_viability_plot.png", plot = NewYorkCity_income_vs_viability, width =10, height =10, dpi =300)NewYorkCity_restaurantes_median_income
NewYorkCity_income_vs_viability
Viabilidad de la ciudad de Nueva York, análisis individual
