Skip to content

Commit ae82f40

Browse files
committed
Tutorial review pass: correctness fixes, notation, ch9 data, ch16 rewrite
Correctness: - ch8: fix the DPP regularization path (loss and penalty were on separate variables, so every path coefficient came out zero); share one variable. - ch13: fix the near-isotonic penalty sign (pos(diff) penalized upward steps and flattened to the mean; pos(-diff) penalizes downward steps). - ch14: clear the fixed start before best_of so the random restarts run (a preset value(x) silently pinned them all to the same result). Math and notation: - ch4: add standard-form definitions for LP, QP, SOCP, SDP. - ch5: normalize the elastic-net loss by 1/(2n) to match the code. - ch6/ch7/ch8/ch10: use n (observations) / p (variables) consistently. - ch15: define monomial and posynomial in the DGP intro. - ch2/ch3: align exercise statements with their code and solutions. Content: - ch9: replace the arbitrary alphabetical stock subset with daily adjusted prices for the 12 largest US companies (Yahoo Finance, 2022-2024); the data is now our own, models cited to Markowitz and Rockafellar-Uryasev; drop all Riskfolio-Lib references and the bib entry. - ch11: present the HiGHS/MIQP solver limitations as current, changeable behavior rather than permanent rules. - ch16: rewrite as a participant-facing open session with honest expectations; drop internal planning notes. - index: clarify the Render button note.
1 parent f3aea90 commit ae82f40

17 files changed

Lines changed: 903 additions & 853 deletions

02-gentle-intro.qmd

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -314,16 +314,16 @@ value(beta)
314314
```
315315
:::
316316

317-
**Exercise 3 (⭐).** Constrain the first four coefficients to sum to at most zero,
318-
$\sum_{i=1}^{4}\beta_i \le 0$, using a matrix row rather than writing out four terms.
317+
**Exercise 3 (⭐).** Constrain the first four coefficients to sum to at least 4,
318+
$\sum_{i=1}^{4}\beta_i \ge 4$, using a matrix row rather than writing out four terms.
319319

320320
::: {.callout-tip collapse="true"}
321321
## Solution
322322
Build a row vector that picks out $\beta_1,\ldots,\beta_4$:
323323
```{r}
324324
#| label: ex-sum
325325
A <- matrix(c(rep(1, 4), rep(0, 6)), nrow = 1)
326-
sum_con <- Problem(obj, constraints = list(A %*% beta <= 0))
326+
sum_con <- Problem(obj, constraints = list(A %*% beta >= 4))
327327
psolve(sum_con)
328328
value(beta)
329329
```

03-dcp.qmd

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -207,9 +207,10 @@ is_dcp(p)
207207
psolve(p)
208208
c(x = value(x), y = value(y))
209209
```
210-
`-square(·)` is concave and the constraint is `affine == affine`, so the problem is
211-
`Maximize(concave)` subject to a legal constraint --- DCP, with optimum
212-
$x = y = 1/2$.
210+
`square` is convex, so `-square(x)` and `-square(y)` are each concave, and a **sum
211+
of concave expressions is concave** --- so the objective `-square(x) - square(y)` is
212+
concave. The constraint is `affine == affine`, so the problem is `Maximize(concave)`
213+
subject to a legal constraint --- DCP, with optimum $x = y = 1/2$.
213214
:::
214215

215216
**Exercise 3 (⭐).** The expression `x^2 / y` (with `y > 0`) is convex --- the

04-problem-types.qmd

Lines changed: 48 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,16 @@ is expressible and which solvers apply.
3333

3434
Linear objective, linear constraints. The workhorse of operations research.
3535

