A provable control of sensitivity of neural networks through a direct parameterization of the overall bi-Lipschitzness
Authors: Yuri Kinoshita
Taro Toyoizumi
Venue: NeurIPS 2024 poster (2024)
Keywords: bi-Lipschitzness
theoretical guarantee
tight control
direct parameterization
inductive bias
convex neural network
Legendre-Fenchel transformation
Primary area: learning_theory

Abstract: While neural networks can enjoy an outstanding flexibility and exhibit unprecedented performance, the mechanism behind their behavior is still not well-understood. To tackle this fundamental challenge, researchers have tried to restrict and manipulate some of their properties in order to gain new insights and better control on them. Especially, throughout the past few years, the concept of *bi-Lipschitzness* has been proved as a beneficial inductive bias in many areas. However, due to its complexity, the design and control of bi-Lipschitz architectures are falling behind, and a model that is precisely designed for bi-Lipschitzness realizing a direct and simple control of the constants along with solid theoretical analysis is lacking. In this work, we investigate and propose a novel framework for bi-Lipschitzness that can achieve such a clear and tight control based on convex neural networks and the Legendre-Fenchel duality. Its desirable properties are illustrated with concrete experiments to illustrate its broad range of applications.

A provable control of sensitivity of neural networks through a direct parameterization of the overall bi-Lipschitzness
Yuri Kinoshita, Taro Toyoizumi
Department of Mathematical Informatics,
Graduate School of Information Science and Technology,
The University of Tokyo. Tokyo, Japan.
Laboratory for Neural Computation and Adaptation,
RIKEN Center for Brain Science, Wako, Japan. puri-kinoshitall. tariiiatl1@@@@@@ ecc. wocccu.okyo... ac.
Abstract
While neural networks can enjoy an outstanding flexibility and exhibit unprece- dented performance, the mechanism behind their behavior is still not well- understood. To tackle this fundamental challenge, researchers have tried to restrict and manipulate some of their properties in order gain new insights and bet- ter control on them. Especially, throughout the past few years, the concept of bi-Lipschitzzess has been proved beneficial inductive bias in many areas. How- ever, due to its complexity, the design and control of bi- bi-Lipschitz architectures are falling behind, and model that precisely designed for bi-Lipschitzness realizing a direct and simple control of the constants along with solid theoretical analysis is lacking. In this work, we investigate and propose novel framework for bi-Lipschitzness that can achieve such a clear and tight control based on convex neural networks and the Segendre-Fenchte duality. Its desirable properties are illustrated with concrete experiments to illustrate its broad range of applications.
Introduction
1.1 Background
Nowadays, neural networks have become indispensable tool in the field of machine learning and artificial intelligence. While they can enjoy an outstanding flexibility and exhibit unprecedented performance, the mechanism behind their behavior is still not well-understood. To tackle this fundamental challenge, researchers have tried to restrict and manipulate some of their properties in order gain new insights and better control them. Especially, throughout the past few years, the concept of sensitivity has been proved as a beneficial inductive bias in many areas.
Sensitivity be translated into the concept of bi- bi-Lipschirzness, which combines two different properties, namely, Lipschitzeess and inverse Lipschitzness. The former describes the maximal, and the latter the minimal sensitivity of a function. Bi-LLipschitzness is attracting interests in various fields not only as it proposes promising solution to avoid unexpected and irregular results caused by a too sensitive or too insensitive behavior of the trained function, but also it achieves approximate isometry that preserves geometries of the input dimension (Li et al., 2020). It plays an essential role in generative models such as normalizing flows to guarantee invertibility (Behrmann al., 2019), in uncertainty estimation to prevent the recurrent problem of feature collapse (Liu et al., 2020a: Van Amersfoort al., 2020) and in inverse problems (Kruse al., 2021) to assure the stability of both the forward and inverse function (Behrmann et al., 2021).
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

Unfortunately, the appropriate design and effective control of bi-Lipschitz neural networks are far from simple, which hinders their application to prospective areas. First of all. the estimation of bi-Lipschitz constants is an NP-hard problem (Scaman and Virmaux, 2018). Second, the design of bi-Lipschitz models is even harder as we cannot straightforwardly extend existing Lipschitz architectures that exploit some unique properties of the concept to inverse Lipschitzness and keep their advantages, and vice-versa. Finally, the control of bi-Lipschitz constants requires particular attention as it is a question of manipulating two distinct concepts in preferably independent and simple manner.
Currently existing bi- Di-Liischitz models still present some issues in terms of design or control. On the one hand, some lack theoretical guarantees because they impose soft constraints by adding regularization terms to the loss function (Van Amersfoort et al., 2020), or because their expressive power not well-understood and may be more limited than expected. On the other hand, others restrict bi- Di-Lipschitzness on layer-wise basis (Behrmann et al., 2019; Liu et al., 2020a). Particularly, this means these approaches can only build in the essence simple bi-Lipschitz function employed as layer of more complex neural network. In practice, this kind of parameterization impacts the generalization ability of the model or leads to loose control (Fazlyab et al., 2019). They can also contain many parameters affecting sensitivity to the same extent that controlling all of them is unrealistic and fixing some may affect the expressive power. See Section for further details.
Therefore, taking into account both the potential and complexity of this inductive bias, we first and foremost need model that is precisely designed for Li-iicschitzeess realizing direct and simple control of the constants of the overall function more complex than single layer neural network equipped with solid theoretical analysis. In this work, we investigate an architecture that can achieve such a clear and tight control and apply it to several problem settings to illustrate its effectiveness.
1.2 Contributions
Our contributions can be summarized as follows. First, we construct model bi- bi-Lipschitz by design based on convex neural networks and the Legendre-Fencher duality, as well as a partially bi-Lipschitz variant. This architecture provides a simple, direct and tight control of the Lipschitz and inverse Lipschitz constants through only two parameters, the ideal minimum, equipped with theoretical guar- antees. These characteristic features are illustrated and supported by several experiments including comparison with prior models. Finally, we show the utility of our model in concrete machine learning applications, namely, uncertainty estimation and monotone problem settings and show that it can improve previous methods.
Organization Section 2, we will first explain in more detail existing architectures around bi- Lipschitzness. In Section 3, we will develop our model followed by theoretical analyses. The next Section 4 will be devoted to experiments and applications of our proposed method.
Notation Throughout this paper, the Euclidean norm is denoted || I| for vectors unless stated otherwise. Similarly, for matrices, || corresponds to the matrix norm induced by the Euclidean norm. For a real-valued function F; RT > R. F is defined as (z) /0xz))/0x........ /0xx,, and V F as its transpose. When the function vector F Rm - R." then its Jacobian is defined as V F (0fi/0xz) The subderivative is denoted as dsub.
Preliminaries
In this section, we explain mathematical backgrounds and existing bi-Lipschitz models clarify the motivation of this work.
2.1 Definition
Let us first start by the definition of bi- bi-Lipschitzness.
Definition 2.1 (bi-Lipschitznees). Ler 0 < L1 S Lz. f R R$ is (L1, Li)- if Ta - vll < (x) - s(y)ll S T21|æ- - yll holds for all x, E TR. The right (left) inequality is the (inverse) Lipschitzness with constant L2 (L1). Li and L2 are called together bi-Lispchitz constants.
2

---

Figure 1: Results of fitting 50r with a Lipschitz model (SN (left) or our model (right)), where the Lipschitz constant is constrained by upper bound L. 50 (red line) is where an L-Lipschitz model with perfect tightness and expressive power should achieve a O loss for the first time. SN achieves this only from around 100 while ours at L 50. See Appendix G.4 for further details.
Several interpretations can be attributed to this definition. Bi-Lipschitzness appears when we want to guarantee a high sensitivity and high smoothness for better generalization such as in uncertainty estimation. This can be also regarded as quasi-isometry with distortion of L2/L1,, which means that the structure of the input is relatively inherited in the output well. Moreover, since bi- Lipschitz function by definition invertible, the inverse Lipschitz can interpreted as the Lipschitz constant of the inverse function. Therefore, controlling the bi-Lipschitz constants is beneficial for the stability of both the forward and inverse.
2.2 Desired Features for a Controllable Inductive Bias
In order gain new insights of a problem and predictable behavior of the bi-Lipschitz model, it is primordial that bi-Lipschitzness works as an inductive bias with high controllability and theoretical foundation. For this goal, the following points are of particular interest: (a) bi-Lipschitzness guaranteed by design, (b) theoretical bounds on the bi-Lipschitz constants, (c) tight bounds, (d) independent control of the Lipschitz and inverse Lipschitz constants, (e) a minimal number of hyper- parameters to tune, (f) theoretical guarantee on the expressive power and (g) direct parameterization of the optimized variables the sense of Definition 2.3 of Wang and Manchester (2023)).
2.3 Related Works and Bi-Lipschitz Models
There currently three methods that are mainly employed achieve and control bi-Lipschitz neural networks. See Appendix A for brief overview on Lipschitz and inverse Lipschitz architectures.
Regularization The first one is based on a regularization method (Gulrajani et al., 2017; Van Amers- foort al., 2020). It adds to the loss a term that incites bi bi-Lipschitzness. While the method is simple, the resulting function may not bi-Lipschitz at the end of the training, we have no theoretical bounds on the bi-Lipschitz constants and cannot achieve real control of them.
i-ResNet The invertible residual network (i-ResNet) of Behrmann et al. (2019) is based on the composition of several T.g.. 1 T.g. bbi..-bii-ipshhhtt layers of the type f((x) = x + (x)), where ge is Lipschitz with constant l.g. < 1. The Lipschitz 9i can be constructed by layer-wise methods such as spectral normalization of weights (SN) (Miyato et al., 2018). However, this kind of layer-wise control of (bis)Lipschitzness is known to be loose in general (Fazlyab et al., 2019). See Figure I for an illustration. The i-ResNet has limited expressive power as it cannot represent the function y -I (Zhang al., 2020), and the construction of gi risks to introduce more restriction on the expressive power in practice (Anil et al., 2019). similar model, the Lipschitz monotone network (LMN) of Nolte al. (2023) combines the GroupSort activation function with SN to create It. These approaches requires to adjust the Lipschitzness of each weight during the training process, which may be sub-optimal (i.c., no direct parameterization).
The concurrent work of Wang et al. (2024) also provides a bi-Lipschitz neural network (BiLipNet), which was mainly used for creating Polyak-E.ojasewiicz function useful in surrogate loss learning. It extends the Lipschitz architecture of Wang and Manchester (2023). The BiLipNet constructed as the composition of monotone Lipschitz layers and orthogonal layers. The former layer is realized based on direct parameterization of the IQC theory of Megretski and Rantzer
BiLipNet

---

