A multiobjective continuation method to compute the regularization path
  of deep neural networks
Sparsity is a highly desired feature in deep neural networks (DNNs) since it
ensures numerical efficiency, improves the interpretability of models (due to
the smaller number of relevant features), and robustness. For linear models, it
is well known that there exists a \emph{regularization path} connecting the
sparsest solution in terms of the $\ell^1$ norm, i.e., zero weights and the
non-regularized solution. Very recently, there was a first attempt to extend
the concept of regularization paths to DNNs by means of treating the empirical
loss and sparsity ($\ell^1$ norm) as two conflicting criteria and solving the
resulting multiobjective optimization problem for low-dimensional DNN. However,
due to the non-smoothness of the $\ell^1$ norm and the high number of
parameters, this approach is not very efficient from a computational
perspective for high-dimensional DNNs. To overcome this limitation, we present
an algorithm that allows for the approximation of the entire Pareto front for
the above-mentioned objectives in a very efficient manner for high-dimensional
DNNs with millions of parameters. We present numerical examples using both
deterministic and stochastic gradients. We furthermore demonstrate that
knowledge of the regularization path allows for a well-generalizing network
parametrization. To the best of our knowledge, this is the first algorithm to
compute the regularization path for non-convex multiobjective optimization
problems (MOPs) with millions of degrees of freedom.
Augustina C. Amakor, Konstantin Sonntag, Sebastian Peitz
# A multiobjective continuation method

###### Abstract

Sparsity is a highly desired feature in deep neural networks (DNNs) since it ensures numerical efficiency, improves the interpretability of models (due to the smaller number of relevant features), and robustness. In machine learning approaches based on linear models, it is well known that there exists a connecting path between the sparsest solution in terms of the \(\ell^{1}\) norm (i.e., zero weights) and the non-regularized solution, which is called the regularization path. Very recently, there was a first attempt to extend the concept of regularization paths to DNNs by means of treating the empirical loss and sparsity (\(\ell^{1}\) norm) as two conflicting criteria and solving the resulting multiobjective optimization problem. However, due to the non-smoothness of the \(\ell^{1}\) norm and the high number of parameters, this approach is not very efficient from a computational perspective. To overcome this limitation, we present an algorithm that allows for the approximation of the entire Pareto front for the above-mentioned objectives in a very efficient manner. We present numerical examples using both deterministic and stochastic gradients. We furthermore demonstrate that knowledge of the regularization path allows for a well-generalizing network parametrization.

## Introduction

Machine Learning (ML) and in particular deep neural networks (DNNs) are nowadays an integral part of numerous applications such as data classification, image recognition, time series prediction, and language processing. Their importance continues to grow at great speed across numerous disciplines and applications, and the increase in available computational capacity allows for the construction of larger models. However, these advances also increase the challenges regarding the construction and training of DNNs, e.g., the required training data, the training efficiency, and the adaptability to changing external factors. This leads to the task of simultaneously fulfilling numerous, sometimes contradictory goals in the best possible way.

Multiobjective optimization (MOO) addresses the problem of optimizing several conflicting objectives. The issue of having to trade between multiple, conflicting criteria is a universal problem, such as the need to have an optimal tradeoff between cost and quality in a production process. In a similar manner, conflicting criteria occur in various ways in ML. The main task of MOO is thus to identify the set of optimal trade-off solutions (the _Pareto set_) between these conflicting criteria. This concerns multitask problems [14], but also the training itself, where we want to minimize the training loss, obtain sparse models and improve generalization.

Interestingly, we found that many papers on multicriteria machine learning do not address the true nature of multiobjective optimization. The reason for this, is that when choosing very large neural network architectures or even considering task-specific layers [14], the different tasks are no longer conflicting. The network is simply too powerful such that both task can be optimally handled simultaneously and there is no tradeoff. From an optimization point of view, the Pareto front collapses into a single point. However, considering the strongly growing carbon footprint of AI [10], there is a growing interest in smaller models that are better tailored to specific tasks. This is why we propose the use of models that are smaller and adapted to a certain task. While this will reduce the general applicability, multicriteria optimization can help us to determine a set of compromising network architectures, such that we can adapt networks to specific situations online (multicriteria decision making).

