ECLipsE: Efficient Compositional Lipschitz Constant Estimation for Deep Neural Networks

Authors: yuezhu  xu
s  sivaranjani
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 23:04:32

Abstract: The Lipschitz constant plays a crucial role in certifying the robustness of neural networks to input perturbations. Since calculating the exact Lipschitz constant is NP-hard, efforts have been made to obtain tight upper bounds on the Lipschitz constant. Typically, this involves solving a large matrix verification problem, the computational cost of which grows significantly for both deeper and wider networks. In this paper, we provide a compositional approach to estimate Lipschitz constants for deep feed-forward neural networks. We first obtain an exact decomposition of the large matrix verification problem into smaller sub-problems. Then, leveraging the underlying cascade structure of the network, we develop two algorithms. The first algorithm explores the geometric features of the problem and enables us to provide Lipschitz estimates that are comparable to existing methods by solving small semidefinite programs (SDPs) that are only as large as the size of each layer. The second algorithm relaxes these sub-problems and provides a closed-form solution to each sub-problem for extremely fast estimation, altogether eliminating the need to solve SDPs. The two algorithms represent different levels of trade-offs between efficiency and accuracy. Finally, we demonstrate that our approach provides a steep reduction in computation time (as much as several thousand times faster, depending on the algorithm for deeper networks) while yielding Lipschitz bounds that are very close to or even better than those achieved by state-of-the-art approaches in a broad range of experiments. In summary, our approach considerably advances the scalability and efficiency of certifying neural network robustness, making it particularly attractive for online learning tasks.

Main paper: 
 introduction: 
 The Lipschitz constant, which quantifies how a neural networks output varies in response to changes in its inputs, is a crucial measure in providing robustness certificates [1, 2] on downstream tasks such as ensuring resilience against adversarial attacks [3, 4] , stability of learning-based models or systems with neural network controllers [5] [6] [7] [8] [9] , enhancing generalizability [10] , improving gradientbased optimization methods and controlling the rate of learning [11] [12] . The problem of calculating the exact Lipschitz constant is NP-hard [13] . Therefore, efforts have been made to estimate tight upper bounds for the Lipschitz constant of feed-forward neural networks (FNNs) [14] [15] [16] [17] [18] and other architectures such as convolutional neural networks (CNNs) [19] [20] [21] . Typical approaches include formulating a polynomial optimization problem [22] or bounding the Lipschitz constant via quadratic constraints and semidefinite programming (SDP) [14] , which in turn requires solving a large-scale matrix verification problem whose computational complexity grows significantly with both the depth and width of the network. These approaches have also motivated the development of methods to design neural networks with certifiable robustness guarantees [19, [23] [24] [25] .Contribution. In this paper, we provide a scalable compositional approach to estimate Lipschitz constants for deep feed-forward neural networks. We demonstrate steep reductions in computation time (as much as several thousand times faster than the state-of-the-art depending on the experiment), while obtaining Lipschitz estimates that are very close to or even better than those achieved by state-of-the-art approaches. Specifically, we develop two algorithms, representing different levels in the trade-off between accuracy and efficiency, allowing for application-specific choices. The first algorithm, ECLipsE, involves estimating the Lipschitz constant through a compositional layer-bylayer solution of small SDPs that are only as large as the weight matrix in each layer. The second algorithm, ECLipsE-Fast, provides a closed-form solution to estimate the Lipschitz constant, completely eliminating the need to solve any matrix inequality SDPs. Both algorithms provably guarantee the existence of solutions at each step to generate tight Lipschitz estimates. In summary, our work significantly advances scalability and efficiency in certifying neural network robustness, making it applicable to a variety of online learning tasks.Theoretical Approach. We begin with the large matrix verification SDP for Lipschitz constant estimation under the well-known framework LipSDP [14] . To avoid handling a large matrix inequality, we employ a sequential Cholesky decomposition technique to obtain an exact decomposition of the large matrix verification problem into a series of smaller, more manageable sub-problems that are only as large as the size of the weight matrix in each layer. Then, observing the cascade structure of the neural network, we develop (i) algorithm ECLipsE, which characterizes the geometric features of the optimization problem and enables us to provide an accurate Lipschitz estimate and (ii) algorithm ECLipsE-Fast, which further relaxes the sub-problems, and yields a closed-form solution for each sub-problem that altogether eliminates the need to solve any SDPs, resulting in extremely fast implementations.Related Work. The simplest way to estimate the Lipschitz constans is to provide a naive upper bound using the product of induced weight norms, which is rather conservative [26] . Another approach is to utilize automatic differentiation to approximate a bound, which is not a strict upper bound, although it is often so in practice [13] . Additionally, compositions of nonexpansive averaged operators and affine operators [16] , Clarke Jacobian based approaches and other methods focusing on local Lipschitz constants [17] [27] have also been studied. Recently, optimization-based approaches such as sparse polynomial optimization [22] and SDP methods such as the canonical LipSDP framework [14] have been successful in providing tighter Lipschitz bounds. SDP-based methods specifically exploit the slope-restrictedness of the activation functions to cast the problem of estimating a Lipschitz constant as a linear matrix verification problem. However, the computational cost of such methods explodes as the number of layers increases. A common strategy to address this is to ignore some coupling constraints among the neurons to reduce the number of decision variables, yielding a more scalable algorithm at the expense of estimation accuracy [14] . Another strategy is to exploit the sparsity of the SDP using graph-theoretic approaches to decompose it into smaller linear matrix inequalities (LMI) [15] [28] . Along similar lines, [21] and [29] employ a dissipativitybased method and dynamic convolutional partition respectively to derive layer-wise LMIs that are applicable to both FNNs and CNNs. Very recent developments also focus on enhancing the scalability of SDP-based implementations through eigenvalue optimization and memory improvement [20] , which are compatible with autodiff frameworks such as PyTorch and TensorFlow. 
problem formulation and background: 
 Notation. We define Z N = {1, . . . , N }, where N is a natural number excluding zero. A symmetric positive-definite matrix P ∈ R n×n is represented as P > 0 (and as P ≥ 0, if it is positive semidefinite). We denote the largest singular value or the spectral norm of matrix A by σ max (A). The set of positive semi-definite diagonal matrices is written as D + .We consider a feedforward neural network (FNN) of l layers with input z ∈ R d0 and output y ∈ R d l defined as y = f (z). The function f is recursively formulated with layers L i , i ∈ Z l , defined asEQUATION)where v (i) = W i z (i-1) +b i with W i and b i representing the weight and bias for layer L i respectively, and ϕ : R di → R di is a nonlinear activation function that acts element-wise on its argument. The last layer L l is termed the output layer. We denote the number of neurons in layerL i by d i , i ∈ Z l . Definition 1. A function f : R d0 → R d l is Lipschitz continuous on Z ⊆ R d0 if there exists a constant L > 0 such that ∥f (z 1 ) -f (z 2 )∥ 2 ≤ L∥z 1 -z 2 ∥ 2 , ∀z 1 , z 2 ∈ Z.The smallest positive L satisfying this inequality is termed the Lipschitz constant of the function f .Without loss of generality, we assume W i ̸ = 0, i ∈ Z l , as any weights being 0 will lead to the trivial case where the output corresponding to any input will remain the same after that layer. Our goal is to provide a scalable approach to give an efficient and accurate upper bound for the Lipschitz constant L > 0. Note that the proofs of all the theoretical results in this paper are included in Appendix A.We begin with a slope-restrictedness property satisfied by most activation functions, which is typically leveraged to to derive SDPs for Lipschitz certificates [14] .Assumption 1 (Slope-restrictedness). For the neural network defined in (1) , the activation functionϕ is slope-restricted in [α, β], α < β in the sense that ∀v 1 , v 2 ∈ R n , we have α(v 1 -v 2 ) ≤ ϕ(v 1 ) -ϕ(v 2 ) ≤ β(v 1 -v 2 ) element-wise. Consequently, we have that for ∀Λ ∈ D + , v 1 -v 2 ϕ(v 1 ) -ϕ(v 2 ) T pΛ -mΛ -mΛ Λ v 1 -v 2 ϕ(v 1 ) -ϕ(v 2 ) ≤ 0, p = αβ, m = (α + β)/2. (2)Now, we can obtain an upper bound for the Lipschitz constant as follows; this result is equivalent to the well-known LipSDP framework [14] . Theorem 1 (LipSDP). For the FNN (1) satisfying Assumption 1, if there exists F > 0 and positive diagonal matricesEQUATIONthen z(l) 2 -z (l) 1 2 ≤ 1/F z (0) 2 -z (0) 1 2, which provides a sufficient condition for the Lipschitz constant L to be upper bounded by 1/F . Remark 1. LipSDP provides three variants that tradeoff accuracy and efficiency, namely, LipSDP-Network, LipSDP-Neuron, and LipSDP-Layer, whose scalability increases sequentially at the expense of decreased accuracy. However, [30] provides a counterexample showing that the Lipschitz estimate from LipSDP-Network is not a strict upper bound; thus, only LipSDP-Neuron, and LipSDP-Layer are valid. Theorem 1 here directly corresponds to LipSDP-Neuron. If all Λ i , i ∈ Z l-1 in (3) are set to multiples of identity matrices, that is, λ i I, i ∈ Z l-1 , then it corresponds to LipSDP-Layer.Assumption 1 holds for all commonly used activation functions; for example, it holds with α = 0, β = 1, that is, p = 0, m = 1/2 for the ReLU, sigmoid, tanh, exponential linear functions. Therefore, we focus on this case in this work. 
