\n",
"\n",
"\n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
"

\n",
"

"
],
"text/plain": [
" probability demand\n",
"sunny skies 0.1 650.0\n",
"good weather 0.6 400.0\n",
"poor weather 0.3 200.0"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Expected demand = 365.0\n"
]
}
],
"source": [
"# price information\n",
"r = 40\n",
"c = 12\n",
"w = 2\n",
"\n",
"# scenario information\n",
"scenarios = {\n",
" \"sunny skies\": {\"probability\": 0.10, \"demand\": 650},\n",
" \"good weather\": {\"probability\": 0.60, \"demand\": 400},\n",
" \"poor weather\": {\"probability\": 0.30, \"demand\": 200},\n",
"}\n",
"\n",
"df = pd.DataFrame.from_dict(scenarios).T\n",
"display(df)\n",
"\n",
"expected_demand = sum(df[\"probability\"] * df[\"demand\"])\n",
"print(f\"Expected demand = {expected_demand}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Subsequent calculations to obtain the EVM can be done directly within the pandas dataframe holding the scenario data."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"probability | demand | |
---|---|---|

sunny skies | 0.1 | 650.0 |

good weather | 0.6 | 400.0 |

poor weather | 0.3 | 200.0 |

\n",
"\n",
"\n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
"

\n",
"

"
],
"text/plain": [
" probability demand order sold salvage profit\n",
"sunny skies 0.1 650.0 365.0 365.0 0.0 10220.0\n",
"good weather 0.6 400.0 365.0 365.0 0.0 10220.0\n",
"poor weather 0.3 200.0 365.0 200.0 165.0 3950.0"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Expected value of the mean demand (EVM) = 8339.0\n"
]
}
],
"source": [
"df[\"order\"] = expected_demand\n",
"df[\"sold\"] = df[[\"demand\", \"order\"]].min(axis=1)\n",
"df[\"salvage\"] = df[\"order\"] - df[\"sold\"]\n",
"df[\"profit\"] = r * df[\"sold\"] + w * df[\"salvage\"] - c * df[\"order\"]\n",
"\n",
"EVM = sum(df[\"probability\"] * df[\"profit\"])\n",
"display(df)\n",
"print(f\"Expected value of the mean demand (EVM) = {EVM}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"No scenario shows a profit loss, which appears to be a satisfactory outcome. However, can we find an order resulting in a higher expected profit?"
]
},
{
"cell_type": "markdown",
"metadata": {
"tags": []
},
"source": [
"## Value of the stochastic solution (VSS)\n",
"\n",
"In order to answer this question, let us formulate the problem in mathematical terms. Let $x$ be a non-negative number representing the number of items that will be ordered, and $y_s$ be the non-negative variable describing the number of items sold in scenario $s$ in the set $S$ comprising all scenarios under consideration. The number $y_s$ of sold items is the lesser of the demand $d_s$ and the order size $x$, that is\n",
"\n",
"$$\n",
"\\begin{align*}\n",
"y_s & = \\min(d_s, x) & \\forall s \\in S.\n",
"\\end{align*}\n",
"$$\n",
"\n",
"Any unsold inventory $x - y_s$ remaining after the event will be sold at the salvage price $w$. Taking into account the revenue from sales $r y_s$, the salvage value of the unsold inventory $w(x - y_s)$, and the cost of the order $c x$, the profit $f_s$ for scenario $s$ is given by\n",
"\n",
"$$\n",
"\\begin{align*}\n",
"f_s & = r y_s + w (x - y_s) - c x & \\forall s \\in S\n",
"\\end{align*}\n",
"$$\n",
"\n",
"Using the constants introduced earlier, the profit $f_s$ for scenario $s \\in S$ can then be written as \n",
"\n",
"$$\n",
" f_s = \\underbrace{r y_s}_\\text{sales revenue} + \\underbrace{w (d_s - y_s)}_\\text{salvage value} - \\underbrace{c x}_\\text{order cost}.\n",
"$$\n",
"\n",
"The expected profit is given by $\\mathbb E(F) = \\sum_s p_s f_s$. Operationally, $y_s$ can be no larger the number of items ordered, $x$, or the demand under scenario $s$, $d_s$. \n",
"The optimization problem is to find the order size $x$ that maximizes expected profit subject to operational constraints on the decision variables. The variables $x$ and $y_s$ are non-negative integers, while $f_s$ is a real number that can take either positive or negative values. Putting these facts together, the optimization problem to be solved is\n",
"\n",
"$$\n",
"\\begin{align*}\n",
" \\text{EV} = \\max \\quad & \\mathbb E(F) = \\sum_{s\\in S} p_s f_s \\\\\n",
" \\text{s.t.} \\quad \n",
" &f_s = r y_s + w(d_s - y_s) - c x & \\forall s \\in S\\\\\n",
" &y_s \\leq x & \\forall s \\in S \\\\\n",
" &y_s \\leq d_s & \\forall s \\in S\\\\\n",
" &y_s \\in \\mathbb{Z}_+ & \\forall s \\in S\\\\\n",
" &x \\in \\mathbb{Z}_+,\n",
"\\end{align*}\n",
"$$\n",
"where $S$ is the set of all scenarios under consideration.\n",
"\n",
"We can implement this problem in AMPL as follows."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Overwriting pop.mod\n"
]
}
],
"source": [
"%%writefile pop.mod\n",
"\n",
"param r;\n",
"param c;\n",
"param w;\n",
"\n",
"# set of scenarios\n",
"set S;\n",
"\n",
"param p{S};\n",
"param d{S};\n",
"\n",
"# decision variables\n",
"var x >= 0;\n",
"var y{S} >= 0;\n",
"var f{S};\n",
"\n",
"# objective\n",
"maximize EV: sum{s in S} p[s] * f[s];\n",
"\n",
"# constraints\n",
"s.t. profit {s in S}: f[s] == r * y[s] + w * (x - y[s]) - c * x;\n",
"s.t. sales_less_than_order {s in S}: y[s] <= x;\n",
"s.t. sales_less_than_demand {s in S}: y[s] <= d[s];"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"cbc 2.10.7: cbc 2.10.7: optimal solution; objective 8920\n",
"0 simplex iterations\n",
"Solver Termination Condition: solved\n",
"\n"
]
},
{
"data": {
"text/html": [
"probability | demand | order | sold | salvage | profit | |
---|---|---|---|---|---|---|

sunny skies | 0.1 | 650.0 | 365.0 | 365.0 | 0.0 | 10220.0 |

good weather | 0.6 | 400.0 | 365.0 | 365.0 | 0.0 | 10220.0 |

poor weather | 0.3 | 200.0 | 365.0 | 200.0 | 165.0 | 3950.0 |

\n",
"\n",
"\n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
"

\n",
"

"
],
"text/plain": [
" demand probability order sold salvage profit\n",
"sunny skies 650.0 0.1 400.0 400 0.0 11200\n",
"good weather 400.0 0.6 400.0 400 0.0 11200\n",
"poor weather 200.0 0.3 400.0 200 200.0 3600"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Expected Profit: 8920.0\n"
]
}
],
"source": [
"# price and scenario information\n",
"r = 40\n",
"c = 12\n",
"w = 2\n",
"\n",
"scenarios = {\n",
" \"sunny skies\": {\"demand\": 650, \"probability\": 0.1},\n",
" \"good weather\": {\"demand\": 400, \"probability\": 0.6},\n",
" \"poor weather\": {\"demand\": 200, \"probability\": 0.3},\n",
"}\n",
"\n",
"# create a data frame with the data and rename the columns to match the model\n",
"scenarios_df = pd.DataFrame.from_dict(scenarios).T.rename(\n",
" columns={\"demand\": \"d\", \"probability\": \"p\"}\n",
")\n",
"\n",
"# Create AMPL instance and load the model\n",
"ampl = AMPL()\n",
"ampl.read(\"pop.mod\")\n",
"\n",
"# load the data\n",
"ampl.param[\"r\"] = r\n",
"ampl.param[\"c\"] = c\n",
"ampl.param[\"w\"] = w\n",
"\n",
"ampl.set_data(scenarios_df, \"S\")\n",
"\n",
"# solve the problem\n",
"ampl.option[\"solver\"] = SOLVER\n",
"ampl.solve()\n",
"\n",
"print(\"Solver Termination Condition:\", ampl.get_value(\"solve_result\"))\n",
"print()\n",
"\n",
"# display solution using Pandas\n",
"df = pd.DataFrame.from_dict(scenarios).T\n",
"df[\"order\"] = ampl.var[\"x\"].value()\n",
"df[\"sold\"] = ampl.var[\"y\"].get_values().toPandas()\n",
"df[\"salvage\"] = df[\"order\"] - df[\"sold\"]\n",
"df[\"profit\"] = ampl.var[\"f\"].get_values().toPandas()\n",
"\n",
"display(df)\n",
"print(\"Expected Profit:\", ampl.obj[\"EV\"].value())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Optimizing over all scenarios provides an expected profit of 8,920 €, an increase of 581 € over the naive strategy of simply ordering the expected number of items sold. The new optimal solution places a larger order, that is $x=400$. In poor weather conditions, there will be more returns and lower profit that is more than compensated by the increased profits in good weather conditions. \n",
"\n",
"The additional value that results from solve of this planning problem is called the **Value of the Stochastic Solution (VSS)**. The value of the stochastic solution is the additional profit compared to ordering to meet the expected demand. In this case,\n",
"\n",
"$$\\text{VSS} = \\text{EV} - \\text{EVM} = 8,920 - 8,339 = 581.$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Expected value with perfect information (EVPI)\n",
"\n",
"Maximizing expected profit requires the size of the order be decided before knowing what scenario will unfold. The decision for $x$ has to be made \"here and now\" with probablistic information about the future, but without specific information on which future will actually transpire.\n",
"\n",
"Nevertheless, we can perform the hypothetical calculation of what profit would be realized if we could know the future. We are still subject to the variability of weather, what is different is we know what the weather will be at the time the order is placed. \n",
"\n",
"The resulting value for the expected profit is called the **Expected Value of Perfect Information (EVPI)**. The difference EVPI - EV is the extra profit due to having perfect knowledge of the future.\n",
"\n",
"To compute the expected profit with perfect information, we let the order variable $x$ be indexed by the subsequent scenario that will unfold. Given decision varaible $x_s$, the model for EVPI becomes\n",
"\n",
"$$\n",
"\\begin{align*}\n",
"\\text{EVPI} = \\max_{x_s, y_s} \\quad & \\mathbb E[f] = \\sum_{s\\in S} p_s f_s \\\\\n",
"\\text{s.t.} \\quad\n",
"& f_s = r y_s + w(x_s - y_s) - c x_s & \\forall s \\in S\\\\\n",
"& y_s \\leq x_s & \\forall s \\in S \\\\\n",
"& y_s \\leq d_s & \\forall s \\in S\n",
"\\end{align*}\n",
"$$\n",
"\n",
"The following implementation is a variation of the prior cell."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Overwriting pop_evpi.mod\n"
]
}
],
"source": [
"%%writefile pop_evpi.mod\n",
"\n",
"param r;\n",
"param c;\n",
"param w;\n",
"\n",
"# set of scenarios\n",
"set S;\n",
"\n",
"param p{S};\n",
"param d{S};\n",
"\n",
"# decision variables\n",
"var x{S} >= 0;\n",
"var y{S} >= 0;\n",
"var f{S};\n",
"\n",
"# objective\n",
"maximize EV: sum{s in S} p[s] * f[s];\n",
"\n",
"# constraints\n",
"s.t. profit {s in S}: f[s] == r * y[s] + w * (x[s] - y[s]) - c * x[s];\n",
"s.t. sales_less_than_order {s in S}: y[s] <= x[s];\n",
"s.t. sales_less_than_demand {s in S}: y[s] <= d[s];"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"cbc 2.10.7: cbc 2.10.7: optimal solution; objective 10220\n",
"0 simplex iterations\n",
"Solver Termination Condition: solved\n",
"\n"
]
},
{
"data": {
"text/html": [
"demand | probability | order | sold | salvage | profit | |
---|---|---|---|---|---|---|

sunny skies | 650.0 | 0.1 | 400.0 | 400 | 0.0 | 11200 |

good weather | 400.0 | 0.6 | 400.0 | 400 | 0.0 | 11200 |

poor weather | 200.0 | 0.3 | 400.0 | 200 | 200.0 | 3600 |

\n",
"\n",
"\n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
" \n",
"

\n",
"

"
],
"text/plain": [
" demand probability order sold salvage profit\n",
"sunny skies 650.0 0.1 650.0 650 0.0 18200\n",
"good weather 400.0 0.6 400.0 400 0.0 11200\n",
"poor weather 200.0 0.3 200.0 200 -0.0 5600"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Expected Profit: 10220.0\n"
]
}
],
"source": [
"# Create AMPL instance and load the model\n",
"ampl = AMPL()\n",
"ampl.read(\"pop_evpi.mod\")\n",
"\n",
"# load the data\n",
"ampl.param[\"r\"] = r\n",
"ampl.param[\"c\"] = c\n",
"ampl.param[\"w\"] = w\n",
"\n",
"ampl.set_data(scenarios_df, \"S\")\n",
"\n",
"# solve the problem\n",
"ampl.option[\"solver\"] = SOLVER\n",
"ampl.solve()\n",
"\n",
"print(\"Solver Termination Condition:\", ampl.get_value(\"solve_result\"))\n",
"print()\n",
"\n",
"# display solution using Pandas\n",
"df = pd.DataFrame.from_dict(scenarios).T\n",
"df[\"order\"] = ampl.var[\"x\"].get_values().toPandas()\n",
"df[\"sold\"] = ampl.var[\"y\"].get_values().toPandas()\n",
"df[\"salvage\"] = df[\"order\"] - df[\"sold\"]\n",
"df[\"salvage\"] = df[\"salvage\"].round(2)\n",
"df[\"profit\"] = ampl.var[\"f\"].get_values().toPandas()\n",
"\n",
"display(df)\n",
"print(\"Expected Profit:\", ampl.obj[\"EV\"].value())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Summary\n",
"\n",
"To summarize, have computed three different solutions to the problem of order size:\n",
"\n",
"* The expected value of the mean solution (EVM) is the expected profit resulting from ordering the number of items expected to sold under all scenarios. \n",
"\n",
"* The expected value of the stochastic solution (EVSS) is the expected profit found by solving an two-state optimization problem where the order size was the \"here and now\" decision without specific knowledge of which future scenario would transpire.\n",
"\n",
"* The expected value of perfect information (EVPI) is the result of a hypotherical case where knowledge of the future scenario was somehow available when then order had to be placed. \n",
"\n",
"For this example we found\n",
"\n",
"| Solution | Value (€) |\n",
"| :------ | ----: |\n",
"| Expected Value of the Mean Solution (EVM) | 8,399.0 | \n",
"| Expected Value of the Stochastic Solution (EVSS) | 8,920.0 |\n",
"| Expected Value of Perfect Information (EVPI) | 10,220.0 |\n",
"\n",
"These results verify our expectation that\n",
"\n",
"$$\n",
"\\begin{align*}\n",
"EVM \\leq EVSS \\leq EVPI\n",
"\\end{align*}\n",
"$$\n",
"\n",
"The value of the stochastic solution \n",
"\n",
"$$\n",
"\\begin{align*}\n",
"VSS = EVSS - EVM = 581\n",
"\\end{align*}\n",
"$$\n",
"\n",
"The value of perfect information\n",
"\n",
"$$\n",
"\\begin{align*}\n",
"VPI = EVPI - EVSS = 1,300\n",
"\\end{align*}\n",
"$$\n",
"\n",
"\n",
"As one might expect, there is a cost that results from lack of knowledge about an uncertain future."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.12"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
demand | probability | order | sold | salvage | profit | |
---|---|---|---|---|---|---|

sunny skies | 650.0 | 0.1 | 650.0 | 650 | 0.0 | 18200 |

good weather | 400.0 | 0.6 | 400.0 | 400 | 0.0 | 11200 |

poor weather | 200.0 | 0.3 | 200.0 | 200 | -0.0 | 5600 |