Experiment design

Страница в процессе перевода.

This tutorial was generated using Literate.jl. Download the source as a .jl file.

This tutorial was originally contributed by Arpit Bhatia and Chris Coey.

This tutorial covers experiment design examples (D-optimal, A-optimal, and E-optimal) from section 7.5 of Boyd2004.

The tutorial uses the following packages

using JuMP
import SCS
import LinearAlgebra
import Random

This tutorial uses sets from MathOptInterface. By default, JuMP exports the MOI symbol as an alias for the MathOptInterface.jl package. We recommend making this more explicit in your code by adding the following lines:

```julia
import MathOptInterface as MOI
```

We set a seed so the random numbers are repeatable:

Random.seed!(1234)

Random.TaskLocalRNG()

The relaxed experiment design problem

The basic experiment design problem is as follows.

Given the menu of possible choices for experiments, , and the total number of experiments to be carried out, choose the numbers of each type of experiment, that is, to make the error covariance small (in some sense).

The variables must, of course, be integers and sum to the given total number of experiments. This leads to the optimization problem:

The basic experiment design problem can be a hard combinatorial problem when , the total number of experiments, is comparable to , since in this case the are all small integers.

In the case when is large compared to , however, a good approximate solution can be found by ignoring, or relaxing, the constraint that the are integers.

Let which is the fraction of the total number of experiments for which or the relative frequency of experiment . We can express the error covariance in terms of as:

The vector satisfies and also, each is an integer multiple of . By ignoring this last constraint, we arrive at the problem:

Several scalarizations have been proposed for the experiment design problem, which is a vector optimization problem over the positive semidefinite cone.

q = 4 # dimension of estimate space
p = 8 # number of experimental vectors
n_max = 3 # upper bound on lambda
n = 12

V = randn(q, p)

eye = Matrix{Float64}(LinearAlgebra.I, q, q);

A-optimal design

In A-optimal experiment design, we minimize tr , the trace of the covariance matrix. This objective is simply the mean of the norm of the error squared:

The A-optimal experiment design problem in SDP form is

aOpt = Model(SCS.Optimizer)
set_silent(aOpt)
@variable(aOpt, np[1:p], lower_bound = 0, upper_bound = n_max)
@variable(aOpt, u[1:q], lower_bound = 0)
@constraint(aOpt, sum(np) <= n)
for i in 1:q
    matrix = [
        V*LinearAlgebra.diagm(0 => np ./ n)*V' eye[:, i]
        eye[i, :]' u[i]
    ]
    @constraint(aOpt, matrix >= 0, PSDCone())
end
@objective(aOpt, Min, sum(u))
optimize!(aOpt)
@assert is_solved_and_feasible(aOpt)
objective_value(aOpt)

5.103223362935165

value.(np)

8-element Vector{Float64}:
 2.949521135024129
 1.7790921964770634
 8.104871349218623e-6
 4.706163791527242e-6
 2.10826698418427
 1.698346379809937
 1.2635010842547454
 2.2012292330175063

E-optimal design

In -optimal design, we minimize the norm of the error covariance matrix, that is, the maximum eigenvalue of .

Since the diameter (twice the longest semi-axis) of the confidence ellipsoid is proportional to , minimizing can be interpreted geometrically as minimizing the diameter of the confidence ellipsoid.

E-optimal design can also be interpreted as minimizing the maximum variance of , over all with . The E-optimal experiment design problem in SDP form is:

eOpt = Model(SCS.Optimizer)
set_silent(eOpt)
@variable(eOpt, 0 <= np[1:p] <= n_max)
@variable(eOpt, t)
@constraint(
    eOpt,
    V * LinearAlgebra.diagm(0 => np ./ n) * V' - (t .* eye) >= 0,
    PSDCone(),
)
@constraint(eOpt, sum(np) <= n)
@objective(eOpt, Max, t)
optimize!(eOpt)
@assert is_solved_and_feasible(eOpt)
objective_value(eOpt)

0.4353863856695074

value.(np)

8-element Vector{Float64}:
  2.999998678509857
  2.7528101637391016
 -3.609558428104128e-7
 -1.487238376393808e-6
  2.1818462791653928
  2.3252635463929763
  0.21896764670154034
  1.521116305115523

D-optimal design

The most widely used scalarization is called -optimal design, in which we minimize the determinant of the error covariance matrix . This corresponds to designing the experiment to minimize the volume of the resulting confidence ellipsoid (for a fixed confidence level). Ignoring the constant factor in , and taking the logarithm of the objective, we can pose this problem as convex optimization problem:

dOpt = Model(SCS.Optimizer)
set_silent(dOpt)
@variable(dOpt, np[1:p], lower_bound = 0, upper_bound = n_max)
@variable(dOpt, t)
@objective(dOpt, Max, t)
@constraint(dOpt, sum(np) <= n)
E = V * LinearAlgebra.diagm(0 => np ./ n) * V'
@constraint(dOpt, [t; 1; triangle_vec(E)] in MOI.LogDetConeTriangle(q))
optimize!(dOpt)
@assert is_solved_and_feasible(dOpt)
objective_value(dOpt)

0.30845393467249654

value.(np)

8-element Vector{Float64}:
 0.4275824266205948
 2.9100163168941116
 2.916872474296307e-6
 2.6291979913959343e-6
 2.9157806299833644
 2.673375664593826
 2.7353950122202764
 0.3378388086254706