Localización plana con orloca

Manuel Munoz-Marquez

2024-02-01

Introducción

En un problema de localización se busca encontrar la ubicación óptima de un servicio, o un conjunto de ellos, de forma que la calidad que dicho servicio presta a un conjunto de puntos de demanda sea, según cierta medida, óptima. Algunos ejemplos de problemas de localización son:

Son numerosos los contextos en los que se plantean problemas de localización, debido a ello la teoría de localización ha sido objeto de gran atención en los últimos años, pudiendo decirse que es un tema de gran actualidad y vigencia. A ello contribuye la aparición de facetas del problema hasta ahora no estudiadas. Por ejemplo, junto a los ya clásicos criterios de minimización de costos, aparecen nuevos criterios: ambientales, sociales, calidad de vida, etc. Estas nuevas vertientes del problema hacen que sea un campo abierto de estudio.

El paquete que se presenta está dedicado a la resolución del problema de localizar un único punto en el plano usando como objetivo la minimización de la suma de las distancias ponderadas a los puntos de demanda. Nuevas versiones del paquete incluirán nuevos modelos de localización.

La clase de objetos loca.p

En un problema de localización plana el conjunto de puntos de demanda viene dado por las coordenadas de dichos puntos. Opcionalmente, se puede asignar a dichos puntos una ponderación, que da mayor importancia a unos puntos que a otros, dado que el objetivo que se considera es minimizar la suma ponderada de las distancias entre el punto de servicio y dicho conjunto demandante. Por ejemplo, si se busca la localización de un hospital comarcal, los puntos de demanda pueden ser las localidades a las que el hospital debe atender y las ponderaciones la población de cada localidad.

Para la resolución de estos problemas se ha definido una clase de objetos designada loca.p, de forma que un objeto loca.p almacena las coordenadas de los puntos de demanda y las ponderaciones de cada uno de los puntos. Cada objeto loca.p tiene tres slots, x e y que almacenan las coordenadas y w que almacena las ponderaciones. Cuando las ponderaciones no se den de forma explícita, se considerará que todos los puntos de demanda tienen igual importancia.

En el resto de esta sección se expondrá la forma de hacer las operaciones básicas con objetos loca.p.

Creación de objetos de clase loca.p

Consideremos un problema de localización en el que el conjunto de puntos de demanda es \((0,0)\), \((4,0)\) y \((2,2)\). Para crear un objeto loca.p que represente a dicho conjunto, se puede hacer llamando a la función constructora usando como argumentos el vector con las coordenadas \(x\) y el vector con las coordenadas \(y\) del conjunto de puntos:

loca.p(c(0, 4, 2), c(0, 0, 2))
#> An object of class "loca.p"
#> Slot "x":
#> [1] 0 4 2
#> 
#> Slot "y":
#> [1] 0 0 2
#> 
#> Slot "w":
#> [1] 1 1 1
#> 
#> Slot "label":
#> [1] ""

o alternativamente:

loca.p(x = c(0, 4, 2), y = c(0, 0, 2))

El constructor tiene dos argumentos opcionales más, el tercero w se usa para especificar un vector de pesos y el cuarto para especificar una etiqueta que se usará para identificar el objeto. Si, usando el mismo conjunto de puntos, se quiere asignar los pesos 1, 1, 3, a dichos puntos y la etiqueta “Problema 1”, se usa:

loca.p(x = c(0, 4, 2), y = c(0, 0, 2), w = c(1, 1, 3), label = "Problema 1")
#> An object of class "loca.p"
#> Slot "x":
#> [1] 0 4 2
#> 
#> Slot "y":
#> [1] 0 0 2
#> 
#> Slot "w":
#> [1] 1 1 3
#> 
#> Slot "label":
#> [1] "Problema 1"

Un objeto loca.p también se puede obtener convirtiendo un objeto data.frame que tenga las columnas x e y, y opcionalmente w. Partiendo del data.frame

d
#>    x y w
#> 1  0 0 1
#> 2 10 0 3
#> 3  2 8 1

se puede construir un objeto loca.p llamando a la función as:

as(d, "loca.p")

o alternativamente:

as.loca.p(d)

Recíprocamente, un objeto loca.p se puede convertir en un objeto data.frame mediante:

p1 <- loca.p(x = c(0, 4, 2), y = c(0, 0, 2), w = c(1, 1, 3), label = "Problema 1")
as(p1, 'data.frame')
#>   x y w
#> 1 0 0 1
#> 2 4 0 1
#> 3 2 2 3

o alternativamente

as.data.frame(p1)

En las conversiones, el slot label del objeto loca.p se almacena como un atributo del objeto data.frame. La etiqueta se puede leer y modificar:

p1@label
#> [1] "Problema 1"
dp1 <- as.data.frame(p1)
attr(dp1, "label")
#> [1] "Problema 1"