(1997). While they can achieve a certified control of the bi-Lipschitz constant, theoretical guarantee on its expressive power is lacking. Moreover, BiLipNet composes many layers which may harm the tightness of the bounds.
Therefore, current models fail to satisfy the desired features of Subsection 2.2. Notably, all these points cannot be satisfied as long as we rely on the simple layer-wise control which offers too conservative bounds and number of tunable hyper-parameters proportional to that of the layers in general. In this work, we will establish bi-Lipschitz architecture that satisfies almost all these expected attributes mentioned above. In terms of mathematical foundations, it takes inspiration from the inverse Lipschitz architecture of Kinoshita et al. (2023).
3 Bi-Lipschitz Neural Network
We introduce our novel bi-Lipschitz model, followed by some theoretical analyses and discussion.
3.1 Additional Definitions
We clarify notions of smoothness, convexity and the Lengendre-Fenchel transformation (LFT) which will be a core process in our construction.
Definition 3.1. Let y > 0. F: Rm + IR is y-soooth if F is differentiable and F is y-Lischhitz Definition 3.2. Let μ > 0. F; R - R is p-strongly convex if F(x) Y/2 $1 is convex.
Definition 3.3. Let F >> IR a convex function over I C CR.. Its Legendre- Fenchel transformation
Ft is defined as F* F*)x) == suPyer ((,, F(y)).
3.2 Construction of Bi-Lipschitz Functions
We extend the approach of Kinoshita et al. (2023) and propose a method that creates bi-Lipschitz functions, where it is sufficient to control two parameters to manipulate the overall bi-Lipschitzness without needing to dissect the neural network and control each layer, a feature that existing methods for (bi- (bschhipssh cannot deliver. The first step is to notice that the gradient of a real-valued u-strongly convex function becomes y-inverse Lipschitz, and that of a 7-smooth function 7-Lipschitz by definition. Therefore, aim to compose function which is both strongly convex and smooth. Interestingly, the LFT of a 1/ -troooyy convex function is 3-soohth. This leads to our main theorem.
Theorem 3.4. Let F be closed 1/8-strongly convex function and o 2 O. Then the following function is a-strongly convex and a + ß-smooth: supyer ((y.x) F(y)} + 3 ||x]|22 Thus, its derivative is (a, a + -i--bi--isschhizzz which equals f (xx == argmax,, (w,xr) F(y)} + az.
See Appendix B.2 for the proof. The term a |/1222| has the effect of turning convex function into an a-strongly convex function by definition of strong convexity. This [=()) constructed as described in the above theorem is precisely the bi-Lipschitz model propose in this paper.
3.3 Implementation of the Forward Pass
Based on Theorem 3.4, we can create a bi-Lipschitz function parameterized by neural networks, and thus applicable to machine learning problems. Here, we clarify the implementation of a strongly convex function and the LFT. The overall explicit formulation of our bi- bi-Lispchitz neural network (BLNN) is already shown in Algorithm I for convenience.
Strongly Convex Neural Networks A justtoogly convex neural network can constructed by adding regularization term to the output of the Input Convex Neural Network (ICNN) from Amos et al. (2017). The resulting structure can be written as Fy) = 2 + |1yll?, where 2++1 =
(: 2, w(") y + b,) 0, ....... - 1). (W( Jv are non-negative, W/(") 0. and all functions 9: are convex and non non-decreasing. This architecture is universal approximator of ju-strongly convex functions defined on a compact domain endowed with the sup norm (Chen et al., 2019). Note that any other choice for the convex architecture is possible. This is only an example, also employed in Huang et al. (2021) and Kinoshita et al. (2023). Indeed, we could opt for a convolutional version by using the convolutional ICNN proposed by Amos al. (2017).

---

Algorithm 1: Forward pass of (a, B--BLNN
Input:
imput data æ, an ICNN Ga with parameters and constants a, 3 2
Output: image of æ by (a, +
function
Step O: Construct / --strongly convex function Fo(y) = Go(y) + 2,
1 llyl
Segendre-Fenchhee transformation Step: Find yÿ (z) argmax, { x) - Fo(
Gradient Step: Compute fā (z) Vx (Fø (z) ||=]1²)) yo (z)
return (x)
Algorithms for LFT By Theorem 3.4, we only need compute the optimal point of the LFT, i... vā (z) == argming { (Fo()) x) which is strongly convex optimization. This computation is thus rather fast, and we can use various algorithms with well-known convergence guarantees such as the steepest gradient descent (Shamir and Zhang, 2013; Bansal and Gupta, 2017). Such convex solvers will generate for a fixed x a sequence of points (v+(+))+=o.....7) and its last point will be an estimate of v (I)) However, this kind of discrete finite time approximation could compromise the Li-ppschitzeess of the whole algorithm. The bi-Lipschitz behavior of ye (x) can be explicitly described for the steepest gradient descent as follows. For other algorithms,, see Appendix C.2.1.
Theorem 3.5. Let the symbols defined as in Algorithm 1. Consider the steepest gradient descent of sup, {y 2² Fe Fy)) generating points {z){ at the t-th iterations and ya (I) is the global maximum. If F'g p-strrogly convex and 7-smooth then the point lim yo is (1/r, 1/px)bii- Lipschitz without any bias. Moreover, with Te = 1/ (/((ttttttttt
(o(t++1) 1)) as step size and (2x4 ) yoo as initial point, then for all y3, Jy [19++ Mt++ (I1) h(t]l=;; - x;I where lim h(t) 1/p.
See Appendix C.2.2 for the concrete formulation of htt) and the proof. The above theorem guarantees that the optimization scheme will ultimately provide a ye (e(x) that is bi bi-Lipschitz with the right constants. This was not trivial as bias may have persisted due to the non-zero step size. More importantly, this statement offers a non-asymptotic bound for Lipschitzness. This greatly useful when we theoretically need to precisely assure certain degree of Lipschitzness such as in robustness certification against adversarial attacks (Szegedy al., 2013). The step size v = 1/(u(t + 1)) can be precisely calculated as = 1/3 for an (o,, BL..N. Experimental results suggest similar behavior for inverse- inverse-Lipschitee Lipschitzness. Note the strong convexity of F' is always satisfied in our setting and smoothness can be fulfilled if Since this simple gradient descent has been thoroughly analysed and assured to properly behave in practice as well, we will employ it as the algorithm of LFT in the remainder of this paper.
an ICNN with softplus activation functions log(1 I + e*).
3.4 Expressive Power
As we mentioned earlier, layer- wise approaches realize a bi-Lipschitz neural network by restricting the sensitivity of each layer. Nevertheless, this kind of construction may be sub-optimal as it does not take into account the interaction between layers, limiting more than expected the expressive power, which was often not considered in the original papers. Layer-wise bi-Lipschitz approaches cannot inherit the proofs and properties of the original network, which makes the understanding of their behavior more difficult. As for our model, the (a,, BB-LLN,, we can guarantee the following universality theorem.
Theorem 3.6. For any proper closed atroooyyy convex and a + 3-smooth function on compact domain, there is BLNN without taking the gradient at the end, i.e., sup,, ( (ly. z) - (Go(v) - llvl2/(233)) + al|x||//2 where Gg is an ICNN with ReLU softplus-type activation function, that approximates within e in terms of the sup norm.
Thus, after taking the gradient, any function which is a + f-Liischitz, a-strongly monotone (Definition B.12) and the derivative of real-valued function. See Appendix B.3 for the proofs and further discussion.
create a sequence of functions that converges point- wise to
3.5 Backward Pass of BLNN
The most straightforward way to compute the gradient of a BLNN is to track the whole computation of the forward pass including the optimization of the LFT and back pack-propagate through it. However,
5

---

this engenders crucial bottleneck since back-pree poock-rpgggating over the for-loop of the optimization involves many computations of the Hessian and the memorization of the whole computational graph. Nevertheless, if the activation function of the ICNN is chosen as softplus, the convergence of the L.F quite fast making this strategy scalable to moderately large data.
Interestingly, when Fa == F(::;) 0) is C2, the backward pass can be computed only with the information of the core ICNN, which means we do not need track the whole forward pass involving many recurrent computations:
Theorem 3.7. Suppose a loss function L + L(z), and the output of the BLNN is f*;;;) == V (z) + ox as defined in Algorithm 1. If F is C2 and F* is differentiable, then the gradient (x; V (xr" 0)) can be expressed as V [(x) (x);))0))
0j F(y (r; ();0))
See Appendix B.4 for the proof and formulation with more complex architectures. In practice, coding is simple since we can directly employ the backward algorithm provided in libraries such as Pytorch (Corollary 19) The advantage of this second method is threefold. First, it can reduce the computational and memory cost of the backward process as it does not depend the iteration number and involves fewer matrix manipulations. Second, this leaves the freedom to choose the optimizer calculate the LFT, while the brute force method requires the gradient be tracked during the optimization process, which is a feature that many current solvers do not provide. Third, this (f)
method approximates the gradient by taking the derivative along the curve point yo (a)) close to Vo (z), which is a fair approximation. In contrast, the brute force method approximates the gradient based on the recurrent computation of Va (t) (z) but we do not know whether this is really a
good estimate.
3.6 Comparison with Deep Equilibrium Models
Interestingly, our BLNN can be regarded as bi-Lipschitz version of broader family of models called deep equilibrium models (DEQs) from Bai et al. (2019). A DEQ an implicit model whose output is defined as the fixed point of neural network ha, i.e., input. Similarly, the output of our algorithm can be rerfouuuldeed as finding the solution & of which corresponds to DEQ with he(2, x) = Fe + The general properties of a DEQ concerning the computational complexity (both time and space) are thus also inherited in model. One major difference is that the iteration to find this fixed point 2 is in our case unique and guaranteed to converge, which addresses one of the main concerns of general DEQs. We believe this interpretation is promising for future work to increase the generality of our model as it enables to converge to the larger flow of work around DEQs and apply the various improvements for DEQs researched so far. However, we will not pursue this direction further since our primary goal is to develop a bi-Lipschitz model with direct and simple control.
ho(2, 2), where x is the
= (z)
3.7 Increased Scalability and Expressive Power: Partially BLNN
Despite all the theoretical guarantees mentioned above and the correspondence of the proposed method to DEQs, its main drawbacks persist in computational efficiency and expressive power. However, these weaknesses can be alleviated by imposing the bi-Lipschitz requirement on limited number of variables. This can be realized by using the partially input convex neural network (PICNN) instead of the ICNN in our architecture (Amos et al., 2017). A PICNN is convex with respect to a pre-defined limited set of variables. Based on this architecture, we can proceed similarly to Algorithm 1 and obtain partially bi- bi-Lipschitz neural network (PBLNN). See Appendix F.2 for further details. As a result, all heavy operations such as the LFT and the gradient are applied on this smaller set of variables, which makes the architecture much more scalable higher dimensions. Moreover, a PICNN with is layers can represent any ICNN with & layers and any purely feedforward network with k layers, as shown in Proposition 2 of Amos et al. (2017). This also enhances the expressive power of our model with larger liberty on the precise construction of the architecture.
3.8 Computational Complexity
Regarding the time and space complexity of our model, it is largely equivalent to that of a DEQ (Bai et al., 2019) except that the core function is the derivative of neural network, which adds some computational burden. In Figure 2, we present a comparison of the computational complexity in

---

siua Tul Dhease
10
Figure 2: Comparison of the time (left) and space (right) complexity for a single iteration between traditional feedforward network and various BLNN variants.
floating point number operations (FI.Ps) and memory for a single iteration between traditional feedforward network and various BLNN variants: BLNN (with brute force backpropagation), BLNN with Theorem 3.7, PBLNN with only one variable constrained to be bi-Lipschitz. Comparing our model with Theorem 3.7 to a traditional feed-forward neural network, we can conclude that their difference is only a factor of order 10. Theorem 3.7 greaty contributes reducing both computational and space requirements (improvement of order of 10 and 102, respectively). This explains the scalability of our model to large datasets. Finally, the PBLNN clearly decreases the complexity of the model by limiting the number of variables we impose bi bi-Lipschitzness. See Appendix H.1 for details on the experimental setup.
3.9 Discussion
Our bi-Lipschitz model BLNN (and PBLNN) possesses interesting properties. First, it provides by design a theoretical bound of the bi-Lipschitz constants. Some methods cannot necessarily afford this. Furthermore, it enables direct parameterization of the overall bi-Lipschitzness through the addition of two regularization terms with respective coefficients a and /3 at different outputs. This translation of bi-Lipschitzness into strong convexity plays indispensable role as we do not disturb the formulation of the core function, do not track the transformation of the input down to the output and only have two hyper-paramettrs to monitor, the strict minimum. As a result, we can create bi-Lipschitz functions with known expressive power and more complex than a single layer which others had to compose many times create a deep bi-Lipschitz neural network. We also avoid the risk of losing tightness in the bounds of the bi-Lipschitz constants caused by this kind of layer-wiss wise parameterization. Therefore, our model satisfies the desired points of Subsection 2.2.
While the PBLNN solves several issues of the BLNN, optimization inside the forward pass is still an expensive procedure. Nevertheless, there are plenty of approaches accelerate the LFT and backpropagation through approximations as well. Since the objective function of the LFT is always strongly concave, its convergence speed does not depend on the dimension, and the only bottleneck is the computation of the gradient information of the objective function (i.c., the ICNN). Approximations by zeroth-order methods can improve the procedure, for example. Moreover, amortizing the LFT as Amos (2023) did may also be a promising simplification of the forward pass. On the other hand, backpropagation through the whole forward pass can also be simplified based on Theorem 3.7. For instance, the Hessian matrix can be replaced by its diagonal or upper triangular elements, making the inversion operation easier. Other improvements, including those for the backward pass, can be adapted from research on DEQs (Bai et al., 2019).
Another limitation, which is overcome by the PBLNN, is that the BLNN cannot represent function that is not the gradient of a convex function. However, such a type of function is the core of machine learning problems related optimal transport (Santambrogio, 2015) and some physical systems (Huang et al., 2021). In the following chapter, we will show that our model can be applied to various problems and its performance is not overshadowed by these limitations. In short, it is a necessary price to pay to gain a bi- bi-Liischitz model with features such as known expressive power and high controllability so that it can outperform other methods with looser constraint but lower controllability and fewer guarantees. Indeed, if we define the unit of a bi-Lipschitz model as the basic architecture that requires the minimal number hyperparameters, i.e., one for Lipschitzness and one for inverse-l Line--ipssssssss....mmm most existing models are limited to constructing simple (bi-) Lipshhitz units with low expressive power (e.g. only linear) and they have to compose those units to achieve higher expressive power, which leads to looser bounds. In that sense, our model still has higher expressive
7

---

Table 1: Tightness of Lipschitz bound when fitting f(z) = z (x < 0), x + 1 (x 2 0)
Mean over
five trials. See Table 5 for further results.
Models
== 5
L = 10
= 50
SN
80.1%
68.9 %
32.5%
Orthogonal
75.1%
48.9%
14.6%
AOL
65.9%
46.5%
15,4%
SLL
77.6%
50.0%
15.7%
Sandwich
84.3%
60.2%
16.4%
LMN
100.0%
98.5%
26.0%
BiLipNet
98.0%
58.9%
6.8%
Ours
100.0%
99.4% %
99.9%
Figure 3: Results of fitting f(x) = x (x <), 0), x + (x: 2 0) with SLL (left), Sandwich (middle) and our method (right) with specified Lipschitzness of 50. See Figure 16 for further details.
power and tighter bounds than other (bis)Lipsshitz units as ours can with only one unit produce complex functions and its parameterization is not layer wise. Now, if we can afford to sacrifice the tightness of the bounds further increase expressive power, we can proceed like other methods by stacking multiple BLNNs or combining them with other architectures to suit the characteristics of the problem at hand (see Appendix F).
Further substantial extensions of our method can be found in Appendix F, such as generalization to different input and output dimensions, non-homogeneous functions and to other norms.
4 Experiments
Experimental setups can be found in Appendix H and codes in https://github. com/yuri-kkin -k11/ teaipschttz-euural---a Network.. Further detailed results are summarized in Appendix
4.1 Bi-Lipschitz Control
In this subsection, the goal to empirically verify that (i) our model achieves tight control of bi-Lipschitz constants, and (ii) it provides beneficial behaviors different from other methods. We focus on simple problems as they effectively convey key ideas.
Tightness of the Bounds Here, we verify the quality of the bi-Lipschitz bounds when the model undergoes training. Especially, focus the Lipschitz bound since it the most affected by the approximation of LFT, and there exist many works to compare with it. Inspired by the experiment
of Wang and Manchester (2023), we aim to learn the function f(x) x
(x: 0) that has a discontinuity at x = 0. We take as comparison the LMN (Nolte et al., 2023), the BiLipNet (Wang al., 2024) and the i-ResNet network represented by its substructures: spectral normalization (SN) (Miyato et al., 2018), AOL (Prach and Lampert, 2022), Orthogonal (Trockman and Kolter, 2021), SLL (Araujo et al., 2023), Sandwich (Wang and Manchester, 2023). The Lipschitzness of each model is constrained by a constant L. A model with tight Lipschitz bound should achieve in order to reproduce the behavior of f. The percentage between the
0),
that upper bound around empirical Lipschitz constant and the imposed upper bound L can be found in Table 1. Interestingly, our method is the only one that achieves an almost perfect score for settings, while for others tightness is decreasing. This can be understood as the result of the direct control of bi- bi-L.isschitzeess without relying on the information of the individual layers and the construction which uses only direct parameterizations. See Figure 3 for visualization and Appendix G.2 for more results.

---

Figure 4: Results of fitting the linear function y/ = x' with (from left) AOL, Sandwich, BiLipNet and method with specified Lipschitzness of 1000. See Figures 17 and 18 for further results.
(a) Deep ensemble
(b) DUQ no reg.
(c) DUQ
(d) (0, 4) BLLNN
(e) (1, BLLN
0 (2, 4) BLNN
Figure 5: Uncertainty estimation with the two moons data set with several models. Blue indicates high
uncertainty, and yellow low uncertainty. (d)-(f) are with DUQ+BLNN, where (ax, 3) are clarified.
Flexibility of the Model If we underestimate the Lipschitz constant of the target function the previous experiment, we may also overestimate In this case, we observe that previous methods are influenced by the imposed Lipschitz bound presenting high fluctuations in the learned representation or slower learning speed. This is illustrated when learning the identity function constraint of L 1000 in Figure 4. Our model learn without any problem and the learning speed is hardly affected by L. This may be due to the way we control the Lipschitz constant of the model. In layer- wise models, the Lipschitz constant is adjusted by scaling the input, while in ours we add a regularization term at the end, resulting in different loss landscapes with seemingly better regularization performance for the latter strategy. See Appendix G.3 for more results.
x with a Lipschitz
4.2 Uncertainty Estimation
We now illustrate that our model can be efficiently used in areas where bi-Lipschitzness already plays essential role inductive bias. In the estimation of uncertainty based single neural network, distibbbbts and in-distributions points may overlap in the feature space making them indistinguishabl. Imposing inverse Lipschitzness the underlying neural network of the architecture thus important as it can avoid this phenomenon called feature collapse. Lipschitzness is also beneficial to improve the generalization performance. In one of the state-of-the-art approaches called deterministic uncertainty quantification (DUQ) from (Van Amersfoort et al., 2020), bi- Lipschitzzesss constrained through regularization method. Therefore, a hard restriction by changing the neural network to a BLNN, resulting in a method we call here DUQ+BLNN, See Appendix for precise formulation of this method and theme.
replace the gradient penalty with
Two moons The first experiment we lead is with the two moons dataset. Here, we compare DUQ+BLNN with the deep ensembles method (Taasshminrrayanan et al., 2017), DUQ and DUQ with no bi-Lipschitz regularization (DUQ with no reg.). Results are plotted in Figure 5, The ideal scenario is for the boundary of the yellow area be close possible the training data. As we can observe, by adding a bi-Lipschitz constraint, the area of certainty is decreased around the training data. DUQ+BLNN with (a, 5) (2, achieves tighter area than DUQ. We chose the value of and / based on grid search. We took the apper-parameters with the highest accuracy and a since higher is expected to create a tighter yellow area as shown in Figure 5. Note this strategy only relies on the training data, and such strategy is possible thanks to our direct and simple parameterization of the bi-Lipschitz constants. This clearly shows the advantage of the unique features of our model. See Appendix G.6.1 for further results and discussion on this experiment.
FaShion-MNIST Next, we use real world data of FashionMNIST (Xiao et al., 2017), MNIST (Le- Cun and Cortes, 2010) and NotMNIST (Bulatov, 2011). The task is to learn to classify FashionMNIST, at the end of training we verify whether the uncertainty of the model significantly increases when other types of data such as MNIST or NotMNIST is given to the model. This task to distinguish FNSH from MNIST datasets is known to be a complicated task (Van Amersfoort al., 2020). We compute the AUROC for the detection performance. The result is shown in Table 2. Our

---

Table 2: detection task of FashionMNIST MNIST and FashionMNIST Vs NotMNIST with DUQ and (0))+(0,,)BBLL),) Means over five trials.
Models
Accuracy
BCE
AUROC MNIST
AUROC NotMNIST
DUQ
0.889
0.064
0.862
0.900
DUQ+BLNN
0.899
0.060
0.896
0.964
Table 3: Comparison of our model with state-of-the thettrrrtttttt monotone models in benchmark datasets.
Means over three trials. Results of LMN and SMNN are from the original papers. BF BlogFeedBack, LD = Leanttfauleer, HD HeartDisease, Acc. accuracy. See Table 9 for complete results.
Models COMPAS (Acc.) BF (RMSE)
LD (Acc.)
HD (Acc.)
AutoMPG (MSE)
LMN
69.3 %
0.160
65.4 %
89.6 %
7.58
SMNN
69.3
0.150
65.0
88.0 %
7.44
Ours
69.4 %
0.157
65.5 %
90.2 %
7.13
model achieves not only higher performance for FashionMNIST but also better detection of MNIST and NotMNIST dataset. See Appendix G.6.2 for further results and discussion on this experiment.
4.3 Partially Monotone Settings
Sometimes, it happens that we have preliminary knowledge some type of monotone behaviors of the dataset (Nolte et al., 2023). For example, COMPAS (Angwin et al., 2016). BlogFeedBack LoanDefaulter (Nolte et al., 2023), HeartDisease (Janosi et al., 1988) and AutoMPG (Quinlan, 1993) are benchmark datasets that possess a monotone inductive bias on some variables. See Nolte et al. (2023) for further details on the dataset. As a result, it is more efficient to tune the architecture of the trained model so that it successfully reflects this inductive bias, and various models have been proposed address this challenge. This is another field where we apply our architecture that can create monotone (or inverse Lipschitz) functions. We can also control the Lipschitzness to improve the generalization performance. We compare our PBLNN with two state-of-the-ar methods: LMN (Nolte et al., 2023) and SMNN (Kim and Lee, 2024) in Table 3, and with other models in Table 9 as well. As we can observe, our model is competitive with the others.
Generalization and Scalability Test Furthermore, we used the dataset provided by Nolte et al. (2023), CIFARIO1, which is a slight augmentation the original CIFAR dataset and designed to exhibit a monotone behavior with respect to one variable. We adopted their training scheme, utilizing the entire dataset for training and intentionally overfitting the data in order to assess both the scalability and expressive power of the model. Successfully, we achieved 0 loss and 100% accuracy for this experiment, and the convergence was faster than that of Nolte et al. (2023). This illustrates the high scalability and expressive power of the PBLNN.
Conclusion
We built a model called BLNN based convex neural networks and the LFT that bi-Lipschitz functions more complex than a single layer can be constructed and its bi-Lipschitz constants manipu- lated through the coefficient of two regularization terms added different outputs. That way, BLNN not only achieves such tight, direct and simple control but also provides rigorous straightforward analysis on the expressive power and on approximations involved in practice. We illustrated with experiments its distinctive advantageous features compared to prior models. While the primary focus of this paper was establish a framework suited for solid analyses and for the practical control of bi- Li-Lipschitzness, it is, of course, essential to complement this effort with more varied machine learning applications future work. We still believe this work on its own contributes to the further exploitation of bi-siischiseeess in prospective fields and to the deeper understanding of neural networks by delivering model with unique features for the control of this central inductive bias.
We also tried combine DUQ with BiLipNet but could not find competitive parameter settings.
10

---

Acknowledgments and Disclosure of Funding
Y.K. was partially supported by Grant-in-Aid for JSPS Fellows Grant Number JP24KJ0862 and JST BOOST Japan Grant Number JPMJBS2418. T.T. was supported by RIKEN Center for Brain Science, IST CREST program JPMJCR23N2 and RIKEN TRIP initiative (RIKEN Quantum). We also thank anonymous reviewers for their valuable feedback.
References
A.-M. Acu, I. Rasa, and A. B. Steopoaie. Strongly convex squared norms. Dolomites Research Notes on Approximation, 16(2), 2023.
B. Amos. On amortizing convex conjugates for optimal transport. International Conference on Learning Representations, 2023.
B. Amos, L. Xu, and Z. Kolter. Input convex neural networks. D. Precup and Y. W. Teh, editors, Proceedings of the 34th International Conference Machine Learning, volume 70 of Proceedings of Machine Learning Research, pages 146-155. PMLR, 06-11 Aug 2017.
J. Angwin, J. Larson, S. Mattu, and L. Kirchner. Machine bias: There's software used across the country predict future criminals. And it's biased against blacks. Propublica, 2016.
C. Anil, J. Lucas, and R. Grosse. Sorting out Lipschitz function approximation. In K. Chaudhuri and R. Salakhutdinov, editors, Proceedings of the 36th International Conference Machine Learning, volume 97 of Proceedings of Machine Learning Research, pages 291-301. PMLR, 09-15 Jun 2019.
A. Araujo, A. Havens, B. Delattre, A. Allauzen, and Hu. A unified algebraic perspective on Lipschitz neural networks. International Conference on Learning Representations, 2023.
M. Arjovsky, S. Chintala, and L. Bottou. Wasserstein generative adversarial networks. In D. Precup and Y. W. Tch, editors, Proceedings of the 34th International Conference on Machine Learning, volume 70 of Proceedings of Machine Learning Research, pages 214-223. PMLR, 06-11 Aug 2017.
S. Bai, J. Z. Kolter, and V. Koltun. Deep equilibrium models. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
N. Bansal and Gupta. arXiv: 111.1112.444111 2017.
Potential-function proofs for first-order methods. arXiv preprint
L. Bartletr, D. J. Foster, and M. Telgarsky. Spectrally-normalized margin bounds for neural networks. In I. Guyon, U. V. Luxburg. S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017.
J. Behrmann, W. Grathwohl, R. T. Q. Chen, D. Duvenaud, and J.-H. Jacobsen. Invertible residual networks. In K. Chaudhuri and R. Salakhutdinov, editors, Proceedings of the 36th International Conference Machine Learning, volume 97 of Proceedings of Machine Learning Research, pages 573-582. PMLR, 09-15 Jun 2019.
J. Behrmann, P. Vicol, K.-C Wang, R. Grosse, and J.-H. Jacobsen. Understanding and mitigating exploding inverses in invertible neural networks. In Banerjee and K. Fukumizu, editors, Proceedings of the 24th International Conference on Artificial Intelligence and Statistics, volume 130 of Proceedings of Machine Learning Research, pages 1792-1800. PMLR, 13-15 Apr 2021.
Y. Bulatov. notMNIST dataset. 2011. URL httts:://yarsslatt cotmmist-dataset. html..
Y. Chen, Y. Shi, and B. Zhang. Optimal control via neural networks: convex approach. In International Conference Learning Representations, 2019,
11

---

M. Cisse, P. Bojanowski, E. Grave, Y. Dauphin, and N. Usunier. Parseval networks: Improving robustness adversarial examples. In D. Precup and Y. W. Teh, editors, Proceedings of the 34th International Conference on Machine Learning, volume 70 of Proceedings of Machine Learning Research, pages 854-863. PMLR, 06-11 Aug 2017.
J. Duchi, E. Hazan, and Y. Singer. Adaptive subgradient methods for online learning and stochastic optimization. Journal of Machine Learning Research, 12(61):2121-211 2011.
F. Farnia, J. M. Zhang, and D. Tse. Generalizable adversarial training via spectral normalization. International Conference on Learning Representations, 2019.
M. Fazlyab, A. Robey, H. Hassani, M. Morari, and G. Pappas. Efficient and accurate estimation of Lipschitz constants for deep neural networks. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d' d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
Y. Gal and Z. Ghahramani. Dropout Bayesian approximation: Representing model uncertainty in deep learning. In M. F. Balcan and K. Q. Weinberger, editors, Proceedings of The 33rd International Conference Machine Learning, volume 48 of Proceedings of Machine Learning Research, pages 1050-1059, New York, New York, USA, 20-22 Jun 2016. PMLR.
X. Glorot and Y. Bengio. Understanding the difficulty of training deep feedforward neural networks. In Y. W. Teh and M. Titterington, editors, Proceedings of the 13th International Conference on Artificial Intelligence and Statisties, volume 9 of Proceedings of Machine Learning Research, pages 249-256, Chia Laguna Resort, Sardinia, Italy, 13-15 May 2010. PMLR.
H. Gouk, E. Frank, B. Pfahringer, and M. J. Cree. Regularisation of neural networks by enforeing Lipschitz continuity. Machine Learning, 110:393-466, 2021.
I. Gulrajani, F. Ahmed, M. Arjovsky, V. Dumoulin, and A. C. Courville. Improved training of Wasserstein GANs. In I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017.
G. Hinton, N. Srivastava, and K. Swersky. Neural networks for machine learning lecture 6a overview of mini-batch gradient descent. Cited on, 14(8):2, 2012.
C.-W. Huang, R. T. Chen, C. Tsirigotis, and A. Courville. Convex potential flows: Universal probability distributions with optimal transport and convex optimization. International Conference on Learning Representations, 2021.
L. Huang, L. Liu, F. Zhu, D. Wan, Z. Yuan, B. Li, and L. Shao. Controllable orthogonalization in training DNNs. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), pages 6429-6438, 2020.
H. Huggins, T. Campbell, M. Kasprzak, and T. Broderick. Practical bounds on the error of Bayesian posterior approximations: A nonasymptotic approach. arXiv preprint arXiv:/809.00 2018.
J. Jacobsen, Behrmann, Zemel, and M. Bethge. Excessive invariance causes adversarial vulnerability. International Conference on Learning Representations, 2019.
A. Janosi, W. Steinbrunn, M. Pfisterer, and Detrano. Heart Disease. UCI Machine Learning Repository, 1988. DOI:
M. Jin and Lavaei. Stability-certifeed reinforcement learning: A control-theoretic perspective. IEEE Access, 8:229086-229100 2020.
H. Kim and J.-S. Lee. Scalable monotonic neural networks. In International Conference on Learning Representations, 2024.
D. Kingma and J. Ba. Adam: A method for stochastic optimization. International Conference on Learning Representationss 2015.
12

---

Y. Kinoshita, K. Oono, K. Fukumizu, Y. Yoshida, and S.-i. Maeda. Controlling posterior collapse by an inverse Lipschitz constraint on the decoder network. In Proceedings of the 40th International Conference on Machine Learning, Proceedings of Machine Learning Research. PMLR, 2023.
J. Kruse, L. Ardizzone, C. Rother, and U. Köthe. Benchmarking invertible architectures on inverse problems. arXiv preprint arXiv:2101. 107663 2021.
B. Lakshminarayanan, A. Pritzel, and C. Blundell. Simple and scalable predictive uncertainty estimation using deep ensembles. In I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017.
Y. LeCun and C. Cortes. MNIST handwritten digit database. 2010. URL http://yaannnn. lecun. com/ eeeeeddbbmiisii
Q. Li, S. Haque, C. Anil, J. Lucas, R. B. Grosse, and J.-H. Jacobsen. Preventing gradient attenuation in Lipschitz constrained convolutional networks. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d' d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
S. Z. Li. Z. Zang, and L. Wu. Markov-Lipschitz deep learning. arXiv preprint arXiv:2006.08256. 2020.
J. Liu, Z. Lin, S. Padhy, D. Tran, T. Bedrax Weiss, and B. Lakshminarayanan. Simple and principled uncertainty estimation with deterministic deep learning via distance awareness. In H. Larochelle, M. Ranzato, R. Hadsell, M. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 7498-7512. Curran Associates, Inc., 2020a.
X. Liu, X. Han, N. Zhang, and Q. Liu. Certified monotonic neural networks. In H. Larochelle, M. Ranzato, R. Hadsell, M. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 15427-15438. Curran Associates, Inc., 2020b.
A. Megretski and A. Rantzer. System analysis via integral quadratic constraints. IEEE Transactions on Automatic Control, (2(6))819-88)) 1997. doi: 10. 109/9.587335.
L. Meunier, J. Delattre, A. Araujo, and A. Allauzen. A dynamical system perspective for Lipschitz neural networks. In K. Chaudhuri, S. Jegelka, L. Song, C. Szepesvari, G. Niu, and S. Sabato, editors, Proceedings of the 39th International Conference on Machine Learning, volume 162 of Proceedings of Machine Learning Research, pages 15484-15500. PMLR, 17-23 Jul 2022.
T. Miyato, T. Kataoka, M. Koyama, and Y. Yoshida. Spectral normalization for generative adversarial networks. arXiv preprint arXiv: I 1802 2018
R. M. Neal. Bayesian learning for neural networks, volume 118. Springer Science & Business Media, 2012.
Y. E. Nesterov. A method of solving convex programming problem with convergence rate o(1/k²). In Doklady Akademii Nauk, volume 269, pages 543-547. Russian Academy Sciences, 1983.
N. Nolte, O. Kitouni, and M. Williams. Expressive monotonic neural networks. In International Conference Learning Representations, 2023.
B. Prach and C. H. Lampert. Almost-orthogonal layers for efficient general-purpose Lipschitz networks. In Computer Vision - ECCV 2022: 17th European Conference, Tel Aviv, Israel, October 23-27, 2022, Proceedings, Part XXI, page 350-365, Berlin, Heidelberg, 2022. Peringer-Verlaa, ISBN 978-3-031-19802-1.
doi: 10. 10.1007/978-3-031-19990
R. Quinlan. Auto MPG. http://dto.org11122442///CC
UCI
Machine
Learning Repository, 1993.
DOI:
M. Raissi, P. Perdikaris, and G. Karniadakis. Physics-informed neural networks: deep learning framework for solving forward and inverse problems involving nonlinear partial differential
equations. Journal of Computational Physics, 378:686-707, 2019. ISSN 0021-9991. doi: https://bb..or/000001//...g2222s
13

---

