Resolvendo o problema das p-Medianas com Python e Pyomo

Exemplo 1

import pyomo.environ as pyo
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.neighbors import DistanceMetric
df_dados = pd.read_csv("lojas_1.csv", delimiter=";")
df_dados
dist = DistanceMetric.get_metric('euclidean')
dij = df_dados[['X', 'Y']].to_numpy()
dij = dist.pairwise(dij)
pd.DataFrame(dij)
plt.figure(figsize=(12, 6))
plot = sns.scatterplot(data=df_dados, x="X", y="Y", s=100)
for index, row in df_dados.iterrows():
plot.text(row['X']+3, row['Y'], index, horizontalalignment='left')

O modelo

Modelagem e resolução através do Pyomo

modelo = pyo.ConcreteModel()
# Índices:
modelo.M = range(len(dij))
modelo.N = range(len(dij))
# Parâmetros
modelo.d = pyo.Param(modelo.M, modelo.N, initialize=lambda modelo, i, j: dij[i][j])
p = 3
modelo.y = pyo.Var(modelo.N, within=pyo.Binary)
modelo.x = pyo.Var(modelo.M, modelo.N, within=pyo.Binary)
def f_obj(modelo):
return sum(modelo.x[i,j] * modelo.d[i,j] for i in modelo.M for j in modelo.N)
modelo.objetivo = pyo.Objective(rule=f_obj, sense=pyo.minimize)
modelo.restricao_a = pyo.Constraint(expr=sum(modelo.y[j] for j in modelo.N) == p)
modelo.restricao_b = pyo.ConstraintList()for i in modelo.M:
modelo.restricao_b.add(sum(modelo.x[i,j] for j in modelo.N) == 1)
modelo.restricao_c = pyo.ConstraintList()for i in modelo.M:
for j in modelo.N:
modelo.restricao_c.add(modelo.x[i,j] <= modelo.y[j])
resultado = pyo.SolverFactory('glpk').solve(modelo)
print(resultado)
modelo.y.pprint()
list_y = list(modelo.y.keys())
[j for j in list_y if modelo.y[j]() == 1]
dados_modelo = df_dados.copy()
dados_modelo['Mediana'] = [modelo.y[i]() for i in list_y]
dados_modelo
modelo.x.pprint()
list_x = list(modelo.x.keys())
alocacoes = [i for i in list_x if modelo.x[i]() == 1]
alocacoes.sort(key=lambda x:x[0])
print(alocacoes)
medianas = [alocacao[1] for alocacao in alocacoes]
dados_modelo['Alocacao'] = medianas
dados_modelo
dados_modelo['Distancia'] = [dij[alocacao[0], alocacao[1]] for alocacao in alocacoes]
dados_modelo
plt.figure(figsize=(12, 6))
plot = sns.scatterplot(data=dados_modelo,
x="X", y="Y",
hue="Alocacao",
alpha=.7,
s=100,
palette="tab10")
highlights = dados_modelo[dados_modelo.index == dados_modelo.Alocacao]
for index, row in highlights.iterrows():
plot.text(row['X']+3, row['Y'], row['Alocacao'], horizontalalignment='left')
dados_modelo_resumo = dados_modelo.copy()
resumo = dados_modelo_resumo.groupby('Alocacao', as_index=False).agg({"Distancia": "sum"})
resumo
resumo.Distancia.sum()

Exemplo 2

df_dados = pd.read_csv("lojas_2.csv", delimiter=";")
print(df_dados)

O modelo

Modelagem e resolução através do Pyomo

# Declaração do modelo:
modelo2 = pyo.ConcreteModel()
# Matriz de distâncias:
dist = DistanceMetric.get_metric('euclidean')
dij = df_dados[['X', 'Y']].to_numpy()
dij = dist.pairwise(dij)