methodology: 
 We now develop two fast compositional algorithms based on LipSDP-Layer and Lipschitz-Neuron respectively. Both algorithms are not only scalable and significantly faster, but also provide comparable estimates for the Lipschitz constant.We circumvent direct solution of the large matrix inequality in (3) , which becomes computationally prohibitive as the FNN (1) grows deeper. Instead, we develop a sequential block Cholesky decom-position method, akin to the technique introduced in [31] , also expanded in [32, 33] . We first restate Lemma 2 of [31] below. Theorem 2 (Restatement of Lemma 2 of [31] ). A symmetric block tri-diagonal matrix defined asEQUATION)is positive definite if and only ifX i > 0, ∀i ∈ {0} ∪ Z l-1, whereX i = P i if i = 0, P i -R T i X -1 i-1 R i if i ∈ Z l-1 .(EQUATIONTheorem 3. Let P l be defined as in (3) with p = 0, m = 1/2. Then, the Lipschitz certificate P l > 0 holds if and only if the following sequence of matrix inequalities is satisfied:EQUATIONwhereEQUATION)Theorem 3 provides an exact decomposition of (3), and allows us to establish necessary and sufficient conditions through small matrix inequalities that scale with the size of the weight matrices of each layer, rather than that of the entire network. To accurately estimate the Lipschitz constant, we need to decide on Λ i , i ∈ Z 1-1 that generate a tight upper bound at the last stage. In other words, we want M l-1 -F W T l W l > 0 to yield the smallest estimate for 1/F . In the following subsection, we provide compositional algorithms to decide the appropriate Λ i , i ∈ Z 1-1 sequentially, so that we only need to solve one small problem corresponding to each layer.We first propose two practical algorithms here. The theory supporting the algorithms and the geometric intuition are deliberately deferred, and will be thoroughly discussed in a the next subsection.The first algorithm, ECLipsE, explores the geometric features that enables us to provide an accurate Lipschitz estimate by solving small semidefinite programs (SDPs), which are of the size of the weight matrices on each layer. The second algorithm, ECLipsE-Fast relaxes the sub-problems at each stage and yields a closed-form solution for each sub-problem that makes it extremely fast. These algorithms represent different trade-offs between efficiency and accuracy; one may choose ECLipsE if pursuing accuracy, and ECLipsE-Fast for applications where time is of the essence.We observe in (7) that M i is obtained in a recursive manner and depends on Λ i and M i-1 , i ∈ Z l-1 . Therefore, we decide Λ i and then calculate M i for i ∈ Z l-1 sequentially. Thus, these two algorithms can be implemented layer-by-layer in a compositional manner.Concretely, for ECLipsE, we obtain Λ i , i ∈ Z l-1 at each stage i using the information from the next layer, i.e. W i+1 , by solving the following small SDP:EQUATIONFor ECLipsE-Fast, Λ i is reduced to λ i I, i ∈ Z l-1 and λ i is calculated in closed-form asEQUATION)Note that this completely eliminates the need to solve matrix inequality SDPs altogether. At last, after all Λ i s, i ∈ Z l-1 are decided, we obtain the smallest 1/F , which yields the smallest Lipschitz estimate L = 1/F , as followsEQUATION)Remark 2. We choose to directly calculate the smallest 1/F rather than first derive the largest F . This is because obtaining the largest F first involves taking the inverse of W T l W l , which can cause numerical issues due to potential singularity of W T l W l . In contrast, directly calculating the smallest 1/F involves taking the inverses of M l-1 , which is already guaranteed to be strictly positive definite at layer l -1 when deciding Λ l-1 .We summarize the algorithms as one in Algorithm 1. Algorithms ECLipsE and ECLipsE-Fast are respectively preferable based on whether the priority is on accuracy or speed. Obtain Λ i from the optimal solution of (8)5:else if ECLipsE-Fast (pursuing speed) then 6:Obtain λ i from (9)7:Λ i ← λ i I 8:end if9:Obtain M i from ( 7) with Λ i and M i-1 10: end for 11: Obtain 1/F from (10) 12: Return L = 1/FNow we dive into the cascade structure of feed-forward neural networks and demonstrate the theory behind the two algorithms. We analyze the compositional algorithms in Section 3.2 in a backward manner, starting with the output layer. After all Λ i , i ∈ Z l-1 are decided, M i > 0, i ∈ Z l-2 hold. From Theorem 3, it remains to guarantee that M l-1 -F W T l W l > 0, and consequently, (10) , for which we state the following result. Proposition 1. For given Λ i , i ∈ Z l-1 that satisfies M i > 0, i ∈ Z l-2 , the tightest upper bound forLipschitz constant is L = σ max W T l W l (M l-1 ) -1 .Now, at stage l -1, when deciding Λ l-1 , Λ i , i ∈ Z l-2 are fixed and thus M l-2 is fixed. According to Proposition 1, we would like to choose Λ l-1 such that σ max W T l W l (M l-1 ) -1 , where M l-1 is a function of Λ l-1 , is as small as possible. We have the following result. Lemma 1.If M i > 0, then W T i+1 W i+1 (M i ) -1 and W i+1 (M i ) -1 (W i+1 ) T share the same non-zero eigenvalues. Note that at stage i, it is guaranteed that M i > 0. Taking i = l -1, Lemma 1 infers that it is equivalent to minimize σ max W l (M l-1 ) -1 W T lwhen deciding on Λ l-1 . Note that M l-1 > 0, and consequently, the existence of M -1 l-1 is already guaranteed when we reach the last stage. For the sake of conciseness, we defineF i ≜ W i (M i-1 ) -1 W T i i ∈ Z l-1 . From (7), M i = Λ i -1 4 Λ i F i Λ i .We further write out the recursive expression for F i asEQUATION)Lemma 2. For any constant γ ∈ (0, 1), anyΛ i ∈ D + that satisfies M i = Λ i -1 4 Λ i F i Λ i > 0 is also a feasible solution for Mi ≜ Λ i -1 4 Λ i (γF i )Λ i > 0.In other words, the feasible region{Λ i : M i > 0, Λ i ∈ D + } ⊆ {Λ i : Mi > 0, Λ i ∈ D + }.Lemma 2 gives us the observation that a contraction F i → γF i , γ ∈ (0, 1) yields a larger feasible space for Λ i ∈ D + to ensure M i > 0. Meanwhile, (11) shows that for any given Λ i , a smaller F i leads to a smaller F i+1 for the next stage. We can characterize how 'small' F i is by its spectral norm σ max (F i ). Then, minimizing σ max (F i ) aligns with our goal of minimizingσ max W T l W l (M l-1 ) -1 = σ max W l (M l-1 ) -1 W T l = σ max (F l )at the last stage. In other words, a smaller F 1 at the start will generally translate to a tighter Lipschitz estimate at output layer if we always choose to minimize the spectral norm σ max (F i ) at each stage. Now we focus on how to specifically optimize Λ i , i ∈ Z l-1 . At stage i, the goal is to seek for the Λ i that minimizes σ max (F i+1 ), where (11) . Note that M i-1 and F i are already fixed and can be regarded as constants at the i-th stage.F i+1 = W i+1 (Λ i -1 4 Λ i F i Λ i ) -1 W T i+1 as inProposition 2. If there exists a singular matrix N ≥ 0 such thatM i = c i W T i+1 W i+1 + N , with constant c i > 0, then σ max (F i+1 ) = 1/c i , ∀i ∈ Z l-1 .In other words, we need to find the largestc i > 0 to minimize σ max (F i+1 ) = 1/c i . Recall that M i = Λ i -1 4 Λ i F i Λ i is a function of Λ i .We state the following proposition that is used to derive the small sub-problems at each stage. Proposition 3. Consider the following optimization problem for ∀i ∈ Z l-1 .max ci c i s.t. Λ i - 1 4 Λ i W i (M i-1 ) -1 W T i Λ i -c i (W T i+1 W i+1 ) > 0, Λ i ∈ D + , c i > 0 (12)Then, the optimal value c i is the largest constant such that M i can be written asM i = c i W T i+1 W i+1 + N ,where N is some singular matrix such that N ≥ 0. Moreover, the feasible region for the optimization problem is always nonempty.We illustrate the process of achieving the largest c i > 0 in Fig. 1 . We geometrically represent a positive semidefinite matrix by the ellipsoid generated by the transformation of a unit ball in the Euclidean space by the matrix. For simplicity of exposition, we refer to this ellipsoid as the 'shape' of the matrix. We plot the shapes of M i and W T i+1 W i+1 in green and blue, respectively, in 2D. The positive definiteness of the constraint in ( 12) is equivalent to the ellipsoid of W T i+1 W i+1 being contained in the ellipsoid corresponding to M i /c i . Specifically, when c i > 1, Fig. 1a demonstrates the maximum contraction of M i , corresponding to the largest c i , such that ellipsoid of W T i+1 W i+1 is still contained in ellipsoid of c i M i . Similarly, for the case where c i < 1, Fig. 1b demonstrates the minimum extent (the smallest 1/c i ) to which M i needs to expand, such that the ellipsoid of W T i+1 W i+1 is contained. Algebraically, in both cases, c i is the ratio of the lengths of the green and pink arrows. By Proposition 2, the resulting ellipsoid (depicted in pink) is M i /c i = W T i+1 W i+1 + N/c i for both cases, and is tangent to the ellipsoid of W T i+1 W i+1 . Moreover, the vector pointing from the origin to the tangency point aligns with the direction of eigenvectors (the grey vector v in the plots) corresponding to the zero eigenvalues of the singular matrix N ≥ 0. Combining Proposition 2 and 3, we can derive an optimization problem to sequentially find the appropriate Λ i , i ∈ Z l-1 . The first constraint in (12) is quadratic in Λ i , which makes it unattractive for practical purposes. Therefore, we apply the Schur Complement to transform it into the linear matrix inequality (LMI) constraint in (8) . Thus, the optimization problem in Proposition 3 becomes equivalent to the SDP in (8), yielding algorithm ECLipsE. Notice that there are several ways to write the Schur complement of the constraint in (12) . We choose this specific structure to avoid singularity of the diagonal entries and ensure positive definiteness.ECLipsE-Fast achieves remarkable speed by further reducing Λ i , i ∈ Z l-1 to a multiple of identity matrix λ i I, where λ i > 0, and by relaxing the sub-problems. While our goal remains to minimizeσ max (F i+1 ) = σ max W i+1 (Λ i -1 4 Λ i F i Λ i ) -1 W T i+1, we intentionally disregard information from W i+1 , and instead focus solely on minimizing the spectral norm of(Λ i -1 4 Λ i F i Λ i ) -1 .Roughly speaking, a smallerσ max (Λ i -1 4 Λ i F i Λ i ) -1 yields a smaller σ max (F i+1). This relaxation allows us to derive a closed-form solution for Λ i , i ∈ Z l-1 as follows. Proposition 4.Choosing λ i = 2 σmax(Fi) > 0 minimizes σ max (Λ i -1 4 Λ i F i Λ i ) -1 where Λ i = λ i I under the constraint that M i = Λ i -1 4 Λ i F i Λ i > 0. Moreover, this closed-form solution for λ i always satisfies M i > 0, i ∈ Z l-1 .By the definition of F i , Proposition 4 matches with (9) , yielding algorithm ECLipsE-Fast. Although this relaxation may result in a loss of tightness, the closed-form solution offers the advantage of significantly increased computational speed.We now demonstrate the geometric analysis behind the development of ECLipsE-Fast and compare it with ECLipsE in the case where c i > 1 (Fig. 2 ). We also include the case c i < 1 in Appendix B. The key idea behind ECLipsE-Fast is that instead of keeping the shape of M i fixed, and contracting the ellipsoid itself, as in ECLipsE, we first find the largest inscribed ball (dark green) for the ellipsoid of M i . Then, we contract this ball to the maximum extent such that it still contains W T i+1 W i+1 . The resulting ball (dark blue) is precisely the smallest circumscribing ball for the ellipsoid of W T i+1 W i+1 . Note that this approach serves as an approximation for the process of contraction depicted in Fig. 7b (corresponding to ECLipsE), thus yielding a smaller c i . We use this approximation to achieve a closed-form solution, which significantly increases the computational speed. Remark 3. In Lemma 2, the analysis initially fixes the shape of F i . However, when optimizing Λ i , the shape of the feasible region depends on F i , which can vary with different Λ i-1 , i ∈ Z l . Thus, this approximation, which allows for a scalable distributed algorithm to solve the centralized problem (3) introduces an unavoidable but minor tradeoff in achieving global optimality. 
