Chapter 13Smoothing by optimization

When a function is non-differentiable (or worse, discontinuous), a reasonable approach is to replace it by a differentiable approximation (or at least, by a continuous relaxation). We refer to the process of transforming a non-differentiable function into a differentiable one as “smoothing” the original function. In this chapter, we begin by reviewing a smoothing technique based on infimal convolution. We then review an equivalent dual approach, based on the Legendre-Fenchel transform. We illustrate how to apply these techniques to compute smoothed ReLUs and smoothed max operators, as well as continuous relaxations of step functions and argmax operators.

13.1 Primal approach

We first review how to smooth functions in the original, primal space of the function, using the infimal convolution and more particularly the Moreau envelope, a.k.a. Moreau-Yoshida regularization. In this chapter, we consider functions taking potentially infinite positive values, that is, functions taking values in the half-extended real line $\R\cup \{\infty\}$. For a function $f:\RR^M \rightarrow \RR \cup \{\infty\}$, we define its domain as $$\mathrm{dom}(f) = \{\u\in \RR^M: f(\u) < \infty\}.$$

13.1.1 Infimal convolution

Sometimes abbreviated inf-conv, the infimal convolution between two functions $f$ and $g$ creates a new function $f \square g$. It is defined as follows.

It is easy to check that the three definitions are indeed equivalent, by using the change of variable $\u \coloneqq \muv + \z$, which is a location-scale transform; see Section 12.4.1.

The infimal convolution can be seen as a counterpart of the usual convolution, in which integration has been replaced by minimization (hence its name). Similarly to the classical convolution, it is commutative, meaning that for all $\muv \in \RR^M$, we have $$(f \square g)(\muv) = (g \square f)(\muv).$$

Computing the infimal convolution involves the resolution of a minimization problem, that may or may not enjoy an analytical solution. Some examples are given in Table 13.1.

Existence

The infimal convolution $(f \square g)(\muv)$ exists if the infimum $\inf_{\u \in \RR^M} f(\u) + g(\muv - \u)$ is finite (Bauschke and Combettes, 2017, Proposition 12.6). A sufficient condition to achieve this is that $\u \mapsto f(\u) + g(\muv - \u)$ is convex for all $\muv \in \RR^M$. However, this is not a necessary condition. For example, the infimum can be finite even if $f$ or $g$ are nonconvex, for example if their domain is a compact set.

Infimal convolution with a regularization function

When a function $f$ is non-differentiable, a commonly-used technique is to replace it by its infimal convolution $f \square R$, with some regularization $R$. The most used regularization is the squared $2$-norm, leading to the Moreau envelope, as we now review.

13.1.2 Moreau envelope

When $R(\z) \coloneqq \frac{1}{2} \|\z\|^2_2$, the infimal convolution $f \square R$ gives the so-called Moreau envelope of $f$, which is also known as Moreau-Yoshida regularization of $f$.

Intuitively, the Moreau envelope is the minimal value over $\u \in \RR^M$ of a trade-off between staying close to the input $\muv$ according to the proximity term $\frac{1}{2}\|\muv-\u\|_2^2$ and minimizing $f(\u)$. Provided that the minimizer exists and is unique, we can define the associated proximal operator of $f$ as $$\mathrm{prox}_{f}(\muv) \coloneqq \argmin_{\u \in \RR^M} \frac{1}{2} \|\muv - \u\|_2^2 + f(\u).$$ In other words, we have for $\mathrm{prox}_{f}(\muv)$ well defined, $$\env_f(\muv) = f(\mathrm{prox}_f(\muv)) + \frac{1}{2} \|\muv - \mathrm{prox}_f(\muv)\|_2^2. \tag{13.1}$$

Properties

A crucial property of the Moreau envelope $\env_f$ is that for any convex function $f$, it is always a smooth function, even when $f$ itself is not smooth. By smooth, we formally mean that the resulting function $\env_f$ is differentiable everywhere with Lipschitz-continuous gradients. We say $L$-smooth, if the gradients are $L$-Lipschitz continuous. Such a property can determine the efficiency of optimization algorithms as reviewed in Section 15.4. We recap below useful properties of the Moreau envelope.

