6  Планиране на диета

6.1 Диета

Първият проблем, решен с използването на симплекс алгоритъма за решаване на линейни оптимизационни задачи, е била задачата за оптимално съставяне на диета, която да отговаря на определени хранителни изисквания и да бъде с минимални разходи. Нека да разгледаме следната задача:

Food Порция kCal Протеин (гр.) Калций (mg) Цена на порция (EUR)
Овесена каша 28 г 110 4 2 1.50
Пилешко месо 100 г 205 32 12 4.25
Яйца две 160 13 54 3.0
Мляко 237 мл 160 8 285 1.20
Черешов пай 170 г 420 4 22 5.8
Свинско месо 260 г 260 14 80 8.70

От характеристиките на различните цени искаме да намерим такава комбинация, че калорийното съдържание на менюто да е поне 2 000 килокалории, да съдържа поне 55 гр. протеини и 0.8 гр. калций.

Формулирайте линейна оптимизационна задача, която да намери най-евтиното меню и намерете оптималното меню с Excel.

В допълнение, въведете ограничение за максимален брой порции от 4 от всяко ястие и намерете оптималното меню отново. Сверете решението си с Таблица 6.1.

Нека x_1, x_2, \ldots x_6 са количествата от всяка храна, които ще съставят менюто. Стойността на менюто е равна на количеството от всяка храна, умножена по цената й в евро.

\min z(x_1,\ldots,x_6) = 1.5x_1 + 4.25 x_2 + \ldots + 8.70 x_6

Ограниченията за менюто са изискванията за минимална енергийна стойност (килокалории), както и съдържанието на протеини и калций. При съставянето на неравенствата обърнете внимание на мерните единици: изискването за съдържание на калций е дадено в грамове, а коефициентите в таблицата са в милиграмове.

\begin{align*} 110 x_1 + 205 x_2 + ... \ldots 260 x_6 & \geq 2000 \quad \text{килокалории}\\ 4 x_1 + 32 x_2 + \ldots 14 x_6 &\geq 55 \quad \text{протеини} \\ 2 x_1 + 12 x_2 + \ldots 80 x_6 & \geq 1000 \cdot 0.8 \quad \text{калций} \\ \end{align*}

Променливите трябва да са неотрицателни:

x_1, x_2, \ldots, x_6 \geq 0

При допълнителното ограничение за максимум четири порции от всяка храна имаме:

x_1 \leq 4, x_2 \leq 4, \ldots x_6 \leq 4

Решение в Excel

Покажи
import gurobipy as gp
from gurobipy import GRB
import pandas as pd

m = gp.Model("Diet")
m.Params.LogToConsole = 0

# Create variables

v_names = ["Oatmeal", "Chicken", "Eggs", "Whole Milk", "Cherry pie", "Pork"]

x = m.addVars(len(v_names), name=v_names, lb=0, ub=4.0)

c = [1.5, 4.25, 3.0, 1.20, 5.8, 8.7]

# Set objective

m.setObjective(gp.quicksum(c[i] * x[i] for i in range(len(v_names))), GRB.MINIMIZE)

# Constraints

# Calories

m.addConstr(110 * x[0] + 205 * x[1] + 160 * x[2] + 160 * x[3] + 420 * x[4] + 260 * x[5] >= 2500, "Calories")

# Protein

m.addConstr(4 * x[0] + 32 * x[1] + 13 * x[2] + 8 * x[3] + 4 * x[4] + 14 * x[5] >= 80, "Protein")

# Calcium

m.addConstr(2 * x[0] + 12 * x[1] + 54 * x[2] + 285 * x[3] + 22 * x[4] + 80 * x[5] >= 800, "Calcium")

m.optimize()

vars_df = pd.DataFrame(
    data=[(var.varName, var.x, var.RC) for var in m.getVars()],
    columns=["Variable", "Value", "RC"]
    )
constr_df = pd.DataFrame(
    data=[(constr.constrName, constr.slack, constr.pi) for constr in m.getConstrs()],
    columns=["Constraint", "Slack", "Dual"]
)

# # Write model to file for inspection
# m.write("diet.lp")

# # Print the model file
# with open("diet.lp", "r") as f:
#     print(f.read())
Restricted license - for non-production use only - expires 2027-11-29
Set parameter LogToConsole to value 0
Покажи
vars_df
Таблица 6.1: Оптимално меню (с ограничение за максимум 4 бройки от свяка храна)
Variable Value RC
0 Oatmeal 4.000000 -0.158479
1 Chicken 0.614897 0.000000
2 Eggs 0.000000 0.248494
3 Whole Milk 4.000000 -1.315372
4 Cherry pie 3.080824 0.000000
5 Pork 0.000000 4.565293

6.2 Диета с най-ниска цена в Нова Зеландия

  • Източник: https://doi.org/10.1093/cdn/nzab132

  • Данни за храни

  • Данни за хранителните изисквания

  • ID: Уникален идентификатор на храната

  • Name: Име на храната

  • Group: Група, към която принадлежи храната

  • Cost: Цена на 100 g от храната

  • Amount: Максимално количество за консумация (в 100 g)

  • Dry matter: Сухо вещество на храната (в g)

  • Energy, total metabolisable (kJ): Енергийно съдържание на храната (в kJ)

  • Fat, total: Общо съдържание на мазнини (в g)

  • Protein, total; calculated from total nitrogen (в g)

  • Total carbohydrates by summation (в g)

  • Други характеристики