The joint consideration of loss and \(\ell^{1}\) regularization is well-studied for linear systems. However, it is much less understood for the nonlinear problems that we face in deep learning. When treating the \(\ell^{1}\) regularization problem as a multiobjective optimization problem (MOP), a popular approach to obtain the entire solution set (the _Pareto set_) is via continuation methods [15, 16]. These continuation methods usually consist of a predictor step (along the tangent direction of the Pareto set) and a corrector step that converges to a new point on the Pareto set close by. However, as the \(\ell^{1}\) norm is non-smooth, classical manifold continuation techniques fail. Due to this fact, a first extension of regularization paths from linear to nonlinear models was presented only very recently in [1], where continuation methods were extended to non-smooth objective functions. Although this extension provides a rigorous treatment of the problem, it results in a computationally expensive algorithm, which renders it impractical for DNNs of realistic dimensions.

As research on regularization paths for nonlinear loss functions has focused on small-scale learning problems until
now, this work is concerned with large-scale machine learning problems. The main contributions of this paper include:

* The extension of regularization paths from linear models to high dimensional nonlinear deep learning problems.
* The demonstration of the usefulness of MOO continuation methods for the generalization properties of DNNs.
* A step towards more resource-efficient ML by beginning with very sparse networks and slowly increasing the number of weights. This is in complete contrast to the standard pruning approaches, where we start very large and then reduce the number of weights.

A comparison of our approach is further made with the standard approach of weighting the individual losses using an additional hyperparameter.

The remainder of the paper is organized as follows. First, we discuss related works, before introducing our continuation method. We then present a detailed discussion of our extensive numerical experiments.

## Related works

### Multiobjective optimization

In the last decades, many approaches have been introduced for solving nonlinear [11], nonsmooth [12], non-convex [11], or stochastic [12] multiobjective optimization problems, to name just a few special problem classes. In recent years, researchers have further begun to consider multiobjective optimization in machine learning [21, 22, 23] and deep learning [24, 25]. We provide an overview of the work that is most pertinent to our own and direct readers to the works done by [24] and [2]. Ehrgott2008 described MOO as the process of optimizing multiple, sometimes conflicting objectives simultaneously by finding the set of points that are _not dominated_1 by any other feasible points. Different methods have been proposed to obtain these Pareto optimal solutions such as _evolutionary algorithms_ which use the evolutionary principles inspired by nature by evolving an entire population of solutions [1]. _Gradient-based techniques_ extend the well-known gradient techniques from single-objective optimization to multiple criteria problems [10]. _Set Oriented Methods_ provide an approximation of the Pareto set through box coverings and often suffer from the curse of dimensionality which makes applications in ML very expensive [13]. Another deterministic approach are _scalarization methods_ which involve the transformation of MOP into (a series of) single objective optimization problems [1]. These can then be solved using standard optimization techniques. Examples include the weighted sum method, epsilon constraint, Chebyshev scalarization, goal programming, min-max and augmented epsilon constraint method [12, 23]. Some drawbacks exist for the latter methods, most notably the incapability of some to handle non-convex problems, as well as the challenge to obtain equidistant coverings. Notably, these drawbacks are most severe in the weighted sum method, which is by far the most widely applied approach in ML when considering multiple criteria at the same time (such as regularization terms). Moreover, the weighting has to be made a priori, which makes selecting useful trade-offs much harder [10, 24].

Footnote 1: A point is said to _dominate_ another point if it performs at least as good as the other point across all objectives and strictly better than it in at least one objective. It is referred to as a _Pareto optimal_ or efficient point.

_Continuation methods_ are another approach for the computation of the solutions of MOPs. These methods were initially developed to solve complex problems starting from simple ones, using homotopy approaches2. Also, they were used in predictor-corrector form to track manifolds [11, 12]. In the multiobjective setting, one can show that the Pareto set is a manifold as well3, if the objectives are sufficiently regular [11]. Hence, the continuation method then becomes a technique to follow along the Pareto set in two steps; **a)** predictor step along the tangent space of the Pareto set which forms in the smooth setting a manifold of dimension _m - I_ (where \(m\) is the number of objective functions) [12] and **b)** a corrector step that obtains the next point on the Pareto set which is often achieved using multiobjective gradient descent. Schutze, Dell'Aere, and Dellnitz (2005) used the predictor-corrector approach to find points that satisfy the Karush-Kuhn-Tucker (KKT) condition and to further identify other KKT points in the neighbourhood.

Footnote 2: These are approaches that involve starting at a simple-to-calculate solution and then continuously vary some parameter to increase the problem difficulty step by step, until finally arriving at the solution of the original problem, which is often very hard to compute directly [12].

Footnote 3: To be more precise, the set of points satisfying the Karush-Kuhn-Tucker (KKT) necessary conditions for optimality along with the corresponding KKT multipliers form a manifold of dimension \(m-1\).