A. Rakhlin, O. Shamir, and K. Sridharan. Making gradient descent optimal for strongly convex stochastic optimization. 2012.
R. Rockafellar. On the maximal monotonicity of subdifferential mappings. Mathematics, 33(1):209-216, 1970.
Pacific Journal of
R. T. Rockafellar. Convex analysis, volume 11. Princeton university press, 1997.
D. Runje and M. Shankaranarayana. Constrained monotonic neural networks. In A. Krause, E. Brunskill, K. Cho, B. Engelhardt, S. Sabato, and J. Scarlett, editors, Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 29338-29353 PMLR, 23-29 Jul 2023.
Santambrogio. Optimal transport for applied mathematicians. Birkäuser, NY, 55(58-63):94, 2015.
K. Scaman and A. Virmaux. Lipschitz regularity of deep neural networks: Analysis and efficient estimation. In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa- Besa-Bianchi, and R. Gamett, editors, Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018.
and Y. Singer. On the equivalence of weak learnability and linear separability: S. New relaxations and efficient boosting algorithms. Machine learning, 10:141-163, 2010.
O. Shamir and T. Zhang. Stochastic gradient descent for non-smooth optimization: Convergence results and optimal averaging schemes. In S. Dasgupta and D. McAllester, editors, Proceedings of the 30th International Conference on Machine Learning. volume 28 of Proceedings of Machine Learning Research, pages 71-79, Atlanta, Georgia, USA, 17-19 Jun 2013. PMLR.
S. Singla and S. Feizi. Skew orthogonal convolutions. In M. Meila and T. Zhang, editors, Proceedings of the 38th International Conference on Machine Learning, volume 139 of Proceedings of Machine Learning Research, pages 7766-9766. PMLR, 18-24 Jul 2021.
C. Szegedy, W. Zaremba, I. Sutskever, Bruna, Erhan, I. Goodfellow. and Fergus. Intriguing properties of neural networks. arXiv preprint 1312.6131 2013.
A. Trockman and J. Kolter. Orthogonalizing convolutional layers with the cayley transform. International Conference on Learning Representations, 2021.
Van Amersfoort, L. Smith, Y. W. Tch, and Gal. Uncertainty estimation using a single deep deter- ministic neural network. In H. D. III and A. Singh, editors, Proceedings of the 37th International Conference on Machine Learning, volume 119 of Proceedings of Machine Learning Research, pages 1690-9700. PMLR, 13-18 Jul 2020.
A. van den Oord, O. Vinyals, and kavukcuoglu. Neural discrete representation learning. In I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett, editors, Advances Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017.
Wang, Y. Chen, R. Chakraborty, and X. Yu. Orthogonal convolutional neural networks. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), June 2020.
R. Wang and I. Manchester. Direct parameterization of Lipschitz-bounded deep networks. In A. Krause, Brunskill, K. Cho, Engelhardt, Sabato, and Scarlett, editors, Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pages 36093-36110, PMLR, 23-29 Jul 2023.
R. Wang, K. D. Dvijotham, and I. Manchester. Monotone, bi- bipsiisc and polyak- Lojyske-ijasiiiizzzz networks. R. Salakhutdinov, Z. Kolter, K. Heller, A. Weller, N. Oliver, J. Scarlett, and E. Berkenkamp. editors, P Proceedings of the 41st International Conference Machine Learning, volume 235 of Proceedings of Machine Learning Research, pages 50379-50399, PMLR, 21-27 Jul 2024.
14

---

Y. Wang, D. Blei, and J. P. Cunningham. Posterior collapse and latent variable non von-identifiability. In M. Ranzato, A. Beygelzimer, Y. Dauphin, P. Liang, and J. W. Vaughan, editors, Advances Neural Information Processing Systems, volume 34, pages 5443-5455. Curran Associates, Inc., 2021.
H. Xiao, K. Rasul, and R. Vollgraf. Fashion-MNIST: a novel image dataset for benchmarking machine learning algorithms. arXiv preprint 1wt../:170807777000 2017.
X. Xu, L. Li, and B. Li. Lot: Layer-wise orthogonal training improving l2 certified robustness. In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, Advances in Neural Information Processing Systems, volume 35, pages 18904-18915. Curran Associates, Inc., 2022.
T. Yu, Li, Y. Cai, and P. Li. Constructing orthogonal convolutions in explicit manner. In International Conference on Learning Representations, 2022.
B. Zhang, D. Jiang, D. He, and L. Wang. Rethinking Lipschitz neural networks and certified robustness: A boolean function perspective. In S. Koyejo, Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, Advances in Neural Information Processing Systems, volume 35, pages 19398-19413. Curran Associates, Inc., 2022.
H. Zhang, X. Gao, J. Unterman, and T. Arodz. Approximation capabilities of neural ODEs and invert- ible residual networks. In H. D. III and A. Singh, editors, Proceedings of the 37th International Conference on Machine Learning, volume 119 of Proceedings of Machine Learning Research, pages 11086-11095.. PMLR, 13-18 Jul 2020.
X. Zhou. On the Fenchel duality between strong convexity and Lipschitz continuous gradient. arXiv preprint 1r03.v18033.00 06557,, 2018.
S. Zlobec. On the flouuoudas convexification of smooth programs. Journal of Global Optimization, 32(3):401-40) 2005.
15

---

A Lipschitz and Inverse Lipschitz models
In this appendix, we review existing architectures that aim to control the Lipschitzness and inverse Lipschitzness separately as existing bi-Lipschitz models borrow a lot ideas from prior works on Lipschitzness and inverse Lipschitzness. Moreover, this appendix may constitute reference for future bi-Lipschitz architectures.
Let us start by the Lipschitz property which has been extensively researched over the past few years.
A.1 Lipschitzness
A.1.1 Definition
Definition A.1 (Lipschitzness), Let L > 0. f -R' is L-Lipschitz (continuous) (f) f(vll & -| - yll holds for all z, R.. L is called the Lipschitz constant and denoted as Lip(
As we can observe, this definition represents the concept of maximal sensitivity with respect to the input as changes in the input upper-bounds changes the output. This is also referred as smoothness since when the function is differentiable the Lipschitz constant upper bounds the Ln-norm of its Jacobian.
convenient property which is often used in Lipschitz constrained neural networks is its closure under composition and addition. Thanks to these properties, the problem of be decomposed into smaller ones for neural networks constructed by the composition and addition of many simple layers.
Lipschitz regulation can
A.1.2 Models
Due to its omnipresence in the field of machine learning, especially in the certification of robustness against adversarial attacks, Lipschitz constrained neural networks have been the focus of many prior works. Interestingly, their realization greatly varies, delivering different solutions for the same ultimate goal.
A.1.3 Gradient Clipping
One of the earliest works trying to control Lipschitzness in the context of machine learning proposes to clip the weights of the function to lie within range of [-c, c| after pack-propagation (Arjovsky al., 2017). Consequently, the set of possible parameters would be limited to a compact set, meaning that the overall function is guaranteed to be Lipschitz with a constant dependent on However, we cannot precisely know this value, and Miyato et al. (2018) pointed out this kind of brutal restriction favors low-rank weights, leading to serious decrease of the expressive power of the neural networks.
Regularization Another simple approach but quite effective is to add regularization term to the loss function so that the behavior of the overall function is encouraged to become Lipschitz in a certain way. On the one hand, gradient of the network as Gulrajani et al. (2017) did in the context of Wasserstein GAN (Arjovsky al., 2017) with the following term:
introduce regularization term which penalizes the whole
E ((IVF - 1)?}.
(1)
On the other hand, we can just focus on the weights. For example, Cisse et al. (2017) proposed a Lipschitz architecture called the Parseval network which uses the following regularization term:
| - -11?.
As result, all the weight matrices are incited to become orthogonal, which means I-Lipschitz.
The main downside of these regularization methods is that they can only impose soft constraints and cannot provide any theoretical guarantees on the Lipschitzness of the neural network. This is crucial in some cases such as adversarial attacks where sensitivity has be exactly controlled. While the implementation is not so complicated as we only need to add term to the loss function, the computational cost heavily depends on the formulation on the whole gradient like that of equation (1) may not be truly efficient since the expectation can only computed for a limited number of points, which means that it unclear whether or not the regularization effect will propagate throughout the whole function.
the regularization term. Moreover, penalty
16

---

Spectral Normalization value since
The Lipschitz constant of a linear layer is equivalent to its largest singular
I x Wyll s [[ ~ <
||ll| - vil-
As result, the problem amounts to normalizing the spectral norm of the weights. This is one of the most famous approaches in this sector. The estimation of the largest singular value is thus crucial for this technique, and power iteration is often used in practice.
Once we evaluated the spectral norm of each layer, we can for example normalize (Miyato et al., 2018 1):
W = =VIIVI
or rescale only weights whose norm is above a certain threshold (Gouk et al., 2021):
W ( V.
This can be applied to most types of layers used so far in practice including convolutional layers which are only essence linear operators (Farnia al., 2019).
During training, since the weights are constantly changing, it is not realistic to exactly compute the spectral norm at each iteration. Nevertheless, it has been observed that weights evolve more slowly than the convergence of the power-iteration, As a consequence, we can only execute a few steps of the power-iteration each step and at the end of training arrive linear layers normalized as desired. This drastically reduces the computational cost.
Spectral normalization enables to simplify the problem of Lipschitz control of the overall function into that of the singular value of the linear units, as the overall Lipschitz constant becomes the product of the spectral norm of each layer. However, it has been repeatedly shown that this approach considerably overestimates the true Lipschitz constant of the whole function (Fazlyab et al., 2019). For example, consider the composition of the following two matrices:
100
)/100)
0
B 1/100
A
1/100
0
100
The spectral norm is determined by the direction of the maximal expansion of the layer but this direction may not be aligned throughout all the layers, especially when there is a nonlinear activation function between each of them. As result, it engenders a gap between the overall maximal sensitivity and the product of that of each layer, and intentional control becomes really difficult.
Similarly, since we are manipulating layers without taking into account the complex interaction between each other, the magnitude of the gradient may vanish (Anil et al., 2019). Indeed, when the weights are restricted to be I-Lipschitz, the norm the gradient can only decrease through the layers during the backward pass, leading potential vanishing gradients (Anil al., 2019). Notably, it was proved that this kind of 1-Lipschitz neural networks with spectral normalization and with ReLU activation functions cannot represent the absolute value [ar| (Anil et al., 2019) which is I-Lipschitz.
In short, layer-wise spectral normalization risks overestimating the general Lipschitz constant, which is problematic if we want to control it, may result in the problem of vanishing gradient and introduce issues into the expressive power of the neural network.
The problem of vanishing gradient can be solved by restricting the matrices to be orthogonal since the eigenvalues an orthogonal matrix always equal 1. The fact that all matrices are isotropic help to stabilize the training (Prach and Lampert, 2022). This leads to large body of work that investigates orthogonal linear layers.
Orthogonalization
There exist several direct realizations of orthogonality. For instance, Xu et al. (2022); Huang et al. (2020) uses the form of
W (VVT)-1/2V.
Singla and Feizi (2021) provides a parameterization with skew matrix as follows:
W exp (V VV)...
The latter method needs to approximate the exponential operator with finite sum. Anil al. (2019) also construct an orthogonal weight from matrix power series:
1+ 2 Qx +...... (-1)² (
-1/2
AR+1
A&
p
17

---

where - Ag. Orthogonal convolutional layers are also studied (Li et al., 2019; Wang al., 2020), and other realizations were considered based on Cayley transform including convolu- tions (Trockman and Kolter, 2021; Yu et al., 2022).
The Altoss-orthogonal layer (AOL) (Prach and Lampert, 2022) reduces the computational cost of creating orthogonal layers by approximating them by the following formulation:
1/2
W = Vddig Σ VVT
Empirically, these weights were nearly orthogonal at the end of the training.
In order to create L-Lipschitz functions, we can simply multiply the output or input by L. Never- theless, while orthogonal layers able to stabilize the training. this regularization destroys the information about the spectrum by setting all the singular values to one. Moreover, the problem of looseness of the Lipschitz bound persists since nonlinear activations are considered apart, and do not know how they interact with the linear layers.
Nonlinear Lipschitz Layer So far, have reviewed methods that deal with linear layers. Recently, some that incorporate nonlinear functions have started to appear.
Meunier et al. (2022) suggest to use residual type layers of the form
f(x) = 2/ oWWx b).
where W (WW is the derivative of s(W + b) and is 1-smooth convex function that s' = a. This type of layer is guaranteed to be I-Lipschitz. This has been generalized by Araujo et al. (2023) who built nonlinear layer of the form
f(z) Hr + Go(Wa + b).
If there exists diagonal matrix A with non-negative scalars that satisfies
1- - H - (a B)AW
T H + 20 AW -a BT G
))= =0.
2A G
where a is slope restricted to Ja, B). then f(x) is L--isschitz. Especially, when 7 1, H = and G -2W7-1 where T' is diagonal with non-negative entries so that T ł W W, the corresponding layer is called SDP-based Lipschitz Layer (SLL). A choice of T is
TN
W W Wlygs//.
where qi 0. See Theorem 3 of their work for further detail.
Fazlyab et al. (2019) proposes LipSDP which estimates the Lipsehitz constant a multi-layer neural network through semi-definite programming (SDP). is based on a SDP derived from integral quadratic constraint. Suppose we reformulate the simple multi-layer neural network as
BX = (AX + +b),
f(x) = CX + bz.
where a is [a, O slope-restrieted, X (o) 2
x
TT
Wo
0
W
0
0
A
WL-1
L-11
0
0
In
0
0
In
B
:
0
Ioz.
18

---

and
C (O..... ................................ 0, WT),
b b= ..... TT
Now, if there is a diagonal matrix with non-negative entries such that
? In
0
( AA B
))C
T
-208A
(3)
)
A B
0
0
3
3)+
-2A
W W7
satisfied, then f is Vi-Lipschitz.
Wang and Manchester (2023) provides a direct parameterization for this SDP. I-layer version is called the Sandwich layer and formulated as follows:
f(z) V2A Bx
where A and B are produced from the Cayley transformation of an arbitrary matrix with correct dimensions, and a diagonal matrix with positive entries.
While these methods provide new interesting possibilities to parameterize Lipschitz functions, they mainly layer-wise, still leading overestimation the overall Lipschitz constant, or use the typical structure of the Lipschitzness, meaning that we cannot extend them to bi-Lispchitznes for our purpose. Notably, the SDP-based approach that comes from Fazlyab al. (2019) cannot handle more general structures such as skip connections.
A.2
Inverse Lipschitzness
A.2.1
Definition
Definition A.2 (inverse Lipschitzness). Let L/ 0. f : IR' is L'-inverse Lipschitz if
Ilfzz) f(v) L'I 1'Il=-y/ yll
holds for all z, E R!. L' is called the inverse Lipschitz constant and denoted as invLip( f))
The above property implies that an inverse Lipschitz function is always injective, which means that has an inverse which is 1/T.-Lipschitz.
The inverse Lipschitzness
mathematical description of the minimal sensitivity of the function. By increasing the inverse Lipschitz constant, we can dilute the function and make it more sensitive to small changes of the input. Unfortunately, the inverse Lipschitzness is not closed under addition as we can understand with the simple example of 0
A.2.2 Models
To the best of our knowledge, the only model that built specially for the control of inverse Lipschitzness is that Kinoshita et al. (2023). The same model was used earlier by Huang et al. (2021) but for the construction normalizing flows compatible with Brenier's theorem in optimal transport.
Kinoshita et al. (2023) observed that the derivative of a-strongly convex function is a-inverse Lipschitz. As a result, they propose a model the type
F (c))
|>:|
VF(x) az.
where F is any convex function parameterized by neural network. This approach has the interesting property that we do not need to know what happens between the layers we only add a term to the output, which is possible thanks to the convexity of That is, layer-wise control not required, resulting in simple and tight control of the inverse Lipschitz constant thanks to the fact that before taking the derivative can also freely reproduce any convex function that not necessarily strongly convex. In our work, we start from this model and extend it to bi- bi-Lipschitzness.
19

---

A.3
Bi-Lipschitzness
A.3.1
Definition
Finally, function which is both Lipschitz and inverse- Liprss-tips is called bi-Lipschitz. Note that by definition the inverse Lipschitz constant cannot exceed the Lipschitz constant.
Definition A.3 (bi-Lipschitznness) Ler 0 < L1 - L2. f R + R is (L1, - Lii-hih-ppschitz if
_ - yl - - f(y)ll - L2||= yl
holds for all z, y E RR. Ln and L2 are called together bi-Lispchitz constants.
A.3.2 Applications
Several interpretations can be attributed to this definition. Here, we provide some examples with clearer explanations of the importance of bi-Lipschitzness which were omitted in the main pa- per. In these applications, both Lipschitzness and inverse Lipschitzness are useful, or sometimes indispensable, and Lipschitzness alone becomes insufficient.
1. Injectivity and Dut-of-Distributtio Detection: Bi-Lipschitz functions are injective thanks to inverse Lipschitzness. This helps distinguish in the feature space, making uncertainty estimation possible (Van Amersfoort et al., 2020). Without inverse Lipschitzness, the problem of feature collapse can occur and compromise the detection of outliers. The inverse Lipschitz constant can control the sensitivity to out-of-distributtin points. Please see Appendix E for further details on this topic.
and in-distribution data
2. Quasi-Isometry and Dimensionality Reduction: Bi-Lipschitz functions can be regarded as a quasi-isometry, meaning that the structure of the input is inherited in the output as well. This property is used for adequate dimensionality reduction or embeddings (Li et al., 2020).
3. Invertibility and Solving Inverse Problems: Bi Bi-Lipschitz functions are invertible, and the inverse Lipschitz constant serves as the Lipschitz constant of the inverse function. In that sense, imposing inverse Lipschitzness helps guarantee good properties of the inverse function, just as Lipschitzness assures them for the forward function (Behrmann et al., 2021).
This aspect is used to create normalizing flows (Behrmann et al., 2019) and solve inverse problems (Kruse al., 2021).
4. Balancing Sensitivity and Robustness: Bi- Bi-Lipschitz functions can avoid overly insensitive behavior with respect to their input by controlling the inverse Lipschitz constant. An overly invariant function is also vulnerable to adversarial attacks, as pointed out by Jacobsen et al. (2019). To the best of our knowledge, the application of bi-Lipschitzness in this direction underexplored, but this concept may provide an effective solution.
A.3.3 Comparison with Prior Models
In this section, we compare our model BLNN with other models in more details.
Convex Potential Flow If set 0, obtain BLNN whose output is V the derivative of the LFT of a 1/S-strongly convex ICNN Fg in terms of Algorithm 1. Since the BLNN is still injective, it has an inverse which only F.. As a consequence, in applications where we need to compute the inverse, this model can provide an interesting solution since we can create invertible functions with known inverse. This model was precisely used by Huang al. (2021) in the context of normalizing flows. It is equivalent to that of Kinoshita et al. (2023), but the motivation is different. Ours can thus be regarded as an extension of their model well.
Residual Network Interestingly, our model ultimately takes the form of ozz + g(z), which can be compared with a residual network. The main differences are that the skip connection is scaled by and that the formulation of g is restricted to derivatives of convex functions. These two features were crucial components for direct parameterization of the bi-Lipschitz constants. Behrmann et al. (2019) composed many layers of the form of residual network to guarantee high expressive power to their invertible residual network. This superposition layers is only sub-optimal as leads to a looseness in the bounds of the bi-Lipschitz constants of the overall function and it is not enough to
20

---

represent some functions such as y =z, In a sense, our work shows that by restricting ourselves to the derivative of convex functions, this kind of heuristics not necessary all and that this condition leads to tighter control.
B Proofs of Statements and Further Discussion
In this appendix, we provide further details and proofs statements in the main paper.
B.1 Additional Definitions
We first remind some definitions.
Definition B.1. Let y > 0. F Rm + y-smooth if F is differentiable and
[VF(() - V S 7|1x yll
holds for all x, y E R. That is, V --iischiitz
In the remainder of this chapter, we will concisely refer smooth function when we do need to specify the smoothness constant. Furthermore, in this work, we will often deal with smooth convex functions. For this type of function, there exist four equivalent characterizations, See Appendix D for the proof.
Theorem B.2. Let 7 > 0 and F Rm + I differentiable convex function convex domain. Then the following are equivalent:
1. F ymsoooth the meaning of Definition B.l:
[V - VF < >|1 yll
2. The following holds for any y domF:
(VF(x) (z) 71 - yll
(x
3. The following holds for any I, E domF:
F(y) F(z)+ VF(x) (y x) + 21ly x112.
4. (co-coercivity) The following holds for any x, y E domF:
1
(VF(=) -F)()FF() v) = (VF(=)-VFFL)
Convexity and strong convexity is defined follows:
Definition B.3. F: m R is convex if
F(tz + (1 - t)y) (1-
holds for all E [0,1] and x, E domF.
Definition B.4. Let 0. F:Rm
R + R is ju-strongly convex if F(z) -/2 || is convex.
We also introduce the Cangendre-Fenchhee transformation which will be a core process in our con- struction.
Definition B.5. Let F: I R convex functtio over T Its Legendre- transformation F' is defined as follows:
F (x)
sup (1,, F(v)). yET
21

---

B.2
Construction of Bi-Lipschitz Functions: Proof of Theorem 3.4
The first step notice that the gradient of a real-valued L-strongly convex function F becomes L-inverse Lipschitz as pointed out by Kinoshita et al. (2023).
Proposition B.6. Let F be an atsoongyl convex differentiable function. Then V is a-inverse Lipschitz.
Proof. Since F is a-strongly convex,
F(y) 2
x)) -x)
2 lly x1?
Similarly,
F(x) 2 F(y) + VF (x y) + 2 +
y12.
By summing both inequalities side by side, obtain
alı y| < (VF(x)) VF(()) (I - y) < VF(2)) VF(v))lly x||,
where we used Cauchy- inequality for the right inequality. As result,
|VF(()) (FF(z))VF(v)) Fz)il 2 a+| yll-
Therefore, an a-strongly convex neural network can be first built, and then its gradient calculated in order to construct a neural network which is guaranteed to be a-inverse Lipschitz. Now, since the derivative of smooth function is Lipschitz by definition, function which is both Lipschitz and inverse Lipschitz. That is, we aim to compose a function which is both strongly convex and smooth.
similarly proceed to construct
Interestingly, smoothness and strong convexity are closely related through the Legendre- Fenchel transformation.
Proposition B.7. If F is closed 1/ B/B-styyygly convex function. Then its Gegendre-Fencher transfor- mation is 8-smooth.
See Zhou (2018) for proof. Importantly, the smoothness of the Legendre-Fenchel transform does not depend on that of as long as it is strongly convex. This means that this statement holds also for strongly convex neural networks with ReLU activation functions, and consequently non-differentiable.
Proposition B.8. The resulting function F of a Fenchel transformation is also convex as long as its domain is convex.
Proof. For all t E [0, 1] and æ1. E domF"
F" (tx1 + (1 t)x²) sup ((y, txj + (1 t)x2) - F(y)}
vET
sup (1 (11,x1) F(w)) + (1- 1) ((1, - F()))}
yel
</sup ((,,m) F(w)) + (1 1) sup (11,02) F(y))
yel
vet
(m1) + (1-t)F (x2).
This leads to the following statement.
Theorem B.9. Let F be a closed 1/ S-stroyyyyy convex function and a 2 0. Then the following function is atronggyy convex and + 3-smooth:
Fdt sup (1,, F(y)} + |0r| yEt
22

