Search approach to solve optimization problem with only a minimum where time series get scaled

Question

Currently, I am working on a relatively simple optimization problem:

There is a set of time series (red) that get summed up to a cumulated time series (blue). The red time series have different forms (simplified here).

During the optimization process, I want to scale the time series up or down so that my optimization solution is a set of $n$ factors for scaled the individual time series.

It applies: $\boldsymbol y = \sum\limits_{i=1}^{n}(\boldsymbol x_i \cdot s_i)$, where the $\boldsymbol x_i$ are the individual time series and $\pmb y$ is the cumulated time series. The factors are not allowed to be negative and together they must add up to $1$.

The function for evaluating the cumulated time series calculates the MSE between the cumulated time series $\boldsymbol y$ and another time series $\boldsymbol z$: $$\text{MSE} = \frac{1}{m}\sum_{i = 1}^{m}{(z_i - y_i)^2}.$$

For my first tests, I used the heuristic Simulated Annealing to set the factors iteratively. This works, however, it seems to be a bit oversized for this problem, because it has only one minimum.

How would you solve this problem?

Answer 1

Your problem actually comes down to a constrained linear regression problem where $z$ is your dependent variable, the $x_j$ for $j=1,\dots,n$ are your independent variables and $s$ is your vector with regression coefficients. Without any constraints, the unique solution optimizing the MSE is the least squares estimator: $$\hat{s}=(X'X)^{-1}X'z, \quad \text{where}\quad X=\begin{bmatrix}x_1&x_2&\cdots& x_n\end{bmatrix}$$

The requirement that all weights add up to one can be formalized as the constraint $\iota's = 1$, with $\iota$ the all-ones vector. There also exists a closed-form solution to the least squares problem with such an equality constraint, see Wikipedia.

Adding the inequality constraints is more complicated. My approach would be to use an Active Set algorithm. Here, we start by solving the problem without the inequality constraints. Then, we check whether there are any time series (regression coefficients) with negative weights. We then remove these time series from the problem and resolve. We continue until all weights are positive.

I'm pretty confident this method would result in finding the global optimum (although it would take me some time to write down a proof).

As an alternative method, if you do not like programming an algorithm yourself, you could also pass your original problem to a commercial solver. For example, both CPLEX and Gurobi can solve convex quadratic programs with linear constraints fairly well.

Answer 2

The bounty convinced me to compete with Rolf's excellent answer, which is exactly how I would approach the problem myself. Next to CPLEX and Gurobi, it also worth noting that MATLAB and Octave provide the function fmincon, which can also be used to solve your problem directly, and SPSS provides Constrained Nonlinear Regression (which also allows for constrained linear regression). Because of all the available software, the method that I propose below is probably not practically useful.

You have the following constraints: \begin{align} \sum_{i=1}^n s_i &= 1\\ s_i &\ge 0 \end{align}

Without these constraints, your problem would be an ordinary least squares regression problem, as pointed out by Rolf. Without the equality constraint, your problem would be a non-negative least squares problem (NNLS). I will assume that you have access to statistics software that can solve NNLS, but is somehow unable to handle the equality constraint.

Let $\lambda \in \mathbb{R}$ be a parameter, with initial value $\lambda = 1$. Add the artificial time period $m+1$ to every time series. At time $m+1$, give $z$ the value $\lambda$, and give each $x_i$ the value 1. Due to this extra time period, we get the additional MSE term $$\frac{1}{m+1} (z_{m+1} - y_{m+1})^2 = \frac{1}{m+1} \left(\lambda - \sum_{i=1}^n s_i\right)^2.$$

This MSE term is an incentive for $\sum_{i=1}^n s_i$ to be close to $\lambda=1$. Otherwise, this MSE term would become large. Now use NNLS to solve the problem with the artificial time period. It is guaranteed that $s_i \ge 0$, and hopefully, $\sum_{i=1}^n s_i$ will be close to 1.

If $\sum_{i=1}^n s_i < 1$, then increase $\lambda$ and try again. If $\sum_{i=1}^n s_i > 1$, then decrease $\lambda$ and try again. It can be shown that $\sum_{i=1}^n s_i$ is non-decreasing in $\lambda$, so you can use the bisection search to find the value of $\lambda$ that nudges the sum to be exactly 1. You could even search for $\lambda$ by hand.

When NNLS gives you an answer with $\sum_{i=1}^n s_i = 1$ for some $\lambda$, then it is guaranteed that you have found the optimal weights!

---

Justification

The above method follows from using the augmented Lagrangian method and solving the problem $$\max_{\mu \in \mathbb{R}} \min_{s_i \ge 0} \mathcal{L}(s,\mu) = \text{MSE} + \mu\left(1-\sum_{i=1}^n s_i\right) + \frac{\rho}{2}\left(1-\sum_{i=1}^n s_i\right)^2.$$ We choose $\lambda= \frac{\mu}{2}+1$ and $\rho = 1$, as we may. After rewriting, the inner minimization is equivalent to the NNLS problem, and the outer maximization is the search for $\lambda$. Guaranteed optimality follows from the optimality of the augmented Lagrangian method.