experiments: 
 We implement our algorithms † on randomly generated neural networks and ones trained on the MNIST dataset. The details of the experimental setup, and training of the neural networks (both randomly generated and trained on the MNIST dataset) are described in Appendix D.Baselines. ‡ For ECLipsE, Λ i , i ∈ Z l-1 can have different diagonal entries, which benchmarks to LipSDP-Neuron. For ECLipsE-Fast, Λ i = λ i I, i ∈ Z l-1 , which benchmarks to LipSDP-Layer. Additionally, we compare our Lipschitz estimates to the naive upper bound [26] , CPLip [16] and LipDiff [20] . The codes for these baselines are available at [34, 35, 20] . Note that LipDiff is accelerated using a node with 2 NVIDIA A100 GPUs (80G) and 512 GB of memory.L naive = l i=1 ∥W i ∥ 2We first consider randomly generated networks, where the number of layers are chosen from {2, 5, 10, 20, 30, 50, 75, 100}, and number of neurons are chosen from {20, 40, 60, 80, 100}, amounting to a total of 40 experiments for each algorithm (including the baselines). We quantify the computation time and tightness of the Lipschitz bounds (raw data in Appendix E). The Lipschitz bounds presented in the following figures are normalized to the trivial upper bound for ease of comparison. † https://github.com/YuezhuXu/ECLipsE ‡ Note that SeqLip [13] is also an often-used benchmark; however, we do not consider it since it does not represent a true upper bound for the Lipschitz constant. We also note that we do not include Chordal-LipSDP [15] as a baseline , since only the case where τ = 0 in that work is valid, and all other cases, are no longer valid in certifying the Lipschitz constant as discussed in Remark 1 as well as [15] .Case 1: Varying network depth (number of layers). We select a network with 80 neurons per layer, and demonstrate the scalability of our algorithm as network depth increases. Note that all baseline approaches fail to provide a Lipschitz estimate within a computational cutoff time of 15 min for networks larger than this size (see results in Appendix E). As the number of layers increases, the computation time for CPLip algorithm explodes (the algorithm does not return a Lipschitz estimate within the cutoff time beyond 20 layers); however, CPLip provides the most accurate estimates in smaller networks. LipDiff provides inadmissible Lipschitz estimates even for moderate networks, returning as much as 10-100 times the trivial bound (see Table 2a , Appendix E for the estimates). Also, while LipDiff has similar computational time for smaller networks, computational time grows for deeper networks as recorded in Appendix E Table 2b . Consequently, we do not include these results in the plots. LipSDP-Neuron and LipSDP-Layer are also scalable to some extent; however, they fail for a networks of 30 and 50 layers respectively. In contrast, the computation time for ECLipsE and ECLipsE-Fast stays low and grows only linearly with respect to the number of layers (Fig. 3b ). Notably, ECLipsE-Fast is significantly faster (thousands of times) than LipSDP-Layer, owing to the closed-form solution at each stage, while ECLipsE is also considerably faster than LipSDP-Neuron. The Lipschitz estimates given by algorithms ECLipsE and ECLipsE-Fast are very close to the ones from LipSDP-Neuron and LipSDP-Layer respectively (Fig. 3a ), and outperform the trivial bound. As the number of layers increases, the normalized Lipschitz estimates are smaller, indicating that our algorithms are well-suited to very deep networks. Case 2: Varying neural network width (number of neurons per layer). We now examine the performance of our algorithms for wider (more hidden neurons per layer), rather than deeper networks (with more layers), and demonstrate the results for networks with 20 and 50 layers respectively (Fig. 4 ). While the complete raw data is presented in Appendix E, we discuss the results for 20 and 50 layer networks here, since they represent the network sizes where different baselines fail to return Lipschitz estimates beyond the computation cutoff time of 15 min. Note that while LipDiff also manages to generate estimates for all network sizes in our 50 layers case, it once again provides inadmissible Lipschitz constants, returning as much as 10 4 -10 6 times the trivial bound. Therefore, we do not include these results in Fig. 4 (see Tables 3a and 3b in Appendix E for the estimates and computation time.) We can observe from Figs. 4b and 4d that the computation time needed for CPLip, LipSDP-Layer, and LipSDP-Neuron significantly increases with the number of neurons, while the computation time of our method still grows linearly. Meanwhile, the Lipschitz estimates from algorithms ECLipsE and ECLipsE-Fast are close to the ones from LipSDP-Neuron and LipSDP-Layer respectively (Figs. 4a and 4c ). Thus, we can conclude that our method significantly improves scalability for wider neural networks.Case 3: Comparison with LipSDP implementations. In order to address the scalability issue as the size of the network grows, LipSDP utilizes a splitting approach, where the network is split into smaller sub-networks and the Lipschitz constants for each sub-network are composed at the end to obtain the final estimate. We benchmark our approach with respect to the performance of LipSDP-Layer and LipSDP-Neuron considering different sub-network sizes. Note that our algorithms do not require any splitting, since they remain scalable to large networks. As the FNNs are larger than the ones in previous cases, we change the cutoff time to 30 minutes. We conduct two sets of experiments Fig. 5 . By plotting the normalized Lipschitz estimates and computation times on the two axes, we illustrate how efficient and accurate an algorithm is by how close the corresponding data point is to the origin. We observe that all the data points for ECLipsE-Fast are at the leftmost extreme of the plot, indicating that it is the most efficient algorithm. Further, ECLipsE-Fast also outperforms the red cluster (LipSDP-layer with the network split into 3) in both tightness and speed. Comparing data points of the same shape, ECLipsE-Fast outperforms LipSDP-Layer for all sub-network splits both in terms of the Lipschitz estimate and the computation time. Finally, the data points corresponding to ECLipsE are clustered at the bottom left, demonstrating that it is relatively more accurate and efficient than all LipSDP methods, no matter how the network is split.In the second set of experiments, we explore even wider networks. Specifically, we choose a fairly deep neural network with 50 layers and vary the width from 150 to 1000. The splitting size for LipSDP-Neuron and LipSDP-Layer is 5. The resulting Lipschitz estimates (normalized with respect to trivial upper bounds) and the computation time are provided in Tables 4a and 4b of Appendix E due to space limitations. From these results, we observe that ECLipsE-Fast is extremely fast even for very wide networks, with a running time of only 15.63 seconds for a network width of 1000, while the computation time for LipSDP-Layer grows significantly. Also, while ECLipsE fails when the width reaches 300, it is comparable to LipSDP-Neuron split into 5 sub-networks in terms of time performance. Remark 4. We notice that when the neural networks are significantly wide, ECLipsE takes more than 30 minutes while ECLipsE-Fast remains efficient. This observation can be explained by examining the computational complexity of these algorithms. Note that we directly state the computational complexity of each algorithm here for brevity; the detailed derivations are included in Appendix C. Suppose a neural network has n hidden layers with m neurons. Then, the computational cost for LipSDP and ECLipsE are O(n 4 m 4 ) and O(nm 4 ) respectively. We can observe that the complexity is significantly decreased in terms of the depth, but is the same in terms of the width, immediately indicating the advantage for deep networks. Nevertheless, as m grows, the difference between O(n 4 m 4 ) and O(nm 4 ) is still drastically enhanced, especially with large n. More importantly, for ECLipsE-Fast, the computational cost drops to O(nm 3 ). This is the fastest one can expect if the weights on each layer are treated as a whole.We now demonstrate our algorithms on four networks trained on the MNIST dataset (see Appendix D for details) to achieve an accuracy of at least 97%. The resulting networks are not very deep (3 layers), with 100, 200, 300, and 400 neurons. We set a computational cutoff time of 30 min to obtain Lipschitz estimates. As described in the note on Baselines earlier in this section, ECLipsE is benchmarked against LipSDP-Neuron and ECLipsE-Fast is benchmarked against the faster LipSDP-Layer due to their mathematical structure. From Fig. 6b , we can see that ECLipsE-Fast is significantly faster than LipSDP-Layer, while ECLipsE is also considerably faster than LipSDP-Neuron. Note that all algorithms provide very similar Lipschitz estimates (Fig. 6a ). Therefore, for networks that are not very deep, such as those in this example, ECLipsE-Fast is the optimal choice, since it significantly outperforms all algorithms in terms of speed, while the approximation error due to the closed-form solution is not too significant compared to the baselines.  