# Índices:
modelo2.M = range(len(dij))
modelo2.N = range(len(dij))
# Parâmetros:
modelo2.d = pyo.Param(modelo2.M, modelo2.N, initialize=lambda modelo, i, j: dij[i][j])
p = 3
# Vetores de demanda e capacidade
Ci = list(df_dados.Demanda)
Kj = list(df_dados.Capacidade)
modelo2.C = pyo.Param(modelo2.M, initialize=lambda modelo2, i: Ci[i])
modelo2.K = pyo.Param(modelo2.N, initialize=lambda modelo2, j: Kj[j])
# Variáveis de decisão:
modelo2.y = pyo.Var(modelo2.N, within=pyo.Binary)
modelo2.x = pyo.Var(modelo2.M, modelo2.N, within=pyo.Binary)
# Função objetivo:
def f_obj(modelo):
return sum(modelo.C[i] * modelo.x[i,j] * modelo.d[i,j] for i in modelo.M for j in modelo.N)
modelo2.obj = pyo.Objective(rule=f_obj, sense=pyo.minimize)# Sujeito a:modelo2.restricao_a = pyo.Constraint(expr=sum(modelo2.y[j] for j in modelo2.N) == p)modelo2.restricao_b = pyo.ConstraintList()
for i in modelo2.M:
modelo2.restricao_b.add(sum(modelo2.x[i,j] for j in modelo2.N) == 1.0)
modelo2.restricao_c = pyo.ConstraintList()
for j in modelo2.N:
modelo2.restricao_c.add(sum(modelo2.C[i] * modelo2.x[i,j] for i in modelo2.M) <= modelo2.K[j] * modelo2.y[j])
resultado2 = pyo.SolverFactory('glpk').solve(modelo2)
print(resultado2)
# Medianas
list_y = list(modelo2.y.keys())
dados_modelo2 = df_dados.copy()
dados_modelo2['Mediana'] = [modelo2.y[i]() for i in list_y]
# Alocações:
list_x = list(modelo2.x.keys())
alocacoes = [i for i in list_x if modelo2.x[i]() == 1]
alocacoes.sort(key=lambda x:x[0])
medianas = [alocacao[1] for alocacao in alocacoes]
dados_modelo2['Alocacao'] = medianas
# Distância
dados_modelo2['Distancia'] = [dij[alocacao[0], alocacao[1]] for alocacao in alocacoes]
dados_modelo2
dados_modelo2['Distancia_total'] = dados_modelo2['Distancia'] * dados_modelo2['Demanda']
print(dados_modelo2)
plt.figure(figsize=(12, 6))
plot = sns.scatterplot(data=dados_modelo2,
x="X", y="Y",
hue="Alocacao",
size='Demanda',
sizes=(50, 250),
alpha=.7,
palette="tab10")
highlights = dados_modelo2[dados_modelo2.index == dados_modelo2.Alocacao]
for index, row in highlights.iterrows():
plot.text(row['X']+3, row['Y'], row['Alocacao'], horizontalalignment='left')
dados_modelo2_resumo = dados_modelo2.copy()
resumo2 = dados_modelo2_resumo.groupby('Alocacao', as_index=False).agg({"Demanda": "sum",
"Distancia": "sum",
"Distancia_total": "sum"})
sum_dist2 = resumo2.Distancia.sum()
sum_dist_total2 = resumo2.Distancia_total.sum()
print("Distância:", sum_dist2)
print("Distância_total:", sum_dist_total2)
resumo2

Exemplo 3

df_dados = pd.read_csv("lojas_3.csv", delimiter=";")
df_dados

O modelo

Modelagem e resolução através do Pyomo

# Declaração do modelo:
modelo3 = pyo.ConcreteModel()
# Matriz de distâncias:
dist = DistanceMetric.get_metric('euclidean')
dij = df_dados[['X', 'Y']].to_numpy()
dij = dist.pairwise(dij)

