{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"toc": "true"
},
"source": [
"# Table of Contents\n",
" "
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Julia Version 1.1.0\n",
"Commit 80516ca202 (2019-01-21 21:24 UTC)\n",
"Platform Info:\n",
" OS: macOS (x86_64-apple-darwin14.5.0)\n",
" CPU: Intel(R) Core(TM) i7-6920HQ CPU @ 2.90GHz\n",
" WORD_SIZE: 64\n",
" LIBM: libopenlibm\n",
" LLVM: libLLVM-6.0.1 (ORCJIT, skylake)\n",
"Environment:\n",
" JULIA_EDITOR = code\n"
]
}
],
"source": [
"versioninfo()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Optimization Examples - Linear Programming\n",
"\n",
"## Linear programming (LP)\n",
"\n",
"* A general linear program takes the form\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{c}^T \\mathbf{x} \\\\\n",
"\t&\\text{subject to}& \\mathbf{A} \\mathbf{x} = \\mathbf{b} \\\\\n",
"\t& & \\mathbf{G} \\mathbf{x} \\preceq \\mathbf{h}.\n",
"\\end{eqnarray*}\n",
"Linear program is a convex optimization problem, why?\n",
"\n",
"![](\"./lp.png\")
\n",
"\n",
"* The **standard form** of an LP is\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{c}^T \\mathbf{x} \\\\\n",
"\t&\\text{subject to}& \\mathbf{A} \\mathbf{x} = \\mathbf{b} \\\\\n",
"\t& & \\mathbf{x} \\succeq \\mathbf{0}.\n",
"\\end{eqnarray*}\n",
"To transform a general linear program into the standard form, we introduce the _slack variables_ $\\mathbf{s} \\succeq \\mathbf{0}$ such that $\\mathbf{G} \\mathbf{x} + \\mathbf{s} = \\mathbf{h}$. Then we write $\\mathbf{x} = \\mathbf{x}^+ - \\mathbf{x}^-$, where $\\mathbf{x}^+ \\succeq \\mathbf{0}$ and $\\mathbf{x}^- \\succeq \\mathbf{0}$. This yields the problem\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{c}^T (\\mathbf{x}^+ - \\mathbf{x}^-) \\\\\n",
"\t&\\text{subject to}& \\mathbf{A} (\\mathbf{x}^+ - \\mathbf{x}^-) = \\mathbf{b} \\\\\n",
"\t& & \\mathbf{G} (\\mathbf{x}^+ - \\mathbf{x}^-) + \\mathbf{s} = \\mathbf{h} \\\\\n",
"\t& & \\mathbf{x}^+ \\succeq \\mathbf{0}, \\mathbf{x}^- \\succeq \\mathbf{0}, \\mathbf{s} \\succeq \\mathbf{0}\n",
"\\end{eqnarray*}\n",
"in $\\mathbf{x}^+$, $\\mathbf{x}^-$, and $\\mathbf{s}$.\n",
"\n",
" Slack variables are often used to transform a complicated inequality constraint to simple non-negativity constraints.\n",
"\n",
"* The **inequality form** of an LP is\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{c}^T \\mathbf{x} \\\\\n",
"\t&\\text{subject to}& \\mathbf{G} \\mathbf{x} \\preceq \\mathbf{h}.\n",
"\\end{eqnarray*}\n",
"\n",
"* Some softwares, e.g., `solveLP` in R, require an LP be written in either standard or inequality form. However a good software should do this for you!\n",
"\n",
"* A _piecewise-linear minimization_ problem\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\max_{i=1,\\ldots,m} (\\mathbf{a}_i^T \\mathbf{x} + b_i)\n",
"\\end{eqnarray*}\n",
"can be transformed to an LP\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& t \\\\\n",
"\t&\\text{subject to}& \\mathbf{a}_i^T \\mathbf{x} + b_i \\le t, \\quad i = 1,\\ldots,m,\n",
"\\end{eqnarray*}\n",
"in $\\mathbf{x}$ and $t$. Apparently \n",
"$$\n",
"\t\\text{minimize} \\max_{i=1,\\ldots,m} |\\mathbf{a}_i^T \\mathbf{x} + b_i|\n",
"$$\n",
"and\n",
"$$\n",
"\t\\text{minimize} \\max_{i=1,\\ldots,m} (\\mathbf{a}_i^T \\mathbf{x} + b_i)_+\n",
"$$\n",
"are also LP.\n",
"\n",
"* Any _convex optimization problem_\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& f_0(\\mathbf{x}) \\\\\n",
"\t&\\text{subject to}& f_i(\\mathbf{x}) \\le 0, \\quad i=1,\\ldots,m \\\\\n",
"\t&& \\mathbf{a}_i^T \\mathbf{x} = b_i, \\quad i=1,\\ldots,p,\n",
"\\end{eqnarray*}\n",
"where $f_0,\\ldots,f_m$ are convex functions, can be transformed to the _epigraph form_\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& t \\\\\n",
"\t&\\text{subject to}& f_0(\\mathbf{x}) - t \\le 0 \\\\\n",
"\t& & f_i(\\mathbf{x}) \\le 0, \\quad i=1,\\ldots,m \\\\\n",
"\t& & \\mathbf{a}_i^T \\mathbf{x} = b_i, \\quad i=1,\\ldots,p\n",
"\\end{eqnarray*}\n",
"in variables $\\mathbf{x}$ and $t$. That is why people often say linear program is universal.\n",
"\n",
"* The _linear fractional programming_\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\frac{\\mathbf{c}^T \\mathbf{x} + d}{\\mathbf{e}^T \\mathbf{x} + f} \\\\\n",
"\t&\\text{subject to}& \\mathbf{A} \\mathbf{x} = \\mathbf{b} \\\\\n",
"\t& & \\mathbf{G} \\mathbf{x} \\preceq \\mathbf{h} \\\\\n",
"\t& & \\mathbf{e}^T \\mathbf{x} + f > 0\n",
"\\end{eqnarray*}\n",
"can be transformed to an LP\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{c}^T \\mathbf{y} + d z \\\\\n",
"\t&\\text{subject to}& \\mathbf{G} \\mathbf{y} - z \\mathbf{h} \\preceq \\mathbf{0} \\\\\n",
"\t& & \\mathbf{A} \\mathbf{y} - z \\mathbf{b} = \\mathbf{0} \\\\\n",
"\t& & \\mathbf{e}^T \\mathbf{y} + f z = 1 \\\\\n",
"\t& & z \\ge 0\n",
"\\end{eqnarray*}\n",
"in $\\mathbf{y}$ and $z$, via transformation of variables\n",
"\\begin{eqnarray*}\n",
"\t\\mathbf{y} = \\frac{\\mathbf{x}}{\\mathbf{e}^T \\mathbf{x} + f}, \\quad z = \\frac{1}{\\mathbf{e}^T \\mathbf{x} + f}.\n",
"\\end{eqnarray*}\n",
"See Section 4.3.2 of Boyd and Vandenberghe (2004) for proof."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## LP example: compressed sensing\n",
"\n",
"* **Compressed sensing** [Candes and Tao (2006)](https://doi.org/10.1109/TIT.2006.885507) and [Donoho (2006)](https://doi.org/10.1109/TIT.2006.871582) tries to address a fundamental question: how to compress and transmit a complex signal (e.g., musical clips, mega-pixel images), which can be decoded to recover the original signal?\n",
"\n",
"![](\"./david-donoho.jpg\")
\n",
"![](\"./emmanuel-candes.png\")
\n",
"![](\"./terrence-tao.png\")
\n",
"\n",
"* Suppose a signal $\\mathbf{x} \\in \\mathbb{R}^n$ is sparse with $s$ non-zeros. We under-sample the signal by multiplying a (flat) measurement matrix $\\mathbf{y} = \\mathbf{A} \\mathbf{x}$, where $\\mathbf{A} \\in \\mathbb{R}^{m\\times n}$ has iid normal entries. [Candes, Romberg and Tao (2006)](https://doi.org/10.1002/cpa.20124) show that the solution to\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\|\\mathbf{x}\\|_1 \\\\\n",
"\t&\\text{subject to}& \\mathbf{A} \\mathbf{x} = \\mathbf{y}\n",
"\\end{eqnarray*}\n",
"exactly recovers the true signal under certain conditions on $\\mathbf{A}$ when $n \\gg s$ and $m \\approx s \\ln(n/s)$. Why sparsity is a reasonable assumption? _Virtually all real-world images have low information content_.\n",
"\n",
"![](\"./movie-scene.png\")
\n",
"\n",
"* The $\\ell_1$ minimization problem apparently is an LP, by writing $\\mathbf{x} = \\mathbf{x}^+ - \\mathbf{x}^-$,\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{1}^T (\\mathbf{x}^+ + \\mathbf{x}^-) \\\\\n",
"\t&\\text{subject to}& \\mathbf{A} (\\mathbf{x}^+ - \\mathbf{x}^-) = \\mathbf{y} \\\\\n",
"\t& & \\mathbf{x}^+ \\succeq \\mathbf{0}, \\mathbf{x}^- \\succeq \\mathbf{0}.\n",
"\\end{eqnarray*}\n",
"\n",
"* Let's try a numerical example.\n",
"\n",
"### Generate a sparse signal and sub-sampling"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
],
"text/html": [
"\n",
"\n"
],
"text/plain": [
"Plot(...)"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Gadfly, Random\n",
"\n",
"# random seed\n",
"Random.seed!(280)\n",
"# Size of signal\n",
"n = 1024\n",
"# Sparsity (# nonzeros) in the signal\n",
"s = 20\n",
"# Number of samples (undersample by a factor of 8) \n",
"m = 128\n",
"\n",
"# Generate and display the signal\n",
"x0 = zeros(n)\n",
"x0[rand(1:n, s)] = randn(s)\n",
"# Generate the random sampling matrix\n",
"A = randn(m, n) / m\n",
"# Subsample by multiplexing\n",
"y = A * x0\n",
"\n",
"# plot the true signal\n",
"plot(x=1:n, y=x0, Geom.line, Guide.title(\"True signal x_0\"))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Solve LP by calling Gurobi directly (not recommended)\n",
"\n",
"Gurobi model formulation: \n",
"\\begin{eqnarray*}\n",
" &\\text{minimize}& (1/2) x' H x + f' x \\\\\n",
" &s.t.& A x \\le b \\\\\n",
" & & Aeq \\cdot x = beq \\\\\n",
" & & lb <= x <= ub\n",
"\\end{eqnarray*}\n",
"Refer to [Gurobi.jl](https://github.com/JuliaOpt/Gurobi.jl) documentation for setting up the model."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Academic license - for non-commercial use only\n",
"Optimize a model with 128 rows, 2048 columns and 262144 nonzeros\n",
"Coefficient statistics:\n",
" Matrix range [1e-08, 3e-02]\n",
" Objective range [1e+00, 1e+00]\n",
" Bounds range [0e+00, 0e+00]\n",
" RHS range [1e-04, 1e-01]\n",
"\n",
"Concurrent LP optimizer: dual simplex and barrier\n",
"Showing barrier log only...\n",
"\n",
"Presolve time: 0.10s\n",
"Presolved: 128 rows, 2048 columns, 262144 nonzeros\n",
"\n",
"Ordering time: 0.00s\n",
"\n",
"Barrier statistics:\n",
" AA' NZ : 8.128e+03\n",
" Factor NZ : 8.256e+03 (roughly 1 MByte of memory)\n",
" Factor Ops : 7.073e+05 (less than 1 second per iteration)\n",
" Threads : 3\n",
"\n",
" Objective Residual\n",
"Iter Primal Dual Primal Dual Compl Time\n",
" 0 2.51327671e+03 0.00000000e+00 1.40e-14 0.00e+00 2.08e+00 0s\n",
" 1 4.36716581e+02 2.93300751e+00 1.90e-14 2.22e-15 2.12e-01 0s\n",
"\n",
"Barrier performed 1 iterations in 0.16 seconds\n",
"Barrier solve interrupted - model solved by another algorithm\n",
"\n",
"\n",
"Solved with dual simplex\n",
"Solved in 161 iterations and 0.16 seconds\n",
"Optimal objective 1.830728192e+01\n"
]
},
{
"data": {
"image/svg+xml": [
"\n",
"\n"
],
"text/html": [
"\n",
"\n"
],
"text/plain": [
"Plot(...)"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Gurobi\n",
"env = Gurobi.Env()\n",
"setparams!(env, OutputFlag=1) # display log\n",
"\n",
"# Construct the model\n",
"model = gurobi_model(env;\n",
" name = \"cs\",\n",
" f = ones(2 * n),\n",
" Aeq = [A -A],\n",
" beq = y,\n",
" lb = zeros(2 * n))\n",
"\n",
"# Run optimization\n",
"optimize(model)\n",
"\n",
"# Show results\n",
"sol = get_solution(model)\n",
"xsol = sol[1:n] - sol[n + 1:end]\n",
"\n",
"plot(x=1:n, y=x0, Geom.point)\n",
"plot(x=1:n, y=xsol, Geom.line, Guide.title(\"Reconstructed signal overlayed with x0\"))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Solve LP by DCP (disciplined convex programming) interface Convex.jl\n",
"\n",
"Check [Convex.jl documentation](https://convexjl.readthedocs.io/en/latest/operations.html) for a list of supported operations."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Problem\n",
" Name : \n",
" Objective sense : min \n",
" Type : LO (linear optimization problem)\n",
" Constraints : 2177 \n",
" Cones : 0 \n",
" Scalar variables : 2049 \n",
" Matrix variables : 0 \n",
" Integer variables : 0 \n",
"\n",
"Optimizer started.\n",
"Presolve started.\n",
"Linear dependency checker started.\n",
"Linear dependency checker terminated.\n",
"Eliminator started.\n",
"Freed constraints in eliminator : 1024\n",
"Eliminator terminated.\n",
"Eliminator started.\n",
"Freed constraints in eliminator : 0\n",
"Eliminator terminated.\n",
"Eliminator - tries : 2 time : 0.00 \n",
"Lin. dep. - tries : 1 time : 0.02 \n",
"Lin. dep. - number : 0 \n",
"Presolve terminated. Time: 0.04 \n",
"Problem\n",
" Name : \n",
" Objective sense : min \n",
" Type : LO (linear optimization problem)\n",
" Constraints : 2177 \n",
" Cones : 0 \n",
" Scalar variables : 2049 \n",
" Matrix variables : 0 \n",
" Integer variables : 0 \n",
"\n",
"Optimizer - threads : 8 \n",
"Optimizer - solved problem : the primal \n",
"Optimizer - Constraints : 1152\n",
"Optimizer - Cones : 0\n",
"Optimizer - Scalar variables : 3073 conic : 0 \n",
"Optimizer - Semi-definite variables: 0 scalarized : 0 \n",
"Factor - setup time : 0.03 dense det. time : 0.00 \n",
"Factor - ML order time : 0.00 GP order time : 0.00 \n",
"Factor - nonzeros before factor : 1.42e+05 after factor : 1.42e+05 \n",
"Factor - dense dim. : 1 flops : 3.59e+07 \n",
"ITE PFEAS DFEAS GFEAS PRSTATUS POBJ DOBJ MU TIME \n",
"0 5.1e+00 1.0e+00 1.2e+03 1.00e+00 1.448154688e+03 0.000000000e+00 1.0e+00 0.09 \n",
"1 2.0e+00 1.0e+00 7.9e+02 0.00e+00 1.868721000e+02 5.017718138e+00 3.7e+00 0.10 \n",
"2 8.9e-01 4.5e-01 3.6e+02 9.38e-01 9.836670682e+01 1.394257795e+01 1.7e+00 0.10 \n",
"3 2.8e-01 1.4e-01 1.1e+02 9.69e-01 4.300059794e+01 1.627916054e+01 5.3e-01 0.11 \n",
"4 9.9e-02 5.0e-02 4.0e+01 9.85e-01 2.693148355e+01 1.743259784e+01 1.9e-01 0.11 \n",
"5 2.0e-02 1.0e-02 8.0e+00 9.94e-01 2.000060021e+01 1.808359996e+01 3.8e-02 0.12 \n",
"6 1.6e-03 8.2e-04 6.5e-01 9.98e-01 1.843978822e+01 1.828480266e+01 3.0e-03 0.13 \n",
"7 1.2e-05 5.9e-06 4.7e-03 1.00e+00 1.830822487e+01 1.830710759e+01 2.2e-05 0.13 \n",
"8 5.9e-08 3.0e-08 2.4e-05 1.00e+00 1.830728673e+01 1.830728103e+01 1.1e-07 0.14 \n",
"9 3.5e-08 1.5e-10 1.7e-07 1.00e+00 1.830728195e+01 1.830728191e+01 5.8e-10 0.14 \n",
"Basis identification started.\n",
"Basis identification terminated. Time: 0.05\n",
"Optimizer terminated. Time: 0.20 \n",
"\n",
"\n",
"Interior-point solution summary\n",
" Problem status : PRIMAL_AND_DUAL_FEASIBLE\n",
" Solution status : OPTIMAL\n",
" Primal. obj: 1.8307281953e+01 nrm: 2e+01 Viol. con: 2e-10 var: 0e+00 \n",
" Dual. obj: 1.8307281912e+01 nrm: 1e+01 Viol. con: 3e-16 var: 6e-12 \n",
"\n",
"Basic solution summary\n",
" Problem status : PRIMAL_AND_DUAL_FEASIBLE\n",
" Solution status : OPTIMAL\n",
" Primal. obj: 1.8307281917e+01 nrm: 2e+01 Viol. con: 2e-14 var: 0e+00 \n",
" Dual. obj: 1.8307281912e+01 nrm: 2e+01 Viol. con: 4e-16 var: 6e-12 \n",
" 0.235640 seconds (95.13 k allocations: 23.659 MiB, 3.46% gc time)\n"
]
},
{
"data": {
"image/svg+xml": [
"\n",
"\n"
],
"text/html": [
"\n",
"\n"
],
"text/plain": [
"Plot(...)"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Convex\n",
"\n",
"# Use Mosek solver\n",
"using Mosek\n",
"solver = MosekSolver(LOG=1)\n",
"\n",
"# # Use Gurobi solver\n",
"# using Gurobi\n",
"# solver = GurobiSolver(OutputFlag=1)\n",
"\n",
"# # Use Cplex solver\n",
"# using CPLEX\n",
"# solver = CplexSolver(CPXPARAM_ScreenOutput=1)\n",
"\n",
"## Use SCS solver\n",
"#using SCS\n",
"#solver = SCSSolver(verbose=1)\n",
"#set_default_solver(solver)\n",
"\n",
"# Set up optimizaiton problem\n",
"x = Variable(n)\n",
"problem = minimize(norm(x, 1))\n",
"problem.constraints += A * x == y\n",
"\n",
"# Solve the problem\n",
"@time solve!(problem, solver)\n",
"\n",
"# Display the solution\n",
"plot(x=1:n, y=x0, Geom.point)\n",
"plot(x=1:n, y=xsol, Geom.line, Guide.title(\"Reconstructed signal overlayed with x0\"))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## LP example: quantile regression\n",
"\n",
"![](\"./growth-2-20-boys.png\")
\n",
"![](\"./quantregbin3.png\")
\n",
"\n",
"* In linear regression, we model the mean of response variable as a function of covariates. In many situations, the error variance is not constant, the distribution of $y$ may be asymmetric, or we simply care about the quantile(s) of response variable. Quantile regression offers a better modeling tool in these applications.\n",
"\n",
"\n",
"* In $\\tau$-quantile regression, we minimize the loss function\n",
"\\begin{eqnarray*}\n",
"\tf(\\beta) = \\sum_{i=1}^n \\rho_\\tau (y_i - \\mathbf{x}_i^T \\beta),\n",
"\\end{eqnarray*}\n",
"where $\\rho_\\tau(z) = z (\\tau - 1_{\\{z < 0\\}})$. Writing $\\mathbf{y} - \\mathbf{X} \\beta = \\mathbf{r}^+ - \\mathbf{r}^-$, this is equivalent to the LP\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\tau \\mathbf{1}^T \\mathbf{r}^+ + (1-\\tau) \\mathbf{1}^T \\mathbf{r}^- \\\\\n",
"\t&\\text{subject to}& \\mathbf{r}^+ - \\mathbf{r}^- = \\mathbf{y} - \\mathbf{X} \\beta \\\\\n",
"\t& & \\mathbf{r}^+ \\succeq \\mathbf{0}, \\mathbf{r}^- \\succeq \\mathbf{0}\n",
"\\end{eqnarray*}\n",
"in $\\mathbf{r}^+$, $\\mathbf{r}^-$, and $\\beta$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## LP Example: $\\ell_1$ regression\n",
"\n",
"* A popular method in robust statistics is the median absolute deviation (MAD) regression that minimizes the $\\ell_1$ norm of the residual vector $\\|\\mathbf{y} - \\mathbf{X} \\beta\\|_1$. This apparently is equivalent to the LP\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\mathbf{1}^T (\\mathbf{r}^+ + \\mathbf{r}^-) \\\\\n",
"\t&\\text{subject to}& \\mathbf{r}^+ - \\mathbf{r}^- = \\mathbf{y} - \\mathbf{X} \\beta \\\\\n",
"\t& & \\mathbf{r}^+ \\succeq \\mathbf{0}, \\mathbf{r}^- \\succeq \\mathbf{0}\n",
"\\end{eqnarray*}\n",
"in $\\mathbf{r}^+$, $\\mathbf{r}^-$, and $\\beta$. \n",
"\n",
" $\\ell_1$ regression = MAD = 1/2-quantile regression."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## LP Example: $\\ell_\\infty$ regression (Chebychev approximation)\n",
"\n",
"* Minimizing the worst possible residual $\\|\\mathbf{y} - \\mathbf{X} \\beta\\|_\\infty$ is equivalent to the LP\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& t \\\\\n",
"\t&\\text{subject to}& -t \\le y_i - \\mathbf{x}_i^T \\beta \\le t, \\quad i = 1,\\dots,n\n",
"\\end{eqnarray*}\n",
"in variables $\\beta$ and $t$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## LP Example: Dantzig selector\n",
"\n",
"* [Candes and Tao (2007)](https://www.doi.org/10.1214/009053606000001523) propose a variable selection method called the Dantzig selector that solves\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\|\\mathbf{X}^T (\\mathbf{y} - \\mathbf{X} \\beta)\\|_\\infty \\\\\n",
"\t&\\text{subject to}& \\sum_{j=2}^p |\\beta_j| \\le t,\n",
"\\end{eqnarray*}\n",
"which can be transformed to an LP. Indeed they name the method after George Dantzig, who invented the simplex method for efficiently solving LP in 50s.\n",
"\n",
"![](\"./george-bernard-dantzig.jpg\")
"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## LP Example: 1-norm SVM\n",
"\n",
"* In two-class classification problems, we are given training data $(\\mathbf{x}_i, y_i)$, $i=1,\\ldots,n$, where $\\mathbf{x}_i \\in \\mathbb{R}^p$ are feature vectors and $y_i \\in \\{-1, 1\\}$ are class labels. [Zhu, Rosset, Tibshirani, and Hastie (2004)](https://www.google.com/url?sa=t&rct=j&q=&esrc=s&source=web&cd=1&cad=rja&uact=8&ved=2ahUKEwi-_6C2sMLiAhUNGDQIHU1RD3MQFjAAegQIBRAC&url=https%3A%2F%2Fpapers.nips.cc%2Fpaper%2F2450-1-norm-support-vector-machines.pdf&usg=AOvVaw1cTOsmPjpxKJHCs93iLUkn) propose the 1-norm support vector machine (svm) that achieves the dual purpose of classification and feature selection. Denote the solution of the optimization problem\n",
"\\begin{eqnarray*}\n",
"\t&\\text{minimize}& \\sum_{i=1}^n \\left[ 1 - y_i \\left( \\beta_0 + \\sum_{j=1}^p x_{ij} \\beta_j \\right) \\right]_+ \\\\\n",
"\t&\\text{subject to}& \\|\\beta\\|_1 = \\sum_{j=1}^p |\\beta_j| \\le t\n",
"\\end{eqnarray*}\n",
"by $\\hat \\beta_0(t)$ and $\\hat \\beta(t)$. 1-norm svm classifies a future feature vector $\\mathbf{x}$ by the sign of fitted model\n",
"\\begin{eqnarray*}\n",
"\t\\hat f(\\mathbf{x}) = \\hat \\beta_0 + \\mathbf{x}^T \\hat \\beta.\n",
"\\end{eqnarray*}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Many more applications of LP: Airport scheduling (Copenhagen airport uses Gurobi), airline flight scheduling, NFL scheduling, match.com, $\\LaTeX$, ...\n",
"\n",
"Apparently any loss/penalty or loss/constraint combinations of form \n",
"$$\n",
"\\{\\ell_1, \\ell_\\infty, \\text{quantile}\\} \\times \\{\\ell_1, \\ell_\\infty, \\text{quantile}\\},\n",
"$$\n",
"possibly with affine (equality and/or inequality) constraints, can be formulated as an LP."
]
}
],
"metadata": {
"@webio": {
"lastCommId": null,
"lastKernelId": null
},
"kernelspec": {
"display_name": "Julia 1.1.0",
"language": "julia",
"name": "julia-1.1"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.1.0"
},
"toc": {
"colors": {
"hover_highlight": "#DAA520",
"running_highlight": "#FF0000",
"selected_highlight": "#FFD700"
},
"moveMenuLeft": true,
"nav_menu": {
"height": "30.66666603088379px",
"width": "252px"
},
"navigate_menu": true,
"number_sections": true,
"sideBar": true,
"skip_h1_title": true,
"threshold": 4,
"toc_cell": true,
"toc_section_display": "block",
"toc_window_display": true,
"widenNotebook": false
}
},
"nbformat": 4,
"nbformat_minor": 2
}