The continuation method has further been extended to regularization paths in machine learning (ML) for linear models [13, 12]. Recently, regularization paths for nonlinear models have been introduced although limited to small dimensions since dealing with non-smoothness is difficult [10].

### Multicriteria machine learning

In the context of multicriteria machine learning, several approaches have been used such as _evolutionary algorithms_[1, 1, 12] or _gradient-based methods_[24, 12]. However, only a few of these methods address high-dimensional deep learning problems or attempts to compute the entire Pareto front. Furthermore, as discussed in the introduction, many researchers have introduced architectures that are so powerful that even with the inclusion of the task-specific parts,
the Pareto front collapses into a single point (e.g., Sener and Koltun (2018)). As we want to pursue a more sustainable path here, we are interested in truly conflicting criteria. The entire regularization path for DNNs (i.e., a MOP with training loss versus \(\ell^{1}\) norm) was computed in (Bieker, Gebken, and Peitz 2022). However, even though the algorithm provably yields the Pareto front, the computation becomes intractable in very high dimensions. Hence, the need to develop an efficient method to find the entire regularization path and Pareto front for DNNs.

For large-scale MOO, gradient-based methods have proven to be the most efficient. Examples are the steepest descent method (Fliege and Svaiter 2000), projected gradient method (Drummond and Iusem 2004), proximal gradient method (Tanabe, Fukuda, and Yamashita 2019) and recently accelerated proximal gradient (Tanabe, Fukuda, and Yamashita 2022). Previous approaches to MOO problems often assume smoothness/differentiability of the objective functions but the \(\ell^{1}\) norm is not differentiable, so we use the multiobjective proximal gradient (MPG) to ensure convergence. MPG has been described by Tanabe, Fukuda, and Yamashita (2019) as a descent method for unconstrained MOPs where each objective function can be written as the sum of a convex and smooth function, and a proper, convex, and lower semicontinuous but not necessarily differentiable one. Simply put, MPG combines the proximal point and the steepest descent method and is a generalization of the iterative shrinkage-thresholding algorithms (ISTA) (Combettes and Wajs 2005; Beck and Teboulle 2009) to multiobjective optimization problems.

## Continuation Method

### Some basics on multiobjective optimization

A multiobjective optimization problem can be mathematically formalized as

\[\min_{\theta\in\mathbb{R}^{n}}\left[\begin{array}{c}F_{1}(\theta)\\ \vdots\\ F_{m}(\theta)\end{array}\right],\] (MOP)

where \(F_{i}:\mathbb{R}^{n}\rightarrow\mathbb{R}\ \forall\ i=1,\ldots,m\) are the generally conflicting objective functions and \(\theta\) the parameters we are optimizing over. In general, there does not exist a single point that minimizes all criteria simultaneously, the solution to MOP is the entire set of optimal compromises, the so-called Pareto set.

**Definition 1** (Miettinen (1998)): _Consider the optimization problem (MOP)._

1. _A point_ \(\theta^{*}\in\mathbb{R}^{n}\) _is_ Pareto optimal _if there does not exist another point_ \(\theta\in\mathbb{R}^{n}\) _such that_ \(f_{i}(\theta)\leq f_{i}(\theta^{*})\) _for all_ \(i=1,\ldots,m\)_, and_ \(f_{j}(\theta)<f_{j}(\theta^{*})\) _for at least one index_ \(j\)_. The set of all Pareto optimal points is the_ Pareto set_, which we denote by_ \(P\)_. The set_ \(F(P)\subset\mathbb{R}^{m}\) _in the image space is called the_ Pareto front_._
2. _A point_ \(\theta^{*}\in\mathbb{R}^{n}\) _is_ weakly Pareto optimal _if there does not exist another point_ \(\theta\in\mathbb{R}^{n}\) _such that_ \(f_{i}(\theta)<f_{i}(\theta^{*})\) _for all_ \(i=1,\ldots,m\)_. The set of all weakly Pareto optimal points is the_ weak Pareto set_, which we denote by_ \(P_{w}\) _and the set_ \(F(P_{w})\) _is the_ weak Pareto front_._

The above definition is not very helpful for gradient-based algorithms. In this case, we need to rely on first order optimality conditions.

**Definition 2** (Hilermeier (2001)): _A point \(\theta^{*}\in\mathbb{R}^{n}\) is called Pareto critical if there exist an \(\alpha\in\mathbb{R}^{m}\) with \(\alpha_{i}\geq 0\) for all \(i=1,\ldots,m\) and \(\sum_{i=1}^{m}\alpha_{i}=1\) satisfying_