Examples

To illustrate smoothing from the Moreau envelope perspective, we show how to smooth the $1$-norm. In this case, we obtain an analytical expression for the Moreau envelope.

We also illustrate in Fig. 13.2 that the Moreau envelope of nonconvex functions can be approximately computed numerically.

13.1.3 Vector-valued functions

The Moreau envelope is defined by $\env_f(\muv) \coloneqq \inf_{\u \in \RR^M} f(\u) + \frac{1}{2} \|\muv - \u\|_2^2$. As such, it is limited to scalar-valued functions $f \colon \RR^M \to \RR$. To extend the Moreau envelope to vector-valued functions $f \colon \RR^M \to \RR^T$, where $f(\u) = (f_1(\u), \dots, f_T(\u))$ and $f_i \colon \RR^M \to \RR$ for $i \in [T]$, we may choose to smooth each $f_j$ separately to define $$\envv_f(\muv) \coloneqq (\env_{f_1}(\muv), \dots, \env_{f_T}(\muv)),$$ where $$\env_{f_i}(\muv) = \inf_{\u_i \in \RR^M} f_i(\u_i) + \frac{1}{2} \|\muv - \u_i\|_2^2.$$ This approach requires solving $T$ separate minimization problems and performs the smoothing of each output coordinate $i \in [T]$ independently. From Proposition 2.9, we then have that the VJP of $\envv_f(\muv)$ with any direction $\d = (d_1, \dots, d_T) \in \RR^T$ is $$\begin{align*} \partial \envv_f(\muv)^*[\d] &= \sum_{i=1}^T \partial \env_{f_i}(\muv)^*[d_i] \\ &= \sum_{i=1}^T d_i \nabla \env_{f_i}(\muv). \end{align*}$$ In the particular case $f(\u) = (f_1(u_1), \dots, f_T(u_T))$, we obtain $$\partial \envv_f(\muv)^*[\d] = \sum_{i=1}^T d_i \nabla \env_{f_i}(\mu_i).$$

An alternative was proposed by Roulet and Harchaoui (2022). For a differentiable function $f \colon \RR^M \to \RR^T$, we recall that the VJP of $f$ with a direction $\d \in \RR^T$ reads $$\partial f(\u)^*[\d] = \nabla \langle f, \d \rangle(\u),$$ where we defined the scalar-valued function $\langle f, \d \rangle(\u) \coloneqq \langle f(\u), \d \rangle$. As a result, if $f$ is non-differentiable, a natural idea is to approximate its VJP $\partial f(\u)^*[\d]$ (had it existed) by the gradient $\nabla \env_{\langle f, \d \rangle}(\muv)$ of the Moreau envelope $$\env_{\langle f, \d \rangle}(\muv) = \inf_{\u \in \RR^M} \langle f(\u), \d \rangle + \frac{1}{2} \|\muv - \u\|^2_2. \tag{13.2}$$ This requires solving a single optimization problem, independently of the number of outputs $T$. Moreover, for $\d = \e_i$, this recovers $\env_{f_i}(\muv)$ as a special case.

This approach allows in principle to perform reverse-mode autodiff (gradient backpropagation) on a neural network whose layers use the Moreau envelope. Indeed, following Proposition 13.1, the approximate VJP of $f$ with a direction $\d$ is given by $$\partial f(\muv)^*[\d] \approx \nabla \env_{\langle f, \d \rangle}(\muv) = \muv - \u^\star,$$ where $\u^\star$ is the solution of the minimization problem in Eq. (13.2). However, we emphasize that this minimization problem could be difficult to solve in general. Indeed, when performing gradient backpropagation, the direction $\d$ is not necessarily non-negative, therefore the function being minimized in Eq. (13.2) could be nonconvex, even if each $f_i$ is convex. Another potential caveat is that the direction $\d$ influences the smoothing strength, while in principle we should be able to smooth a function independently of whether we compute its VJP or not. To see that, for example in the particular case $f(\u) = (f_1(u_1), \dots, f_T(u_T))$, one easily checks that for $\d = (d_1, \dots, d_T)$, we get $$\env_{\langle f, \d \rangle}(\muv) = \sum_{i=1}^T \env_{d_i f_i}(\mu_i).$$ Smoothing vector-valued functions by Moreau envelope (or more generally, by infimal convolution) remains an open area of research. We will see in Chapter 14 that smoothing by convolution more naturally supports vector-valued functions.

