Instalamos y cargamos la librería
install.packages('RgoogleMaps')
library(RgoogleMaps)
Posterior a ello, definimos algunos parámetros como el zoom de nuestro mapa, la latitud y la longitud centro, y adicional a ello una distancia (Euclidiana) que será el radio de círculos concéntricos a partir de nuestro centro.
z0 <- 8 # zoom
y0 <- 19.040442 # Latitud zero
x0 <- -98.250160 # longitud zero
k <- 0.4 # Define el salto del circulo
# Saltos de 11 Km por cada décima
Colocaremos en el centro del mapa (o cualquier coordenada), un marcador que nos indique alguna referencia adicional.
markers = paste0("&markers=color:blue|label:I|19.040442,-98.250160")
Acto seguido, obtener nuestro mapa base:
PuebMap <- GetMap(center= c(y0,x0), zoom=z0, GRAYSCALE = TRUE,
size=c(640,640),markers=markers,
destfile = file.path(tempdir(),"interm.png"),
maptype = 'mapmaker-hybrid')
Enseguida, importamos nuestra base de datos de las unidades económicas, la cual deberá tener una estructura parecida a esto:
> head(DF)
id lat lon x1 x2
1 3535207 18.87717 -97.72655 87.5724 3
2 3398690 18.49201 -97.40435 1167.6320 40
3 3263793 19.04774 -98.04557 87.5724 3
4 3467483 19.12102 -98.24399 1167.6320 40
5 3316191 18.95277 -97.91986 87.5724 3
6 3457112 19.06409 -98.23185 87.5724 3
x1 y x2 representan variables cuantitativas a representar mediante un gráfico de burbujas. Lo siguiente es elaborar el gráfico sobre nuestro mapa base.
tpm <- bubbleMap(DF,coords=c("lon","lat"),
crs=sp::CRS("+proj=longlat +datum=WGS84"),
max.radius=500,legendLoc="bottomright",
filename = file.path(tempdir(),"interm.png"),
map=PuebMap,zcol='zinc',
key.entries = c(100,250,500,750,2500,6000));
Lucirá de la siguiente manera:
Lo siguiente que haremos será graficar n círculos concentricos de radio n*k en nuestro mapa, con las siguientes líneas de comando:
r1 <- k*1
r2 <- k*2
r3 <- k*3
n <- seq(0,359,0.05)
a <- seq(1,360,0.05)
theta <- (360*n)/a
x1 <- x0 + r1*cos(theta)
y1 <- y0 + r1*sin(theta)
x2 <- x0 + r2*cos(theta)
y2 <- y0 + r2*sin(theta)
x3 <- x0 + r3*cos(theta)
y3 <- y0 + r3*sin(theta)
circle <- data.frame(c(x1,x2,x3),c(y1,y2,y3))
N <-length(x1) # Tamaño del vector
color <- c() # Inicializo vector
color[1:N] <- 'gold'
color[(N+1):(2*N)] <- 'darkorange'
color[(2*N+1):(3*N)] <- 'red'
{PlotOnStaticMap(PuebMap, lat = circle[,2],
lon = circle[,1], cex=0.5, pch=19,
destfile = "interm.png",
col= color ,add=TRUE)};
Quedará algo como esto:
Podemos definir un nuevo nivel de zoom:
Por último y con fines de análisis geoespacial, para definir en que radio se encuentra cada una de las Unidades Económicas (UE), empleamos las siguientes líneas:
D <- dim(DF)[1]; D
# Este vector me indica si la Unidad Económca se encuentra en el área del radio 1 o no
DF$PotR1 <- NULL
for(j in 1:D){
if (((DF[j,2] - y0)^2 + (DF[j,3] - x0)^2) <= r1^2){
DF[j,'PotR1'] <- 'dentro'
} else DF[j,'PotR1'] <- 'fuera'
}
# Este vector me indica si la Unidad Económca se encuentra en el área del radio 2 o no
DF$PotR2 <- NULL
for(j in 1:D){
if (((DF[j,2] - y0)^2 + (DF[j,3] - x0)^2) <= r2^2){
DF[j,'PotR2'] <- 'dentro'
} else DF[j,'PotR2'] <- 'fuera'
}
# Este vector me indica si la Unidad Económca se encuentra en el área del radio 3 o no
DF$PotR3 <- NULL
for(j in 1:D){
if (((DF[j,2] - y0)^2 + (DF[j,3] - x0)^2) <= r3^2){
DF[j,'PotR3'] <- 'dentro'
} else DF[j,'PotR3'] <- 'fuera'
}
Finalmente, en este ejemplo, se analiza a las variables cuantitativas en función de su pertenencia o no al área definida por cada radio, tal como el número de posibles empleados en cada radio:
sumER_1 <- 0
for (j in 1:D){
if (DF[j,'PotR1']=='dentro'){
sumER_1 <- DF[j,'pers']+ sumER_1
}
}
sumER_1 # Potencial Radio 1
sumER_2 <- 0
for (j in 1:D){
if (DF[j,'PotR2']=='dentro'){
sumER_2 <- DF[j,'pers']+ sumER_2
}
}
sumER_2 <- sumER_2 - sumER_1
sumER_2 # Potencial Radio 2
sumER_3 <- 0
for (j in 1:D){
if (DF[j,'PotR3']=='dentro'){
sumER_3 <- DF[j,'pers']+ sumER_3
}
}
sumER_3 <- sumER_3 - sumER_2 - sumER_1
sumER_3 # Potencial Radio 3
sumER_x <- sum(DF['pers']) - sumER_1 - sumER_2 - sumER_3
sumER_x # Potencial fuera de cualquier radio
Información que finalmente representamos en una tabla:
> potencial <- data.frame(rbind(c(sumER_1, sumER_2, sumER_3, sumER_x)
+ )
+ )
> colnames(potencial) <- c(' Radio 1', ' Radio 2', ' Radio 3', ' Radio 4')
> rownames(potencial) <- 'Empleados'
> round(potencial)
Radio 1 Radio 2 Radio 3 Radio 4
Empleados 1539 372 398 128
Lo que nos indica que existen aproximadamente 1539 empleados en el radio 1, que enfatizemos es de 44 Km a partir de la coordenada central que hemos fijado.