36+
$$
37+
\begin{array}{ll}
38+
\text{minimize} & c^\top x\\
39+
\text{subject to} & Gx \preceq h,\quad Ax = b,
40+
\end{array}
41+
$$
42+
43+
with variable $x$ (length $n$); the data $c, G, h, A, b$ are constant and every
44+
expression is affine in $x$. ($\preceq$ is the componentwise $\le$ for vectors.)
45+
3646
```{r}
3747
#| label: lp
3848
x <- Variable(2)
@@ -47,6 +57,16 @@ value(x)
4757
A convex **quadratic** objective with linear constraints. Least squares (and ridge,
4858
and the Markowitz portfolio of @sec-portfolio) live here.
4959

60+
$$
61+
\begin{array}{ll}
62+
\text{minimize} & \frac{1}{2}\,x^\top P x + q^\top x\\
63+
\text{subject to} & Gx \preceq h,\quad Ax = b,
64+
\end{array}
65+
$$
66+
67+
with $P \succeq 0$ (positive semidefinite) --- that condition is exactly what makes
68+
the objective convex. Taking $P = 0$ recovers an LP.
69+
5070
```{r}
5171
#| label: qp
5272
set.seed(1)
@@ -60,8 +80,20 @@ value(z)
6080
## Second-order cone program (SOCP)
6181

6282
Adds constraints involving the Euclidean norm --- "second-order cones." Anything with
63-
a `p_norm(., 2)`, `huber`, or robust-constraint flavor is typically an SOCP. Minimizing
64-
a norm subject to a budget is the canonical toy:
83+
a `p_norm(., 2)`, `huber`, or robust-constraint flavor is typically an SOCP.
84+
85+
$$
86+
\begin{array}{ll}
87+
\text{minimize} & f^\top x\\
88+
\text{subject to} & \| A_i x + b_i \|_2 \le c_i^\top x + d_i,\quad i = 1,\dots,m,\\
89+
& Fx = g.
90+
\end{array}
91+
$$
92+
93+
Each "norm $\le$ affine" row is a *second-order cone* constraint --- it forces the
94+
pair $(A_i x + b_i,\; c_i^\top x + d_i)$ into the cone
95+
$\{(y, t) : \|y\|_2 \le t\}$. Minimizing a norm subject to a budget is the canonical
96+
toy:
6597

6698
```{r}
6799
#| label: socp
@@ -74,8 +106,20 @@ round(value(xs), 4) # 1/3 each -> the min-norm point on the simplex
74106

75107
The most general of the four: constraints that a **matrix is positive semidefinite**
76108
($\succeq 0$), written with `%>>%`. Eigenvalue optimization, covariance/correlation
77-
problems, and relaxations of hard combinatorial problems are SDPs. Here we compute the
78-
largest eigenvalue of a matrix $M$ as $\min\{t : tI - M \succeq 0\}$:
109+
problems, and relaxations of hard combinatorial problems are SDPs.
110+
111+
$$
112+
\begin{array}{ll}
113+
\text{minimize} & c^\top x\\
114+
\text{subject to} & x_1 F_1 + \cdots + x_n F_n + G \;\succeq\; 0,\quad Ax = b,
115+
\end{array}
116+
$$
117+
118+
where $F_1,\dots,F_n, G$ are symmetric matrices and $\succeq 0$ means positive
119+
semidefinite; that matrix inequality is a *linear matrix inequality* (LMI).
120+
Equivalently, the variable can itself be a symmetric matrix $X \succeq 0$. Here we
121+
compute the largest eigenvalue of a matrix $M$ as $\min\{t : tI - M \succeq 0\}$ (an
122+
LMI in the scalar $t$):
79123

80124
```{r}
81125
#| label: sdp

05-elastic-net.qmd

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -32,7 +32,7 @@ combines them:
3232

3333
$$
3434
\underset{\beta}{\mbox{minimize}}\quad
35-
\frac{1}{2m}\|y - X\beta\|_2^2
35+
\frac{1}{2n}\|y - X\beta\|_2^2
3636
\;+\;\lambda\Big(\tfrac{1-\alpha}{2}\|\beta\|_2^2 + \alpha\|\beta\|_1\Big).
3737
$$
3838

06-change-the-loss.qmd