13.2 Legendre–Fenchel transforms, convex conjugates

The Legendre-Fenchel transform, a.k.a. convex conjugate, is a way to turn a function $f$ into a new function, denoted $f^*$. We now review it in detail, as it plays a major role for the dual approach to smoothing.

13.2.1 Definition

Consider the class of affine functions of the form $$\u \mapsto \langle \u, \v \rangle - b.$$ These functions are parametrized by their slope $\v \in \RR^M$ and their intercept $-b \in \RR$. Now, suppose we fix $\v$. Given a function $f(\u)$, affine lower bounds of $f(\u)$ are all the functions of $\u$ such that $b$ satisfies for all $\u \in \RR^M$, $$\langle \u, \v \rangle - b \le f(\u) \iff \langle \u, \v \rangle - f(\u) \le b.$$ The tightest lower bound is then defined by $b$ such that $$b \coloneqq \sup_{\u \in \dom(f)} \langle \u, \v \rangle - f(\u),$$ where we recall that the domain of $f$ is defined by $$\dom(f) \coloneqq \{\u \in \RR^M \colon f(\u) < \infty\}.$$ This leads to the definition of Legendre-Fenchel transform, a.k.a. convex conjugate.

We use a $\sup$ rather than a $\max$ to indicate that $f^*(\v)$ is potentially $\infty$. Following the previous discussion, $-f^*(\v)$ is the intercept of the tightest affine lower bound with slope $\v$ of $f(\u)$. This is illustrated in Fig. 13.3.

The Legendre-Fenchel transform is a function transformation, as it produces a new function $f^*$. It can be seen as a dual representation of a function: instead of representing a convex function $f$ by its graph $(\u, f(\u))$ for $\u \in \dom(f)$, we can represent it by the set of tangents with slope $\v$ and intercept $-f^*(\v)$ for $\v \in \dom(f^*)$, as illustrated in Fig. 13.4. As the name “convex conjugate” indicates, it is convex, even if the original function is not.

13.2.2 Closed-form examples

Computing $f^*(\v)$ involves the resolution of a maximization problem, which could be difficult in general without assumption on $f$. In some cases, however, we can compute an analytical expression, as we now illustrate.

See for instance Boyd and Vandenberghe (2004) or Beck (2017) for many more examples.

Constraining the domain

We can incorporate constraints using an indicator function with values in the extended real line $\RR \cup \{\infty\}$, $$\indic_\cC(\u) \coloneqq \begin{cases} 0 &\mbox{ if } \u \in \cC \\ +\infty &\mbox{ otherwise } \end{cases}.$$

13.2.3 Properties

The conjugate enjoys several useful properties, that we now summarize.

13.2.4 Conjugate calculus

While deriving a convex conjugate expression can be difficult in general, in some cases, it is possible to use simple rules to derive conjugates in terms of other conjugates.

13.2.5 Fast Legendre transform

When an analytical expression is not available, we can resort to numerical schemes to approximately compute the transform / conjugate. When $f$ is convex, because $-f$ is concave, the maximization in Definition 13.3 is that of a concave function. Therefore, the conjugate can be computed to arbitrary precision in polynomial time using classical iterative algorithms for constrained optimization such as projected gradient descent (Section 16.3) or conditional gradient a.k.a. Frank-Wolfe (Jaggi, 2013). Without convexity assumption on $f$, $f^*(\v)$ can be approximated by $$f^*(\v) \approx \sup_{\u \in \cU} \langle \u, \v \rangle - f(\u),$$ where $\cU \subseteq \dom(f)$ is a discrete grid of values. We can then compute $f^*(\v)$ for several inputs $\v \in \cV$ using the linear-time Legendre transform algorithm (Lucet, 1997), where $\cV \subseteq \dom(f^*)$ is another discrete grid. The complexity is $O(|\cU| \cdot |\cV|)$, which is linear in the grid sizes. However, the grid sizes are typically $|\cU| = |\cV| = O(N^M)$, for $N$ equally-distributed points in each of the $M$ dimensions. Therefore, this approach is limited to small-dimensional settings, e.g., $M \in \{1,2,3\}$.