conclusion: 
 We propose a scalable approach to estimate Lipschitz constants for deep neural networks by developing a new matrix decomposition that yields two fast algorithms. Our experiments demonstrate that our algorithms significantly outperform the state-of-the-art in terms of computation speed, while providing comparable Lipschitz estimates. We envision that further computational speedup can be achieved through sparse matrix multiplication and eigenvalue estimation techniques, and leveraging autodiff frameworks, along the lines of [20] . While we can unroll the convolutional layers in CNN structure to a large fully connected neural network layer to apply ECLipsE and ECLipsE-Fast to estimate Lipschitz constant, better compositional methods that are tailored to feature the convolutional layers are expected for future work. Similarly, other architectures, such as residual networks, present additional challenges due to their unique structures and will be considered in future research. 
 
 Appendix: Proof of Relation in Assumption 1. We show how the slope-restrictedness of the activation function implies (2) . The inequalityα(v 1 -v 2 ) ≤ ϕ(v 1 ) -ϕ(v 2 ) ≤ β(v 1 -v 2 ) that holds elementwise yields λ i [(ϕ(v 1 ) -ϕ(v 2 )) -α(v 1 -v 2 )] i × [(ϕ(v 1 ) -ϕ(v 2 )) -β(v 1 -v 2 )] i ≤ 0, ∀i ∈ R n , ∀λ i ≥ 0,where subscript i indexes the i th element of the vector. Summing all the inequalities for i ∈ R n and letting Λ = diag(λ 1 , λ 2 , ..., λ n ), we have[(ϕ(v 1 )-ϕ(v 2 ))-α(v 1 -v 2 )] T Λ [(ϕ(v 1 )-ϕ(v 2 ))-β(v 1 -v 2 )] ≤0.In other words, we have the following quadratic inequality(ϕ(v 1 )-ϕ(v 2 )) T Λ(ϕ(v 1 )-ϕ(v 2 ))- α + β 2 (ϕ(v 1 )-ϕ(v 2 )) T Λ(v 1 -v 2 ) - α + β 2 (v 1 -v 2 ) T Λ(ϕ(v 1 )-ϕ(v 2 ))+αβ(v 1 -v 2 ) T Λ(v 1 -v 2 ) ≤ 0,which can be directly rewritten as (2) with p = αβ and m = α+β 2 . Proof of Theorem 1. For any two inputs z(0) 1 and z (0) 2 , let z (i) 1 , z (i) 2 , v (i) 1 , v (i) 2 , i ∈ Z l be computed as in (1). Define ∆z (i) = z (i) 1 -z (i) 2 , i ∈ {0} ∪ Z l and ∆v (j) = v (j) 1 -v (j)2 , j ∈ Z l . By both left and right multiplying the left hand side in (3) by the vector [∆z (0) , ∆z (1) , ..., ∆z (l-1) ] T , we haveEQUATIONNow, we show that every summand in the above inequality ( 13) is negative semidefinite. In fact, with Assumption 1, using notation∆v (i) = v (i) 1 -v (i)2 and ∆ϕ(v(i) ) = ϕ(v (i) 1 ) -ϕ(v (i)2 ), i ∈ Z l , and taking Λ = Λ i , we can writeEQUATION)Note that ∆v (i) = (W i z 1) and ∆ϕ(v i ) = ∆z i . We can express these relationships in matrix form asEQUATION+ b i ) = W i ∆z (i-EQUATION)Substituting ( 15) into ( 14), we haveEQUATION)Combining ( 13) and ( 16), we haveEQUATION)Similarly, as the last layer is a linear layer, ∆z (l) = ∆v (l) = W l ∆z (l-1) . Then ( 17) is exactly(∆z l ) T ∆z l ≤ 1 F (∆z (0) ) T ∆z (0) ,yielding a upper bound 1/F for the Lipschitz constant.Proof of Theorem 2. Applying Lemma 2 in [31] , we defineEQUATION)andEQUATIONThen, with P i and R i defined above, we directly haveEQUATIONIn other words, the sufficient and necessary conditionEQUATIONwhich is the same as (21) .Proof of Proposition 1. As M i > 0, i ∈ Z l-2 has been guaranteed, it remains to ensure thatM l-1 -F W T l W l > 0 by Theorem 3. This is equivalent to M l-1 /F > W T l W l . Therefore, the smallest possible 1/F is σ max (W T l W l (M l-1 ) -1). Then by Theorem 1, the upper bound for the Lipschitz constant isEQUATIONmeaning that all eigenvalues of W i+1 (M i ) -1 W T i+1 are non-negative. As W i+1 ̸ = 0, we know that the largest eigenvalue is positive and should be the same as σ max W i+1 (M i ) -1 W T i+1 > 0. Now consider any non-zero eigenvalue λ a of matrix W i+1 (M i ) -1 W T i+1 and let v a ̸ = 0 be its corresponding eigenvector. Then,EQUATION)Left multiplying both sides with W T i+1 , we haveEQUATION)As λ a ̸ = 0, we know that W T i+1 v a ̸ = 0 from (23) . (Otherwise, λ a v a = 0 with v a ̸ = 0 will lead to λ a = 0). With W T i+1 v a ̸ = 0, (24) implies that λ a is also an eigenvalue ofEQUATION) Let v c = 1 λ b W i+1 (M i ) -1 v b . We have v b = W T i+1 v c .Then we substitute v b on the both sides in ( 25) and obtainW T i+1 W i+1 (M i ) -1 W T i+1 v c = λ b W T i+1 v c . Left multiplying both side with W i+1 (M i ) -1 , we further have W i+1 (M i ) -1 W T i+1 W i+1 (M i ) -1 W T i+1 v c = λ b W i+1 (M i ) -1 W T i+1 v c . Let v d = W i+1 (M i ) -1 W T i+1 v c , we finally have W i+1 (M i ) -1 W T i+1 v d = λ b v d .Similarly, we can conclude that λ b ̸ = 0 is the eigenvalue of W i+1 (M i ) -1 W T i+1 and v d is the corresponding eigenvector.EQUATION)At stage i, M i-1 > 0 is guaranteed. Then, similar to ( 22), we know thatF i = W i (M i-1 ) -1 W T i ≥ 0 Since Λ i ∈ D + and 0 < γ < 1, implying 1-γ 4 > 0, we have 1-γ 4 Λ i F i Λ i ≥ 0. Proof of Proposition 2. Recall that by definition, F i+1 = W i+1 (M i ) -1 W T i+1 . Applying Lemma 1, we have σ max (F i+1 ) = σ max (W T i+1 W i+1 (M i ) -1 ). Meanwhile, with M i = c i W T i+1 W i+1 + N , W T i+1 W i+1 (M i ) -1 = (W T i+1 W i+1 + N c i )(c i W T i+1 W i+1 + N ) -1 - N c i (c i W T i+1 W i+1 + N ) -1 = 1 c i I - N c i (M i ) -1 .(27) WIth N ≥ 0 and M i > 0 after Λ i is decided, we show that N ci (M i ) -1 only has non-negative eigenvalues. As M i > 0 is guaranteed to symmetric according to (7) , there exists a symmetric square root for (M i ) -1 and we denote it to be(M i ) -1 2 . Then N (M i ) -1 is similar to (M i ) -1 2 N (M i ) -1 2, thus sharing the same eigenvalues. Furthermore, for ∀x ̸ = 0,x T (M i ) -1 2 N (M i ) -1 2 x = (M i ) -1 2 x T N (M i ) -1 2 x ≥ 0. It indicates (M i ) -1 2 N (M i ) -1 2only has non-negative eigenvalues. The first equation holds by the symmetry of (M i ) -1 2 and the second inequality is because of the definition of positive semidefiniteness of N . Therefore, withc i > 0, the eigenvalues of N ci (M i ) -1 = 1 ci N (M i ) -1are all non-negative. and all eigenvalues of W T i+1 W i+1 (M i ) -1 should be less than or equal to 1/c i . On the other hand, as N is singular, N (M i ) -1 is also singular, thus having eigenvalue 0. Therefore,σ max (F i+1 ) = σ max (W T i+1 W i+1 (M i ) -1 ) = 1/c i .Proof of Proposition 3. We first show that if M i-1 > 0, the feasible region is always non-empty for ∀i ∈ Z l-1 . We use σ to denoteEQUATIONTherefore, the feasible region at least includes the Λ i and c i we specified, and is thus not empty.To make the feasibility complete, we prove thatM i > 0, ∀i ∈ {0} Z l-1 by induction. When i = 0, M i = M 0 = I > 0.When it comes to stage i, we have by induction that M i-1 > 0 is true. As Λ i are obtained satisfying (12) with c i > 0 and recall the recursive relation forM i , i ∈ Z l-1 to be M i = Λ i -1 4 Λ i W i (M i-1 ) -1 W T i Λ i , we have M i > c i W T i+1 W i+1 ≥ 0.We now prove by contradiction that the optimal value c i is the largest constant such that M i can be written asM i = c i W T i+1 W i+1 + N , where N is some singular matrix that N ≥ 0. M i = c i W T i+1 W i+1 + N .Suppose there exists a ĉi > c i such that it satisfies all the constraints. Then, from the first constraint in (12), we haveEQUATIONLet v ̸ = 0 be the eigenvector of N corresponding to eigenvalue 0. Observe thatEQUATIONThe last step holds because c i -ĉi < 0 by definition of ĉi . It contradicts (29) .Proof of Proposition 4.When Λ i = λ i I, Λ i - 1 4 Λ i F i Λ i -1 = λ i I - 1 4 λ 2 i F i -1.EQUATIONMinimizing the above spectrum is equivalent to maximizing the denominator λ i -λ 2 i σ max (F i )/4 in (31) , which is quadratic in λ i . To find the optimal λ i , we set the derivative of the denominator with respect to λ i to be 0, and obtain the closed-form solution λ i = 2 σmax(Fi) . Moreover, with Λ i = 2 σmax(Fi) I and F i > 0 guaranteed at stage i, we haveM i = Λ i - 1 4 Λ i F i Λ i = 1 σ max (F i ) I > 0.The geometric analysis for algorithm ECLipsE-Fast analogous to Fig. 2 , comparing the case where c i > 1 and in other cases are shown in Fig. 7 . We derive the time complexity for both ECLipsE and ECLipsE-Fast in detail here. Suppose a neural network has n hidden layers with m neurons. Then, the large matrix in Theorem 1 has dimension nm + O(1) and the decision variable is of size nm + O(1). Note that the computational complexity for solving an LMI with the size of the matrix constraint being size A and the number of decision variables being B is O(A 3 + A 2 B 2 ). Therefore, the computational cost for LipSDP isO((nm + O(1)) 3 + (nm + O(1)) 2 (nm + O(1)) 2 ) = O(n 4 m 4). Contrarily, ECLipsE solves n sub-problems as in Eq. ( 8), each involving a matrix of size O(m) and m decision variables. The corresponding total computational cost is n× (O(m 3 + m 2 m 2 )) = O(nm 4). This directly indicates the advantage of ECLipsE for deep networks. Also, as m grows, the difference between O(n 4 m 4 ) and O(nm 4 ) is still significantly enhanced, especially with large n. Regarding ECLipsE-Fast, we note that we do not need to solve any SDPs and the computational cost drops to n × O(m 3 ) = O(nm 3 ). This is the fastest one can expect if the weights on each layer are treated as a whole.Experimental Setup. All experiments are implemented on a Windows laptop with a 12-core CPU with 16GB of RAM.Randomly Generated Neural Networks. We set the input dimension to be 4 and the output dimension to be 1. The activation functionsare chosen to be ReLU, and the number of neurons in each hidden layer is set to be the same. We randomly generate weights for each layer to follow the normal distribution. Also, in order to avoid the case where the Lipschitz constant is too large or too small and may cause numerical issues, we scale the weights on each layer such that the is norm randomly chosen in [0.4, 1.8], following a uniform distribution.For training on the dataset MNIST, the input dimension is 784 and output dimension is 10, which is compatible with the dataset. he activation functionsare chosen to be ReLU, and the number of neurons in each hidden layer is set to be the same. We train neural networks using the SGD optimizer with a learning rate of 0.01 and momentum of 0.9 until they achieve at least 97% accuracy on test data.The Lipschitz constant estimates and computation times for randomly generated neural networks with the number of layers chosen from {2, 5, 10, 20, 30, 50, 75, 100}, and number of neurons are chosen from {20, 40, 60, 80, 100}, are provided below. This work is primarily theoretical and pertains to obtaining upper bounds on the Lipschitz constant, which can serve as a measure of the robustness of deep neural networks, and does not have any direct societal impact.(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.In the case of closed-source models, it 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 to have some path to reproducing or verifying the results. • 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., if there is a special consideration due to laws or regulations in their jurisdiction). 10. Broader Impacts Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?• 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 a 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, paperswithcode.com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a 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 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: All the code and data are submitted along with the paper, and will be released publicly with detailed documentation upon publication. Guidelines:• The answer NA means that the paper does not release new assets.• Researchers should communicate the details of the dataset/code/model as 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.  