---

Fu(r): 1/B-strongly convex network
Legendre Transformation
Fe(K): fismooth neural network
Fgin) Elltt flarmooth ir-strongly neural network
Derivative
(ar, a (D-bbi-iipschitzz neural network
Figure 6: Construction flow of bi-Lipschitz neural network through Legendre-Fencher transformation and Brenier map.
Proof. Since F is a closed 1/ strrongly convex function, we know that F" convex and 3-smooth. Now clearly F + 2 || is strongly convex by definition of strong convexity. Moreover, since
V F (x) - -1=1²) v
F
at
|b
F* (x)
(VF (v) az)ll
$ VF (x) VF ||a: azll
$ (3 + +>)l== - yll.
+ + | 2 is a + 3-smooth.
Corollary B.10. Suppose F is constructed Theorem B.9, then V is (a, a Bi- and f = V F (z) = argmax, {(y F(y)}
} + axc.
See Zhou (2018) for proof of (F) = argmax, (4, F()). See Figure 6 for a summary of the construction flow.
Remark B.11. This construction can only handle functions with the same input and output dimensions. However, this common problem in bi-Lipschitz architectures such as that of Behrmann al. (2019). It can be addressed by using composition of bi-Lipschitz functions which will be explained Appendix F.
B.3 Expressive Power of BLNN
In this subsection, we discuss the expressive power of BLNN. It is crucial for this type of construction elarify this. As we mentioned earlier, layer-wise constraints realize bi-Lipschitz neural network by restricting the sensitivity of each layer. Nevertheless, this kind of construction is only sub-optimal it limits the expressive power of the function more than expected. For instance, neural networks with spectral normalization and ReLU activation functions cannot represent = |æ| (Anil et al., 2019), and invertible residual network cannot express the linear function y = (Zhang al., 2020), which were not mentioned in the original paper (Behrmann et al., 2019). Layer-wise bi-Lipschitz approaches cannot directly inherit the proofs and properties of the original network, which makes the understanding of their behavior more difficult.
As for model, the (a, BLLNN,, we can prove that
1. Before taking the gradient, our model can approximate any a-strongly convex and a + 8- smooth functions on a compact domain endowed with the sup norm (Theorem 3.6).
2. After taking the gradient, any function which a + 3-Lisschit,, a-strongly monotone (see Definition B. 12) and the derivative of real-valuee function.
create sequence of functions that converges point-wise to
The proof is almost straightforward thanks to the fact that our model only makes modification at the output of the core neural network, and we can build the proofs of previous works (Chen al., 2019; Huggins et al., 2018).
Definition B.12. monotone function / R - RC. defined as (f(x) - f(v), x v) 2 Furthermore, f is a-strongly monotone if (f((x) - f(v), y) all - 211
23

---

B.3.1 First Part: Proof of Theorem 3.6
Theorem B.13. For any proper closed a-strongly convex and a + B-smooth function on compact domain, there BLNN without taking the gradient at the end, i.e., sup, (y, (Go(y) llyl2/(233)) + a||x>|/// /2 where Gg is an ICNN with ReLU or softplus-tye typp activation finction,,,, that approximates it within in terms of the sup norm.
Proof. From Proposition of Huggins et al. (2018), we already know that for any convex function on compact domain, there is an ICNN with ReLU softplus-type activation function that approximates it within in terms of the sup norm." 2 The key point in this proof is to show that the Legendre-Fenne Fenchel transformation does not deteriorate the approximation quality and that introducing strongly convexity and smoothness constraints does not cause any unexpected limitations on the expressive power.
Let f be proper closed a-strongly convex and a + f-smooth function. Since f(x)) is a-strongly convex and a B-smooth, f(())
Tul |ar|| becomes a convex 3-soooth function. Indeed, by definition of strong convexity, [ax - 2 ||r]] is convex. Moreover,
(<x()) our (V az)) y) S B y/1²
As a result, by definition of smoothness (Theorem B.2), f(x $ |,||² is 3smooth. Now, since fa)) : (x) ||~|| is proper closed convex function, f f (Theorem 12.2 of Rockafellar (1997)). That is, the Legendre- transforms is an involution. Taking f", we obtain a 1/3- strongly convex function. Similarly, (w) z |lyll is convex. Therefore, there exists an ICNN f so that
sup f(e) (f'x)= 111111
<e,
or
1 111112 - (w) <c 2B
sup f(v) +
(4)
This precision is preserved after Legendre-Fencher transformation. Indeed, if we define f3 as the Legendre-Fencher transformation of f(y) + 전 || then
sup [f - (x)) <c
(5)
which is equivalent
sup (f))f F(x)) <e
This can be derived as follows:
(w,x) f(v) + //1
23 |11133 }
1
fz( = sup
y/
(w,x)) (((x) + (x)) -(oo) + " |lvll }
sup
(w. (x)) + sup ((())) ((()) - + 23 1u12)
1
}
sup
<x)
where we used the definition of (**( and equation (4) for the last inequality. Since we can change the role of f and f and the above inequality holds for all x, we obtain inequality (5). Finally, by adding =|| //////2/ to (x), we obtain function that approximates of within
sup fz()) + 1|2 - f(x) =sup (z) - (f(x) a a||x||//2))
ssup fz] f(x))
Se.
'eere softplus-typ activation function is function that satisfies the following conditions: s 2 ReLU, convex and (<x)) - ReLU(x)] 0 as |x| 30 (Chen et al., 2019).
24

---

In order words, ff is the ICNN that creates a BLNN (before derivation) that approximates f within
As we can conclude, our model with parameters and / (before taking the gradient) can approximate any a-strongly convex and a + B-smooth functions on a compact domain endowed with the sup norm. Importantly, the approximation quality of the BLNN equals that of the model we use to create the core convex function.
B.3.2 Second Part
The second point of this section is proved by the following theorem.
Theorem B.14 (Huggins et al. (2018), Theorem 2 adapted). Suppose G .R + R a proper convex function and almost everywhere differentiable. there is sequence F, B R of BLNN before
laking the derivative so that Fn - G. Then, for almost every x E Rd,, (x) + VG(x).
If a monotone function can be represented the derivative of the real-valued function, then it is the gradient of a convex function. Hence, we can use the above theorem.
Remark B.15. The existence of such sequence F, is assured by Theorem B.13 since we can let Fy approximate G with a uniform error of 1/n on the compact domain [-n,,]]
Therefore, an (ax, B)-BLNN can represent any a + 3-Lipschitz a-strongly monotone (i.e., (ar, or +o-+000 B)- bi-Lipschitz) function that is the derivative of another function almost everywhere.
B.3.3 A Brief Discussion on the Expressive Power BLNN
In this section, we briefly discuss the difference between the whole set of bi bi-Lipschitz functions and the expressive power of the BLNN and the possibility of designing an architecture that can approximate any bi-Lipschitz functions based on our model.
There three relevant classes of functions: (1) bi-Lipschitz functions, (2) monotone bi-Lipschitz functions, and (3) cyclically monotone bi-Lipschitz functions. Our model can represent any function of class (3) since cyclically monotone bi- Lipsschiit functions are equivalent to the class of derivatives of (strongly) convex functions (Rockafellar, 1970). In dimensions larger than 1, these three classes different: [(x, y) (-y, 2r) is in (1) but not in (2), and f(x y) (a: + 2y, y) is in (2) but not in (3). Current bi-L.ipschizz models are supposedly in class (1) or (2), like that of Behrmann et al. (2019) and Nolte et al. (2023). Here, we will only focus on bi-Lipschitz functions with the same input and output dimensions, as most of the bi-Lipschitz models fall into this category. Let us now discuss how class (3) can be used to produce functions of class (1) (and (2)).
First, if we suppose the function can be represented by the gradient of another real-valued function, class (1) and (3) are equivalent. Interestingly, under this condition, we can provide even stronger statement based on Theorem 2.1 of Zlobec (2005). It says that If f : R" + R is continuously differentiable and its derivative is Lipschitz on convex
K with some Lipschitz constant L, then there are a convex function C(x) on K and a > L so that f(x) = C(x) ax x/2. In other words, be written
under the condition that function class (3) is equivalent to all Lipschitz functions up to a quadratic concave term. Therefore, our model can still express quite large family of functions.
the gradient of another real-valued function,
Next, the condition of being the gradient of a function can be interpreted as being rotation-free. ie., V x F O. It is known that any vector field can be decomposed into rotation-fres and divergence free component (Hodge decomposition theorem) under some regularity conditions. In R33 this means that function f in class (1) can be decomposed follows:
-=Vx
where is a vector field, and scalar function. From the other direction, we can take as function of class (3), and A as an arbitrary function with bounded gradient. That way, we can reproduce large variety of functions so that f is bi-Lipsehitz without necessarily being representable as the gradient of another function. The choice of A is still ambiguous that f is effectively
25

---

bi-Lipsehitz with explicit bi bi-Lipschitz parameter control, but we believe this approach is promising for the generalization of our method.
In conclusion, while our model may have theoretically restrictive expressive power, it still leaves a large liherty to express more general bi-Lipschitz functions both in theory and practice.
B.4 Backward Pass of BLNN
In this subsection, we provide the explicit formulation of the gradient of the BLNN with respect to the parameters. The first half considers the simple BLNN and the second a more complete case analysis with a complex architecture where BLNN only one component of it.
B.4.1 Proof of Theorem 3.7
First of all, the gradient of the loss with respect to the parameters can be more concretely formulated follows:
Lemma B.16. Suppose loss function L L:za - L(z), and the output of the BLNN is f" 0) =
V (z) + ozz as defined in Algorithm 1. Then the gradient can be expressed as follows:
V T (x: 0)) = L(2 Vo (x:
where: = (x;
This is only application of the chain rule. Now, interestingly, VoF*(x; 0) can be written function of the core ICNN as the following statements shows.
Theorem B.17. Suppose F and F are both differentiable, then the following holds: (x; 0) -Ôn -0, F(y" (x; 0);
where On is the partial derivative with respect to 0.
Proof. Since, y" (a;; 0) argmax wer (ly x)) Fly: 0)), 3² satisfies the following stationary point condition:
v, (, F,, =
(8)
Now, taking the derivative of ()) ) with respect to leads to:
V&F* (x: (x;
Vo (z: (()))=;)):
x a,F(y"(x:0):0) (V (z)))(=)))
= (V (z; 0)) x 0g F(
(x:
FOyF (z; 0);
where we used equation (8) for the last equality.
This results in the following representation of the gradient (Theorem 3.7).
Theorem B.18. Suppose a loss function L 2 + L(2), and the outpur of the BLNN fxx;;) == V F. (x) + or defined in Algorithm 1. If F and F" are both differentiable, then the gradient can be expressed as follows:
(TT()) - -V. M22)
(9)
Proof. From Lemma B.16 and Theorem B.17, we know that
V [x: 0)) - L(x) VoF'
Continuing the procedure, we obtain:
V (x; = v. L(2)0 VzF()'xx
L[ F,F(y"
= (1
F,)"((((())))
where the last equality, we used the fact that y/ V. F and V Fa are inverse with respect to each other (Zhou, 2018).
26

---

INPUT
GLAN
Figure 7: A generalization of an architecture including our model.
Therefore, compute the gradient without back-propagating through the whole optimization process the Legendre-Ferne Fenchel transformation but only with the information the ICNN Fg. In practice, the following observation facilitates even more coding since backward algorithm provided libraries such as Pytorch.
directly employ the
Corollary B.19. Under the assumptions of Thoerem B.18, the gradient of the loss L is equivalent to the following real-valued function:
- IJ (1;; rreuires_grad_(PP V F(G f.. (eqqurre_grad_[(F) ((F(((());
where O. reqqures grad_ indicates that this is considered constant and not a parameter be differentiated with respect to.
B.4.2 A more complete case analysis
In practice, a model may possess more complex architecture where parameters to be optimized do not only come from the ICNN but also from other components. Consider the schema of Figure 7. Mathematically, this can be formulated as follows. Suppose 0, á and v are parameters, fa; 0) BLNN, (dd;; Φ) is a parameterized function like a neural network, and L(z, vb) is a real- valued function returning the loss. In this problem setting, the overall loss can be expressed as (h((; o) 0), hdd: Φ): v). This means, for a model that incorporates a BLNN, we can classify parameters into three groups:
1. those that define the ICNN used in BLNN,
2. those that detine the transformation whose output is transferred to the input of the BLNN,
3. and those that unrelated the BLNN.
For instance, this formulation includes VAEs with a BLNN for the underlying neural network of the decoder. The derivative the first group of parameters was already discussed above, and that last one is straightforward. As for the second family, we can proceed in similar way Theorem 1.
Theorem B.20. Suppose model with loss included L(f*(h(d Φ); (); The output of the BLNN is f= 0) == V (x) + az defined in Algorithm 1. 1f F and F* are both differentiable, then the gradient of 6 can be expressed as follows:
Vol (h(d: (); 0); h(d: 0): v! ((,
[ ( (x;)00) +ai] vJh(d:
V
where: = = ("(hdd; 0) and h(d;
Proof. We can proceed similarly as the previous theorems. By chain rule,
V. V ()) Liz ( 0): w, v.)
27

---

where 2 (h(d: 6)) and w = h(d; 6). The second term does not require further computation. We will focus the first term.
V (1(( V +
(dd + [V/F)((=;)))))) tal (od(d))).
This provides the desired result.
To summarize, gradients used for the update each type of parameters can be written as follows: Corollary B.21. Let I. TA (h(d: Φ): 0). h(d; vb). Then under assumptions of Theorem B. B.20,
V V. v) (1) (x:)0)) 0 v F(;
V L V_ v) [V (x; al V (Td(d: Φ) v,
V V LZ, b:)
where x = h(d: の) and (h(d:)));)
C Algorithms for Legendre-Fenchel Transformation
In this appendix, we discuss the implementation of the LFT as an optimization problem and properties derived from the choice the optimization algorithm. As a reminder, we suppose that convex solver generates for fixed a sequence of points {y)()))t=o...) based on the objective function (y. .) - Fyy), converging to the true optimizer y" (x) == argmax, (v, F(y)).
C.1 Influence of Approximate Optimization on Bi-Lipschitz Constants: Experiments
We first run experiments with the following common convex solvers: steepest gradient descent (GD), Nesterov's accelerated gradient (AGD) (Nesterov, 1983), Adagrad (Duchi et al., 201 1), RM- Sprop (Hinton et al., 2012), Adam (Kingma and Ba, 2015) and the Newton method. We calculated the bi-Lipschitz constants of the generated curve (a:) for each optimization scheme at each iteration F was chosen as a two-dimensional ICNN with two hidden layers. The activation was set as the softplus function log(1 + exp(x)). The derivative of the softplus function is the sigmoid function, which means the overall convex neural network becomes smooth. At the output of the ICNN, we added regularization term ||>|| 2 (2 10) so that the overall function becomes 10-strongly convex as well. As a result, the Legendre-Fente Fenchel transformation is also smooth and strongly convex with respective constants c' and that estimated beforehand. In other words, the Lipschitz and inverse Lipschitz constants of the function yr (x) with respect to should converge to é' and c, respectively. As for the step size, chose convex objective function, it is known that GD converges with decreasing step size of p(t-Tī- Thus, we chose this for GD but also for Adam and RMSprop since it helped the algorithm to converge. For Adagrad, we set 1. for the Newton method as and for AGD as 1/c, the smoothness constant of the objective function. Results are shown in Figures 8 and See Appendix H for further details on the experimental setup.
that the corresponding convex solver converged. For p-strongly
28

---

Edirated ipestr
Folmmatd ipestr
Yean lponte
The inaere ipereer
16x 10
10"
0
200
1006
2200
2000
2200
300
3000
1002
2000
2100
lleration
teration
(a) GD
21:00
Litimated Jbic1tz Lpactit
Lubmette lnad UIIE
14x
307
teration
(b) AGD
10x
33"
16x 12
7x.10*
653
1060
2000
2500
1008
1992
2060
iteration
teration
(c) Newton
Figure 8: Evolution of bi- bi-Lipschitzness (Lipschitz: left, inverse Lipschitz: right) through the iteration of several optimization algorithms: GD (top row), AGD (middle row) and Newton (bottom row).
29

---

Edirated ipestr
Yean lponte
Eelimated iive UIic
107
0
200
1006
2200
2000
2200
300
3000
1202
2000
2100
lleration
teration
(a) Adagrad
Latimated Jpic1tz Lpactit
invoss
307
3000
1943
tteration
teration
(b) Adam
Jovcntz
10"
30a
30"
2000
2500
1008
1992
2000
iteration
teration
(c) RMSprop
Figure 9: Evolution of bi- bi-Lipschitzness (Lipschitz: left, inverse Lipschitz: right) through the iteration of several optimization algorithms: Adagrad (top row), Adam (middle row) and RMSprop (bottom row).
As can observe, the simple gradient descent, AGD, Adagrad and the Newton method perform well in the perspective of conservation of bi-Lipschitzness while RMSprop and Adam are not able to provide such feature, especially for the inverse Lipschitz constant. This incites to GD for simplicity or Adagrad as convex solver of the Legendre-Foud Fenchel transformation.
C.2 Influence of Approximate Optimization on Bi-Lipschitz Constants: Theoretical Analysis
Now that we have an idea of the behavior well-known convex solvers, we proceed theoretical analysis of the approximation quality of optimization schemes concerning the bi-Lipschitz constants predefined in our model. The aim is to clarify when bi-L.ipschszzness preserved and evaluate the effective bi-Lipschitz constants through the iterations.
30

---

C.2.1
Proof for General Algorithms
Here, we focus on deriving general non-asymptotic bounds, i.e., the evaluation and evolution of the bi-Lipschitzness of (x) for finite period of timer for general optimization algorithms. This is crude bound.
Theorem C.1. Let the symbols defined as in Algorithm 1. Consider an optimization scheme of sup, (y. - = generating points {ve(r)) that achieves an error of (lue(()--"" (x))] & c(t, x) after I iterations, where y (x) is the global maximum. Then, for all wi, 2i such that ||2 xz 2
for ó > 0,
11 (z.) 1/0 (x,)) 1 - 1
a (e(((1,)))) + et, 2;)1/6
sa+ B (e(l, + 2,1/,,
where (x)= = (x) + oz is the finite time approximation of fö (z).
Proof. Concerning the upper bound,
11 (xt) (x;) 1 (z,))] (2,) + at; (Ve +
Sa [x - xyره + (x))(x;)) ,r(x1))
<a |z æyll + (2,) 5² (z+)) + (zy) y" (22))
+ Illy" (xrz) - (x,)]]
<a æzıl + <(t, + v E(t,)
<a - æj]l + <(1,,,)) + (VF](,)) (F) + xl,
(
Sa [2; æyı + =((,,;)) + B||2; 2yyl (z)
=a [2 - 1, + B||2; - x,ا + <(t, v.) + +e(t,zg),
where Fÿ is the Legendre- Fenchel transformation of Fg. As for the lower bound, we begin by evaluating (we Ya (1; ti xz).
(x;) - Me (xz)), mタ} = (vx( - y (zi), H2- mÿ) + (y" (x.) (
(x,)) (x(x)), - xz))
2 (xxe - V² (xx)))]] 1ij;;;--)] 2y|l
(VF (I:)) - (1,(),
(1;)) - (x - 2z
2 =((,,x:)]<:; - x-]l+ (VF (1,)), xz))
<(t,,;)]lı; 10:--00 =,,1
2 c(t, - æ (t,x,,)); -2yl1-
where for the last inequality we used that for any convex function (VF (x,) xz)) 0.
Now since f! (x) = (z)x)) + ax
(f( (I,) - (1() (1())) - kz) ((y((x)) + az, - (u + + - 23}
=(v)(xr)) ye (e(r))) x4 - 2z)) + a|| - 11
2a -2,112 - ((t,,;))12; 1yا - -e(t, xz, |z, -xyll- By Stuchy-Schwarz inequality, (F/²) (I:) f/10 (zy), 2,) s 11 (2.) 1/²0 (2,))l|ı, x,
Thus, we finally obtain
11 -
2 a||=; -0,--,11--44 - - (tt, m.y). (I,,) (
In order to derive the inequalities of the statement, it suffices to divide all sides by || - 1351 and use ||2, æ 2 5 on the two obtained bounds.
L
Since (l., xx) - 0 as + oo for any optimization scheme that converges, we can recover an arbitrarily close approximation of a and + 8, the inverse Lipschitz and Lipschitz constants, with enough number of iterations. This is a general bound that can be applied to any algorithm with known behavior. Below, we provide two concrete evaluations applied to the gradient descent.
31

---

Corollary C.2. Suppose Fø is p-strongly convex and z - Fo S G(z). If we employ gradient descent with tre 1/(p(t + 1)) as a step size and a fixed initial point yo = 0, 4t, (G(z)/(uVB) Therefore.
a ---22 -[(((x)) + ((T))) S 115 (0) (2:) S/0 (1x)))))
Sa+ 8 2((()) (((5))/6wV)
12; - -
Proof. This follows from Lemma 1 of Rakhlin et al. (2012).
similar bound can be provided when x - Folvell S (z)) is satisfied almost surely. See Lemma 2 of Rakhlin et al. (2012).
Corollary C.3. Suppose F6 u-sttrongly convex and y-soooth. If we employ gradient descent with
7) 1/7 as step size and fixed initial point yo, we can choose c(t, x) = (1 y" (2)
a²
t/2
llvo(z))
72
Proof. Since F is 7-smooth, h(y) = (u. x) - Fly) is smsoohhh and ju-strongly concave. As result, || - y" 1| = |1/l + 7/ v)((())) - y²11²
= |Iy/e y² 112 - 2n/ (Vh(ve) ) V11 -
=I +II/= n 112 - 2n7 Vh(yt) Vh(y"),,, 10).) + + - Vh(y"))
S SIIye n 112 20 [h(ve)) - Vh('') -WW)"\)
1
=Ilye - y 8²1
[Vh(v)) - Vhyy 1
42
SlIue
1ive 3² 11?
722
1121
11ve - 1²11
?2
++1
2²
11:/0 - y 11
where we used that Vh(y") = 0 for the third equality, the co-coercivity of smooth convex functions for the first inequality, and inverse Lipschitzness for the second inequality.
This is a crude estimation of the bi-Lipschitz constants. As long as we assume that all values can only take discrete ones with intervals bigger than 6, then this theorem provides meaningful evaluation. For example, during training this does not pose a problem since the number of training data is usually finite. Technically speaking, will not have any problem in practice either as the computer also deals with discrete values.
The limitation of the above theorem appears once we authorize infinite precision we could take two values arbitrarily close to each other, i.e., 8 - 0, leading the lower and upper bound to explode. This makes the model prone to adversarial attacks. One way to remedy this is to evaluate the Legendre- Fenchel transformation for some discrete values and interpolate with some bi-Lipschitz functions (e.g., linear) so that the whole function satisfies the expected bounds.
Ideally, we would like to derive bound of the type
h (1 xzllst) < 115 (t)
(x:) (0) (15))l
< h2(11 - :))
(10)
I1 x3
where h and h2 converge as ||x;; - 1 + 0 and converge to 0 as - 30. We could not achieve this in Theorem C. which is predictable as we did not take into account the similarity of the optimization paths of z; and x, and derived the bounds by isolating each point. The next step is thus to take into account this closeness. Nevertheless, this kind of analysis is more involved and we do not know if it applicable to all optimization schemes. This could explain the existence of some algorithms where the bi- bi-Lipschitzness do not evolve desired as we have pointed out in Figure 9.
32