Покажи
import pandas as pd
# Прочитане на данните

df = pd.read_csv('https://raw.githubusercontent.com/febse/data/refs/heads/main/opt/NZFoods.csv').set_index('ID')
df.head()
Source Group SubGroup Cost Amount Chapter Name Description Alanine Alpha-carotene ... Vitamin A, retinol equivalents Vitamin B12 Vitamin B6 Vitamin C Vitamin D; calculated by summation Vitamin E, alpha-tocopherol equivalents Vitamin K Water Zeaxanthin Zinc
ID
M1170 Animal Beef Beef blade 2.042333 1.0 M Beef, forequarter bolar, separable lean & fat,... Beef, forequarter bolar, separable lean & fat,... 1505.04 0.0 ... 15.95 1.53 0.18 0.0 0.23 0.83 0.0 57.57 0.0 5.93
M1154 Animal Beef Beef blade 2.042333 1.0 M Beef, forequarter bolar, separable lean & fat,... Beef, forequarter bolar, separable lean & fat,... 1036.61 0.0 ... 8.88 1.71 0.32 0.0 0.14 0.45 0.0 70.26 0.0 3.46
M1069 Animal Beef Beef blade 2.042333 1.0 M Beef, forequarter bolar, separable lean, braised Beef, forequarter bolar, separable lean, braised 1627.87 0.0 ... 13.49 1.65 0.19 0.0 0.21 0.76 0.0 60.46 0.0 6.38
M1044 Animal Beef Beef blade 2.042333 1.0 M Beef, forequarter bolar, separable lean, raw Beef, forequarter bolar, separable lean, raw 1072.58 0.0 ... 6.25 1.80 0.33 0.0 0.13 0.40 0.0 73.11 0.0 3.63
M1171 Animal Beef Beef blade 2.042333 1.0 M Beef, forequarter brisket navel end, separable... Beef, forequarter brisket navel end, separable... 870.80 0.0 ... 27.91 0.84 0.09 0.0 0.31 1.11 0.0 42.28 0.0 3.70

5 rows × 145 columns

Покажи
# Einlesen der Ernährungsrichtlinien

reqs = pd.read_csv('https://raw.githubusercontent.com/febse/data/refs/heads/main/opt/NZPersonReq.csv')
reqs.head()
Alanine Alpha-carotene Alpha-tocopherol Arginine Asparagine Available carbohydrate by difference Beta-carotene Beta-carotene equivalents Beta-tocopherol Beta-tocopherol + Gamma-tocopherol ... Vitamin A, retinol equivalents Vitamin B12 Vitamin B6 Vitamin C Vitamin D; calculated by summation Vitamin E, alpha-tocopherol equivalents Vitamin K Water Zeaxanthin Zinc
0 0 0 0 0 0 0 0 0 0 0 ... 800 2.400000e+00 1.3 45 5 8.5 65 0 0 11
1 99999999999 99999999999 99999999999 99999999999 99999999999 99999999999 99999999999 99999999999 99999999999 99999999999 ... 3000 1.000000e+11 50.0 1000 80 300.0 99999999999 99999999999 99999999999 40

2 rows × 134 columns

Покажи
# Данни за съдържанието на определена храна

df.loc["M1170"]
Source                                          Animal
Group                                             Beef
SubGroup                                   Beef blade 
Cost                                          2.042333
Amount                                             1.0
                                              ...     
Vitamin E, alpha-tocopherol equivalents           0.83
Vitamin K                                          0.0
Water                                            57.57
Zeaxanthin                                         0.0
Zinc                                              5.93
Name: M1170, Length: 145, dtype: object
Покажи
import gurobipy as gp
from gurobipy import GRB

m = gp.Model("Diet")
m.Params.LogToConsole = 0

# Променливи

x = m.addVars(df.index, name="x")

# Целева функция

m.setObjective(x.prod(df['Cost'].to_dict()), GRB.MINIMIZE)

# Ограничения

for col in reqs.columns:
    m.addConstr(x.prod(df[col].to_dict()) >= reqs.loc[0, col], col)
    m.addConstr(x.prod(df[col].to_dict()) <= reqs.loc[1, col], col)

m.optimize()

print(f"Total Cost: {m.objVal}")
Set parameter LogToConsole to value 0
Total Cost: 2.804167442632024
Покажи
df["Optimal"] = [x[i].x for i in df.index]
df[df["Optimal"] > 0][["Name", "Group", "Optimal", "Amount", "Cost"]]
/tmp/ipykernel_3812/968555577.py:1: PerformanceWarning:

DataFrame is highly fragmented.  This is usually the result of calling `frame.insert` many times, which has poor performance.  Consider joining all columns at once using pd.concat(axis=1) instead. To get a de-fragmented frame, use `newframe = frame.copy()`
Name Group Optimal Amount Cost
ID
T1026 Mussel, green, meat, boiled Seafood 0.148860 1.0000 0.475000
F1026 Milk, cow, standard 3.3% fat, fluid, Christchu... Milk 3.500815 2.5750 0.173746
X23 Cabbage, leafy vegetable, red Vegetable 0.446254 0.6000 0.209083
X1039 Cabbage, leafy vegetable, raw, tat soi Vegetable 1.307076 0.6000 0.209083
X1137 Potato, root vegetable,  flesh & skin,  baked ... Vegetable 0.366369 1.3500 0.223250
X75 Chickpea, seed vegetable, dried Legumes 1.071478 1.3500 0.523571
A21 Bread, white, toasted Grain 0.346834 0.2600 0.218750
E17 Flour, wheat, white, standard, upper North Island Grain 4.672461 0.5000 0.135167
J1035 Oil, soya bean Fats and oils 0.137830 0.1380 0.360326
J1028 Margarine, original, Flora, fortified vitamin D Fats and oils 0.111599 0.1395 0.454167
J1006 Margarine, polyunsaturated, 70% fat, reduced s... Fats and oils 0.030950 0.1395 0.454167
Q9 Coconut, desiccated Nuts 0.023885 0.3000 0.899091
Q43 Seed, sunflower, dry roasted, no salt added Seeds 0.221816 0.3000 1.228571
Покажи
# Analyse der Lösung

df["EnergyOptimal"] = df["Energy, total metabolisable (kJ)"] * df["Optimal"]

# Die Fettmenge in der optimalen Lösung berechnen

df["FatOptimal"] = df["Fat, total"] * df["Optimal"]

# Die Proteinmenge in der optimalen Lösung berechnen

df["ProteinOptimal"] = df["Protein, total; calculated from total nitrogen"] * df["Optimal"]

# Die Kohlenhydratmenge in der optimalen Lösung berechnen

df["CarbohydrateOptimal"] = df["Total carbohydrates by summation"] * df["Optimal"]

df.aggregate({"EnergyOptimal": "sum", "FatOptimal": "sum", "ProteinOptimal": "sum", "CarbohydrateOptimal": "sum"})
/tmp/ipykernel_3812/1387002755.py:3: PerformanceWarning:

DataFrame is highly fragmented.  This is usually the result of calling `frame.insert` many times, which has poor performance.  Consider joining all columns at once using pd.concat(axis=1) instead. To get a de-fragmented frame, use `newframe = frame.copy()`

/tmp/ipykernel_3812/1387002755.py:7: PerformanceWarning:

DataFrame is highly fragmented.  This is usually the result of calling `frame.insert` many times, which has poor performance.  Consider joining all columns at once using pd.concat(axis=1) instead. To get a de-fragmented frame, use `newframe = frame.copy()`

/tmp/ipykernel_3812/1387002755.py:11: PerformanceWarning:

DataFrame is highly fragmented.  This is usually the result of calling `frame.insert` many times, which has poor performance.  Consider joining all columns at once using pd.concat(axis=1) instead. To get a de-fragmented frame, use `newframe = frame.copy()`

/tmp/ipykernel_3812/1387002755.py:15: PerformanceWarning:

DataFrame is highly fragmented.  This is usually the result of calling `frame.insert` many times, which has poor performance.  Consider joining all columns at once using pd.concat(axis=1) instead. To get a de-fragmented frame, use `newframe = frame.copy()`
EnergyOptimal          11150.000000
FatOptimal                58.491066
ProteinOptimal            92.515608
CarbohydrateOptimal      475.393903
dtype: float64
Покажи
# Активни ограничения

constr_df = pd.DataFrame(
    [(c.ConstrName, c.Sense, c.Slack, c.Pi, c.SARHSLow, c.SARHSUp) for c in m.getConstrs() if c.pi != 0],
    columns=["Name", "Sense", "Slack", "Shadow", "Lower", "Upper"]
    )
constr_df
Name Sense Slack Shadow Lower Upper
0 Biotin > 0.0 0.006518 25.646329 61.741632
1 Calcium > 0.0 0.000927 967.712384 1179.545899
2 Dietary folate equivalents < 0.0 -0.000277 789.784567 1060.486688
3 Energy, total metabolisable (kJ) > 0.0 0.000006 10102.074880 11150.000000
4 Fatty acid 18:3 omega-3 > 0.0 0.026536 0.972971 1.362049
5 Fatty acid cis,cis 18:2 omega-6 > 0.0 0.004134 7.506915 11.256698
6 Pantothenic acid > 0.0 0.109966 4.125095 5.216972
7 Potassium > 0.0 0.000258 3051.327614 3931.002695
8 Selenium > 0.0 0.004293 63.649627 72.645262
9 Sodium > 0.0 0.000071 380.892200 1611.752799
10 Vitamin B12 > 0.0 0.000896 1.807678 2.512669
11 Vitamin C > 0.0 0.000879 21.228732 108.700819
12 Vitamin D; calculated by summation > 0.0 0.012653 3.793052 10.318384