Decision: This paper proposes a novel, efficient algorithm to estimate the Lipschitz constant of neural networks. The algorithm is based on SDP based methods, developed to estimate the Lipschitz constant. The main idea is breaking the SDP into sequential subproblems which involve matrix inequality on smaller matrices.  Using this decomposition, the paper proposes two algorithms with different computational cost and accuracy, allowing to establish cost-accuracy trade-offs. In particular, the proposed method significantly reduces the computational complexity with depth. Experiments show that this proposed method can achieve almost the same estimate while being scalable to deeper networks. The paper is well-motivated, well-written with a solid contribution. While the proposed method is limited to feedforward networks, this limitation is included in the discussion. Thus, I recommend to accept this paper for a spotlight presentation.

Reviews: review: paper_summary: paper_summary: The authors are able to decompose a particular case of LipSDP-neuron exactly into a series of sub-problems leading to the proposed algorithm ECLipSE. In the case of the relaxed LipSDP-layer, it can be shown that each sub-problem can be solved analytically and eliminate the need for solving an SDP. The proposed algorithm ECLipSE-fast can provide Lipschitz estimates of deep-NN quite fast at the expense of increased conservativeness.
strength_weakness: strength_weakness: -It seems that the approach doesn’t yet apply to residual networks or CNN commonly found in state-of-the-art vision models. For this reason, it seems difficult to show improved certified robustness on common image classification benchmarks which would benefit the most from increased scalability (CIFAR10-100, Tiny-imagenet, etc).