---

C.2.2 Proof of Theorem 3.5
The following theorem shows the first half of the statement of Theorem 3.5. That is, show that in the limit there is no bias for GD terms of bi-Lipschitz constants.
Theorem C.4. Let the symbols defined as in Algorithm 1. Consider an optimization scheme of sup, Uy x) F(y) generating points {y(())} the 1-th iterations. If F is p-strongly convex, 7-smooth and twice differentiable and we employ gradient descent with a step size Tp so that the discrete optimization converges to the global maximum y (2), then for all
(o) P
Proof. The recurrence relation of GD is as follows:
(x) = (x)) + fs - V, (()))))
Since F is twice differentiable, we can take the Jacobian of both sides:
V {)) V (e( + ne {I V.
By replacing ye+1 (x) and (x) by Booo (z:) as the sequence is converging by hypothesis, we obtain
V Voc (x) =V =V y
= V (x)) V (z) m
F2F()) since Therefore, since is invertible as v (iy),
V Uox (x) = (+)))
By using μ = V2FF()) I y, we arrive at the desired result.
Remark C.5. Note that (x) and y² (x) are different in the sense that the former is defined by the limit (t + oo) of the recurrence relation of the optimizer, while the latter is defined as the global maximum of sup, (yy, .²) - Fly)} which also satisfies V Fly" (x)) = x. short, they differ the formulation and explicit dependence
Since F' is p-strongly convex and 7-moohhh its Segendre-Fenchhe transform is 1 s/n-strongly convex and 1 //u-smooth, leading to (1/7, 1 1/pj)bii-Lisschizz function by taking the derivative of the latter. This behavior is inherited by the optimization scheme in the limit of f. - oo as the above theorem shows. The point is that there is no bias, which may have occurred the influence of the step size persisted in the final result. A similar result can be obtained for AGD. For the others, the analysis is more complicated as the gradient accumulated throughout the iterations. Since our main interest non-asympiotic behavior, we will not further develop this.
We will now prove the second half of the theorem and show that equation (10) can be actually established for some parameter settings at least for Lipschitzness.
Theorem C.6. Let the symbols defined as in Algorithm 1. Consider an optimization scheme of sup, Ky x) Fly)} generating points (y+( the t-th iterations and (x:) is the global maximum. If F is j-sttongly convex and y-smooth and we employ gradient descent with ne 1/1(t+1) as a step size and yo(xz) = Yo as initial point, then for all His Zz
(1b)+1(z)) - (+2((1))1 S V1 2mp4 + (2?2)1 - ye yx(2))) + melp - xzy--
As a result,
11:4+1 - (1+1(1;))1 h(tt]ll= Zzy1.
where
lim ((t))
+1000 J
33

---

Proof. Let us first prove the recurrence relation of the statement. Throughout the proof, we use the following notations:
Ay == ye - ye yr())))
Ar 2, I,j.
DD; == VF(()(I;)) (F(()(2)))
Since
|| || A 17c AD; n. S IA 77 AD + Te |
(11)
it suffices to evaluate
|Ay; 1/: Now,
| - ",,,, 112 - =1 12 21 (Ay. AD,) + 7/ |
= Ay 1 - 2 1| + n/ A D, 1
1 A 12 - -2MM A 1| + nh > |1 1
(1 2n,/1 + 7² /
where for the first inequality we used that F is u-strongly convex, i.e., (VF(y)- - VF()), 1/ ally 1/1 and for the second inequality we used that / is y-smohth.
As μ S 7 always holds,
1 1 2nd + np?² =1 2014 + 7/ + (1? (2) = =(1 m0 + 2 - 12) 20. 0.
As result
| - 7,D S /1 22 + 1
By substituting this to equation (11), we obtain
| $ / 2n1 + n121211yee11 + 17. |Ax],
which was the desired inequality.
Next, we assume ||rl|| > 0 and divide both sides by || Denoting as == /ll///tal//AAzll we arrive at
a S / 20 + 3 a + Teft.
We want to show that lim4-+00 0g 1.
By the above-mentioned argument when / 2 1, we can further develop this as follows:
0++1 S
24/ + 2 + Te
(1 1//4122 + nス (7, 222 ) + nel
1 n/ (72 112)
(1 1111)) +
Or + Bel
2 77/44
m/ (7² 122)
(1 m(11 +
a2t + 7/W..
1 1/
Since n/ 1/(L+ + 1), we get
04+1 = { t
22/12 -1
2 t(t + 1) ) a4 + 1
t +++1*
++1 1
By multiplying both sides by f + 1, defining A, == toy and k: := 72 /22p222 1(2 0), we can conclude that
k 2
A411 S
1
A + 1.
(12)
Let us prove that from equation 12, we have
K 1
A.
exp
2
(13)
keist-1
34

---

for all t = 1,2,.... We will proceed by mathematical induction.
When / 1, A1 on. Since | || 7 Do A and the initial point was constant
independent of r;, ||y || = 7/ | A = - Thus, A = 01 = 1. This means, equation (13) is satisfied.
Now suppose equation (13) holds for 1 m. Then by equation 12,
K 1
Am+1 5
1 +
Am + 1
2 m
K 1
K 1
<
1 +
2 m2
+
2
exp
2
+1
k=
Icism-1
K 2
K
K
R:
Sm +
+
exp
+1
m²
2
k
2
k= Rcism-1
K 1
K 1
A
K 2.
Sm + 1+
2 m
+ exp
2 m2
2
exp
k=
k<ism-1
Sm + 1+
K
+
II R<ism
exp
K: 2
m
k=1
K
m
1
5m+1+
2
exp
2
k<ism
Therefore, by mathematical induction equation (13) holds for all / 1, .....
Finally, since
n²
1,2,
1
6
t=1
and
Ollgg j),
can conclude that
A
0g
1+ +/4 Et=l F I exp(#1)
t+ $ E exp E Exeist (#1))
t 5
- E b exp 5
+ 0(log) 1)
t +1 (+- 00).
This theorem proves that the generated curve (xx) at fixed time is Lipschitz with a constant that converges to the true value with a convergence speed of O(log l/t). This theorem is interesting it not only guarantees that ye (I)) is bi-L.ipshitz through the whole iteration and converges to the desired values but also is equipped with non-asymptotic bound offering concrete convergence speed and
35

---

value at each iteration. The advantage of this proof is that setting the step size as 1) = 1/(((t 1)) realistic, since our setting we know the strong convexity constant which corresponds 1/3. In practice, the convergence of the bi-L Di-Lipschitz constants seems to be faster, which means that there may exist tighter bounds but we let this investigation for future work since our goal was to just assure a convergence with good properties.
As for the lower bound, we could not prove a similar property for GD. Nevertheless, we can guarantee that it will not diverge to -0 when I - xyYl + 0 as
(va( (x)), I;; mn} S - (xx) ve (z,) -xy1l
and (1:) (x;))) was just proved to be upper-bounded by constant. We let further theoretical investigation for future work. Since GD shows and assures rather good performance, we will mainly use this algorithm this paper.
C.2.3 Analogue of Theorem 3.5 with optimal step size
An analogue of Theorem 3.5 can be shown with a different step size, i.e., 1) = 1/7, where ? is the smoothness constant of the objective function. This step size is optimal as it assures the fastest convergence. However, in practice, it is actually useless as estimating the smoothness constant y for the step size is computationally demanding and during training this constant constantly changing. Only after training. it suffices to estimate 7 once, and all forward passes can be executed with this step size and convergence speed. The theorem is as follows:
Theorem C.7. Let the symbols defined as in Algorithm 1. Consider an optimization scheme of sup, fyy x) F(y)} generating points {ye(r)} at the t-th iterations and y" (x) the global maximum. If p-sttrongy convex and y-smooth and we employ gradient descent with d) = 1/7 as step size and yo (2xa) = yo as initial point, then for all xı, 2;3
[1:++(x:)) - 1 - (1 (1- [1ba(a.) - ye (x,)))] + Th 12 - xz|l-
As a result,
11:)+ 1141 (x,))1 S h(t)]l=;; -
where
lim h(t) 1
1-000
Proof.
Let us first prove the following inequality:
- VFW)) 1/)
11y 2112 1+
- - (14)
| y
If 7/, then by strong convexity of F'
112
(VF()-
v²) 2 plly 3/ 112 M4, 11y y /12
1y 3/ 11.
11 + 7
/ +
However, since μ = Y, we have م fly 3/1 = [[ - FF Thus, the above incquality leads
to
(VF()) VF()).v v/) 2 pl1p v
17 1
bly 1111 + リ -FF
When > Jr, since F is ju-stronlyy convex and smoohhhh F(y) ~ーj ||,| 2 is convex and 7 م smooth. This is straightforward from Definition 2 of Theorem B.2. Now, by applying Definition 4 of Theorem B.2 to F(y) - _ ||
we obtain
1
FV) - py (VF())) - pn)) w/>
- py - -
(VF()) n3 11.
7
By developing both sides, we arrive
(VF()) - VF(f). w/) - plly - 2212
5/) + --1) 7- (IVV(() -VFG))) - VF()).
36

---

Rearranging the terms, we have
Y + (VF())
VF()), 1/)
Nn7 |1
7 7- 1
7 p
7- 1
which leads to
M7 1y
y 11
(VF())
VF(
p -
1+7
Let us now prove the recurrence relation of the statement. Throughout the proof, we use the following notations:
Ay == yı - 30(x)),
Ar 2; Ijj:
AD; == VF((((z,)) VF(U)(Z))).
Since
|| || Ayy MAD nA $ I4 n. + me /
(15)
it suffices to evaluate
|Ay 71 D,l. Now,
|| Ayr - 1/1 A D,
= Ayell' 11? 2n (AUr. AD, + th ||
< Ayel? - -2m ( |A +
MY
IADD1I)) + IADA?
11 14 7
27047
|Awelll + "1 (Th xt+) ||AD.IP.
2
where we used inequality (14) for the inequality. Since we set 1/+ = 1/7,
2p
/1
2
| 2D
1
|
+
-
AADi
#1 + 7
? "7
# ut?
2p
/1
<
1
||Ayell
13110aa11
1
+
?
7
1 ++7
/
1
where for the second inequality we used that aty
2 s 0 (p S 7) and [|V - VF(33] 2
ally- y'll (inverse Lipschitzness). Substituting this to equation (15), we conclude that
A $ (1 (1-4) || + n |
as p/7 1. Since | = ||ااا To Do + 7o 7orr| and the initial point was a constant independent of x,,
| nol ||r| = + Axll, and
- #4 7) ||
A <
Finally,
lim Avel
1
Σ
x
1
11-) ||r|||
1/7 1- (1
|
(1/1))
37

---

Number of needed GD
Number of tterations naaded for GD
home valipor.
Training
Traningg lteration
Figure 10: Evolution of the number of iterations needed for the gradient descent of the Legendre- Fenchel transformation to satisfy the stopping condition during training, with two regimes for the choice of initial points: a fixed or adaptive initial points during training on the FashionMNIST dataset from Appendix G.6.2.
C.3 On the Choice of the Initial Point
In practice, as the iteration proceeds and evolves, the objective function changes. That is, the value of argmax, {y, x - Fo(y)} } for the same training data z differs between each training epoch. However, is often observed that the evolution of the weights relatively slow. This means that the objective function does not change too much compared to the previous time. As result, in order to accelerate the computation, we can use the final estimation of argmax, {y,, - Fo(y)) of the previous epoch the initial point of the next epoch. This may be really simple trick, and there could exist more involved methods, but Figure 10 illustrates that this approach already effective. See Appendix H for further details on the experimental setup. We also leave this research of better approaches for future work.
D Equivalent Definitions of Smoothness under Convexity
In this Appendix, we prove Theorem B.2 by showing general theorem extended to dual norms.
D.1 Theorem and Proof
Dual norms are defined follows:
Definition D.1. The dual norm of norm
is defined
|lyl- = sup {xx 9)
Let us now prove the following theorem. When the norm in question is L2, dual is also the same and arrive at the initial statement.
Theorem D.2. Let 7 > O, || and ||. pair of dual norms, and F R > Re a differentiable convex function on convex domain. Then the following are equivalent:
DI. The following holds for any x, E domF:
[VF((z) VF())]+ 3 <>11= vll
D2. The following holds for any y, E domF:
(VF(x) v y) S </10 - y112.
D3. The following holds for any x, E domF:
F F(z) VF(x) (y 2) + 2f1ly x11
D4. (co-coercivity) The following holds for any x. y E domF:
_ -
(VF(x) VF())) (x
38

---

Proof. D1.= D2.
When x y. the inequality trivially holds. Schwarz inequality, since
When x. # y, this is straightforward by the Cauchy-
VF(w). 2 -ー <| اl sup ((, -VF()) -|=|-= 13
| -z-
-|
<>||> y||
where we used D1. in the last inequality. D2.= D3.
Consider the function G(t) := F(s Mly - x)), which is well-defined as the domain of convex. Since F is differentiable and F(u) = G(1).
d G(t)dt [
dt
-G(0)+ + + - x)). y - d)
=F(x) + (VF()),) v - -x) + (VF(x + My- VF(x), x)
f" ty/1= yll2tt
<F(x) + (VF()), - x) +
-F(x) + (V F (x)) v) + 2 1 y|13
where we used D2. for the inequality. D3.= D4.
Consider the following function:
which is y-smooth function. Moreover, since F is convex, is the minimizer of Fz. Consequently, the following PL inequality generalized to dual norms holds:
Fa inf (w)
< inf { (2(()) + (VF(())) 2)+ w| 21²²
inf inf I &
+ 2 [1
inf inf [2(2))+ ((t))+(VF() (z). 2 2
|oll=1
inf Fx(()) 1 ((VE.()),))
|+||=1 27
=Fa(2) 24
where we used D3. for the inequality. Since this holds for all z, by substituting y. we obtain
2 (VF()) - v - IIF < F.(v) (x) F()) - F(x) - (VF()),))))))) z).
2-
Since we can interchange x and y, we also have
- V
[(F())
2
39

---

Summing side by side these two inequalities, we arrive at
VF v=()) & (VF(x) - VF(U). x v)-
D4.= D1.
When x = y. the inequality trivially holds. When x + y. by Snwpy--Sshaarzzzuuulllll
|VF()) V (VF(>) (FF())-F(().2) x 3) S (VF()> VF y/l-
E
Uncertainty Estimation
E.1
Background
Deep neural networks are nowadays used in many applications such as self-driving cars and large language models. However, in the real world, they are constantly subjected to a large amount of ambiguous scenarios, and assuring fail-safes is a top priority. One approach to realize this is to design the agent so that quantifies the reliability of its own outputs and behavior, while achieving high performance. However, it is well-known that usual deep neural networks over-confidently extrapolate unknown data or poorly perform in uncertainty tasks, As result, the quantification of uncertainty has become the subject of many research, leading to Bayesian neural networks (Neal, 2012), Monte Carlo dropouts (Gal and Ghahramani. 2016) and deep ensembles (Lakshminarayanan al., 2017). Unfortunately, these methods quickly become computationally expensive as they require to process multiple forward passes or to retain distribution samples over a large number of parameters. Furthermore, such ensemble methods cannot in the essence avoid that all constitutive members make the same mistake or place high confidence on the same diitt---iitiitiiinnn data.
Therefore, a line of work that concentrates on the uncertainty estimation using a single neural network has recently started to draw attention (Liu al., 2020a; Van Amersfoort et al., 2020). In these works, they are interested in quantifying the uncertainty through single forward pass and analyzing its behavior. On their own, they can not only provide new methods with low computational cost but also address the challenges of ensemble methods as successful individual agents can again be integrated those models.
Thanks to these works, some properties necessary for an accurate uncertainty quantification have been discovered. Notably, bi- Li-Liischitzeess is an indispensable inductive bias. Intuitively, bi-Lipschitzness guarantees the neural network to be distance aware, i.c., distance is moderately preserved between the input and the feature space, resulting in correct detection of out-of-distribution points. Without this property, it is known that the problem of feature collapse occurs, which refers to the phenomenon that nut-of-distribution and in-distribution points overlap the feature space, making ouitt--diittibutioo detection impossible in principle (Van Amersfoort et al., 2020).
E.2 DUQ
We will focus on a recent model proposed by Van Amersfoort et al. (2020) as is a model that has been shown perform well as deep ensemble methods that are the state-of-the-art in this area. It also performs well on the FashionMNIST dataset.
Their model is called Deep Uncertainty Quantification (DUQ) and is mainly used for classification tasks. The idea to create class e represented by a centroid e, and to calculate the distance K. between data point and all centroids. The data will be classified as the label of the closest centroid to that point. They draw inspiration from the radial basis function kernel. The mathematical formulation of Ke is as follows:
Ko(So(z),e)) = exp W (Wff((x) - -
20²
where fo is general neural network also called feature extractor, W. is weight matrix, and a is hyper parameter called length scale.
40

---

The loss function is the sum of the binary cross entropy over all classes. During the training, the centroids are updated by the following rule:
Ne. -NN (1 Me,t,
Ma,t -, (1 1), WeJoeee))))
Cct =m,,t/Na,t
where Re,t is the number of data points assigned to class c, ae,t,á the i-th element of the minibatch corresponding to class c, and is a hyper parameter (van den Oord et al., 2017).
Bi-LLipschitzness incorporated through gradient penalty term follows:
V
Σ
Kell
The discussion about this regularization can be found Section By deleting this regularization term, and replacing the neural network fo by our bi-Lipschitz neural network, we obtain model called here DUQ+BLNN.
F
Extensions
F11
Extension 1: Different Input and Output Dimensions
With only one BLNN, we can only provide an output that has the same dimension as the input. This is common problem in bi-Lipschitz models. In this extension, we provide an architecture that extends to a more realistic case with different input and output dimensions. While we will lose many theoretical guarantees we have provided so far, some will be retained. Furthermore, throughout this paper, we will show that this method works well in practice.
The idea follows Wang et al. (2021) and Kinoshita et al. (2023). Let f+ and få (2;02) be two BLNN with input (and output) dimension dj and d2, respectively. By composing them as
(s
(16)
where D is a d2 x dN matrix with all diagonal components to we obtain a function with input dimension dj and output dimension d2. Note that if d1 < d2 then the whole bi- bi-Lipschitz property is preserved.
On the other hand, if d1 2 d2 then the Lipschitzness remains at least. Nevertheless, this parameteri- zation (16) can represent any bi-Lipschitz function 1: Rd - & with di 2 d2. It suffices fī the identity map and consider the extension / RRd Rd whose image coincides with that of / in the da dimenssionl subspace. This remains bi- bi-Lipschitz like 1. Therefore, this representation is reasonable some extent.
The gradient computation can also be explicitly formulated as follows:
Theorem F.1. Suppose model with loss included
L: (Mdd 01; h(d: o): v),
where fi (x1; 81) + 01 and f == (Fj 02) + 0x2322 are two C2 BLNN defined in Algorithm 1 with input dimensions dj and d2, respectively. D is d2 x d1 diagonal matrix with all diagonal elements set to If and F are both differentiable, then the gradients with respect to each parameter can be expressed as follows:
V -
V (E (
02) 02) -,,0)))_____ (2 Vol. -V. &,),,,,, V x) + L. h)
V - ((2, 2;)),
where E, (x 01);0)) + + 5, = (v + #22l =
([D) w Db(((h(d: 01:0) and x == h)d;
41

---

F.2 Extension 2: Non-Homogeneous F Functions
As pointed out by Jin and Lavaei (2020), in some cases, we may have knowledge of specific physical constraints or requirements to design con-homogeneous control of the sensitivity like in prssics-informed neural networks (Raissi et al., 2019). However, the current model imposes the same bi-L.isschitz constraint for each dimension. In the linear case, this corresponds to bounding all the eigenvalues of the matrix with the same constants from above and below respectively. For example, consider the following function:
f(1 x2) (2x1, 100r2).
we can only choose parameter for Lipschitzness, then would be 100. However, this ignores the non-homogeneity of each dimension and their sparseness. In the above example, the Lipschitzness with respect to 21 is only 2. Therefore, we would like to be able to control the bi-Lipschitzness the function with respect to each dimension of the input.
This can be executed by introducing a diagonal matrix and B instead of a and | in our model. That is, instead of adding a: |>/a:1/// /2 and ||v|1 (2 (see Algorithm 1), we add x: A²x /2 and y B B-2y//2 respectively, which makes possible more flexible control of the sensitivity of the overall function. The drawback is that we obtain more parameters to tune, but this can be avoided by designing these matrices as parameters to learn throughout the training too.
In the same idea, it may also be interesting to impose the convexity requirement on limited number of variables. The regendre-Fenche' transformation is then executed those variables solely, As result, we obtain partially bi bi-Lipschitz functions. This idea of incomplete convexity can be realized by the partially input convex neural network (PICNN) proposed by Amos et al. (2017). PICNN fæ, y; f), convex only with respect to with L layers, can be formulated as follows:
Us+1 = Waas + bi),
21+1 =9. (w!+*)
(z) (W: us + +b(.211 6,²)1+)
+ wr(w)
(wlww) us + b. b.00)) + W. us b.) (1 (1= O, 1),
f(x,, y: =2,, Ugg = x;
where @ denotes the Hadamar product, only Wx(*) are non-negative, di are convex non-decreasing and Wx(=) O.
F33 Extension General Norms
In this paper, we mainly discussed notions of Lispchitzness, inverse Lispchitzness, strong convexity and smoothness in terms of the L However, in some cases, it may be more interesting to work in other norms such as the Loc norms (Zhang et al., 2022). In this extension, we briefly discuss how other norms can be introduced in our theoretical framework and its advantages. The previous extension can be regarded as such an example where the L2 norm was changed to weighted L2 norm. Here, we treat the more general case. In this section, || will denote general norm.
The dual norm can be defined as follows. Definition F.2. Let ||- be a norm. Its dual || |l. is defined as follows:
|ll|+ == sup { (zx y) ||x||| = 1}
For example, the dual of the norm for any p 2 1 is the T.g-norm such that 1/p + 1/9 = 1. Definitions of strong convexity and smoothness accordingly modified as follows:
Definition F.3. Let 0. F Rm + IR is p-strongly convex with respect to a norm || for all t [0, 1]
F(Lx- + (1 ())) < (z))) + (1 - ()F t1
1) vll2.
Definition F.4. Lety > 0. F: R + R is 7-mmoohh with respect to norm || 1111 if F is differentiable
and
F(y) S F(x) + (VF()). z))
? | /
42

