The Diet Problem

Summary

The goal of the Diet Problem is to select foods that satisfy daily nutritional requirements at minimum cost. This problem can be formulated as a linear program, for which constraints limit the number of calories and the amount of vitamins, minerals, fats, sodium, and cholesterol in the diet. Danzig (1990) notes that the diet problem was motivated by the US Army's desire to minimize the cost of feeding GIs in the field while still providing a healthy diet.

Problem Statement

The Diet Problem can be formulated mathematically as a linear programming problem using the following model.

Sets

$F$ = set of foods
$N$ = set of nutrients

Parameters

$c_i$ = cost per serving of food $i$, $\forall i \in F$
$a_{ij}$ = amount of nutrient $j$ in food $i$, $\forall i \in F, \forall j \in N$
$Nmin_j$ = minimum level of nutrient $j$, $\forall j \in N$
$Nmax_j$ = maximum level of nutrient $j$, $\forall j \in N$
$V_i$ = the volume per serving of food $i$, $\forall i \in F$
$Vmax$ = maximum volume of food consumed

Variables

$x_i$ = number of servings of food $i$ to consume

Objective

Minimize the total cost of the food
$\min \sum_{i \in F} c_i x_i$

Constraints

Limit nutrient consumption for each nutrient $j \in N$.
$Nmin_j \leq \sum_{i \in F} a_{ij} x_i \leq Nmax_j$, $\forall j \in N$

Limit the volume of food consumed
$\sum_{i \in F} V_i x_i \leq Vmax$

Consumption lower bound
$x_i \geq 0$, $\forall i \in F$

In [1]:
import matplotlib.pyplot as plt
import numpy as np
from optlang import Model, Variable, Constraint, Objective
In [2]:
alimentos = ["Cheeseburger",
             "Ham Sandwich",
             "Hamburger",
             "Fish Sandwich",
             "Chicken Sandwich",
             "Fries",
             "Sausage Biscuit",
             "Lowfat Milk",
             "Orange Juice"]
cs = [1.84,2.19,1.82,1.44,2.29,.77,1.29,.60,.72]
Vs = [4.0 , 7.5, 3.5, 5.0, 7.3,2.6, 4.1,8.0,12.0]
nutrientes = ["Cal", "Carbo", "Protein", "VitA", "VitC", "Calc", "Iron"]
A = np.array([[510,34,28,15,6,30,20],
              [370,35,24,15,10,20,20],
              [500,42,25,6,2,25,20],
              [370,38,14,2,0,15,10],
              [400,42,31,8,15,15,8],
              [220,26,3,0,15,0,2],
              [345,27,15,4,0,20,15],
              [110,12,9,10,4,30,0],
              [80,20,1,2,120,2,2]
])
Vmax = 75
Nmin = [2000,350,55,100,100,100,100]
#Los 0 quieren decir que no hay cantidad maxima
Nmax = [0,375,0,0,0,0,0]

Nalimentos = len(alimentos)
Nnutrientes = len(nutrientes)

Plantear el problema

In [3]:
x = [Variable('N_'+alimento.replace(' ','_'), lb=0, type='integer') 
     for alimento in alimentos]

obj = Objective(sum(x[j]*cs[j] for j in range(Nalimentos) ), direction='min')

constraints = [
    #Volumen maximo de la nevera
    Constraint(sum(x[j]*Vs[j] for j in range(Nalimentos) ), ub=Vmax),
]
for k in range(Nnutrientes):
    #Cantidad minima de cada nutriente
    constraints.append(
        Constraint(sum(x[j]*A[j,k] for j in range(Nalimentos) ), lb=Nmin[k])
    )
    #Cantidad maxima de cada nutriente
    if Nmax[k] > Nmin[k]:
        constraints.append(
        Constraint(sum(x[j]*A[j,k] for j in range(Nalimentos) ), ub=Nmax[k])
    )

model = Model(name='Diet problem')
model.objective = obj
model.add(constraints)

status = model.optimize()

print("status:", model.status)
print("objective value:", model.objective.value)

for var_name, var in model.variables.iteritems():
    print(var_name, "=", var.primal)
status: optimal
objective value: 15.05
N_Cheeseburger = 4.0
N_Chicken_Sandwich = 0.0
N_Fish_Sandwich = 1.0
N_Fries = 5.0
N_Ham_Sandwich = 0.0
N_Hamburger = 0.0
N_Lowfat_Milk = 4.0
N_Orange_Juice = 0.0
N_Sausage_Biscuit = 0.0
In [4]:
#Escribimos el coste, volumen total, 
print("cost:",sum(x[j].primal*cs[j] for j in range(Nalimentos)))
print("volume:",sum(x[j].primal*Vs[j] for j in range(Nalimentos)))
for k in range(Nnutrientes):
    print(nutrientes[k])
    print(sum(x[j].primal*A[j,k] for j in range(Nalimentos) ))
cost: 15.05
volume: 66.0
Cal
3950.0
Carbo
352.0
Protein
177.0
VitA
102.0
VitC
115.0
Calc
255.0
Iron
100.0
In [ ]:
 
In [ ]:
 

Ejercicio 1

¿Cuál es el resultado si eliminamos la restricción de que las cantidades tienen que ser números enteros? Es decir, si nos permiten comprar el 70% de un sandwich, etc.

In [ ]:
 
In [ ]:
 
In [ ]:
 
In [ ]:
 

Ejercicio 2

¿Cuál es el coste de la dieta si ponemos una cantidad máxima de calorías de 3500?

In [ ]:
 
In [ ]:
 

Ejercicio 4

Visto que la dieta óptima no satisface la restricción de volumen máximo, nos planteamos ahorrar reduciendo el tamaño de la nevera:

¿Hasta qué tamaño podemos reducir la nevera si estamos dispuestos a asumir un coste de 16 euros? Es suficiente con una respuesta aproximada.

In [ ]:
 
In [ ]:
 

Ejercicio 5

Podemos conseguir patatas a 10 céntimos el kg. Lee las propiedades nutricionales de las patatas y estudia cómo afectan al coste total de la dieta:

In [ ]:
 
In [ ]:
 

Ejercicio 6

Busca en páginas como:

Añade otros alimentos, investiga cómo afectan al resultado: ¿puedes conseguir una dieta que cumpla las restricciones originales y que además tenga menos de 3000 calorías?

In [ ]:
 
In [ ]:
 

Ejercicio 7

Reflexiona sobre qué tipo de nutriente es el más crítico: Relaja un 10% la cantidad mínima de cada nutriente y observa cómo varía el coste total de la dieta.

In [ ]:
 
In [ ]:
 

pyomo

Hay una librería en python, algo más difícil de usar, que mencionaron en una charla en la Escuela en Septiembre de 2017. Se llama pyomo y la usan en Airbus, por ejemplo. En realidad es un lenguaje de modelado, que puede llamar a varios solvers.

Si te interesa, compara el mismo problema en el lenguaje de modelado pyomo:

http://nbviewer.jupyter.org/github/Pyomo/PyomoGallery/blob/master/diet/DietProblem.ipynb

Hemos instalado el software necesario para que puedas ejecutar ese archivo en este mismo servidor.

cvxpy

Otra librería interesante, que también he oído mencionar en esta escuela, y que he usado con buenos resultados, es cvxpy que permite plantear problemas de optimización convexa, y tiene una sintaxis potente y elegante, pero no sirve para problemas MILP como éste.