\[\sum_{i=1}^{m}\alpha_{i}\nabla F_{i}(\theta^{*})=0. \tag{1}\]

_Equation (1) is generally referred to as the Karush-Kuhn-Tucker (KKT) condition._

**Remark 3**: _All Pareto optimal points are Pareto critical but the reverse is in general, not the case. If all objective functions are convex then each Pareto critical point is weakly Pareto optimal._

In this work, we consider two objective functions, namely the empirical loss and the \(\ell^{1}\) norm of the neural network weights. The Pareto set connecting the individual minima (at least locally), is also known as the regularization path. In terms of a MOO, we are looking for the Pareto set of

\[\min_{\theta\in\mathbb{R}^{n}}\left[\begin{array}{c}\mathbb{E}_{(x,y) \sim\mathcal{D}}[\mathcal{L}(f(\theta,x),y)]\\ \frac{1}{n}\|\theta\|_{1}\end{array}\right],\] (MoDNN)

where \((x,y)\in\mathcal{X}\times\mathcal{Y}\) is labeled data following a distribution \(\mathcal{D}\), the function \(f:\mathbb{R}^{n}\times\mathcal{X}\rightarrow\mathcal{Y}\) is a parameterized model and \(\mathcal{L}(\cdot,\cdot)\) denotes a loss function. The second objective is the weighted \(\ell^{1}\) norm \(\frac{1}{n}\|\theta\|_{1}\) to ensure sparsity as in the well known LASSO. Our goal is to solve (MoDNN), by which we are able to obtain the regularization path. However, this problem is challenging due to the fact that the \(\ell^{1}\) norm is not differentiable. In order to solve this problem, we explicitly take the non-smoothness of the second objective into account, which is done by using a proximal gradient method (Tanabe, Fukuda, and Yamashita 2019).

**Remark 4**: _A common approach to solve the problem (MOP) is the use of the weighted sum method. Which includes all types of regularization when considering multiple criteria in DNN training terms. E.g.,_

\[F(\theta)=\sum_{i=1}^{m}\lambda_{i}F_{i}(\theta), \tag{2}\]

_with \(\lambda_{i}\geq 0\) and \(\sum_{i=1}^{m}\lambda_{i}=1\), where \(\lambda\) is an additional hyperparameter and we solve \(\min_{\theta\in\mathbb{R}^{n}}F(\theta)\)._

### Proximal gradient method

Given functions of the form \(F_{i}=f_{i}+g_{i}\), such that \(f_{i}:\mathbb{R}^{n}\rightarrow\mathbb{R}\) is convex and smooth, and \(g_{i}\) is convex and non-smooth with computable proximal operator \(\text{prox}_{g}\),

**Definition 5**: _proximal operator;_

_Given a convex function \(g:\mathbb{R}^{n}\rightarrow\mathbb{R}\), the proximal operator is expressed as_

\[\text{prox}_{g}(\theta)=\operatorname*{arg\,min}_{\phi\in\mathbb{R}^{n}}\left\{g (\phi)+\frac{1}{2}\|\phi-\theta\|_{2}^{2}\right\}\]
Tanabe, Fukuda, and Yamashita (2019) proposed an algorithm that aims at efficiently solving the resulting MOP using _Algorithm 1_.

**Theorem 6** (Tanabe, Fukuda, and Yamashita (2019)): _Let \(f_{i}\) be convex with \(L\)-Lipschitz continuous gradient and let \(g_{i}\) be proper, convex and lower semi continuous for all \(i=1,\ldots,m\). Then, every accumulation point of the sequence \(\{\theta^{k}\}\) computed by Algorithm 1 is Pareto critical. In addition, if the level set of one \(F_{i}\) is bounded, the sequence \(\{\theta^{k}\}\) has at least one accumulation point._

**Remark 7**: _Proximal gradient iterations are "forward-backward" iterations, "forward" referring to the gradient step and "backward" referring to the proximal step. The proximal gradient method is designed for problems where the objective functions include a smooth and a nonsmooth component, which is suitable for optimization with a sparsity-promoting regularization term._