---

Then, the following theorem holds:
Theorem F.5 (Shalev-Shwards and Singer (2010), Lemma 18). If is a closed and p-strongly convex function with respect to a norm ||, then its Legendre-Ferne Fenchel transformation F is 1/ 1/p-smooth with respect to the dual norm 1+
Following Theorem D.2, the Legendre-Fencen Fendhel transformation of a closed and a-strongly convex function with respect to norm satisfies:
[VF VF & - yl
(17)
This means that by taking the derivative, we obtain a monotone Lipschitz function with respect to the dual norm. If] || || || [| then || ||. = || || and arrive at Lipschitzness in terms of the Tox norm which is preferred in some situations (Zhang et al., 2022). We can thus create convex smooth functions with relatively wide choice the norm.
In many applications, since all norms are equivalent a finite dimensional vector space, it may not be interesting to be able to deal with different norms." In other words, controlling the sensitivity of function f in terms of a norm can be executed simply through the calculation of the L2-norm as
C y s || yll2 S C2lla vll-
where C and C2 are constant independent of However, it turns that usually C₁ and C2 dependent on the dimension. For example, between the Lp and Lg norms where 2 the following holds:
| S ||||| S n2//0//P/|<///
As result, in high dimensional spaces, the equivalence of norms may become impractical as we will need to deal with extremely large or small values when translating the requirements into Lę norm, making the training unstable. The main advantage of our approach is that we can avoid this as we can directly characterize the function in any desired norm without any dependency on the dimension.
On the other hand, the control of the lower bound (inverse Lipschitzness) becomes more difficult. Indeed, when p or p > 2, the Lp norm is not strongly convex with respect its own norm (Acu et al., 2023). As a result, we cannot similarly proceed like our approach for the L2 norm and introduce the generalized inverse Lipschitzness.
Note that adding the squared norm |l|= to a function F satisfying equation (17), will lead for the lower bound to
ali - (allll/)))))))) -y)
-
but for the upper bound to
(x) ((F))
+alu/l/l))))) v) S011x - y11² + alIx x - v112
<(C2a + +B]l 1 yll?,
where the constant C22 of the left hand side is highly dependent on the dimension. In high dimension, this may attenuate the importance of 3, making control of the Lipschitzness difficult and unclear in practice.
The above discussion is rather theoretical. We let conerete applications, as well as deeper investiga- tions for future work. However, this may be a promising avenue since we can use other norms and the dependency on the dimension completely dropped.
F.4 Extension 4: Improved Expressive Power through Superposition and Combination
Superposition As mentioned in the main paper, while the expressive power of our bi-Lipschitz unit constrained to a + L-piishhht,,, a-strongly monotone functions which are themselves the derivative of real-valued function, this can be alleviated by superposing several BLNNs. For example, we have tested to the sign function, and the composition of 5 BLNNs could achieve an accuracy around O.
3 Any two norms || I|y and || (|( on tinite dimensional space X over C are equivalent. That is, there exist constants C1 and C2 that for all E X. Ci]|=111) (| S Ca][æ]](a)-
43

---

Table 4: Out-of-distribution detection task of CIFARIO SVHN with DUQ and BLNNconv.
Models
Accuracy
Loss
AUROC SVHN
DUQ BLNNconv then DUQ
0.929
0.04
0.83
0.930
0.04
0.89
Combination Another solution improve the expressive power of our BLNN is combine it with other architectures, While there may be various approaches, our model could be used as pre-processing module for difficult tasks since the BLNN is by definition bi-Lipschitz which means that the geometric properties of the input data is relatively preserved in the output as well. Therefore, for some problems. we could implement model where we first process the data through a BLNN and then transfer it to other networks. For example, we conducted an experiment using the convolutional version of BLNN (BLNNconv) on the problem of uncertainty estimation (as Subsection 4.2) with the CIFAR-10 10 vs. SVHN dataset to illustrate the scalability of our method (Van Amersfoort al., 2020). For this problem, we implemented the model so that we first process the data through BLNNconv and then transfer it to the DUQ. The result compared to DUQ can be found in Table 4. Our model not only scalable large-scale networks but also improves -fttef-detection performance (the AUROC of SVHN). Using BLNNconv instead of the fully connected BLNN also improved the computation time (c.g., 2.5 times faster for the 5 first iterations).
F5 A Simpler Architecture
In this work, the Legendre-Fenccel transformation used in order to provide direct and single parameterization the Lipschitzness in contrast to prior methods that controlled it on layer-wise basis. At the expense losing this benefit, is also possible to characterize Lipschitzness layer. This may be more practical and faster to apply in some contexts. In this case, we will have to control the spectral norm of all the weights. Since we calculate the derivative, the relation between the overall Lipschitz constant and that of each layer is complex though. We provide it in the following proposition as reference.
each
Proposition F.6. Define the ICNN Gg as follows
2 (((W(()) 24 + (v) y y+b) (i= (i=0,...,k 1),
Go = 2k.
where gs is 1-Lipschitz, fanction, with bounded and Lipschitz derivative. Then the Lipschitz constant of Go can be obtained by solving the following set of inequalities, where og is the Lipschitz constant:
81 -|| l|2
0++1 < V/ 12:;
[allw: + [ 2Lip(V,.)\Ix/"" "al [12] B
+ +Lpp(Vg.)[
Bi11 <l| + || 11= 1).
Proof. Let denote, 2;( the i-th layer of the ICNN with input y. By construction,
(11+1((11)0 2+ (12)
Vg: (w:(3) 26(1) + w y1 + b.) (w! Vyz((er) + w!!')
-Vg (W!( =((12) (1(())2((2)) w. u2 + b,) (W! V + wi
44

---

By defining h4(1) Vgs (w. (x) + w/ y + b.). we can further develop as follows:
V (1)))((11) - v (32) - | W (uz((nn) h2((2) vy
+ +1[A:01))w (2(02) (,0
The second term of the right hand side becomes:
(1,(w1) wi (v2)
<Lip(Vo:) ( (011) + Win0 y1 + - (w/(') 2122) + w 1/2 + +b,)) ||
<Lip(VV:)1I<( (1w - 22222) + || |||y - 3221
As for the first term,
(14(11) W V , -
- -|AA:Oo)W (, (v)())))))))(())))))))))))))))))
+h+(22) (1)) 14 (82
[1h,(11) - h
IIV (22(
<Tip(Vgs) -|]]
|ll 3all) |
+ V> [ V - (32)l1
where in the last inequality we used the definition of (v) On the other hand, we have
1o S V + w
By replacing | by B,, we obtain
8 < | |W/
Since 2y) is differentiable,
||=(v])) (x2))l S B,.
Finally, by setting oq := V (v1) 2, (1244112221111)) // //ll/llll// y2l and substituting the results into the first incquality, we obtain the desired result.
Remark F.7. Note that in the above theorem gi has to be not only convex non dend-deereiinng but also 1 I-Lipschitz with bounded and Lipschitz derivative. The layer-wise formulation requires thus further conditions on gi that were unnecessary in the BLNN (Algorithm 1). This means that our original model has more freedom and flexibility than the above simple approach.
G Additional Experiments
See Appendix for further details on the setup each experiment.
G.1 Simple Estimation of Bi-Lipschitz Constants at Initialization
First of all, we verify that BLNN behaves expected with different values of and 3. We set all ga as ReLU or softplus functions. The Legendre-Fenchel transformation step executed with GD (since it is simple and has good theoretical guarantees) until 1000 iterations are reached or the stopping condition, i.e., ||V, ((y, x)) F(y)) < 10-3 is satisfied. The BLNN has input dimension 2, and a is set to 4 as varying it does not severely impact the fundamental quality of the theoretical bounds.
45

---