Los objetos loca.p también pueden transformase en o construirse desde objetos de tipo matrix.

Generación aleatoria de objetos de clase loca.p

Se pueden crear objetos aleatorios de clase loca.p usando la función rloca.p. El primer argumento, n indica el número de puntos a generar. Por defecto, dichos puntos se generan en el cuadrado unidad \([0,1] \times [0, 1]\). Así, para generar un objeto loca.p con 5 puntos en el cuadrado unidad, se usa:

set.seed(161236)
rloca.p(5)
#> An object of class "loca.p"
#> Slot "x":
#> [1] 0.12405385 0.50263749 0.15097058 0.02051109 0.57587429
#> 
#> Slot "y":
#> [1] 0.2435999 0.5788917 0.7562588 0.1144796 0.2721178
#> 
#> Slot "w":
#> [1] 1 1 1 1 1
#> 
#> Slot "label":
#> [1] ""

Los argumentos xmin, xmax, ymin e ymax permiten especificar el rectángulo en el que se generarán los puntos. Además, la función rloca.p permite especificar la etiqueta para el nuevo objeto. Por ejemplo, para generar los puntos en el rectángulo \([-1, 1] \times [-5,5]\) con etiqueta “Rectángulo” se usa:

rloca.p(5, xmin = -1, xmax = 1, ymin = -5, ymax = 5, label = "Rectángulo")
#> An object of class "loca.p"
#> Slot "x":
#> [1] -0.9452062  0.5769137  0.8819392 -0.6214939 -0.1391079
#> 
#> Slot "y":
#> [1] -3.325587  3.931197 -4.497122 -1.304617 -2.815923
#> 
#> Slot "w":
#> [1] 1 1 1 1 1
#> 
#> Slot "label":
#> [1] "Rectángulo"

Los puntos generados por la función rloca.p se pueden generar en grupos repartidos espacialmente. El argumento groups permite especificar el número de grupos mediante un número o el número de puntos en cada grupo a través de un vector. En este segundo caso, el valor dado al argumento n se ignora. Para generar aleatoriamente un conjunto de demanda con tres grupos de igual tamaño:

rloca.p(9, groups = 3, label = "Tres tamaños iguales")
#> An object of class "loca.p"
#> Slot "x":
#> [1] 1.0435201 1.1954117 1.6470744 0.7868431 0.7648655 0.8508291 0.7502594
#> [8] 1.0788618 1.1116052
#> 
#> Slot "y":
#> [1] 1.1400092 0.6076905 0.8163944 1.3954599 1.3186699 0.9262516 0.9044921
#> [8] 0.9464590 0.7384576
#> 
#> Slot "w":
#> [1] 1 1 1 1 1 1 1 1 1
#> 
#> Slot "label":
#> [1] "Tres tamaños iguales"

para tres grupos de tamaños desiguales:

rloca.p(groups = c(2, 2, 5), label = "Tres tamaños desiguales")
#> An object of class "loca.p"
#> Slot "x":
#> [1] 1.6689680 0.9978501 0.5592493 0.4330953 0.7297446 1.1752629 0.9298473
#> [8] 0.4017964 0.8527443
#> 
#> Slot "y":
#> [1] 0.9097321 0.4319344 1.4223225 1.4137017 0.5544933 1.2903588 0.4965480
#> [8] 1.0362726 0.4947184
#> 
#> Slot "w":
#> [1] 1 1 1 1 1 1 1 1 1
#> 
#> Slot "label":
#> [1] "Tres tamaños desiguales"

Para generar los datos en grupos se genera en primer lugar un desplazamiento del centro de cada grupo y luego se generan los puntos sumando a cada punto el desplazamiento que corresponda a su grupo. Por tal motivo, groups = 1 no es equivalente a no especificar dicho parámetro. El desplazamiento de los centros se puede especificar mediante los argumentos xgmin, xgmax, ygmin e ygmax. Para ilustrar mejor el funcionamiento de la función se puede pintar el resultado:

rl <- rloca.p(60, groups = 3, xmin = -1, xmax = 1, ymin = -1, ymax = 1, xgmin = -10, xgmax = 10, ygmin = -10, ygmax = 10, label = "Tres grupos")
plot(rl)

Resumiendo los datos

Para obtener un resumen numérico de un objeto loca.p se puede usar la función summary:

summary(rl)
#>        label  n      xmin   xwmean     xmax      ymin   ywmean     ymax
#>  Tres grupos 60 -1.909368 3.334343 7.010933 -2.148852 1.072656 4.603418

En el resumen se muestran los valores mínimo, máximo y medio de ambas coordenadas, además de la medias ponderadas de las coordenadas de los puntos para cada componente.

Distancia media ponderada