-Certified robustness is mentioned as an application at several points in the paper. While the scalable and tighter Lipschitz estimates on random networks and MNIST certainly suggest some improvements, no practical measures of certified robustness are presented (e.g. robust accuracy).  For certified robust accuracy on MNIST, it is common to use 1-Lipschitz parameterization (SLL, AOL, other direct parameterizations) which eliminates the need for Lipschitz estimation. Would applying ECLipSE to composed 1-Lipschitz networks provide significantly tighter estimates in this case?
questions: questions: -How does the full row rank assumption limit the applications of ECLipSE’s? Is it possible to relax this assumption? It seems like networks considered in this paper are only of constant width or decreasing in output dimension in the case of MNIST which will usually satisfy the full row-rank assumption.

-Can you comment on the slight gap between the Lipschitz constant estimates shown in Appendix E between ECLipsE and LipSDP-Neuron? Is it that LipSDP is not as accurate in larger settings or is ECLipsE somehow slightly conservative?
limitations: limitations: Limitations are adequately addressed.
paper_summary: paper_summary: The paper proposes two novel Lipschitz constant estimation algorithms ECLipsE and ECLipsE-Fast. They are supported by a new decomposition theory developed for the Lip-SDP framework, derived by applying an existing theory (Lemma 2 of [31]). Experiments demonstrate the estimation accuracy and acceleration using toy radom networks and networks trained on MNIST data, by comparing with classical Lipschitz constant estimation method.
strength_weakness: strength_weakness: (1) The proposed algorithm is efficient at addressing network depth, but does not look at the width.  The theory behind self-explains the success of its capability of handling depth.  However, layers with very high numbers of neurons still pose challenges for the sub-problems, e.g., solving for Eq. (6) and Eq. (7). There is a lack of mentioning of this. Also the experiments only studied some modest widths up to 100 neurons. It would be good to see more empirical results with higher neuron numbers in each layer,  to understand the limit of the proposed algorithms on network width.