**Remark 8**: _In the problem (MoDNN) we are in the situation \(F_{1}(\theta)=f_{1}(\theta)=\mathbb{E}_{(x,y)\sim\mathcal{D}}[\mathcal{L}(f( \theta,x),y)]\) and \(F_{2}(\theta)=g_{2}(\theta)=\frac{1}{n}\|\theta\|_{1}\). The proximal operator \(\text{prox}_{\frac{1}{n}\|\cdot\|_{1}}(\theta)\) has a simple closed form. This allows for an efficient implementation of Algorithm 1 as described in (Tanabe, Fukuda, and Yamashita 2022)._

### A predictor-corrector method

This section describes our continuation approach to compute the entire Pareto front of problem (MoDNN). _Figure 1_ shows an exemplary illustration of the Pareto front approximated by a finite set of points that are computed by consecutive predictor and corrector steps. After finding an initial point \(\theta^{0}\) on the Pareto front, we proceed along the front. This is done by first performing a predictor step in a suitable direction. As this will take us away from the front (Hilermeier, 2001; Bieker, Gebken, and Peitz, 2022), we need to perform a consecutive corrector step that takes us back to the front.

Predictor stepDepending on the direction we want to proceed in, we perform a predictor step simply by performing a gradient step or proximal point operator step. Note that this is a deviation from the classical concept of continuation methods as described previously, where we compute the tangent space of a manifold. However, due to the nonsmoothness, the Pareto set of our problem does not possess such a manifold structure, which means that we cannot rely on the tangent space. Nevertheless, this is not necessarily an issue, as the standard approach would require Hessian information, which is too expensive in high-dimensional problems, anyway. The approach presented here is significantly more efficient, even though it may come at the cost that the predictor step is sub-optimal. Despite this fact, we found in our numerical experiments that our predictor nevertheless leads to points close enough for the corrector to converge in a small number of iterations. We formulate in equations (3) and (4) the two different steps that can be computed during the predictor step, depending on the direction we want to proceed in. (3) illustrates the gradient step taken for the loss objective function (i.e., "move left" in _Figure 1_) and (4) represents the shrinkage performed on the \(\ell^{1}\) norm objective function (i.e., "move down" in _Figure 1_).

\[\theta^{k+1} =\theta^{k}-\eta\nabla_{\theta}\mathbb{E}_{(x,y)\sim\mathcal{D} }[\mathcal{L}(f(\theta^{k},x),y)] \tag{3}\] \[\theta^{k+1} =\text{prox}_{\eta\|\cdot\|_{1}}(\theta^{k}) \tag{4}\]

Corrector stepFor the corrector step, we simply perform a multiobjective gradient descent step using the multiobjec

Figure 1: Continuation method; The predictor steps are shown in black and blue for Eq. (3) and Eq. (4) respectively.The red is the corrector step.
tive proximal gradient method (MPG). As our goal is to converge to a point on the Pareto front, this step is identical for both predictor directions.

Both the predictor and corrector steps for both directions of movement along the front are illustrated in _Algorithm 2_ where \(f_{1}\) is the loss and \(g_{2}\) is the \(\ell^{1}\) norm function.

## Numerical experiments

In this section we present numerical experiments for our algorithm and the resulting improvements in comparison to the much more widely used weighted sum approach mentioned in _equation_ (2), when considering multiple criteria in DNN training.

In real-world problems, uncertainties and unknowns exist. Our previous approach of using the full gradient on every optimization step taken, oftentimes makes it computationally expensive, infeasible in some situations and unable to handle uncertainties [17]. A stochastic approach is included in our algorithm to take care of these limitations. The stochasticity is applied by computing the gradient on mini-batches of the data, hence randomness is introduced.

### Experimental settings

To evaluate our algorithm, we first perform experiments on the Iris dataset [13]. Although this dataset is by now quite outdated, it allows us to study the deterministic case in detail. We then extend our experiments to the well-known MNIST dataset [1] using the stochastic gradient approach. The Iris dataset contains 150 instances, each of which has 4 features and belongs to one of 3 classes. The MNIST dataset contains 70000 images, each having 784 features and belonging to one of 10 classes. We split the datasets into training and testing sets in a 80-20 ratio.

We use a dense linear neural network architecture with two hidden layers, where each layer has 4 neurons for the Iris dataset and 20 neurons for the MNIST data and the Rectified Linear Unit (ReLU) activation function for both layers. Cross-Entropy is used as the loss function.

We compare the results from our algorithm against the weighted sum algorithm in _equation_ (2), where \(\lambda\) is the weighting parameter. We choose 44 different values for \(\lambda\), chosen equidistantly on the interval \([0,1]\) and solve the resulting problems using the ADAM optimizer [14].