# Índices:
modelo3.M = range(len(dij))
modelo3.N = range(len(dij))
# Parâmetros:# Vetores de demanda e capacidade
Ci = list(df_dados.Demanda)
Kj = list(df_dados.Capacidade)
modelo3.d = pyo.Param(modelo3.M, modelo3.N, initialize=lambda modelo, i, j: dij[i][j])
modelo3.C = pyo.Param(modelo3.M, initialize=lambda modelo3, i: Ci[i])
modelo3.K = pyo.Param(modelo3.N, initialize=lambda modelo3, j: Kj[j])
# Vetor de custos:
Sj = list(df_dados.Custo)
modelo3.S = pyo.Param(modelo3.N, initialize=lambda modelo, j: Sj[j])
B = 50000
# Variáveis de decisão:
modelo3.y = pyo.Var(modelo3.N, within=pyo.Binary)
modelo3.x = pyo.Var(modelo3.M, modelo3.N, within=pyo.Binary)
# Função objetivo:
def f_obj(modelo):
return sum(modelo.C[i] * modelo.x[i,j] * modelo.d[i,j] for i in modelo.M for j in modelo.N)
modelo3.obj = pyo.Objective(rule=f_obj, sense=pyo.minimize)# Sujeito a:modelo3.restricao_b = pyo.ConstraintList()
for i in modelo3.M:
modelo3.restricao_b.add(sum(modelo3.x[i,j] for j in modelo3.N) == 1.0)

modelo3.restricao_c = pyo.ConstraintList()
for j in modelo3.N:
modelo3.restricao_c.add(sum(modelo3.C[i] * modelo3.x[i,j] for i in modelo3.M) <= modelo3.K[j] * modelo3.y[j])
modelo3.restricao_a = pyo.Constraint(expr=sum(modelo3.S[j] * modelo3.y[j] for j in modelo3.N) <= B)
resultado3 = pyo.SolverFactory('glpk').solve(modelo3)
print(resultado3)
# Medianas
list_y = list(modelo3.y.keys())
dados_modelo3 = df_dados.copy()
dados_modelo3['Mediana'] = [modelo3.y[i]() for i in list_y]
# Alocações:
list_x = list(modelo3.x.keys())
alocacoes = [i for i in list_x if modelo3.x[i]() == 1]
alocacoes.sort(key=lambda x:x[0])
medianas = [alocacao[1] for alocacao in alocacoes]
dados_modelo3['Alocacao'] = medianas
# Distâncias
dados_modelo3['Distancia'] = [dij[alocacao[0], alocacao[1]] for alocacao in alocacoes]
dados_modelo3['Distancia_total'] = dados_modelo3['Distancia'] * dados_modelo3['Demanda']
dados_modelo3
plt.figure(figsize=(12, 6))
plot = sns.scatterplot(data=dados_modelo3,
x="X", y="Y",
hue="Alocacao",
size='Demanda',
sizes=(50, 250),
alpha=.7,
palette="tab10")
highlights = dados_modelo3[dados_modelo3.index == dados_modelo3.Alocacao]
for index, row in highlights.iterrows():
plot.text(row['X']+3, row['Y'], row['Alocacao'], horizontalalignment='left')
dados_modelo3_resumo = dados_modelo3.copy()
dados_modelo3_resumo['Custo'] = dados_modelo3_resumo['Custo'] * dados_modelo3_resumo['Mediana']
resumo3 = dados_modelo3_resumo.groupby('Alocacao', as_index=False).agg({"Demanda": "sum",
"Distancia": "sum",
"Distancia_total": "sum",
"Custo": "sum"})
sum_dist3 = resumo3.Distancia.sum()
sum_dist_total3 = resumo3.Distancia_total.sum()
sum_custo3 = resumo3.Custo.sum()
print("Distância:", sum_dist3)
print("Distância_total:", sum_dist_total3)
print("Custo:", sum_custo3)
resumo3

Exemplo 4

O modelo

Modelagem e resolução através do Pyomo

# Declaração do modelo:
modelo4 = pyo.ConcreteModel()
# Matriz de distâncias:
dist = DistanceMetric.get_metric('euclidean')
dij = df_dados[['X', 'Y']].to_numpy()
dij = dist.pairwise(dij)