Lines changed: 5 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -45,12 +45,12 @@ responses** --- nasty outliers by construction.
4545
```{r}
4646
#| label: huber-data
4747
set.seed(1289)
48-
m <- 450; M <- 1; pflip <- 0.10
48+
n <- 450; M <- 1; pflip <- 0.10
4949
beta_true <- 5
50-
X <- matrix(rnorm(m), ncol = 1)
50+
X <- matrix(rnorm(n), ncol = 1)
5151
y_true <- X * beta_true
52-
flip <- 2 * rbinom(m, 1, 1 - pflip) - 1 # -1 flips the sign
53-
y <- flip * y_true + rnorm(m)
52+
flip <- 2 * rbinom(n, 1, 1 - pflip) - 1 # -1 flips the sign
53+
y <- flip * y_true + rnorm(n)
5454
```
5555

5656
Now fit the *same* problem with two losses --- squared and Huber:
@@ -151,7 +151,7 @@ regression is for.
151151

152152
**Exercise (⭐).** The Huber threshold $M$ interpolates between two losses. What does
153153
Huber regression become as $M \to 0$? As $M \to \infty$? Check the small-$M$ case
154-
against the LAD fit from @sec-gentle-intro.
154+
against LAD --- the loss you met in @sec-gentle-intro.
155155

156156
::: {.callout-tip collapse="true"}
157157
## Solution
@@ -175,8 +175,6 @@ The two nearly coincide.
175175
quantiles, absolute → median.
176176
- Each is a **one-line change** to the objective --- the modeling vocabulary you are
177177
building transfers directly.
178-
- A kinked loss can trip a particular solver; naming a robust solver
179-
(`solver = "CLARABEL"`) is sometimes the easiest fix.
180178