13.3 Dual approach

Previously, we presented how to smooth a function by performing its infimal convolution with a primal-space regularization $R$. We now present how to smooth a function by regularizing its Legendre-Fenchel transform (convex conjugate) instead. This dual, equivalent approach, is often mathematically more convenient.

13.3.1 Duality between strong convexity and smoothness

We begin by stating a well-known result that will underpin this whole section: smoothness and strong convexity are dual to each other (Beck, 2017; Hiriart-Urruty and Lemaréchal, 1993; Kakade et al., 2009; Zhou, 2018).

For a review of the notions of smoothness and strong convexity, see Section 15.4. We give two examples of strongly-convex and smooth conjugate pairs in Table 13.2.

13.3.2 Smoothing by dual regularization

The duality between smoothness and strong convexity suggests a generic approach in order to smooth a function $f \colon \RR^M \to \RR$, by going through the dual space.

1. Compute the conjugate $f^*$: $$f^*(\v) \coloneqq \sup_{\u \in \dom(f)} \langle \u, \v \rangle - f(\u).$$

2. Add strongly-convex regularization $\Omega$ to the conjugate: $$f_\Omega^*(\v) \coloneqq f^*(\v) + \Omega(\v). \tag{13.3}$$

3. Go back to the primal space, by computing the conjugate of $f^*_\Omega$: $$f_\Omega(\u) \coloneqq f^{**}_\Omega(\u) = \max_{\v \in \RR^M} \langle \u, \v \rangle - f_\Omega^*(\v).$$

Note that $\u$ and $\v$ belong to different spaces, i.e., $\u \in \dom(f)$ and $\v \in \dom(f^*)$. Following Proposition 13.4, if $\Omega$ is $\mu$-strongly convex, then $f_\Omega(\u)$ is $\frac{1}{\mu}$-smooth. Furthermore, following Proposition 13.2, $f_\Omega(\u)$ is convex, even if $f$ is nonconvex. Therefore, $f_\Omega(\u)$ is a smooth and convex relaxation of $f(\u)$.

Steps $1$ and $3$ are the most challenging, as they both require the derivation of a conjugate. While an analytical solution may not exist in general, in some simple cases, there is, as we now illustrate.

ReLU functions can be smoothed out in a similar way, as we see in more details in Section 13.4.

The dual approach allows us to easily bound the smoothed function in terms of the original function.

13.3.3 Equivalence between primal and dual regularizations

So far, we saw two approaches to obtain a smooth approximation of a function $f$. The first approach is based on the infimal convolution $f \square R$, where $R \colon \dom(f) \to \RR$ denotes primal regularization. The second approach is based on regularizing the Legendre-Fenchel transform (convex conjugate) $f^*$ of $f$ with some dual regularization $\Omega$, to define $f_\Omega = (f^* + \Omega)^*$. It turns out that both approaches are equivalent.

In particular, with $\Omega = \frac{1}{2}\|\cdot\|^2_2 = \Omega^*$, this shows that the Moreau envelope can equivalently be written as $$\env_f = f_\Omega = f_{\Omega^*}.$$ Given the equivalence between the primal and dual approaches, using one approach or the other is mainly a matter of mathematical or algorithmic convenience, depending on the case.

In this book, we focus on applications of smoothing techniques to differentiable programming. For applications to non-smooth optimization, see (Beck and Teboulle, 2012; Nesterov; 2005).

13.3.4 Regularization scaling

Dual approach