# Índices:
modelo4.M = range(len(dij))
modelo4.N = range(len(dij))
# Parâmetros:# Vetores de demanda, capacidade e custo
Ci = list(df_dados.Demanda)
Kj = list(df_dados.Capacidade)
Sj = list(df_dados.Custo)
modelo4.d = pyo.Param(modelo4.M, modelo4.N, initialize=lambda modelo, i, j: dij[i][j])
modelo4.C = pyo.Param(modelo4.M, initialize=lambda modelo4, i: Ci[i])
modelo4.K = pyo.Param(modelo4.N, initialize=lambda modelo4, j: Kj[j])
modelo4.S = pyo.Param(modelo4.N, initialize=lambda modelo4, j: Sj[j])
B = 70000
T = [3,5,7]
# Variáveis de decisão:
modelo4.y = pyo.Var(modelo4.N, within=pyo.Binary)
modelo4.x = pyo.Var(modelo4.M, modelo4.N, within=pyo.Binary)
# Função objetivo:
def f_obj(modelo):
return sum(modelo.C[i] * modelo.x[i,j] * modelo.d[i,j] for i in modelo.M for j in modelo.N)
modelo4.obj = pyo.Objective(rule=f_obj, sense=pyo.minimize)# Sujeito a:modelo4.restricao_a = pyo.Constraint(expr=sum(modelo4.S[j] * modelo4.y[j] for j in modelo4.N) <= B)modelo4.restricao_b = pyo.ConstraintList()
for i in modelo4.M:
modelo4.restricao_b.add(sum(modelo4.x[i,j] for j in modelo4.N) == 1.0)

modelo4.restricao_c = pyo.ConstraintList()
for j in modelo4.N:
modelo4.restricao_c.add(sum(modelo4.C[i] * modelo4.x[i,j] for i in modelo4.M) <= modelo4.K[j] * modelo4.y[j])
modelo4.restricao_d = pyo.Constraint(expr=sum(modelo4.y[j] for j in modelo4.N if j in T) <= 1)
resultado4 = pyo.SolverFactory('glpk').solve(modelo4)
print(resultado4)
# Medianas
list_y = list(modelo4.y.keys())
dados_modelo4 = df_dados.copy()
dados_modelo4['Mediana'] = [modelo4.y[i]() for i in list_y]
# Alocações:
list_x = list(modelo4.x.keys())
alocacoes = [i for i in list_x if modelo4.x[i]() == 1]
alocacoes.sort(key=lambda x:x[0])
medianas = [alocacao[1] for alocacao in alocacoes]
dados_modelo4['Alocacao'] = medianas
# Distâncias
dados_modelo4['Distancia'] = [dij[alocacao[0], alocacao[1]] for alocacao in alocacoes]
dados_modelo4['Distancia_total'] = dados_modelo4['Distancia'] * dados_modelo4['Demanda']
dados_modelo4
plt.figure(figsize=(12, 6))
plot = sns.scatterplot(data=dados_modelo4,
x="X", y="Y",
hue="Alocacao",
size='Demanda',
sizes=(50, 250),
alpha=.7,
palette="tab10")
highlights = dados_modelo4[dados_modelo4.index == dados_modelo4.Alocacao]
for index, row in highlights.iterrows():
plot.text(row['X']+3, row['Y'], row['Alocacao'], horizontalalignment='left')
dados_modelo4_resumo = dados_modelo4.copy()
dados_modelo4_resumo['Custo'] = dados_modelo4_resumo['Custo'] * dados_modelo4_resumo['Mediana']
resumo4 = dados_modelo4_resumo.groupby('Alocacao', as_index=False).agg({"Demanda": "sum",
"Distancia": "sum",
"Distancia_total": "sum",
"Custo": "sum"})
sum_dist4 = resumo4.Distancia.sum()
sum_dist_total4 = resumo4.Distancia_total.sum()
sum_custo4 = resumo4.Custo.sum()
print("Distância:", sum_dist4)
print("Distância_total:", sum_dist_total4)
print("Custo:", sum_custo4)
resumo4

Conclusões

Referências

--

--

--

Cientista de Dados — https://acsjunior.com

Love podcasts or audiobooks? Learn on the go with our new app.

Get the Medium app

A button that says 'Download on the App Store', and if clicked it will lead you to the iOS App store
A button that says 'Get it on, Google Play', and if clicked it will lead you to the Google Play store
António C. da Silva Júnior

António C. da Silva Júnior

Cientista de Dados — https://acsjunior.com

More from Medium

How To Create Desktop GUI Using Python

Day 4.1 Python — Randomness

[FRAMED]: Sami Chakhachiro, Production Coordinator, MTL

On-Balance Volume Indicator Algorithm in Python