4
$\begingroup$

I have a matrix $M$ with non-negative real entries, and I would like to minimize the objective function $$\Phi(v) = \|Mv\|_\infty,$$ where $v$ is constrained to be a probability vector, i.e., $v_1+\dots+v_n=1$ and $v_i\ge 0$. Is there an efficient algorithm for this?

Motivation. This comes up in solving a two-player game, where the first player chooses $v$ (which represents a probability distribution on the columns of $M$), then the second player sees $v$ and chooses a row $i$, and the first player loses $(Mv)_i$ and the second player gains $(Mv)_i$. I want to find an efficient algorithm to find the optimal strategy for the first player.

Best I can do. I can see how to solve this with linear programming. We introduce a variable $\ell$ to represent the final loss of the first player, then we add the linear inequalities $M_i v \le \ell$ (and the constraints on $v$) and we minimize $\ell$. Is there a more efficient algorithm than this?

$\endgroup$
4
  • $\begingroup$ Do you wish to do this inside of a mathematical programming tool, or can you use arbitrary algorithms and libraries? $\endgroup$ – Richard Jul 26 '20 at 14:24
  • $\begingroup$ Looks like you are trying to solve $\min_v\max_i M_iv$ and that is a tough problem to solve, especially when $M$ is nonsingular. I don't know if you can do any better than linear programming $\endgroup$ – Abdullah Ali Sivas Jul 26 '20 at 17:20
  • $\begingroup$ @Richard, I'm happy to use arbitrary algorithms. Abdullah Ali Sivas, yes, that's exactly the problem. Thank you for the feedback. $\endgroup$ – D.W. Jul 26 '20 at 17:51
  • $\begingroup$ Keep in mind that linear programming solves many problems quickly in practice. What kind of input sizes are you looking at? Do you know LP is too slow? $\endgroup$ – Richard Jul 26 '20 at 19:56
2
$\begingroup$

We can do some transformations to your problem to show that it's easily solvable via a linear program:

Given a matrix $M$ with non-negative real entries and a vector $v$ you wish to solve the problem: $$ \begin{align} \min_v \quad & \lVert Mv \rVert_\infty \\ s.t. \quad & v_i\ge0 \\ & \sum_i v_i = 1 \end{align} $$ Now, note that $\lVert Mv \rVert_\infty = max_i |(Mv)_i|$. With this in hand we get to the problem: $$ \begin{align} \min_v \quad & \max_i |(Mv)_i| \\ s.t. \quad & v_i\ge0 \\ & \sum_i v_i = 1 \end{align} $$ We can treat an absolute value $|x|$ in a min-LP by replacing $|x|$ with a variable $y$ and adding the constraints $x\le y$ and $-x\le y$.

We can replace a maximum function in a min-LP by replacing $\max_i (x)_i$ with a variable $y$ and adding the constraints $x_1\le y, x_2\le y, \ldots, x_n \le y$.

We can therefore rewrite the problem as $$ \begin{align} \min_v \quad & y \\ s.t. \quad & v_i\ge0 \\ & \sum_i v_i = 1 \\ & (M_{1,*} v) \le y \\ & -(M_{1,*} v) \le y \\ & (M_{i,*} v) \le y \\ & -(M_{i,*} v) \le y \\ \end{align} $$ Where $M_{i,*}$ is the $i$-th row of the matrix $M$.

Since this is a convex problem, you can solve it using cvxpy like so:

import cvxpy as cp
import numpy as np

M = np.random.rand(10,10)
v = cp.Variable(10)

objective = cp.Minimize(cp.norm(M*v, 'inf'))
constraints = [sum(v)==1, v>=0]
problem = cp.Problem(objective, constraints)
objval = problem.solve()

print("Objective value = ", objval)
print("v values = ", v.value)

Notice that CVXPY has automagically performed all of the transformations we used above.

Now, efficiency. We can judge this by multiple metrics.

  1. Your time. It's hopefully obvious that CVXPY and similar tools offer an extremely efficient and flexible way of solving this kind of problem. If you need to add constraints later you can do so quickly.
  2. Compute time. LP solvers are often highly optimized. You should expect them to work quickly even with large numbers of variables and constraints.

Let's look at this second point by timing the above:

import cvxpy as cp
import numpy as np
import timeit

M = np.random.rand(1000,1000)
v = cp.Variable(1000)

objective = cp.Minimize(cp.norm(M*v, 'inf'))
constraints = [sum(v)==1, v>=0]
problem = cp.Problem(objective, constraints)
timeit.timeit(lambda: problem.solve(), number=4)

This gives:

     Size | Time
  10x10   | 0.39s
 100x100  | 3.37s
1000x1000 | 345s

A lot of this is Python overhead. If we instead use Julia, we get much better timing:

using Convex
using ECOS

M = rand(10,10);
v = Variable(10);
problem=minimize(norm_inf(M*v), [v>=0, sum(v)==1])
@time solve!(problem, ECOS.Optimizer)

Timing results:

     Size | Time
  10x10   | 0.0033s
 100x100  | 0.07s
1000x1000 | 96s

Much better! Note that we're also using the ECOS solver. Other options, especially commercial ones, might be much faster.

I'm skeptical that other approaches would improve much on the times for smaller (10x10, 100x100) problems, or that you'd be able to make meaningful use of those improvements (outside of some HPC context).

Dynamic programming, as another answer suggested, might also be tricky to implement here. DP alone is slow because the game tree expands exponentially for each additional level of recursion. You make DP fast by memoizing states, but that's impractical if your states are continuous (your problem) or don't overlap (chess, Go).

EDIT:

Brian Borchers comments:

Note that since M has nonnegative entries and v≥0, you don't actually need to handle the absolute values

I'd avoided making use of this information initially in order to provide a fully general answer, but if we do leverage it in Julia:

using Convex
using ECOS

N = 1000
M = rand(N,N);
v = Variable(N);
problem=minimize(maximum(M*v), [v>=0, sum(v)==1])
@time solve!(problem, ECOS.Optimizer)

With this simplification of the constraints, the 1000x1000 problem takes only 19s!

$\endgroup$
4
  • 1
    $\begingroup$ Thank you for the details and analysis! This is basically the same as the approach I outline in "Best I can do" (your $y$ is my $\ell$) except I forgot about the absolute value, and you've developed this in considerably more detail. Thank you! $\endgroup$ – D.W. Dec 27 '20 at 21:29
  • $\begingroup$ Note that since $M$ has nonnegative entries and $v \geq 0$, you don't actually need to handle the absolute values- this greatly simplifies the constraints. $\endgroup$ – Brian Borchers Dec 29 '20 at 22:32
  • $\begingroup$ @BrianBorchers: I'd initially ignored that information to provide a fully general answer, but you're right that it was worth exploring: it gives a 5x speed-up. $\endgroup$ – Richard Dec 30 '20 at 4:53
  • $\begingroup$ @D.W.: See updated answer for a 5x speed-up. $\endgroup$ – Richard Dec 30 '20 at 4:53
1
$\begingroup$

You can and should solve this problem without linear programming and apply the Bellman equation instead.

Actually, the minmax theorem -- handled numerically via LP -- is only required to solve the problem where both players simultaneously choose an action.

In contrast, your game consists of a two-step process, and the mathematical model should incorporate this structure. This can be realized by a Markov decision process that is optimized via the Bellman equation. Basically, you there solve two "max" problems instead of one "minmax" problem, which is way easier from both the mathematical and the computational perspective.

$\endgroup$
2
  • 2
    $\begingroup$ This sounds very interesting -- would you mind elaborating on this a little bit? Can you show me what the Bellman equation would be in this context? I've never heard that term before, but Wikipedia seems to indicate it is the recursive formula used in dynamic programming. I'm familiar with dynamic programming but I'm not sure how we'd obtain a suitable recursive formula in this setting, or how we'd set up the problem so that dynamic programming can be applied. Would you mind expanding on your answer? $\endgroup$ – D.W. Jul 28 '20 at 0:54
  • $\begingroup$ The other answer gave reason to look again at your problem, and I must admit I got in wrong in the first place (and associated this case to a problem I once solved, where DP was much more efficient). Reconsidering it again, what @Richard wrote is completely correct: DP would probably work here, but it will be much more cumbersome than LP. $\endgroup$ – davidhigh Dec 30 '20 at 17:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.