{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "title: Experiment Design\n",
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Originally Contributed by**: Arpit Bhatia, Chris Coey"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This tutorial covers experiment design examples (D-optimal, A-optimal, and E-optimal)\n",
    "from section 7.5 of the book Convex Optimization by Boyd and Vandenberghe[[1]](#c1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Relaxed Experiment Design Problem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The basic experiment design problem is as follows.\n",
    "Given the menu of possible choices for experiments, $v_{1}, \\ldots, v_{p}$,\n",
    "and the total number $m$ of experiments to be carried out, choose the numbers of each type of experiment,\n",
    "$i . e ., m_{1}, \\ldots, m_{p}$ to make the error covariance $E$ small (in some sense).\n",
    "The variables $m_{1}, \\ldots, m_{p}$ must, of course, be integers and sum to $m,$ the given total number of experiments.\n",
    "This leads to the optimization problem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{cl}{\\operatorname{minimize}\\left(\\mathrm{w.r.t.} \\mathbf{S}_{+}^{n}\\right)} & {E=\\left(\\sum_{j=1}^{p} m_{j} v_{j} v_{j}^{T}\\right)^{-1}} \\\\ {\\text { subject to }} & {m_{i} \\geq 0, \\quad m_{1}+\\cdots+m_{p}=m} \\\\ {} & {m_{i} \\in \\mathbf{Z}}\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The basic experiment design problem can be a hard combinatorial problem when $m,$ the total number of experiments,\n",
    "is comparable to $n$ , since in this case the $m_{i}$ are all small integers.\n",
    "In the case when $m$ is large compared to $n$ , however, a good approximate solution can be found by ignoring,\n",
    "or relaxing, the constraint that the $m_{i}$ are integers.\n",
    "Let $\\lambda_{i}=m_{i} / m,$ which is the fraction of the total number of experiments for which\n",
    "$a_{j}=v_{i},$ or the relative frequency of experiment $i$.\n",
    "We can express the error covariance in terms of $\\lambda_{i}$ as"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "E=\\frac{1}{m}\\left(\\sum_{i=1}^{p} \\lambda_{i} v_{i} v_{i}^{T}\\right)^{-1}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The vector $\\lambda \\in \\mathbf{R}^{p}$ satisfies $\\lambda \\succeq 0,\n",
    "\\mathbf{1}^{T} \\lambda=1,$ and also, each $\\lambda_{i}$ is an integer multiple of $1 / m$.\n",
    "By ignoring this last constraint, we arrive at the problem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{ll}{\\operatorname{minimize}\\left(\\mathrm{w.r.t.} \\mathbf{S}_{+}^{n}\\right)} & {E=(1 / m)\\left(\\sum_{i=1}^{p} \\lambda_{i} v_{i} v_{i}^{T}\\right)^{-1}} \\\\ {\\text { subject to }} & {\\lambda \\succeq 0, \\quad \\mathbf{1}^{T} \\lambda=1}\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Types of Experiment Design Problems"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Several scalarizations have been proposed for the experiment design problem,\n",
    "which is a vector optimization problem over the positive semidefinite cone."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "using JuMP\n",
    "using SCS\n",
    "using LinearAlgebra\n",
    "using Random\n",
    "\n",
    "Random.seed!(1234);\n",
    "\n",
    "q = 4 # dimension of estimate space\n",
    "p = 8 # number of experimental vectors\n",
    "nmax = 3 # upper bound on lambda\n",
    "n = 12\n",
    "\n",
    "V = randn(q, p)\n",
    "\n",
    "eye = Matrix{Float64}(I, q, q);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### A-optimal design"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In A-optimal experiment design, we minimize tr $E$, the trace of the covariance matrix.\n",
    "This objective is simply the mean of the norm of the error squared:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\mathbf{E}\\|e\\|_{2}^{2}=\\mathbf{E} \\operatorname{tr}\\left(e e^{T}\\right)=\\operatorname{tr} E\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The A-optimal experiment design problem in SDP form is"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{ll}{\\operatorname{minimize}} & {\\mathbf{1}^{T} u} \\\\ {\\text { subject to }} & {\\left[\\begin{array}{cc}{\\sum_{i=1}^{p} \\lambda_{i} v_{i} v_{i}^{T}} & {e_{k}} \\\\ {e_{k}^{T}} & {u_{k}}\\end{array}\\right] \\succeq 0, \\quad k=1, \\ldots, n} \\\\ {} & {\\lambda \\succeq 0, \\quad \\mathbf{1}^{T} \\lambda=1}\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "objective_value(aOpt) = 5.041247589148935\n",
      "value.(np) = [1.7479360444207457, 1.1153135402236287, 1.8896080115566673e-6, 1.6619566409940143, 2.9999969400033906, 0.8414161134946255, 1.3825673956267661, 2.25080400472545]\n"
     ]
    }
   ],
   "source": [
    "aOpt = Model(optimizer_with_attributes(SCS.Optimizer, \"verbose\" => 0))\n",
    "@variable(aOpt, np[1:p], lower_bound = 0, upper_bound = nmax)\n",
    "@variable(aOpt, u[1:q], lower_bound = 0)\n",
    "\n",
    "@constraint(aOpt, sum(np) <= n)\n",
    "for i = 1:q\n",
    "    @SDconstraint(aOpt, [V * diagm(0 => np ./ n) * V' eye[:, i]; eye[i, :]' u[i]] >= 0)\n",
    "end\n",
    "\n",
    "@objective(aOpt, Min, sum(u))\n",
    "\n",
    "optimize!(aOpt)\n",
    "\n",
    "@show objective_value(aOpt);\n",
    "@show value.(np);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### E-optimal design"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In $E$ -optimal design, we minimize the norm of the error covariance matrix, i.e. the maximum eigenvalue of $E$.\n",
    "Since the diameter (twice the longest semi-axis) of the confidence ellipsoid $\\mathcal{E}$\n",
    "is proportional to $\\|E\\|_{2}^{1 / 2}$,\n",
    "minimizing $\\|E\\|_{2}$ can be interpreted geometrically as minimizing the diameter of the confidence ellipsoid.\n",
    "E-optimal design can also be interpreted as minimizing the maximum variance of $q^{T} e$,\n",
    "over all $q$ with $\\|q\\|_{2}=1$.\n",
    "The E-optimal experiment design problem in SDP form is"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{cl}{\\operatorname{maximize}} & {t} \\\\ {\\text { subject to }} & {\\sum_{i=1}^{p} \\lambda_{i} v_{i} v_{i}^{T} \\succeq t I} \\\\ {} & {\\lambda \\succeq 0, \\quad \\mathbf{1}^{T} \\lambda=1}\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "objective_value(eOpt) = 0.4489430774609742\n",
      "value.(np) = [3.0000033649730002, 0.6752098200096417, -2.1874856564735654e-6, 1.0458573723670574, 2.9999992875514727, 1.7869721672722332, 0.30150709122580693, 2.190460469380653]\n"
     ]
    }
   ],
   "source": [
    "eOpt = Model(optimizer_with_attributes(SCS.Optimizer, \"verbose\" => 0))\n",
    "@variable(eOpt, np[1:p], lower_bound = 0, upper_bound = nmax)\n",
    "@variable(eOpt, t)\n",
    "\n",
    "@SDconstraint(eOpt, V * diagm(0 => np ./ n) * V' - (t .* eye) >= 0)\n",
    "@constraint(eOpt, sum(np) <= n)\n",
    "\n",
    "@objective(eOpt, Max, t)\n",
    "\n",
    "optimize!(eOpt)\n",
    "\n",
    "@show objective_value(eOpt);\n",
    "@show value.(np);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### D-optimal design\n",
    "The most widely used scalarization is called $D$ -optimal design,\n",
    "in which we minimize the determinant of the error covariance matrix $E$.\n",
    "This corresponds to designing the experiment to minimize the volume of the resulting confidence ellipsoid\n",
    "(for a fixed confidence level).\n",
    "Ignoring the constant factor 1$/ m$ in $E$, and taking the logarithm of the objective,\n",
    "we can pose this problem as convex optimization problem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{ll}{\\operatorname{minimize}} & {\\log \\operatorname{det}\\left(\\sum_{i=1}^{p} \\lambda_{i} v_{i} v_{i}^{T}\\right)^{-1}} \\\\ {\\text { subject to }} & {\\lambda \\succeq 0, \\quad \\mathbf{1}^{T} \\lambda=1}\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "objective_value(dOpt) = 0.19015174239253715\n",
      "value.(np) = [-8.79402147330283e-7, 2.567393113706525, 4.4145237864313906e-7, 0.2627580211834385, 2.9434800539941457, 2.3925267125050307, 2.8369691609526035, 0.9968787186836513]\n"
     ]
    }
   ],
   "source": [
    "dOpt = Model(optimizer_with_attributes(SCS.Optimizer, \"verbose\" => 0))\n",
    "@variable(dOpt, np[1:p], lower_bound = 0, upper_bound = nmax)\n",
    "@variable(dOpt, t)\n",
    "@objective(dOpt, Max, t)\n",
    "@constraint(dOpt, sum(np) <= n)\n",
    "E = V * diagm(0 => np ./ n) * V'\n",
    "@constraint(dOpt, [t, 1, (E[i, j] for i in 1:q for j in 1:i)...] in MOI.LogDetConeTriangle(q))\n",
    "\n",
    "optimize!(dOpt)\n",
    "\n",
    "@show objective_value(dOpt);\n",
    "@show value.(np);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### References\n",
    "<a id='c1'></a>\n",
    "1. Boyd, S., & Vandenberghe, L. (2004). Convex Optimization. Cambridge: Cambridge University Press. doi:10.1017/CBO9780511804441"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.5.1",
   "language": "julia",
   "name": "julia-1.5"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.5.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}