If $\Omega$ is $1$-strongly convex, then $f_\Omega$ is a $1$-smooth approximation of the original function $f$. To control the smoothness of the approximation, we regularize with $\varepsilon \Omega$ for $\varepsilon>0$, defining $$f_{\varepsilon \Omega}(\u) \coloneqq (f^* + \varepsilon \Omega)^*(\u) = \sup_{\v \in \dom(f^*)} \langle \u, \v \rangle - f^*(\v) - \varepsilon \Omega(\v).$$ Consequently, the gradient is given by $$\nabla f_{\varepsilon \Omega}(\u) = \argmax_{\v \in \dom(f^*)} \langle \u, \v \rangle - f^*(\v) - \varepsilon \Omega(\v).$$ This formulation works for any closed convex $f$. The approximation error induced by the smoothing can be quantified using Proposition 13.5. Provided that $\cL_\Omega \le \Omega(\v) \le \cU_\Omega$ for all $\v \in \dom(\Omega)$, we have $$f(\u) - \varepsilon \cU_\Omega \le f_{\varepsilon \Omega}(\u) \le f(\u) - \varepsilon \cL_\Omega.$$

Primal approach

Following Definition 13.2, if we use dual regularization $\varepsilon \Omega$, where $\varepsilon > 0$ controls the regularization strength, assuming $f$ is convex, the corresponding primal regularization is $R = \varepsilon \Omega^*(\cdot / \varepsilon)$ and we have $$f_{\varepsilon \Omega} = f \square \varepsilon \Omega^*(\cdot / \varepsilon).$$ In the particular case $\Omega(\v) = \frac{1}{2} \|\v\|^2_2$, we have $$R(\u) = \frac{\varepsilon}{2} \|\u / \varepsilon\|^2_2 = \frac{1}{2 \varepsilon} \|\u\|^2_2 = \frac{1}{\varepsilon} \Omega(\u).$$ We therefore get $$f_{\varepsilon \Omega} = f \square \frac{1}{\varepsilon} \Omega = \frac{1}{\varepsilon} (\varepsilon f \square \Omega) = \frac{1}{\varepsilon} \env_{\varepsilon f}.$$

13.3.5 Generalized entropies

A natural choice of dual regularization $\Omega(\piv)$, when $\piv \in \triangle^M$ is a discrete probability distribution, is a negative entropy function, also known as negentropy. Since negentropies play a major role in smoothed max operators, we discuss them in detail here.

Information content and entropy

An entropy function measures the amount of “surprise” of a random variable or equivalently of a distribution. To define an entropy, we must first define the information content $I(E)$ of an event $E$. The value returned by such a function should be $0$ if the probability of the event is $1$, as there is no surprise. Conversely, information content should attain its maximal value if the probability of the event is $0$, as it is maximally surprising. Furthermore, the more probable an event $E$ is, the less surprising it is. Therefore, when $p(E)$ increases, $I(E)$ should decrease. Overloading the notation, we also write the information content of the outcome $y$ of a random variable $Y$ as the information content of the event $\{Y = y\}$, $$I(y) \coloneqq I(\{Y = y\}).$$

Given an information content function, we can then define the entropy $H(Y)$ of a random variable $Y \in \cY$ as the expected information content, $$H(Y) \coloneqq \EE[I(Y)].$$ Different definitions of information content lead to different definitions of entropy.

Shannon’s entropy

A definition of information content satisfying the criteria above is $$I(E) \coloneqq \log\left(\frac{1}{p(E)}\right ) = - \log p(E).$$ Indeed, $-\log 1 = 0$, $-\log 0 = \infty$ and $-\log$ is a decreasing function over $(0,1]$. Using this information content definition leads to Shannon’s entropy (Shannon, 1948) $$H(Y) = \EE[I(Y)] = -\sum_{y \in \cY} p(y) \log p(y).$$ We can therefore define the Shannon entropy of a discrete probability distribution $\piv \in \triangle^M$ as $$H(\piv) = -\sum_{i=1}^M \pi_i \log \pi_i = -\langle \piv, \log \piv \rangle$$ and use the corresponding negentropy as regularization $$\Omega(\piv) = -H(\piv) = \langle \piv, \log \piv \rangle.$$ The function is strongly convex w.r.t. $\|\cdot\|_1$ over $\triangle^M$. However, it is not strongly convex over $\RR_+^M$, since this is not a bounded set; see for instance (Blondel, 2019, Proposition 2). Since $\Omega$ is added to $f^*$ in Eq. (13.3), we can therefore use this choice of $\Omega$ to smooth out a function $f$ if $\dom(f^*) \subseteq \triangle^M$.