Dado un objeto loca.p se puede evaluar la distancia ponderada desde un punto dado. Así mismo, se puede evaluar el gradiente de dicha función y se puede resolver el problema de minimizar dicho objetivo.

Evaluación

La función distancia media ponderada se denomina en el paquete distsum. Dado un punto, por ejemplo: \((3, 1)\) se puede evaluar la distancia media ponderada a un objeto loca.p:

pt3 <- loca.p(x = c(0, 4, 2), y = c(0, 0, 2), label = "Tres puntos")
distsum(o = pt3, x = 3, y = 1)
#> [1] 5.990705

También se puede calcular el gradiente de distsum llamado distsumgra:

distsumgra(o = pt3, x = 3, y = 1)
#> [1] 0.9486833 0.3162278

Resolución

Para encontrar la solución óptima al problema de localización anterior se usa la función distsummin:

s <- distsummin(pt3)
s
#> [1] 2.00000 1.15332

Evaluando la función y el gradiente en el punto obtenido

distsum(o = pt3, x = s[1], y = s[2])
#> [1] 5.464102
distsumgra(o = pt3, x = s[1], y = s[2])
#> [1]  3.110246e-07 -8.970172e-04

Como se puede comprobar por el valor del gradiente, la solución encontrada es un óptimo local y al ser la función objetivo convexa un óptimo global.

Estas tres funciones admiten un argumento opcional lp, si se omite este argumento, se utiliza la norma euclídea, es decir, la norma \(l_2\), si se especifica un valor para lp se utilizará la norma \(l_p\) para dicho valor de \(p\).

Obsérvese que si se especifica lp = 2 se utiliza el algoritmo genérico para la norma l_p con \(p\) igual a 2. La utilización del algoritmo genérico requiere un mayor esfuerzo computacional para la resolución del problema, por lo que no es recomendable especificar dicho argumento para usar la norma euclídea.

Dibujando

Tanto los objetos loca.p como la función objetivo pueden representarse en un gráfico. Para la función objetivo se proporciona una representación basada en curvas de nivel y otra en un gráfico 3D.

Dibujar un objeto loca.p

La gráfica de un objeto loca.p consiste en representar en el plano el diagrama de dispersión del conjunto de puntos de demanda usando la función plot:

plot(pt3)

Gráfico de curvas de nivel

El gráfico de curvas de nivel se realiza con la función contour:

contour(pt3)

En el gráfico se puede observar cómo la función alcanza el mínimo en el punto calculado anterioremente. Ampliando:

contour(pt3, xlim = c(1.9, 2.1), ylim = c(1, 1.2), levels = c(5.465, 5.47, 5.475))

Las funciones plot y contour admiten un argumento opcional img que permite especificar un gráfico raster que se usará como fondo del gráfico.

Gráfico en 3D

Análogamente se puede realizar una representación en tres dimensiones usando la función persp:

persp(pt3)

persp(pt3, col = "lightblue", theta = 45, ltheta = 120, shade = 0.75, ticktype = "detailed")

Las tres funciones de representación pasarán los restantes argumentos opcionales a la función genérica plot.

Ejemplo de localización en Andalucía

Se cargan los datos de las capitales andaluzas y se convierte en un objeto de clase loca.p:

data(andalusia)
o <- loca.p(x=andalusia$x[1:8], y=andalusia$y[1:8])

Se calculan los valores límite para el gráfico:

xmin <- min(andalusia$x)
ymin <- min(andalusia$y)
xmax <- max(andalusia$x)
ymax <- max(andalusia$y)

Se carga el mapa de Andalucía y se representan los puntos con el mapa de fondo

file = system.file('img', 'andalusian_provinces.png', package='orloca')
img = readPNG(file)
plot(o, img=img, main='Andalucía', xleft=xmin, ybottom=ymin, xright=xmax, ytop=ymax)

El gráfico de curvas de nivel es:

contour(o, img=img, main='Andalucía', xleft=xmin, ybottom=ymin, xright=xmax, ytop=ymax)

La solución óptima del problema de localización con las 8 capitales, ocho primeras filas, se obtiene:

andalusia.loca.p <- loca.p(andalusia$x[1:8], andalusia$y[1:8])
sol <- distsummin(andalusia.loca.p)
sol
#> [1] -4.610679 37.248691

La solución óptima que proporciona el algoritmo está localizada a unos 35 Km al norte de Antequera. Recuérdese que usualmente se considera a Antequera como el centro geográfico de Andalucía. El gráfico presenta la solución como un punto de color rojo:

contour(o, img=img, main='Andalucía', xleft=xmin, ybottom=ymin, xright=xmax, ytop=ymax)
points(sol[1], sol[2], type='p', col='red')

Por simplicidad en el ejemplo, no se ha tenido en cuenta la curvatura terrestre.