181179
See [Huber](https://cvxr.rbind.io/examples/regression/huber-regression.html) and
182180
[Quantile](https://cvxr.rbind.io/examples/regression/quantile-regression.html).

07-logistic.qmd

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -27,7 +27,7 @@ log-likelihood
2727

2828
$$
2929
\underset{\beta}{\mbox{maximize}}\quad
30-
\sum_{i=1}^m \Big[\, y_i\, x_i^\top\beta - \log\big(1 + e^{x_i^\top\beta}\big) \Big].
30+
\sum_{i=1}^n \Big[\, y_i\, x_i^\top\beta - \log\big(1 + e^{x_i^\top\beta}\big) \Big].
3131
$$
3232

3333
The first term is linear; the second is a sum of $\log(1 + e^{z})$ terms.
@@ -58,12 +58,12 @@ $\sum_i \log(1+e^{x_i^\top\beta})$ as `sum(logistic(X %*% beta))`:
5858
```{r}
5959
#| label: fit
6060
set.seed(42)
61-
m <- 500; n <- 8
62-
beta_true <- rnorm(n)
63-
X <- matrix(rnorm(m * n), m, n)
64-
y <- as.numeric(1 / (1 + exp(-X %*% beta_true)) > runif(m)) # 0/1 labels
61+
n <- 500; p <- 8
62+
beta_true <- rnorm(p)
63+
X <- matrix(rnorm(n * p), n, p)
64+
y <- as.numeric(1 / (1 + exp(-X %*% beta_true)) > runif(n)) # 0/1 labels
6565
66-
beta <- Variable(n)
66+
beta <- Variable(p)
6767
loglik <- sum(y * (X %*% beta)) - sum(logistic(X %*% beta))
6868
prob <- Problem(Maximize(loglik))
6969
psolve(prob)

08-dpp.qmd

Lines changed: 11 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -68,7 +68,7 @@ good <- Problem(Minimize(g * w), list(w == g * z, z >= 1))
6868
cat("Is DPP?", is_dpp(good), "\n")
6969
```
7070

71-
## The payoff: watch the compile time vanish
71+
## The payoff: watch the compile time drop
7272

7373
Recall from @sec-gentle-intro that `verbose = TRUE` splits each solve into a
7474
**compile time** and a **solver time** in its summary. Parameters target the
@@ -79,9 +79,9 @@ a `Parameter`, and solve it the first time:
7979
#| label: verbose-first
8080
#| message: true
8181
set.seed(1)
82-
n <- 120; m <- 40
83-
A <- matrix(rnorm(n * m), n, m); b <- rnorm(n)
84-
x <- Variable(m)
82+
n <- 120; p <- 40
83+
A <- matrix(rnorm(n * p), n, p); b <- rnorm(n)
84+
x <- Variable(p)
8585
error <- sum_squares(A %*% x - b)
8686
8787
gamma <- Parameter(nonneg = TRUE)
@@ -101,7 +101,7 @@ value(gamma) <- 0.7
101101
psolve(param_prob, verbose = TRUE) # re-solve: compilation is cached
102102
```
103103

104-
The compile time has **collapsed** --- CVXR reused the cached parameter-to-data
104+
The compile time has **dropped** --- CVXR reused the cached parameter-to-data
105105
mapping and went straight to the solver. A freshly built problem, by contrast, pays
106106
the full compile *every* time.
107107

@@ -122,7 +122,7 @@ bench::mark(
122122
```
123123

124124
On our machine the DPP re-solve is roughly **2–3× faster per solve** and allocates
125-
less memory --- precisely the compile time we just watched disappear. Over a sweep of
125+
less memory --- precisely the compile time we just watched drop. Over a sweep of
126126
$K$ values you save that difference $K$ times.
127127

128128
::: callout-note
@@ -144,11 +144,11 @@ re-solve --- the path looks identical, but only one compilation happens:
144144
#| label: dpp-path
145145
lambda_vals <- 10^seq(-3, 1, length.out = 30)
146146
lambda <- Parameter(nonneg = TRUE)
147-
beta <- Variable(m)
148-
path_prob <- Problem(Minimize(error + lambda * p_norm(beta, 1)))
147+
beta <- Variable(p)
148+
path_prob <- Problem(Minimize(sum_squares(A %*% beta - b) + lambda * p_norm(beta, 1)))
149149
150150
betas <- sapply(lambda_vals, function(l) { value(lambda) <- l; psolve(path_prob); value(beta) })
151-
dim(betas) # m coefficients x number of lambda values
151+
dim(betas) # p coefficients x number of lambda values
152152
```
153153

154154
Whenever you plan to solve the *same shape* of problem many times --- a path, a
@@ -157,11 +157,11 @@ grid search, a simulation, a bootstrap --- reach for a `Parameter`.
157157
## The other half of the solve: warm starts
158158

159159
Recall that every solve splits into a **compile time** and a **solver time**. DPP
160-
attacks the *compile* time --- it reuses the canonicalized problem structure so we
160+
mitigates the *compile* time --- it reuses the canonicalized problem structure so we
161161
never re-stuff the matrices. But the solver still starts **cold** each time, refining
162162
an internal guess from scratch until it converges.
163163

164-
A **warm start** attacks the *other* half. A numerical solver is iterative; if you
164+
A **warm start** addresses the *other* half. A numerical solver is iterative; if you
165165
hand it a starting point already near the answer, it converges in fewer iterations.
166166
In a sweep over a fine grid the solution at one $\lambda$ is an excellent guess for
167167
the next, so the previous answer is exactly the start point we want. In CVXR that is

09-portfolio.qmd

Lines changed: 14 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,8 @@
11
# Portfolio Optimization {#sec-portfolio}
22

3-
<!-- LIVE. Mean-variance (Markowitz) efficient frontier + mean-CVaR, on a real
4-
12-asset price subset bundled in data/stock_prices.csv. Example and data
5-
adapted from Riskfolio-Lib Tutorial 1 (Dany Cajas). Originals: Markowitz
3+
<!-- LIVE. Mean-variance (Markowitz) efficient frontier + mean-CVaR, on daily
4+
adjusted prices for the 12 largest U.S. companies by market cap (Yahoo
5+
Finance, 2022-2024), bundled in data/stock_prices.csv. Models: Markowitz
66
(1952) mean-variance; Rockafellar-Uryasev (2000) CVaR. Verified vs CVXR 1.9.1. -->
77

88
::: callout-note
@@ -13,12 +13,11 @@ Value-at-Risk (CVaR)** --- all from real return data.
1313
:::
1414

1515
::: callout-tip
16-
## Credits
17-
The example and the bundled price data are adapted from
18-
[**Riskfolio-Lib**](https://github.com/dcajasn/Riskfolio-Lib) (Dany Cajas) ---
19-
"Tutorial 1: Classic Mean-Risk Optimization" [@riskfolio], itself built on CVXPY. The
20-
underlying models are due to @Markowitz:1952 (mean-variance) and
21-
@RockafellarUryasev:2000 (CVaR).
16+
## Data and models
17+
The models are classics: mean-variance optimization is due to @Markowitz:1952 and the
18+
CVaR reformulation to @RockafellarUryasev:2000. The bundled price data are daily
19+
adjusted closing prices for the 12 largest U.S. companies by market capitalization
20+
(source: Yahoo Finance, 2022--2024).
2221
:::
2322

2423
```{r}
@@ -47,8 +46,9 @@ return; large $\gamma$ avoids risk.
4746

4847
## Real data
4948

50-
We use ~3 years of daily prices for 12 well-known stocks (a subset of the Riskfolio-Lib
51-
dataset), and turn prices into simple returns, a mean vector, and a covariance matrix:
49+
We use ~3 years (2022--2024) of daily adjusted prices for the **12 largest U.S.
50+
companies by market capitalization**, and turn prices into simple returns, a mean
51+
vector, and a covariance matrix:
5252

5353
```{r}
5454
#| label: data
@@ -133,8 +133,8 @@ ggplot(comp, aes(gamma, weight, fill = asset)) +
133133
```
134134

135135
The left-most bar (greedy, low $\gamma$) leans on one or two names; as $\gamma$ grows the
136-
bars fragment into the full palette. That shift *is* the efficient frontier, seen from
137-
the portfolio's side.
136+
weight spreads across a broader, more defensive set of names. That shift *is* the
137+
efficient frontier, seen from the portfolio's side.
138138

139139
## A different risk measure: CVaR
140140

@@ -238,6 +238,5 @@ The cap spreads weight off the assets that the unconstrained solution concentrat
238238
- Swapping variance for **CVaR** is just a new objective --- CVXR makes alternative risk
239239
measures, position limits, and policy constraints one-liners.
240240

241-
For many more portfolio models (risk parity, Black-Litterman, worst-case, drawdown
242-
measures), see [Riskfolio-Lib](https://github.com/dcajasn/Riskfolio-Lib) and the
241+
For more portfolio models and risk measures, see the
243242
[CVXR portfolio example](https://cvxr.rbind.io/examples/finance/portfolio-optimization.html).

10-invcov.qmd

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -55,14 +55,14 @@ invert it to get the covariance $R$, and sample from $N(0, R)$:
5555
```{r}
5656
#| label: data
5757
set.seed(1)
58-
n <- 10 # number of variables
59-
m <- 1000 # number of samples
58+
p <- 10 # number of variables
59+
n <- 1000 # number of samples
6060
61-
A <- rsparsematrix(n, n, 0.15, rand.x = rnorm)
62-
Strue <- as.matrix(A %*% t(A) + 0.05 * diag(n)) # sparse, symmetric, positive definite
61+
A <- rsparsematrix(p, p, 0.15, rand.x = rnorm)
62+
Strue <- as.matrix(A %*% t(A) + 0.05 * diag(p)) # sparse, symmetric, positive definite
6363
R <- solve(Strue) # the covariance to sample from
6464
65-
x <- matrix(rnorm(m * n), m, n) %*% chol(R) # rows are i.i.d. N(0, R)
65+
x <- matrix(rnorm(n * p), n, p) %*% chol(R) # rows are i.i.d. N(0, R)
6666
Q <- cov(x) # sample covariance
6767
6868
sum(abs(Strue) > 1e-4) # how many true nonzeros (of 100)
@@ -76,7 +76,7 @@ each solve, zero out numerically negligible entries to read off the recovered pa
7676
```{r}
7777
#| label: fit
7878
alphas <- c(10, 6, 4, 1)
79-
S <- Variable(c(n, n), PSD = TRUE)
79+
S <- Variable(c(p, p), PSD = TRUE)
8080
obj <- Maximize(log_det(S) - matrix_trace(S %*% Q))
8181
8282
S_est <- lapply(alphas, function(a) {
@@ -106,7 +106,7 @@ mask of the truth alongside each fit:
106106
#| fig-height: 2.4
107107
#| fig-cap: "Nonzero pattern (dark = nonzero) of the true precision matrix and the estimates. Smaller alpha -> sparser. Around alpha = 4 the estimate matches the truth's sparsity."
108108
mask_df <- function(M, lab) {
109-
data.frame(i = rep(seq_len(n), n), j = rep(seq_len(n), each = n),
109+
data.frame(i = rep(seq_len(p), p), j = rep(seq_len(p), each = p),
110110
nz = factor(as.integer(as.vector(M) != 0)), panel = lab)
111111
}
112112
panels <- rbind(

11-mip.qmd

Lines changed: 17 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -114,14 +114,21 @@ a CRAN solver simply by writing the objective in a form SCIP accepts. SOC plus
114114
integers is SCIP's home turf; the second-order cone is doing the work the quadratic
115115
objective could not.
116116

117-
::: callout-warning
118-
## Two traps worth remembering
119-
- **Never send an SOC-containing mixed-integer problem to HiGHS.** It does not error
120-
--- it silently returns a *wrong* `infeasible`. SOC + integers needs **SCIP**.
121-
- **A quadratic objective with integers (MIQP)** routes through CVXR's QP path, which
122-
only commercial solvers serve. On CRAN, **reformulate the quadratic as an SOC**
123-
(above) and use SCIP --- same answer. (Commercial solvers can take the squared form
124-
directly, if you have one.)
117+
::: callout-note
118+
## Two current limitations (both may change)
119+
These reflect the state of the tooling *today*, not permanent rules --- worth knowing
120+
now, but re-test rather than assume, since solvers and CVXR's routing keep evolving.
121+
122+
- **HiGHS with an SOC inside a mixed-integer problem.** As of this writing, HiGHS does
123+
not support a second-order cone in a mixed-integer model, and rather than erroring it
124+
can silently return a *wrong* `infeasible` --- so for now, send SOC + integers to
125+
**SCIP**. HiGHS is under active development, so this is the likelier of the two to
126+
change; check before relying on it.
127+
- **MIQP routing.** A quadratic objective with integers currently travels through
128+
CVXR's QP path, where only commercial solvers are registered. On CRAN today,
129+
**reformulate the quadratic as an SOC** (above) and use SCIP --- same answer.
130+
(Commercial solvers take the squared form directly.) This routing is a CVXR
131+
implementation choice that could also change.
125132
:::
126133

127134
## Exercise
@@ -245,8 +252,8 @@ rules written as constraints. (Swap in any puzzle by editing the `0`-padded matr
245252
- **MI-SOCP** (a cone plus integers) runs on **SCIP**. An integer problem with a
246253
*quadratic* objective (MIQP) goes to a commercial solver through CVXR's QP path ---
247254
but **reformulate the quadratic as an SOC** and SCIP solves it on CRAN, same answer.
248-
- Match the solver to the problem --- and never trust an SOC sent to a pure-MILP
249-
solver like HiGHS.
255+
- Match the solver to the problem --- today a pure-MILP solver like HiGHS can't handle
256+
an SOC (and may not tell you), so send those to SCIP for now.
250257
- Some integer problems have **no objective at all**: logical puzzles like Sudoku are
251258
pure *feasibility* programs, where every rule is an "exactly one" sum of booleans.
252259

0 commit comments

Comments
 (0)