(2) In experiments, there seems a missing comparison with the “parallel implementation by splitting” version of Lip-SDP as reported in their paper [14], which was proposed to address the depth issue.

(3) The  method description can be improved, e.g., being more organised. The paper can present information that is more important and helpful to practitioners and general readers’ understanding in main paper, while leave some analysis and supporting lemmas to appendix. Personally, I find it helpful to see in the main paper a description of the existing Lemma 2 in [31] and the pseudo code of the proposed algorithm.
questions: questions: The authors are invited to address my comments (1) and (2) in weakness section. If possible, it would be good to see some added results to help understand more the width capacity of the proposed algorithms.
limitations: limitations: Discussion on discussion is pretty limited. For instance, it  can be improved around the network width issue.
paper_summary: paper_summary: The paper proposes two algorithms, ECLipsE and ECLipsE-Fast, to estimate the Lipschitz constant of a feed-forward neural network. The estimation of the Lipschitzness plays a crucial role in certifying the robustness of neural networks and is known to be an NP-hard problem. The proposed algorithms are based on the LipSDP of Fazlyab et al. (2019), which describes the semidefinite program (SDP) that an upper bound of such a Lipschitz constant should generally satisfy. The authors decompose the original large SDP into smaller layer-wise SDPs to improve the scalability of the original approach. The validity of the resulting methods, ECLipsE and ECLipsE-Fast, is supported by theoretical analyses. Experiments show a steep reduction in computation time while maintaining competitive accuracy compared to the LipSDP.
strength_weakness: strength_weakness: Overall, in my view, the main weakness of this paper is that there is no explicit comparison with other works trying to improve the scalability of LipSDP (such as [20]). As a result, while the paper indeed improves the original LipSDP in a new way, it is unclear how significant this work is taking into account existing literature. In addition, the limitations should be more carefully discussed. See below for further detailed comments and advice.

### **Quality**
 
(W-Q1) There are some typos: l.60 (constans), p.3 eqation (3) (index of the bottom right element is i+1 but should be l), l.172 (in the matrix WMW there are 2 unnecessary “L”), l.264 (computatinoal), between l.420 and l.421 (i\in \mathbb{R}^n), p.13 equation (18) (if i=0 should be if i=1?), l.477 (functionsare), l.480 (the is norm), 

(W-Q2) I feel that some use of words is misleading.  

1. l.60 “We develop a sequential Cholesky decomposition technique to obtain […]“: If this is about Theorem 2, it directly uses the result of Agarwal et al., so maybe “use, employ” would be better than “develop”.

2. l.10 “The first algorithm [...] enables us to provide a *tight* Lipschitz constant”, l.56 “ algorithm [...] enables us to provide an *accurate* Lipschitz”... What do you mean by “tight” and “accurate”? Since the proposed algorithms are using some simplifications, I think that those adjectives should be relative, i.e., only used in comparison with something else. Notably, the experiments (e.g., Figure 3) show that CPLip is far more precise than the proposed algorithms so the authors should clarify the meaning of “tight” and “accurate” (or delete them if there is no justification) when describing their own algorithms.

3. Between l.462 and l.463, you write “$N/ci(M_i)^{-1}\ge0$”. Semi-positive definiteness was only defined for symmetric matrices but this one is not necessarily symmetric. There should be an easy fix, but it is ambiguous what you mean by this sign. 

(W-Q3) In Proposition 3, the authors simplify Theorem 2 by setting $M_i$ to $c_iW_{i+1}^\top  W_{i+1}+N$ without any proper discussion about this choice (Q2).

(W-Q4) The proof of positive definiteness of $M_i$ is missing in Proposition 3 (Q3). 

(W-Q5) (minor) Some references should be adjusted: [2] was accepted at ICLR2018, some capital letters are missing (L of lipschitz)...

### **Clarity**

(W-C1) Some parts may require clarifications (See also Questions):

1. (minor) l.111 “[30] provides a counterexample to the most accurate approach”: The concrete property of the most accurate approach disproved by [30] could be explained in a few words.
2. (minor) A mathematical comparison of the computational complexity between your approach and LipSDP will largely help the reader to quickly understand the difference in scalability.

(W-C2) While limitations of proposed algorithms are all stated, they are dispersed throughout the paper. A short subsection summarizing them would be useful. (See also Limitations)

(W-C3) It was a bit difficult to understand the idea of the proposed algorithms from Subsection 3.2. Perhaps showing Algorithm 1 in the main text would be better.

### **Significance**
 
(W-S1) (major) One of the main contributions of this work is the scalability of the proposed algorithms. However, the authors do not compare their method with other accelerations of LipSDP. Therefore, it is difficult to situate their work within the broader context of efforts to improve the scalability of LipSDP. At least, comparing the algorithms with that of Wang et al. (2024) [20] is important to clarify these points. 

(W-S2) The algorithms were run on medium-scale neural networks, and it is difficult to imagine the scalability of ECLipsE and ECLipsE-Fast. Experiments with even larger architectures (for example, those for training on CIFAR or ImageNet) would be more convincing.  

(W-S3) I feel experiments of Subsection 4.2 are a little bit redundant (Q4). 

(W-S4) Limitations of the algorithms should be discussed in more detail and more explicitly.
questions: questions: Q1: Is Theorem 1 *necessary* and sufficient? If not, this should be added to the limitations of the work.

