[PYTHON] Let's decide the position of the fire station by combinatorial optimization

problem

You are the new mayor of City A. I decided to build three new fire stations. Where should I build it? (This story is fiction)

Way of thinking 1

Think of it as Facility placement problem without capacity constraints. It is assumed that the location of the citizen's household is known. Decide to minimize the total distance from each household to the location of the nearest fire station.

Solve with Python

First, randomly determine the household position.

python

%matplotlib inline
import numpy as np, pandas as pd, matplotlib.pyplot as plt
from ortoolpy import facility_location_without_capacity
np.random.seed(1)

Number of households,Number of fire stations= 100,3
Household position= np.random.rand(Number of households,2)
plt.scatter(*Household position.T);

I will try to solve it.

python

Result 1= facility_location_without_capacity(Number of fire stations,Household position)
np.unique(Result 1)
>>>
array([ 9, 62, 77])

We found that it would be better to build the 9th, 62nd, and 77th household positions. Visualize.

python

Fire station location 1=Household position[np.unique(Result 1)]
plt.scatter(*Household position.T)
plt.scatter(*Fire station location 1.T);

Let's look at the distribution of distances for each household. The middle line is the average.

python

Distance 1= np.linalg.norm(Household position-Household position[Result 1], axis=1)
plt.hist(Distance 1, range=(0,0.6), bins=20, alpha=0.5)
plt.vlines(Distance 1.mean(), 0, 20);

Way of thinking 2

If it is the sum of the distances, the people who are far and the people who are close will be offset. Try to reduce the number of people far away as much as possible. Here we minimize the sum of the squares of the distances.

python

Result 2= facility_location_without_capacity(Number of fire stations,Household position,
    func = lambda i,j: (Household position[i][0]-Household position[j][0])**2+(Household position[i][1]-Household position[j][1])**2)
np.unique(Result 2)
>>>
array([37, 39, 73])

It's a different place. If you make a mathematical model with PuLP, it will be as follows. (Same result)

python

from pulp import *
from ortoolpy import addbinvars, addvars
r = range(len(Household position))
m = LpProblem()
x = addvars(len(Household position), len(Household position))
y = addbinvars(len(Household position))
m += lpSum(((Household position[i][0]-Household position[j][0])**2+(Household position[i][1]-Household position[j][1])**2) * x[i][j]
           for i in r for j in r)
m += lpSum(y) <=Number of fire stations
for i in r:
    m += lpSum(x[i]) == 1
    for j in r:
        m += x[i][j] <= y[j]
m.solve()
Result 2= [int(value(lpDot(r, x[i]))) for i in r]

Visualize.

python

Fire station position 2=Household position[np.unique(Result 2)]
plt.scatter(*Household position.T)
plt.scatter(*Fire station location 1.T)
plt.scatter(*Fire station position 2.T);

It feels like it's getting closer to the center. Let's compare the distributions.

python

Distance 2= np.linalg.norm(Household position-Household position[Result 2], axis=1)
plt.hist(Distance 1, range=(0,0.6), bins=20, alpha=0.5)
plt.vlines(Distance 1.mean(), 0, 20)
plt.hist(Distance 2, range=(0,0.6), bins=20, alpha=0.5)
plt.vlines(Distance 2.mean(), 0, 20)
print('average',Distance 1.mean(),Distance 2.mean())
print('Distributed',Distance 1.var(),Distance 2.var())
>>>
Average 0.235140814776 0.237069972634
Variance 0.0138310436529 0.0100843497562

You can see that the average has increased slightly, but the variability has decreased and the number of people farther away has decreased.

Way of thinking 3

Try to minimize your expectations, with a 90% chance of getting to the closest place and a 10% chance of getting to the second closest place. It is OK if you prepare variables for the first and second and do not select the first and second at the same time.

python

m = LpProblem()
x1 = np.array(addvars(len(Household position), len(Household position))) #Closest fire station
x2 = np.array(addvars(len(Household position), len(Household position))) #Second closest fire station
y  = addbinvars(len(Household position))
m += lpSum(((Household position[i][0]-Household position[j][0])**2+(Household position[i][1]-Household position[j][1])**2) * x1[i,j] * 0.9
          +((Household position[i][0]-Household position[j][0])**2+(Household position[i][1]-Household position[j][1])**2) * x2[i,j] * 0.1
           for i in r for j in r)
m += lpSum(y) <=Number of fire stations
for i in r:
    m += lpSum(x1[i]) == 1
    m += lpSum(x2[i]) == 1
    for j in r:
        m += x1[i,j] + x2[i,j] <= y[j]
m.solve()
Result 3= [int(value(lpDot(r, x1[i]))) for i in r]
np.unique(Result 3)
>>>
array([37, 39, 93])

Let's visualize it. It changed a little.

python

Fire station position 3=Household position[np.unique(Result 3)]
for i in [Household position,Fire station location 1,Fire station position 2,Fire station position 3]:
    plt.scatter(*i.T)

that's all