Ivenn Unvchiz sound
Jeer
Figure 11: Estimated Lipschitz and inverse Lipschitz constants of BLNN with gradient descent with softplus activation functions. The left figure is with 3 hidden layers and the right with 10. The axis corresponds to j with 3 = 0.05 + 99.95 ://100.110000000000 /100..
Joher bould
Level
Figure 12: Estimated Lipschitz and inverse Lipschitz constants of BLNN with gradient descent with ReLU activation functions. The left figure is with 3 hidden layers and the right with 10. The axis corresponds to j with 3 = 0.05 + 99.95 j/100
On the other hand, 3 was changed following the sequence f0.05 99.95 5/100) See Figures 11 and 12 for results.
As we can observe, our model combined with gradient descent can generate functions with bi- Lipschitz constants that are inside the predefined bounds as expected for both softplus and ReLU activation functions. That is why, it seems reasonable to use the simple gradient descent in order to exccute the Cegendre-Fenchh' transformation.
Theoretically, since the initialization of weights are executed at random, the Lipschitz and inverse Lipschitz constants can be any value between a and a + B. Curiously though, the estimated Lipschitzness is close to the upper bound for almost every model. Note that we are generating the network and then directly estimating the bi-Lipschitz constants, which means this phenomenon is possibly closely related to the initialization of weights. To better understand this phenomenon, consider the linear network x >> Wx where W E Rxxa our architecture, this can regarded as simplification of V Fg. where Fg is the core ICNN. If we add 3 to this network and compute the inverse operation, corresponding to V F7 we obtain y ~ (W + 1/3) 11.8) By letting amin the smallest singular value of W, the Lipschitz constant becomes (amin + 1/ As a result, the initialization provides W with small min close to 0, we inevitably obtain function with Lipschitz constant close the upper bound 3, independently of the distribution of the other eigenvalues. This seems to happen in Figures 1I and 12. In other words, the construction employing the default initialization is biased, more concretely, the derivative of ICNN F, with initialized weights is with high probability e-inverse Lipschitz with << 1. This is supported by Figure 13 which caleulated the bi-Lhipchitz constants of a (4, 60)-BINN initialized over 100 different trials. While the inverse Lipschitz constant is moderately distributed, the Lipschitz constant is concentrated on the theoretical maximum, i.c., 64.
We also modified the initialization procedure of the non-negative weights W(()) to confirm our hypothesis. The default approach was to draw each element from uniform distribution proposed by Glorot and Bengio (2010) and clamp all negative elements O. Figure 14 replaces this uniform distribution interval of [0, 1), and Figure 15 to [1.0,
46

---

Figure 13: Histograms of Lipschitz and inverse Lipschitz constants of 100 (4, 60)-BBNNs with non- negative weights Wi initialized following the uniform distribution of Glorot and Bengio (2010). Negative elements were clamped to O.
Figure 14: Histograms of Lipschitz and inverse Lipschitz constants of 100 (4, 6B)-NN with non-negative weights w(()) initialized following a uniform distribution over (0, 1.0).
29
Figure 15: Histograms of Lipschitz and inverse Lipschitz constants of 100 (4, 60)-BLLNN with non-negative weights w((*) initialized following a uniform distribution over (1.0, 1.1).
As we can conclude, the choice of the initialization impacts the distribution constants for (a, 31)BBNBBNNN while theoretically any values are possible. In the default setting, we have a strong bias for the Lipschitz constant. Nevertheless, the fact that the bi-Lipschitz constants are close to the bounds can be regarded as the representation of a function with rich features. For example, a matrix has the greatest cigenvalue and the smallest eigenvalue close to the limits, implies that the other eigenvalues can assume any values in-between, allowing for considerable flexibility. This is a direct consequence of our meticulous construction of bi-Lipschitz neural networks respecting the desired constraints using convex neural networks and Legendre- Fenccel transformation.
starting bi-Lipschitz
The adequate initialization of the weights of an ICNN in order to obtain BLNNs with a richer distribution of bi-L.ipschitz constants is a whole new problem and outside the scope of our work. We will not investigate further but seems be an interesting research topic since an appropriate initialization may facilitate training. In this paper, we used the original setup.
47

---

Table Tightness Lipschitz bound with different methods. Each model was built with an upper- bound constraint on the Lipschitzness L, and the true Lipschitzness of the model after training was evaluated. The percentage between this value and L reported in the table, with mean and standard deviation over five different seeds.
Models
L = 2
L =
L = 10
= 50
(Nbr. of param.)
SN (4.3K)
+1.
80.1%+5.0
68.9 4.44
32.5%
AOL (4.3K)
96.7%
65.9%+4.7
46.5%
15.4%+2..
Orthogonal (4.3K)
98.7%:+0.4
75.1%+12..0
48.9%-+11.7
14.6%++3.7
SLL (4.2K)
97.5 %.44.4
77.6%:16.5
50.0%18.5
15.7%=
Sandwich (4.5K)
97.5%:10.7
84.3%:1..0
60.2%:14.1
16.4%.3.6
LMN (4.4K)
100.0%110.0
100.0% 1.000
98.5%.10.3
26.0%=3..6
BiLipNet (5.0K)
100.0% 0.0
98.0%:10.8
58.9%
6.8%
Ours (4.3K)
100
100.0 #0.0
99.4%+1.2
99.9%=0.1
G.2 Tightness of the Bounds: Underestimation Case
This subsection provides additional results experiment 4. of the main paper. As a reminder, this experiment verifies the tightness of the proposed architecture by fitting the following function that has a discontinuity at the point x = 0:
(x)
x x x++1
x
Therefore, when constraining the Lipschitzness of the model by positive finite L, a model with tight Lipschitz bounds should achieve that upper bound around = in order to reproduce the behavior of f. We compare our model with other Lipschitz constrained architectures that were introduced in Appendix A. We take comparison the spectral normalization (SN) (Miyato et al., 2018), AOL (Prach and Lampert, 2022), Orthogonal (Trockman and Kolter, 2021), SLL (Araujo al., 2023) and Sandwich (Wang and Manchester, 2023). They can all regulate the Lipschitzness by an upper bound, but it is executed on a layer-wise basis.
Models are built so that I, is set to 2, 5, 10 and 50. For the compared prior works, this is realized by making all layers I-Lipschitz, and scaling the input by L, which is the strategy often employed in practice. After learning was completed, the empirical Lipschitz constant was computed. The percentage between the empirical Lipschitz constant and the true upper bound L can be found in Table 5. In order to be able to legitimately compare results, the number of parameters for each model were all of the same order. The difference tightness is even clearer when setting the upper bound to 50 and plotting the result shown in Figure 16. Only our method can capture the presence of discontinuous jump around We can conclude that our method clearly achieves the tightest bound for Lipschitzness well the most coherent fit of the function.
G.3 Flexibility of the Model: Overestimation Case
This subsection provides additional results of experiment 4. 1 of the main paper. As a reminder, this experiment verifies the flexibility of our model. While we already provided a theoretical guarantee of expressive power of our model, is important to show that it can learn without problem functions with smaller Lipschitz and larger inverse Lipschitz constants than those we set parameters. At the same time, this addresses one concern of the previous experiment since the tightness of our model may be interpreted result of the strong influence the regularization terms add, effectively generating only a + 3-Lipschitz functions. Therefore, we lead this experiment, where we overestimated the Lipschitz constant of the target function. The Lipschitz constant was set to 50 but the function want to learn a linear function with slope 1. Experimental setting was the same as the previous subsection. Results are shown in Figures 17 and 18.
Our model can without any problem fit the curve y = 50x and y = r even in the overestimation scenario. Interestingly, this is not the case for prior work. On the one hand, SN and AOL cannot learn the linear function with slope 50. the other hand, others cannot learn the simple identity function after convergence. While the overall trend is along y I, can observe small fluctuations.
48

---

lAnmeed
Ieareed
14
14-
0.00
-14
-430
-1.00
-079
-000
-925
002
025
630
0.7
100
-0.
-0.2
82
SLL
2.8
The
130
Thet
learted
0.25
0.1
-01
-12
-850
635
(b) Sandwich
03
-420
-12
-1.00
-0.75
-653
525
650
2.0
(c) Ours
Figure 16: Results of fitting the curve with SLL (first row), Sandwich (sccond row) and our method (third row) with a specified Lispchitzness of 50. The right column is a zoom of the right figure.
This may be due to the large scaling factor L we multiply at the input, leading to overfitting or interpolation which usually depends on smoothness. That is, with high Lipschitzness, the function is likely to dynamically fluctuate around unseen points during training.
In order to examine the influence of this scaling factor on the training we run the experiment again with even larger L. namely, 100 and 1000. The evolution of the loss and the fitted function for 1000 can be found in Figures 19 and 20, respectively. As we can observe, only our method is hardly affected by the size of L and achieves a good generalization performance. While we do not have rigorous mathematical explanation for this phenomenon, we believe our approach provides a better regularization than the previous methods. In other words, these layer-wise methods seem to provide loss landscape that has many spurious minima where function is driven. This phenomenon is aggravated with higher L. The difference in the smooth decrease of the loss function also supports this. The presence of these fluctuations due to the large scaling factor thought to be related with the way we incorporated into the model. For previous work, had choice but
49

---

202
Ieareed
-30
-4.0
-208
-100
-10
-15
-20
-03
50
03
11
15
24
-1.00
-379
-1.0
-025
000
022
050
6.7
100
(a) SN
202
The
1.00
Tiet
-25
-.22
-72
-4.7
-1.00
825
(b) AOL
0.00
-4.0
-7
-6.7
-15
-100
055
1.
15
-1.50
-250
060
625
(c) Orthogonal
50x for left column and for right column) with SN (first row), AOL (second row) and Orthogonal (third row) with specified Lipschitzness of 50.
Figure 17: Results of fitting linear function (y
scale the function by multiplying with L. However, our method is adding regularization term without touching the core function. This is an interesting avenue for future research.
In short, our model can perfectly lear both the function with slope 50 and even though greatly overestimated the Lipsehitz constant.
G.4 Summary of the Two Previous Experiments
This experiment summarizes the two previous experiments. We propose to run regression task of the function y = 50r with SN (representing the layer-wise methods) and our model, where the Lipschitz bound L is changed from 25 to 125. If the model provides tight bounds and perfect expressive power, then the loss should equal for T. 2 50 as theoretically y = 50x can be reconstructed, and should increase once L becomes smaller than 50 as the maximal slope
reproduce is limited
50

---

202
leareed
-30
-4.0
-208
-100
-10
-15
-20
-03
50
03
11
15
24
-1.00
-379
-1.0
-025
000
022
050
6.7
100
SLL
202
Tce
1.00
-25
-.22
-72
-6.7
-3.00
-050
825
(b) Sandwich
0.00
-.50
-7
-6.25
-20
-15
-10
65
15
15
-1.50
-250
060
625
(c) Ours
Figure 18: Results of fitting linear function (y 50x for left column and for right column) with SLL (first), Sandwich (second row) and our method (third row) with a specified Lipschitzness of 50.
to L. Moreover, expect that the optimization should proceed faster when L is near 50 as the initialization of the function should also be close to y = 50x. Results are shown in Figure 21.
As can observe, SN, layer-wise model known to provide conservative bounds, only achieves 0 loss from around L = 100, while ours exactly exhibits an increase of the loss from = 50. This clearly shows that layer-wise model has expressive power that is largely lower than expected in the initial works. Notably, the convergence speed drastically decreases around = 50 for our approach. This suggests choosing bi-Lipschitz constants tight as possible tailored to the specific problem to ensure faster convergence.
Interestingly, this difference of behavior around the true Lipschitz constant of our model could provide new method to estimate bi-Lipschitz constants. Indeed, by starting with a high Lipschitz constant a + B and decreasing it, we can find the largest point where the loss starts to increase. This value can become a good estimation of the Lipschitzness of the true function. We can similarly proceed for inverse Lipschitzness. We leave further investigation for future work.
51

---

1000
109
Upper
Upper
Upper
10
Upper
10p
10P
Loss
Loss
10"1
1011
10"4
10-2
30
20
25
10
15 Epochs
20
Epochs
(a) AOL
(b) Sandwich
Uppccitz Upper
Lipschtt Upper
100
Lipachite Upper Bound
Loss
10-1
15
Epochs
(c) Ours
Figure 19: The evolution the loss with different upper bound constraints on the Lipschitz constant
with AOL (top left), Sandwich (top right) and our method (bottom).
G.5 Better Control for Annealing
Since model provide tight bounds and good regularization performance shown in the precedent experiments, the annealing its parameters could be effective in some suitable settings. Let us consider the case where want learn the exponential function. It is not Lipschitz in general since its slope is monotonically diverging to oc as x - oo and not inverse Lipschitz either since e" 0 2 -00 We will start by a low Lipschitz parameter of the model and increase it by evaluating its effective Lipschitzness during training. The amnealing will proceed in an elastic way. If the estimation is close the upper bound then we will relax the bound since means that bound too strong. That way, if we have tight bound, the annealing process should progress faster and the overall function should rapidly converge to (x). The inverse Lipschitz constant is set to 0 and kept fixed. We compared our model with the Sandwich model. See Figure 22. We can verify that the annealing process is faster than the layer-wise model. Ours constantly wants to increase the Lipschitz constant while the Sandwich layer abandons at some point and cannot keep up the increase of the Lipschitzness.
G.6
Uncertainty Estimation
G.6.1 Two Moons
This subsection provides additional results of experiment 4.2 of the main paper. As a reminder, we compare the uncertainty when learning the two moons dataset with not only DUQ and DUQ+BLNN but also the Deep Ensembles method and DUQ with no bi-Lipschitz regularization. Results are plotted in Figure 23. Uncertainty is calculated as the distance to the closest centroid (see Appendix E for the mathematical formulation). Blue indicates high uncertainty, and yellow low uncertainty.
As mentioned in the main paper, DUQ+BLNN achieves tighter bound which is possible because we are able tune the bi-Lipschitz constant unlike DUQ. The experimental setup for DUQ was chosen like the original paper (Van Amersfoort et al., 2020).
52

---

#stted
IANTHI
000
50-
-0.0
-1.0
-10
1.00
-0.5
-250
-0.23
6.00
0.2
050
0.73
1.0
-1.00
-573
-0.30
-4.00
0.00
023
033
235
182
(a) AOL
(b) Sandwich
1.00
Tul
325
-025
-1.0
-035
(c) Ours
Figure 20: Results of fitting the linear function y = x with AOL (top left), Sandwich (top right) and our method (bottom) with specified Lipschitzness of 1000.
Table Uncertainty quantification of two moons dataset with DUQ+BLNN with different a and B. Mean and standard deviation over five different seeds.
0.1
1.0
2.0
3.0
4.0
5.0
0.0
0.716#.111
0.87411006
0.99011002 0.9041.045
0.972+.043
0.95711.048
1.0
0.8791.004
0.9
0.9941.004
0.992:10044
0.895#.195
0.993+.007
2.0
0.878+.006
0.990=.003
0.993+.004
0.936:.087
0.99511.003
0.5181.047
5.0
0.50
0.500=.000
0.541+.080
0.503+.004
0.501#.002
0.503+.004
Table 6 provides details on the accuracy of the learned model with different and 83 used for the grid search. As for the accuracy, the value of a does not matter as long as 3 is chosen large enough. Indeed, BLNN with (a, 3) (0.0,3.0) theoretically includes that with (a, 3) (1.0, 2.0). That is why they exhibit similar performance. Nevertheless, too loose setting of these parameters, especially for 8, still seems affect the training, and tight possible choice is beneficial as can see for () = (2.0, 5.0). sense, decreasing & which mainly in charge of the Lipschitz constant helps better generalization pointed out in prior work well (Bartlett al., 2017). the same direction. a too large value of a leads to nonsensical result as we can understand when a 5. In this case, the model is randomly guessing which means that unsure everywhere. This is also translated in the uncertainty estimation as shown in Figure 24. The area where the model is unsure increases as is increased, meaning that becomes more sensitive to unknown data (see Figure 25). Comparatively. setting a = 0 does not affect so much since the injectivity constraint is already beneficial in this type of task. Finally, too low value of 3 means that the inverse Lipschitz and Lipschitz constants are close to each other and the function behaves more like linear function. Consequently, learning becomes difficult as we can conclude from the table with = 0.1. Similar trends could be also observed in more complex tasks such as FashionMNIST MNIST and FashionMNIST NotMNIST.
53

---

3000
Mm
400
100
120
100 120 120
(a)
Eoraton
LO35
800
100
120
110
(b) Ours
Figure 21: Results of fitting the linear function y 50x with different L for SN and model. The left figure is the loss in function of L, and the right figure is the first iteration that the loss was below 0.5 for each L (if there no point, means that loss below 0.5 never reached for this value of L). The red line emphasizes the point of L = 50.
Table Out-of-distribution detection task with down-sampled data. FashionMNIST vs MNIST and FashionMNIST NotMNIST dataset with DUQ and DUQ+BLNN, For the BLNN, 0.2 and 0.4. Mean and standard deviation over five trials.
AUROC MNIST
AUROC NotMNIST
Model
Accuracy
BCE (loss function)
DUQ
0.85641.0066
0.078=0004 0.078-.004
0.870+L.040
0.852:10110
DUQ+BLNN 0.8595 1.8595-.0032
0.876:1.008
0.960:.000
G.6.2 FaSHionMNIST
This subsection provides additional results of experiment 4.2 the main paper. In this experiment, we use real world data of FashionMNIST (Xiao et al., 2017), MNIST (LeCun and Cortes, 2010) and NotMNIST (Bulatov, 2011). The task learn to classify FashionMNIST but at the same time to detect out-of out-of-distribution points from MNIST and NotMNIST. This task to distinguish FAST from MNIST datasets is known to be a complicated task (Van Amersfoort al., 2020). We compute the AUROC for the detection performance. Results for downsampled dataset and full size dataset are shown in Table and respectively. visualization the ROC curve can found in Figures 26 and 27.
54

---

Interated ppccermees
2
20
00
50
300
(a) Sandwich
The
insurct Ubscc tanca Exttr Lpedeereee.
66-
(b) Ours
Figure 22: Results of fitting y exp(x) with Sandwich (top row) and our method (bottom row). The right figures are the evolution of the upper bound and the estimated Lipschitz constant over the iteration, and the left column the learned function. The x-axis of the right pictures represents the iteration number.
Table 8: Out-of-distribution detection task with full size data. FashionMNIST MNIST and FashionMNIST vs NotMNIST dataset with DUQ and DUQ+BLNN. For the BLNN, a 0 and 3.0. Mean and standard deviation over five trials.
AUROC MNIST
AUROC NotMNIST
Model
Accuracy
BCE (loss function)
DUQ 188933+00037 0.8985=.0040
0.064=.005 0.060=.000
0.862+.035
0.900
DUQ+BLNN
0.896:2008
0.964+.005
3.0
1.-
TPR: Thue positive rate
TPh: True positive
62
50
positive
FPR: positive
(a) FashionMNIST vs MNIST
(b) FashionMNIST vs NorMNIST
Figure 26: ROC between dwonsampled FashionMNIST vs MNIST (left) and FashionMNIST vs
NotMNIST (right).
55

---

(a) Deep ensemble
(b) DUQ reg.
-2
3
(c) DUQ
(d) DUQ+BLNN
Figure 23: Uncertainty estimation with the two moons data set. Left figure is with a simple neural network without any constraints, and the right is with our model constraining bi-Lipschitzness.
2
Figure 24: Uncertainty quantification of two moons dataset with DUQ+BLNN with high 5.0
and B = 3.0. We do not show the points so that the highly certain area is visible.
True positive rate
[P8: True positive
0.6
0.
0.4
3UO JMCH-ALAN
04
[PR: Talse positive rate
FPR: False positive rate
(a) FashionMNIST MNIST
(b) FashionMNIST vs NotMNIST
Figure 27: ROC between full size FashionMNIST vs MNIST (left) and FashionMNIST vs NotMNIST (right).
56

---

(a) (UQ+(0, 4)-BLNN
(b) DUQ+(L, 4)-N
-2
(c) DUQ+(2, 4)-BLNN
Figure 25:
Uncertainty quantification of two moons dataset with DUQ+BLNN with different a.
Table 9: Comparison of our model with state-of-the-art monotone models. Means and standard deviations over three trials. Results of prior models are from the original papers. COMPAS, BF BlogFeeeBBcck. LD LoanDefaulter, HD = HcartDisease, AM AutoMPG, Acc. = accuracy.
Models C (Acc.)
BF (RMSE)
LD (Acc.)
HD (Acc.)
AM (MSE)
Certified
68.8 + 0.2
0.158 + 0.001
65.2 + 0.1
Constrained
69.2 = 0.2
0.154 0.001
65.3 0.0
0.0
8.37
+ 0.1
LMN
69.3 + 0.1
0.160 + 0.001
0.0
1.6
7.58
- 1.2
SMNN
69.3 = 0.9
0.150
0.001
6.10
8.0
7.44 1.2
Ours
69.4 + 0.6
0.157 00.001
6..5
1.6
7.13 + 1.2
G.7 Partially Monotone Settings
Table 9 provides the experimental results with full data. We also compare with the certified monotone neural network of Liu et al. (2020b) (Certified) and the constrained monotone neural network of Runje and Shankaranarayanz (2023) (Constrained). See Tables 2 and 3 of Kim and Lee (2024) for a comprehensive comparison with other general models.
H Experimental Setups
In this appendix, we provide further experimental details, including values of hyper parameter, specific composition of used architectures and optimization schemes, which were omitted in main text as concern of clarity. All experiments were executed with Pytorch. We mainly used the GPU NVIDIA A 1000 with a memory 80GB for the computation. Additional information for each experiments can be found in each subsection below and in the supplementary material.
57

---

Table 10: General details on the architectures of Figures and 9.
Model
Hidden Dimension
Number of Layers
Ours
10
2
H.1 Computational Complexity Comparison
In this experiment, corresponding to Figure 2, presented comparison of the computational complexity for a single iteration between traditional feedforward network and various BLNN variants.: BLNN (with brute force backpropagation). BLNN with Theorem 3.7, PBLNN with only one variable constrained to be bi-Lipschitz. Respective parameter sizes are 1.11M, 1.42M, 1.42M and 1.38M. The input size dataset, and the batch size was varied the set 1, 5, 25, 75, 100. For this experiment, BLNN was implemented according to the template of the DEQ library FLOPs was computed following another github library 5. The feedforward neural network had ReLU activation function and 4 layers with hidden dimension 210, and the BLNN variants had softplus activation function and 3 layers with hidden dimension 150. a and 3 were both set to
randomly generated data of size 3 x 32 x 32, simulating CIFAR
H.2
Algorithms for LFT
H.2.1 Influence of Approximate Optimization on Bi-I-ipschit Constants: Experiments
In this experiment, were interested in the evolution of the bi-Lipschitz constants under different optimization schemes.
Corresponding Figures and Tables
Figures 8 and 9.
Data We generated 5000 2-dimensional points uniformly sampled at random in the interval [-1,9] [-1.9].
Architecture Our BLNN was created with / 10. The value of a does not matter since we focused on the output of the optimization process which occurs before adding the regularization term a/ // We used the softplus function as the activation function. The effective Lipschitz and inverse Lipschitz constants were estimated by simple sampling with sample size 5000 chosen from the data.
Optimization Schemes We used several optimization schemes, namely, steepest gradient descent (GD), Nesterov's 's accelerated gradient (AGD) (Nesterov, 1983), Adagrad (Duchi et al., 201 1), RM- Sprop (Hinton et al., 2012), Adam (Kingma and Ba, 2015) and the Newton method For Adagrad, R.Mpropp and Adam, all parameters were default except the learning rate.
H.2.2 On the Choice of the Initial Point
In this experiment, we were interested the influence of the initial point on the convergence speed of the Legendre-Fennchll transformation. We used as starting point either the final point of the previous epoch for each training data or the point (I....,1) independently of the history. This experiment was run with the architecture explained in Subsection H.4.2.
Corresponding Figures and Tables
Figure 10.
H.3 Bi-Lipschitz Control
The codes were implemented in Python 3.1 with PyTorch 2.0.1
httppe//ggichbb.ccmmoooussaaaa
httt:::://gitabb..o//////////// flops pytorch
58

---

H.3.1 Simple Estimation of Bi-Lipschitz Constants at Initialization
In this experiment, were interested in verifying in more detail whether an (a: B)-BLNN respects the pre-defined bounds of bi-Lipschitzness, namely, a and a + 3. We used GD for the computation of the Legendre-Fenchel transformation since it showed good performance and practical usefulness in the previous sections.
Corresponding Figures and Tables Figures 11, 12, 13, 14 and 15,
Data We randomly created 200 2-dimensional points, and estimated the effective bi-Lipschitz constants with them.
Architecture For a fixed activation function (ReLU or softplus) and number of hidden layers (3 or 10), we created 100 architectures with different and a was to 4 for all 100 models, and / as 0.05 + 1909 where i E
The value of the hidden dimension was set to 10.
H.3.2 Tightness of the Bounds: Underestimation Case
In this experiment, were interested in the tightness of the Lipschitz bound provided by each framework. The experiment was inspired from Figure of Wang and Manchester (2023). Codes also imported and modified from their github combined with others?
Corresponding Figures and Tables
Tables I and 5 and Figures 3 and 16.
Data We created 300 training points (z, f()),, where x was sampled at random in the interval [_ The test dataset consisted of 2000 points, with in the interval [-1, 1].
Optimization Based Wang and Manchester (2023), used the Adam optimizer with learning rate 0.01, and other parameters were set as default. The objective function was the mean squared error.
Architecture In addition our model, we used seven layer-wise Lipschitz architectures, namely, spectral normalization (Miyato et al., 2018). AOL (Prach and Lampert, 2022), Orthogonal (Trockman and Kolter, 2021). SLL (Araujo et al., 2023), Sandwich (Wang and Manchester, 2023). LMN (Nolte et al., 2023) and BiLipNet (Wang et al., 2024). See Appendix A for the mathematical formulation of each Lipschitz approach. Each layer was designed to be I-Lipschitz. and we employed ReLU as activation function for all of them, including ours, except that for LMN. Since the Lipschitz constant of the overall network input by L. As result, we obtain T. T--Lipschitz neural networks. In the experiments, T. was set to 2, 5, 10 and 50. For our model, corresponds to a + B. We set a = 1 and changed 3 accordingly. See Table 11 for further details the architecture. The effective Lipschitz constant simple sampling with sample size 1000 within the range [-1, 1]. As for BiLipNet, we only used their monotone Lipschitz layer and orthogonal layer in order to avoid losing tightness by composing many simple layers.
we included a scaling factor before the first layer which multiplies the
estimated by
H.3.3 Flexibility of the Model: Overestimation Case
In this experiment, we were interested in the influence of the overestimation of the Lipschitz constant when building the model. The dataset, architecture and optimization scheme were mostly the same as H.3.2. We just changed the function to be fitted to linear one with slope and 50. The architectures were designed to have Lipschitz bound of 50, 100 or 1000.
Corresponding Figures and Tables Figures 4, 17, 18, 19, 20.
"attp:///gttabb.cc//cttttt htttp:://ggghuu..oc///aaaaaaat/s
59

---

Table 11: General details on the architectures of Table 5 and Figure 16.
Model
Hidden Dimension
Number of Layers
Size of Model
SN
45
4.3K
AOL
45
4
4.3K
Orthogonal
45
4.3K
SLL
45
2
4.2K
Sandwich
65
4.5K
LMN
64
3
4.4K
BiLipNet
5.0K
Ours
64
2
4.3K
H.3.4 Summary of the Two Previous Experiments
In this experiment, we were interested the relation between the imposed Lipschitz upper bound and the loss function when learning the function y 50x. The dataset, architecture and optimization scheme were mostly the same as the previous section. We just changed the activation function the BLNN from ReLU to softplus. We also reported the first time when the loss reached value below 0.50.
Corresponding Figures and Tables
Figures and 21.
H.3.5 Better Control for Annealing
In this experiment, we were interested in the behavior of our model under an annealing scheme. The dataset, architecture and optimization scheme mostly the same as H.3.2. We just changed the function to be fitted to = e². Moreover, we set & O for our model. The training was executed for 200 epochs.
Corresponding Figures and Tables Figure 22
Annealing We employed simple annealing scheme. We started with Lipschitz constant of y 2. For the Sandwich model, it is the scaling factor at the input, and for ours, it is a + 3 since a = 0. We estimated the Lipschitz constant each iteration, and if the effective Lipschitzness was close to the imposed one, i.e., 7, with an accuracy of 0.05, we relaxed the condition by multiplying by 1.5. This update was executed every five epochs, if the condition was satisfied, in order provide the model the time to learn a new appropriate representation.
H.4 Uncertainty Estimation
In this chapter, the experiments were inspired from Figure of Van Amersfoort al. (2020). Codes were also imported and modified from their github. The codes were implemented in Python 3.10.13 with PyTorch 2.0.1 1ccul
H.4.1
Two Moons
In this experiment, we were interested in the uncertainty estimation when learning the two moons dataset. We compared the performance of our model with other existing models.
Corresponding Figures and Tables Figures 5, 23, 24, 25 and Table 6.
Data The moon dataset was generated using the sklearn toolkit with noise of 0.1. We used 1500 points for training and 200 for tests. Batch size was 64.
Optimization We used SGD with learning rate of 0.01, momentum of 0.9 and weight decay of 10-4, We used the binary cross entropy as loss function.
httt::///gttabl..oo///////eemmmp
60

---

Table 12: General details on the architectures of Figure 23. "reg." " stands for regularization.
Hidden Dimension
Num. of Layers
Centroid Size
Num. of Parameters
Output Dim.
DUQ DUQ no reg.
40
3
40
10
4.2K
40
3
40
10
4.2K
DUQ+BLNN
20
2/2
40
10
3,4K
Deep Ensembles
20
4
0.9K x 5
Table 13: General details on the architectures of Table and Figure 26. For DUQ+BLNN, since the input and output dimensions were not equal we used the architecture explained in Extension 1 of Appendix F. While these parts share almost the
structure, we clarify the differences by the following notation "info. of first BLNN/info. of second BLNN".
Hidden
Num. of Layers
Num. of Parameters
Output Dim.
Centroid Size
Dimension
DUQ DUQ+BLNN
49
3 3/1
256
100
280.9K
40
256
100
293.7K
Architecture We used four models, namely, deep ensemble method (Lakshminarayanan et al., 2017), the original DUQ, DUQ with regularization and DUQ with method (DUQ+BLNN), The latter means that we used the DUQ framework but we deleted the gradient penalty of the loss function and directly replaced the vanilla neural network with BLNN. We used the softplus function for all architectures. See Table 12 for further details on the architecture. As for our model, the best values of and 3 were found through a grid search: and 3 are shown in Table 6.
and 4. Results with other
H.4.2 Fashion-MNIST (downsampled)
In this experiment, we were interested in applying BLNN to the task of out-of-distribution detection of FashionMNIST MNIST and FASHIONMNIST NotMNIST datasets, and comparing its performance with existing models. This is a downsampled version of the next experiment.
Corresponding Figures and Tables
Tables 2, 7 and Figure 26.
Data The FashionMNIST dataset was generated using the torchvision toolkit, with 60000 training data and 2000 test data. Batch size was set to 128. The dataset was downsamlped by pooling data from 28 28 to 14 14.
Optimization
We used SGD with learning rate of 0.05, momentum of 0.9 and weight decay of 10-4, We used the binary cross entropy as loss function.
Architecture We compared two models DUQ and DUQ with our method (DUQ+BLNN), The latter means that we used the DUQ framework but we deleted the gradient penalty of the loss function and directly replaced the vanilla neural network with our BLNN. We used the softplus function for DUQ+ours and ReLU for the other architectures. See Table 13 for further details on the architecture. As for our model, the best values of a and / were found through a grid search: a = 0.2 and 3 = 0.4. For both models, data was down-sampled by max-pooling layer from 28 28 to 14 x 14.
H.4.3 FaShion-MNIST (full size)
In this experiment, were interested in applying BLNN to the task of out-of-distribution detection of FashionMNIST vs MNIST and FASHIONMNIST vs NotMNIST datasets, and comparing its performance with existing models. This is the full size version of the previous experiment.
Corresponding Figures and Tables
Table 8 and Figure 27.
61

---

Table 14: General details on the architectures of Table 8 and Figure 27.
Hidden Dimensions
Num. of Layers
Num. of Parameters
Model
Ouiput Dim.
Centroid Size
DUQ
196/49
256
100
432.3K
DUQ+BLNN
150
50
100
449.2K
Data The FaSHionMNIST dataset was generated using the torchvision toolkit, with 60000 training data and 2000 test data. Batch size was set to 128. We used the original size of the dataset, namely 28 28.
We used SGD with learning rate of 0.05, momentum of 0.9 and weight decay of 10 We used the binary cross entropy loss function.
Optimization
Architecture We compared two models DUQ and DUQ with our method (DUQ+BLNN), The latter means that we used the DUQ framework but we deleted the gradient penalty of the loss function and directly replaced the vanilla neural network with our BLNN. We used the softplus function for DUQ+ours and ReLU for the other architectures. See Table 14 for further details on the architecture. As for our model, the best values of and were found through a grid search: a 0.0 and 3 = 3.0. Particularly, only used one BLNN in this experiment compared to the previous one, and projected it to a smaller sub-space of size 50 by diagonal rectangular matrix with all components set to 1.
H.S Partially Monotone Settings
In this experiment, we were interested in applying PBLNN to several tasks where data exhibited mono- tone behaviors with respect to some variables. The codes were imported and modified from (Nolte et al., 2023). 10 They were implemented in Python 3.10.14 with PyTorch 2c0.1.cu117.
Corresponding Figures and Tables
Table 3.
Data COMPAS (Angwin et al., 2016) is classification task with 13 features, where of them have monotone inductive bias. important to highlight that this dataset known be racially biased (Angwin et al., 2016). BlogFeedBack (Nolte et al., 2023) is prediction task with 276 features, where of them have monotone bias. LoanDefaulter (Nolte et al., 2023) a classification task with 28 features, where 5 of them have monotone inductive bias. HeartDisease (Janosi et al., 1988) is also classification task with features, including with monotone inductive bias. AutoMPG (Quinlan, 1993) contains 7 features where 2 of them shows monotone inductive bias. CIFAR101 is augmented dataset used by Nolte et al. (2023).
Optimization We used the same training scheme Nolte et al. (2023).
Architecture We compared our model, the PBLNN, with that of Nolte et al. (2023) and that of Kim and Lee (2024). As for the former, we constructed the PBLNN based the inductive bias. For HeartDisease and AutoMPG, we took as the final output the average of the output of one PBLNN with bi-Lipschitzness imposed on the monotone variables. For COMPAS and LoanDefaulter, created independent 1-dimensional PBLNNs whose bi- bi-Lipschitzness is imposed on each variable that monotone and then took the average. For CIFAR101, used 101 PBLNNs with respect the last variable. See Table 15 and the codes for further details.
62

---

Table 15: General details on the architectures of PBLNN used for Table 3.
Hidden Dimensions
Num. of Layers
Inverse Lipschitzness
Dataset
Lipschitznsss
COMPAS
45 x4
4
0.3
0
BlogFeedBack
3x 8
1
LoanDefaulter
10x 5
2
0.5
1
HeartDisease
10
10
2
AutoMPG
20
3
1
0
CIFAR101
10x101
NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification:
See Abstract and Introduction.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
is fine to include aspirational goals as motivation as long as is clear that these goals not attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification:
See end of Section 3.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (c.g.. independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g. if the approach was only tested few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might not be used reliably provide closed captions for online lectures because it fails handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimess..
63

---

While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?
Answer:
[Yes]
Justification:
See appendices and main paper.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proof's provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should properly referenced.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided not)?
Answer: [Yes]
Justification: See Appendix H.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily a new model architecture, the paper should describe the architecture clearly and fully.
64

---

reproducible, they (c) If the contribution is a new model (e.g.. a large language model), then there should either way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility. the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described supplemental material?
Answer: [Yes]
Justification:
See supplemental material.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips.. cc/ public/ mublic/guides wetgggg//cdeuubbeee (rodeSubmissionP. oooicyyyy for more details.
While encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https: / //niss. wwwp//pllcgggpiiieeeeeeeeeoeeeee /ttec/gguuddeeoeeuuuuubmiiiiii Cooeeuumiissseeeeee wwriiisssooPoiooooooo for more details.
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code is permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary understand the results?
Answer: [Yes]
Justification:
See supplemental material.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment Statistical Significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [Yes]
Justification:
See experiments in appendix.
65

---

Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by enror bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall
with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
should be clear whether the eπror bar is the standard deviation or the standard error of the mean.
It is OK to report I-sigma error bars, but should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not show in tables figures symmetric error bars that would yield results that are out of range (c.g. negative error rates).
error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
Question:
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Answer: [Yes]
Justification: See Appendix H.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't 't make it into the paper).
9. Code Of Ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics https: htt//enurii...e. http:p/p///pppppp.cc///////////s
Answer: [Yes]
Justification:
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g.. there special consid- eration due to laws or regulations in their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive socictal impacts and negative societal impacts of the work performed?
Answer: [NA]
66

---

Justification: This research is foundational and not tied to particular applications, let alone deployments. This work is primarily proof of concept for a novel paradigm control the sensitivity of neural networks, which could serve to increase robustness against adversarial attacks and uncertainty quantification performance. However, there are no direct societal consequences which we feel must be specifically highlighted here.
Guidelines:
The answer NA means that there is no socictal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there negative societal impacts, the authors could also discuss possible mitigation strategies (e.g.. gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of MI...
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data or models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification:
The paper poses no such risks
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (c.g., code, data, models), used in the paper, properly credited and the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
Justification:
See supplemental material.
Guidelines:
The answer NA means that the paper does not existing assets.
67

---

The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset. For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification:
See supplemental material.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/model part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsoureing nor research with human subjects
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Bthies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsoureing nor research with human subjects.
68

---

Guidelines:
The answer NA means that the paper does not involve crowdsoureing nor research with
human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent)
may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such the institution conducting the review.
69


Reviews:
[{"review_id": "wDIICu5jAD", "reviewer": "Reviewer_u3yR", "summary": "This paper proposes to control the bi-Lipschitzness of a neural-network by parameterizing the output by the Legendre-Fenchel-Dual. This involves parameterizing a strongly convex function and computing the minimum of that function in the forward pass. Several benchmarks are studied in simple regression tasks and uncertainty quantification.", "questions": "-Is it also necessary that BLNN is a strongly monotone function? It seems that many of the regression experiments involve monotone target functions (figure 2 and 3), but I'm not sure if that is because BLNN is not capable of representing monotone functions or just not a great representer due to your approximation theorem. If it can represent non-monotone functions, it would be interesting to see a simple regression comparison to SLL, AOL, etc. Answering this question will greatly help my evaluation.\n\n-The Lipschitz parameterization of SLL, AOL, and Sandwich layers commonly uses compositions of 1-Lipschitz layers for the application of certified robustness. How would the BLNN parameterization compare to existing 1-Lipschitz layer networks in the certified robustness setting? I’d imagine BLNN might be much more expressive than composed 1-Lipschitz layers which could have a big impact.\n\n-Have you considered amortizing the LFT computation as done in this paper? https://arxiv.org/abs/2210.12153 \n\n-I'm curious if there is any possibility of extending BLNN to convolutional layers? These settings are interesting for larger image classification problems like CIFAR10 and Imagenet.", "limitations": "Limitation are adequately addressed.", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "-The framework is interesting because it parameterizes bi-Lipschitz networks in a way that is not layer-wise and instead takes advantage of the Legendre-Fenchel transform LFT / convex conjugate of parameterized strongly-convex functions (ICNN), which only modifies the output of the output.\n\n-Computing the LFT of a given function can be costly, however the paper offers a non-asymptotic bound for the Lipschitz constant and tractable gradient.\n\n-The experimental results show a considerable improvement in tightness and regularity over other Lipschitz controlled networks like spectral normalization, AOL and Sandwich layers on small regression tasks. In particular BiLipNet behaves a lot better when the Lipschitz constant is overestimated in existing parameterizations.", "weaknesses": "-Computing the LFT seems to be quite expensive, hence why the experiments are only on simple 2d problems and fashion-MNIST. For this reason I'm doubtful that it will be used for any large-scale network training pipelines where tight Lipschitz control and estimation is challenging.\n\n-The provable approximation class is limited to alpha-strongly monotone functions and is the derivative of some function almost everywhere. Lipschitz layers like AOL, SLL and Sandwich layer are all solutions to the LipSDP framework which only requires the activations themselves to be alpha-strongly monotone for alpha >= 0 (Fazlyab et al., 2019).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you very much for spending your time on carefully reviewing our paper. We really appreciate all the questions highlighting the significant potential and future directions of our work. Please find below answers to your questions and some clarifications to other important points of your review.\n \n(W=Weakness, Q=Question)\n\n>Q4\n \nYes, we can extend our BLNN to convolutional layers. We just have to change the core neural network of the BLNN, which is the ICNN, to the convolutional ICNN proposed in the original work of Amos et al. (2017). We will definitely clarify this in the revised version. In the next answer, we provide experimental results with this architecture. \n \n> W1\n \nWe conducted an additional experiment using the convolutional version of BLNN (BLNNconv) on the problem of uncertainty estimation (as in Subsection 4.2) with the CIFAR-10 vs. SVHN dataset to illustrate the scalability of our method. For this problem, we implemented the model so that we first process the data through BLNNconv and then transfer it to the DUQ. The result compared to DUQ is as below. Our model is not only scalable to large-scale networks but also improves out-of-detection performance (the AUROC of SVHN). Using BLNNconv instead of the fully connected BLNN also improved the computation time (e.g., 2.5 times faster for the 5 first iterations). Please also refer to the global rebuttal for this topic. \n \n|      | Accuracy  | Loss | AUROC SVHN |\n|-|-|-|-|\n| DUQ   | 0.929    | 0.04 | 0.83|\n| BLNNconv | 0.930    | 0.04 | **0.89**   |\n \n> W2\n \nThis is indeed a limitation of our algorithm but also an improvement compared to prior works in some aspects. From the theoretical perspective, several (bi-)Lipschitz layers do not provide a provable approximation class, and, from the practical perspective, some bi-Lipschitz models are also restricted to strongly monotone functions (e.g., BiLipNet). We would like to add that most existing layer-wise approaches have difficulty in creating partially bi-Lipschitz networks (e.g., controlling the spectral norm of linear matrices cannot handle this case) but ours can easily realize this with the PBLNN. Finally, while out of the scope of this work, there are some straightforward ways to increase the expressive power of BLNN. Please see the next answer.\n \n> Q1\n \nThe following answer is also related to W1 and W2, and is therefore a bit longer. \n \nYes, the vanilla BLNN is necessarily a strongly monotone function. However, our model still has a higher expressive power and tighter bounds than other (bi-)Lipschitz units (= the most basic architecture guaranteed to be bi-Lipschitz without using composition of several bi-Lipschitz entities). Indeed, most existing models can only create simple (bi-)Lipschitz units with low expressive power (e.g., only linear) and they have to compose those units to achieve higher expressive power, but this leads to looser bounds. In comparison, our model can **with only one unit** achieve high expressive power and keep tight bounds. This is an improvement realized thanks to a parameterization not layer wise. Experiments of figures 2 and 3 were designed to illustrate this fact. \n \nNow, we can of course, like any previous (bi-)Lipschitz methods, start to employ our BLNN unit as a component of a deeper network. For example, we can successively compose several BLNNs to improve its expressive power. With this approach, it **can indeed represent non-monotone functions**. We have tested to fit the sign function and the composition of 5 BLNNs achieves an accuracy around 0.  \n\nSome other extensions to increase the flexibility were already mentioned such as the PBLNN, and the BLNN used as a pre-processing stage for the next networks as we did in the answer to W1. Future research should work on how to use this novel unit in various fields inside and outside bi-Lipschitz problems.\n\nNevertheless, this direction is slightly out of the scope for this work. Indeed, our work is primarily foundational, providing a novel paradigm for the creation of bi-Lipschitz architectures and addressing some non-negligeable issues of the field. Therefore, this paper focuses on the fundamental design and properties of our model related to the control of bi-Lipschitzness and its theoretical analysis. One of the main contributions of this work resides thus in the construction of such a general framework for improved bi-Lipschitz control equipped with theoretical guarantees. While we ran several experiments related to bi-Lipschitz control to illustrate the benefits of our model, we acknowledge that there are some applications outside the bi-Lipschitz control framework that we could not investigate. This is also why we have tried to provide an extensive discussion and analyses in the appendix and kept the formulation as general as possible to promote further extensions and encourage other interpretations for future work. \n \n> Q2\n\nThank you for pointing out this interesting possibility. We also believe that BLNN might perform better in this setting by composing many BLNNs. However, this is out of the scope for this paper as explained in the above answer. This is still one of the first avenues we will investigate for future works.\n \n>Q3\n \nAmortizing the LFT computation as the cited paper may be a fundamental solution to the computational complexity of our model (which can already at least scale to CIFAR10). However, this approach is currently out of the scope of this work as the main goal of our paper is to provide a bi-Lipschitz architecture with guaranteed bi-Lipschitz control and build its foundation starting with well-known basic algorithms such as the gradient descent. This amortization technique is more involved as we have to carefully analyse how such amortization influences the provable bi-Lipschitz bounds of the architecture, which is not trivial.\n \nWe hope we addressed all your questions and concerns. We will be glad to provide further explanation and clarification if necessary.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "5WRSfjcRzm", "reviewer": "Reviewer_XcLq", "summary": "This paper proposes a novel neural network architecture called BLNN (Bi-Lipschitz Neural Network) that allows direct control and parameterization of the overall bi-Lipschitzness of the network. The main contributions include: i) a framework that allows tight control of Lipschitz and inverse Lipschitz constants of networks via using convex neural networks and the Legendre-Fenchel transformation, ii) comprehensive theoretical analysis, iii) empirical evaluation showing the nice performance of BLNN on tasks like function fitting, out-of-distribution detection, and monotone regression.", "questions": "1. How does the proposed method perform on larger, more complex datasets like TinyImageNet or ImageNet?\n\n2. Can the authors clarify the computational complexity of their approach?\n\n3. Can the authors provide a more comprehensive study on hyperparameter sensitivity?\n\n4. Can the authors comment on other network structures (e.g. implicit models, DEQs, etc)?", "limitations": "It seems that an improved discussion on potential negative societal impacts or broader ethical considerations is still missing.", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "**Originality:**\n\nThe paper presents a novel approach to constructing bi-Lipschitz neural networks that is distinctly different from existing methods. The use of convex neural networks and Legendre-Fenchel transformation to directly parameterize overall bi-Lipschitzness is quite novel. The extension (e.g. partially bi-Lipschitz networks, etc) is also new.\n\n\n**Quality:**\n\nThe quality of the paper is good. The authors provide detailed proofs and analyses for their key claims, including the bi-Lipschitz properties of their construction and the expressive power of the resulting networks. The experiments cover various scenarios, from simple function fitting to uncertainty estimation and monotone regression. The results are quite competitive.\n\n\n**Clarity:**\n\nThe paper is generally well-structured and clearly written. However, given the technical nature and the length of the paper, understanding the paper fully is still a tough task.\n\n**Significance:**\n\nThe paper's contributions are significant in its solid theoretical developments. The significance is further underscored by the improved performance on tasks like out-of-distribution detection and monotone function learning. In conclusion, this paper presents a novel approach to an important problem in deep learning theory and practice.", "weaknesses": "1. Computational Complexity: A detailed analysis of time and space complexity compared to traditional networks can be helpful.\n\n2. Scalability and Practical Implications: There's insufficient exploration of how the method scales to very large networks or complex datasets (e.g. TinyImageNet).\n\n3. Hyperparameter Sensitivity: More discussions on this issue will be beneficial.\n\n4. The paper could be more explicit about scenarios where the theoretical guarantees might not hold, and could explore potential extensions to other network architectures beyond feedforward networks.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you very much for spending your time on carefully reviewing our paper. We really appreciate all the advice you provide to improve our paper. Please find below answers to your questions. We also summarized the comparison of the time and space complexity of our models in an independent thread.\n\n(W=Weakness, Q=Question)\n\n>W1+Q2\n\nWe apologize for omitting this important point. Please find in the global rebuttal both theoretical and experimental discussions on this topic. This will be added in the updated paper.\n\n> W2+Q1\n\nThe main goal of our work was to create a new paradigm for the bi-Lipschitz control and to analyse its behavior. Nevertheless, we understand this is an important direction of exploration. We conducted an additional experiment using the convolutional version of BLNN (BLNNconv) on the problem of uncertainty estimation (as in Subsection 4.2) with the CIFAR-10 vs. SVHN dataset to illustrate the scalability of our method. We implemented the model so that we first process the data through BLNNconv and then transfer it to the DUQ. The result compared to DUQ is as below. Our model is not only scalable to large-scale networks but also improves out-of-detection performance (the AUROC of SVHN). Using BLNNconv instead of the fully connected BLNN also improved the computation time (e.g., 2.5 times faster for the 5 first iterations). Moreover, the amortization method of Reviewer u3yR may be a fundamental solution to the computational complexity of our model (which can already at least scale to CIFAR10) but out of the scope for this current work.\n\n| | Accuracy   | Loss | AUROC SVHN |\n|-|-|-|-|\n| DUQ   | 0.929    | 0.04 | 0.83|\n| BLNNconv | 0.930    | 0.04 | **0.89**   |\n \n \n> W3+Q3\n \nThank you very much for this advice. Here, we would like to focus on the inverse Lipschitz and Lipschitz hyperparameters ($\\alpha$ and $\\beta$). Note that other hyperparameters such as batchsize, depth…, follow the basic properties of the core neural network ICNN, and therefore, not directly related to the goal of this work (see e.g., Schalbetter, Adrian. *Input convex neural networks for energy optimization in an occupied apartment*. MS thesis. ETH Zurich, 2020.).\n\nWe would like to start the discussion by reminding that there is already an analysis with different $\\alpha$ and $\\beta$ for uncertainty estimation in the Appendix Table 5 p.50. We discuss in detail how the value of $\\alpha$ and $\\beta$ influence the performance of the model. The influence is rather intuitive, and we would like to verify this with additional experiments with a downscaled FashionMNIST for further clarity. The observations are as follows:\n1. Changes in performance due to hyperparameter changes are continuous (no hyper-sensitivity).\n|   $\\alpha$   | $\\beta$   | Accuracy | Loss| AUROC MNIST| AUROC NotMNIST|\n|-|-|-|-|-|-|\n| 0.2  | 0.99    | 0.8600|0.08|0.84|0.95|\n| 0.2|1.0| 0.8630|0.08|0.84|0.95|\n\n2. An increase of the domain of search leads to better performance due to higher expressivity.\n|   $\\alpha$   | $\\beta$   | Accuracy | Loss|\n|-|-|-|-|\n|0.0|0.1| 0.5680|0.24|\n|0.0|0.2 | 0.7485|0.14|\n|0.0|0.3 |0.7825|0.11|\n\n3. But too loose smoothing leads to worse performance.\n|   $\\alpha$   | $\\beta$   | Accuracy | Loss| AUROC MNIST| AUROC NotMNIST|\n|-|-|-|-|-|-|\n|1.0|4.0| 0.1070|0.10|0.50|0.50|\n|1.0|1.0| 0.8455|0.11|0.81|0.98|\n\n4. Too low sensitivity leads to worse out-of-distribution detection.\n|   $\\alpha$   | $\\beta$   | AUROC NotMNIST|\n|-|-|-|\n|0.0|0.2 |0.42|\n|0.3|0.2 |0.83|\n\nMoreover, in Appendix G.3 we also investigated the influence of Lipschitz constant on convergence speed in many different Lipschitz architectures (Figure 18 p.50). While increasing the Lipschitz constant (i.e., decreasing the smoothness of the network) leads to slower convergence, this decrease in speed is the smallest for ours among the compared models. In that sense, our method is more stable with respect to these hyperparameters. \n\n> W4+Q4\n\nConcerning the first part of W4, our theoretical guarantee is almost valid for usual practical scenarios. For example, the setting of Theorem 3.5 agrees with our experiments.  \n\nAs for the potential extensions to other network architectures, we would like to first remind that our work is primarily foundational, providing a novel paradigm for the creation of bi-Lipschitz architectures and addressing some non-negligeable issues of the field. Therefore, this paper focuses on the fundamental design and properties of our model related to the control of bi-Lipschitzness and its theoretical analysis. We acknowledge that there are some applications outside the bi-Lipschitz control framework that we could not investigate. That is also why we have tried to provide an extensive discussion and analyses in the appendix and kept the formulation as general as possible to promote further extensions and other interpretations for future work. \n\nNow, in the context of DEQ and implicit models, the LFT can be indeed re-formulated as finding the solution $z$ of $z=x-\\nabla F(z) +z$ which corresponds to a DEQ. In that sense, our model can be regarded as a *bi-Lipschitz* DEQ. We believe this novel interpretation is interesting for future work to increase the generality of our model.\n\nFurthermore, as mentioned in the first answer, we can extend our BLNN to convolutional layers as well. We just have to change the core neural network of the BLNN, which is the ICNN, to the convolutional ICNN proposed in the original work of Amos et al. (2017).\n\n> It seems that […] discussion on potential negative societal impacts […] is still missing.\n\nWe will ensure that a more detailed discussion is included in the revised version.\n\nWe hope we addressed all your questions and concerns. We will be glad to provide further explanation and clarification if necessary. \nMoreover, clarity is one of our major concerns. If you have any recommendations or some parts of the paper that were difficult to understand, we will be glad to improve them for the next version.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "QCWTSBhfci", "reviewer": "Reviewer_Hqwz", "summary": "This paper investigates and proposes a novel bi-Lipschitz neural network architecture. This architecture provides a simple, direct and tight control of the Lipschitz and inverse Lipschitz constants through the use of two parameters, the ideal minimum, equipped with theoretical guarantees. To devise their architecture the authors exploit convex neural networks and the Legendre-Fenchel duality. The authors also propose a variant of their bi-Lipschitz architecture that is more scalable by exploiting partially input convex neural networks. Finally, the authors propose a set of experiments to showcase the utility of our model in concrete machine learning applications, namely, uncertainty estimation and monotone problem settings and show that it can improve previous methods.", "questions": "NA.", "limitations": "Lack of experiments wrt to computational cost and maybe an experiment with a more larger dataset than the current ones use in the paper.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- The paper is well written. After writing a clearly structured related work (with an extensive background and related work proposed in Appendix A), the authors propose their new design and explicitly explain how the forward pass of their network is computed as well as the expressivity and how the backpropagation can be done. \n- The authors acknowledge that the computational cost of their approach can pose serious limitation and propose to overcome this problem with partially input convex neural networks.", "weaknesses": "- It would be interesting of the authors could provide experiments with both their architectures with respect to computational cost, and highlight time of training etc.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you very much for spending your time on carefully reviewing our paper. We really appreciate all the advice you provide to improve our paper. Please find below answers to your questions. We also summarized the comparison of the time and space complexity of our models in an independent thread.\n\n> It would be interesting of the authors could provide experiments with both their architectures with respect to computational cost, and highlight time of training etc.\n \nThank you very much for this advice. Since a few reviewers were also interested in this topic, both theoretical and experimental discussions are summarized in the global rebuttal. We would appreciate it if you could verify this. This will be added in the updated version of our paper. Please note that we did not directly compare training times in seconds but instead compared the number of floating-point operations (FLOPs) for each architecture, as training time is heavily dependent on factors such as the machine, code, and libraries used.\n \n> maybe an experiment with a more larger dataset than the current ones use in the paper.\n \nWe conducted an additional experiment using the convolutional version of our model (BLNNconv) on the problem of uncertainty estimation (as in Subsection 4.2) with the CIFAR-10 vs. SVHN dataset to illustrate the scalability of our method. For this problem, we implemented the model so that we first process the data through BLNNconv and then transfer it to the DUQ. The result compared to DUQ is as below. Our model is not only scalable to large-scale networks but also improves out-of-detection performance (the AUROC of SVHN). Using BLNNconv instead of the fully connected BLNN also improved the computation time (e.g., 2.5 times faster for the 5 first iterations).\n \n|      | Accuracy  | Loss | AUROC SVHN |\n|-------|-----------|----------|-----------|\n| DUQ   | 0.929    | 0.04 | 0.83|\n| BLNNconv | 0.930    | 0.04 | **0.89**   |\n\nWe hope we addressed all your questions and concerns. We will be glad to provide further explanation and clarification if necessary.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
# Time and Space Complexity of the BLNN and its Variants
 
This global rebuttal discusses the time and space complexity of the BLNN and its variants. Figures can be found in the attached PDF. This discussion will be added to the updated version of the paper.
 
## Theoretical Discussion
 
Concerning the time complexity, a forward pass of BLNN with $T$ iterations for the LFT has a time complexity of $O(T)$, where we supposed that the computation of a neural network is $O(1)$. Based on Theorem 3.7, the backward pass will require $O(d_x^3+d_\theta)$ as we have to solve one linear system and compute one Hessian vector product consecutively. $d_x$ is the dimension of the input and $d_\theta$ is the number of parameters. Concerning the space complexity, using Theorem 3.7, the storage capacity is $O(d_x+d_\theta)$ independent of $T$. If we are in an over-parameterized regime with $d_\theta>>d_x^3$ then we can neglect the dependence on $d_x$. With the brute force method (i.e., without Theorem 3.7), the forward is at least of $\Omega(T(d_x+d_\theta))$ and the storage capacity is at least proportional to $T$. Therefore, if we are in an over-parameterized regime with large iteration steps $T$, then the approach of Theorem 3.7 becomes more scalable in both time and memory.
 
## Experiments
 
We also computed with experiments the time and space complexity for one iteration of BLNN and its variants: BLNN (with brute force backpropagation), BLNN with Theorem 3.7, PBLNN with only one variable constrained to be bi-Lipschitz, and a traditional feed forward network as a control. Respective parameter sizes are 1.42M, 1.42M, 1.11M and 1.38M. The input size was 3x32x32, simulating a CIFAR dataset, and the batch size was varied in the set {1, 5, 25, 75, 100}.  Results are shown in the attached PDF.
 
As we can observe in Figures 1 and 2 of the PDF, while both our BLNN and BLNN with Theorem 3.7 present higher time and storage complexity than a traditional feed-forward neural network, Theorem 3.7 greatly contributes to reducing both computational and space requirements (improvements of order of $10$ and $10^2$, respectively). Comparing our model with Theorem 3.7 to a traditional feed-forward neural network, we can conclude that their difference in complexities is only a constant factor of order 10. This explains the scalability of our model to large datasets. Finally, PBLNN provides evidence that we can considerably decrease the complexity of the model by limiting the number of variables we impose bi-Lipschitzness.

Meta-review: This paper presents a novel bi-Lipschitz network approach. It is heavily theoretical with proof of concept ML experiments. Reviewers agreed that it is original and insightful. Productive rebuttal period yielded meaningful improvements and added essential information (e.g. time complexity). If accepted, please make sure to incorporate them.
Decision: Accept (poster)