Q2: Why did the authors set $M_i$ to $c_iW_{i+1}^\top  W_{i+1}+N$ in Proposition 3?

Q3: How is the positive definiteness of $M_i$ guaranteed in Equation (6) and Proposition 3?

Q4: What was the motivation to run experiments of Subsection 4.2?
limitations: limitations: There are several limitations of the algorithms that are worth mentioning:

1. The algorithms (Lemma 1) need that the last weight matrix is full row rank. So, we cannot blindly apply them to any feed-forward neural network.

2. There is a simplification when transforming Theorem 2 into Proposition 3 by limiting the expression of  $M_i $ to $c_iW_{i+1}^\top  W_{i+1}+N$. This may lead to looser bounds than the original LipSDP.

3. The algorithms cannot be applied to CNNs and residual networks.
paper_summary: paper_summary: This paper tackles the problem of computing the Lipschitz constant of a neural network. Since computing the exact Lipschitz constant is NP-hard, efforts have been made to obtain tight upper bounds on the Lipschitz constant. This paper builds on the work of LipSDP [1], which involves solving a large matrix verification problem. Since that the large matrix verification problem grows significantly for both deeper and wider networks, this paper proposes a compositional approach to estimate the Lipschitz constants of deep feed-forward neural networks more efficiently. First, the authors obtain an exact decomposition of the large matrix verification problem into smaller sub-problems and, then, exploiting the underlying cascade structure of the network the authors develop two algorithms to compute a bound on the Lipschitz constant:
- The first algorithm explores the geometric features of the problem and provides a tight estimate of the Lipschitz constant by solving small semidefinite programs (SDPs) that are only as large as the size of each layer.
- The second algorithm relaxes these subproblems and provides a closed-form solution to each subproblem for extremely fast estimation, altogether eliminating the need to solve SDPs altogether.

Finally, the authors provide extensive experiments to show the different levels of tradeoffs between efficiency and accuracy of the two algorithms. They show that their approach provides a steep reduction in computation time while yielding Lipschitz bounds that are very close to, or even better than, those achieved by state-of-the-art approaches.
strength_weakness: strength_weakness: - it looks like the approach is restricted to a very limited set of neural networks (feedforward neural networks), can the approach be used for convolutional neural networks?
- In the experiments, the authors use randomly generated neural networks for their first set of experiments, in my experience it is usually easier to compute SDP on random weight matrices than on trained weight matrices due to conditioning, could the authors provide results of these experiments with trained networks? 
- Can the authors clarify if ECLipsE has to compute all subproblems at the same time or if a sequential approach is possible? 
- I assume that the ECLipsE algorithm uses Matlab for SPD optimization, have the authors tried using a deep learning framework (e.g. PyTorch) for ECLipsE-Fast? Could ECLipsE-Fast be used during training, e.g. for regularization?
questions: questions: See Weaknesses
limitations: limitations: NA
paper_summary: paper_summary: This paper presents two novel algorithms for computing the Lipschitz constant of feedforward neural networks (NN). The starting point is a previously-known semi-definite programming (SDP) problem which enables to compute the Lipschitz constant. The paper proposes a decomposition of this SDP in sequential subproblems over layers, then relax the subproblems to enable iterative computations across layers, instead of solving a joint problem on all layers. This approach scales much better with both width and depth, as demonstrated by experiments on neural networks at initialization and after training on MNIST.
strength_weakness: strength_weakness: [EDIT (Aug.7): the rebuttal answered my questions adequately. In particular, the method does provide provably-correct upper-bounds.]

- Theoretical results in Section 3.3 are a bit hard to follow, because the section gives the story behind the proposed relaxation, as well as geometric interpretation, but does not provide a main result summarizing the theoretical guarantees of the proposed approach. This is a key point, because provably correct upper-bounds on the Lipschitz constant are of course much preferable. Although it is suggested that the proposed algorithms are indeed provably correct, it is not clearly stated in the paper. So the paper would highly benefit from a clear statement on this fact, as well as a summary of the theoretical results into a theorem (see also Questions).
- The comparison with methods in the literature is limited to SDP methods, which is OK given that the main contribution is to provide a clever relaxation of these methods, but still a broader comparison would have been interesting. 
- The approach only applies to feedforward NN.
questions: questions: - I believe that Proposition 4 shows that ECLipsE-Fast provably gives an upper-bound on the true (unknown) Lipschitz constant of the neural network. Is this correct?
- Does ECLipsE also always give an upper-bound on the true Lipschitz constant? If so, are there assumptions for this to hold? If not, at which steps are approximations made? I guess this should more or less follow from Lemmas 1 and 2, and Propositions 2 and 3, but it is not clearly stated in the paper.

Minor remarks:
- Authors could consider typesetting their algorithms as “Eclipse” and “Eclipse-Fast” to improve readability.
- line 199: “such” missing.
limitations: limitations: The limitations are adequately addressed.

Comments: title: Thank you for the rebuttal
comment: Thank you for your positive evaluation of our work, and for raising your score! We appreciate your feedback on the paper!
Thank you for addressing my concerns about the row-rank assumption and the gap of the estimates of Appendix E. I think ECLipsE is an interesting result with high impact in the neural-network robustness community. I have raised my score.
We thank you for your feedback on our manuscript, particularly on wider vs deeper networks, and for raising your score! We will add the experiments on width and a discussion section on limitations.
I thank the authors for their very clear explanation on their algorithm complexity with respect to neural network depth and width, and providing experiments with higher width to demonstrate algorithm capacity, while acknowledging the limit/boundary. I am happy to see the paper to be accepted, therefore will increase my score to 7. Meanwhile, I recommend the authors to discuss around "width" in their discussion/limitation section.
Thank you for considering our response, and raising your evaluation! We also appreciate your feedback on improving our manuscript!
We agree that unrolling the convolutional layer and applying FNN based methods may not be the most practical solution. Further study is necessary to develop scalable algorithms for other network architectures, which will be the subject of future work.
We thank you for considering our response, and for raising your evaluation. We are truly grateful for the thorough and detailed suggestions on improving our manuscript! Yes, the product of two symmetric matrices is not necessarily symmetric. We can instead focus on the discussing eigenvalues of $N(M_i)^{-1}$ and prove that it can only have non-negative eigenvalues. We will edit the proof accordingly in our final version.
Thank you for the clarifications. I agree that the unrolling approach could be used in the context of this paper - although it might not be the most scalable approach (unrolling a large convolution leads to extremely large matrices). I agree that the exploration of other types of architectures could be left to future work, and I am in favor of accepting this paper and will raise my score accordingly.
Thank you for these clarifications and additional experiments. The mathematical justification of the algorithms now seems valid to me. As you suggested, a clearer explanation of the points in the proof I asked for clarification on may be beneficial for the updated version of your paper. I am in favor of accepting this paper and will raise my score accordingly.

Still just one detail:  

>(b) By the definitions of $N$ and $M_i$, they are indeed guaranteed to be symmetric.

I agree that $N$ and $M_i$ are both symmetric, but the product of symmetric matrices (e.g., $N/c_i (M_i)^{-1}$) is not necessarily symmetric. That is why I pointed out that $N/c_i (M_i)^{-1}\geq 0$ was not well-defined as $\geq$ was only defined for symmetric matrices in the Notation. Shifting the focus of the discussion in l.462-463 on the singular values without using $\geq 0$ should solve the problem anyway.
Thank you for considering our response and for increasing your evaluation. We appreciate your feedback on the paper!
I thank the authors for their precise rebuttal, which answers adequately my questions. I raised my score accordingly.