The experiments are carried out on a machine with 2.10 GHz 12th Gen Intel(R) Core(TM) i7-1260P CPU and 32 GB memory, using Python 3.8.8. The source code is available at [https://github.com/aamakor/continuation-method](https://github.com/aamakor/continuation-method).

To illustrate the behaviour of our algorithm, we first study the Iris dataset in a deterministic setting. To obtain a baseline, we have executed _Algorithm 2_ using very small step sizes. Interestingly, the Pareto set and front consist of multiple components, which we were only able to find by repeated application of the continuation method with random initial conditions (multi-start). The resulting solution is shown in blue in _Figure 2_, where three connected components of the Pareto critical points are clearly visible. As this approach is much too expensive for a realistic learning setting (the calculations took close to a day for this simple problem), we compare this to a more realistic setting in terms of step sizes. The result is shown via the red symbols. Motivated by our initial statement on more sustainable network architectures, we have initialized our network with all weights being very close to zero (the black "**"** in _Figure 2_) and then proceed along the front towards a less sparse but more accurate architecture.4

Footnote 4: Due to symmetries, an initialization with all zeros poses problems in terms of which weights to activate, etc., see [1] for details.

As we do not need to compute every neural network parametrization from scratch, but use our predictor-corrector scheme, the calculation of each individual point on the front is much less time-consuming than classical DNN training. Moreover, computing the slope of the front from consecutive points allows for the online detection of relevant regions. Very small or very large values for the slope indicate that small changes in one objective lead to a large improvement in the other one, which is usually not of great interest in applications. Moreover, this can be used as an early stopping indicator to avoid overfitting, as we will see in the next example.

## Results

Motivated by these promising results, we now study the MNIST data set in a stochastic setting (i.e., using mini batches for the loss function). To obtain the initial point on the Pareto front using _Algorithm 1_, we perform 500 iterations. For the subsequent steps, 7 iterations are used for a predictor step and 20 for a corrector step repeated 43 times i.e., \(N_{cont}=44\). _Figure 3 (a)_ and _(b)_ show the Pareto front and accuracy versus \(\ell^{1}\) norm, respectively. In this setting, we have started with a point in the middle of the front in blue and then apply the continuation method twice (once in each direction). As indicated above, we observe overfitting for non-sparse architectures, which indicates that we do not necessarily have to pursue the regularization path until the end, but we can stop once the slope of the Pareto front becomes too steep. This provides an alternative training proce

Figure 2: Pareto front approximation for the Iris dataset using _Algorithm 2_ (red symbols) versus the reference Pareto front in “blue” (computed using the same algorithm with very small step sizes and many different initial conditions) with unscaled \(\ell^{1}\) norm.
dure for DNNs where in contrast to pruning we start sparse and then get less sparse only as long as we don't run into overfitting.

_Figure 4_ shows the comparison of our results to the weighted sum (WS) approach in _equation_ (2). For the WS approach, we compute 200 epochs per point on the front repeated 44 times i.e., \(N_{\lambda}=44\) equidistantly distributed weights \(\lambda\). In total the WS approach needs 5 times as many epochs as the continuation method (_Algorithm 2_) and therefore 5 times as many forward and backward passes. We see also, that the WS approach results in clustering around the unregularized and very sparse solutions, respectively. In contrast, we obtain a similar performance with a much sparser architecture using continuation. This shows the superiority of the continuation method and also highlights the randomness associated with simple regularization approaches when not performing appropriate hyperparameter tuning.

Finally, we note that once an initial point on the Pareto front for predictor-corrector method (CM) is identified, fewer iterations are required when compared with the WS method.

## Conclusion

We have presented an extension of regularization paths from linear models to high-dimensional nonlinear deep learning models. This was achieved by extending well-known continuation methods from multiobjective optimization to non-smooth problems and by introducing more efficient predictor and corrector steps. The resulting algorithm shows a performance that is suitable for high-dimensional learning problems.

Moreover, we have demonstrated that starting with sparse models can help to avoid overfitting and significantly reduce the size of neural network models. Due to the small training effort of consecutive points on the Pareto front, this presents an alternative, structured way to deep neural network training, which pursues the opposite direction than pruning methods do. Starting sparse, we increase the number of weights only as long as we do not obtain a too-steep Pareto front, as this suggests overfitting.

For future work, we will extend our results to even higher-dimensional problems (e.g., CIFAR-10), as well as consider more than two objectives.

## Acknowledgements

This work was supported by the German Federal Ministry of Education and Research (BMBF) funded AI junior research group "Multicriteria Machine Learning".
