Determinant Estimation under Memory Constraints and Neural Scaling Laws
Authors: Siavash Ameli
Chris van der Heide
Liam Hodgkinson
Fred Roosta
Michael W. Mahoney
Venue: ICML 2025 poster (2025)
Keywords: log-determinants
neural scaling laws
numerical linear algebra
Primary area: probabilistic_methods->gaussian_processes
TL;DR: We compute determinants under memory constraints, and derive neural scaling laws for fast, accurate approximation.

Abstract: Calculating or accurately estimating log-determinants of large positive definite matrices is of fundamental importance in many machine learning tasks. While its cubic computational complexity can already be prohibitive, in modern applications, even storing the matrices themselves can pose a memory bottleneck. To address this, we derive a novel hierarchical algorithm based on block-wise computation of the LDL decomposition for large-scale log-determinant calculation in memory-constrained settings. In extreme cases where matrices are highly ill-conditioned, accurately computing the full matrix itself may be infeasible. This is particularly relevant when considering kernel matrices at scale, including the empirical Neural Tangent Kernel (NTK) of neural networks trained on large datasets. Under the assumption of neural scaling laws in the test error, we show that the ratio of pseudo-determinants satisfies a power-law relationship, allowing us to derive corresponding scaling laws. This enables accurate estimation of NTK log-determinants from a tiny fraction of the full dataset; in our experiments, this results in a $\sim$100,000$\times$ speedup with improved accuracy over competing approximations. Using these techniques, we successfully estimate log-determinants for dense matrices of extreme sizes, which were previously deemed intractable and inaccessible due to their enormous scale and computational demands.

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Siavash Ameli 12 Chris van der Heide" Liam Hodgkinson" 4
Fred Roosta"
Michael W. Mahoney 126
Abstract
1. Introduction
Calculating or accurately estimating log- determinants of large positive definite matrices is of fundamental importance in many machine learning tasks. While its cubic computational complexity can already be prohibitive, in modern applications,
Many quantities of interest in machine learning require ac- curate estimation of the log-(pseudo- jedeeudnnddttemminatt of large dense positive some number of datapoints. These quantities arise in a num- ber of tasks, including training Gaussian processes (GPs) and other kernel-based methods (Wang et al., 2019), graph- ical models (Rue & Held, 2005), determinantal point pro- cesses (Kulesza et al., 2012), and model comparison and selection techniques (Hodgkinson et al., 2023b). Problems where the calculation cannot be avoided can often be re- duced to computing a volume form, which is the case for tasks in statistical mechanics (Mézard & Montanari, 2009) or Bayesian computation (Gelman al., 2013), applica- tions of the Karlin-MeGreee theorem (Karlin & McGre- gor, 1959) including determinantal point processes. Simi- larly, when training GPs via empirical Bayes (Rasmussen & Williams, 2006), the log-determinant term is the most difficult to compute. In many applications, these matrices are not only dense but highly ill-conditioned. This renders methods that leverage sparsity inappropriate and makes ac- curate estimation of small eigenvalues essential, since errors in these small values are magnified at log-scale. Much of the research focus has concentrated on approximation tech- niques that ameliorate the cubic computational complex- ity of these log-determinant calculations. Methods lever- aging stochastic expansions that rely upon matrix-vects multiplications-aa such as Lanczos-basse based methods (Ubaru et al., 2017)-have have been particularly successful, yielding approximate methods that scale linearly (Dong et al., 2017), and methods leveraging Taylor's expansions (Fitzsimons et al., 2017) have also been proposed. However, accuracy can suffer when they encounter common ill-conditioning pathologies. Crucially, time complexity isn't tt the only bot- tleneck play: large-scale matrices found modern ap- plications often hit the memory wall before computational cost becomes issue (Gholami et al., 2024).
matrices, often indexed by
storing the matrices them- selves can pose a memory bottleneck. To address this, we derive novel hierarchical algorithm based on block-wise computation of the LDL decomposition for large-scale log-determinant calculation memory-constrrinee settings. In extreme cases where matrices are highly ill- conditioned, accurately computing the full matrix itself may be infeasible. This is particularly relevant when considering kernel matrices at scale, including the empirical Neural Tangent Kernel (NTK) of neural networks trained on large datasets. Under the assumption neural scaling laws in the test error, show that the ratio of pseudo-determinnnts satisfies power-law relationship, allowing us derive corresponding scaling laws. This enables accurate estimation of NTK log-determinants from a tiny fraction of the full dataset; in our experiments, this results 100,,000 speedup with improved accuracy over competing approximations. Using these techniques, we successfully estimate log-determinants for dense matrices of extreme sizes, which were previously deemed intractable and inaccessible due to their enormous scale and computational demands.
Department Statistics, University of California, Berkeley CA. USA 'International Computer Science Institute, Berkeley CA. USA 'Dept. of Electrical and Electronic Engineering, Univer- sity of Melhourne, Australia "School of Mathematies and Statis- tics, University of Melbourne, Australia 'CIRES and School of Mathematics and Physics, University of Queensland, Australia 'Lawrence Berkeley National Laboratory, Berkeley CA, USA. Correspondence to: Chris van der Heide <chris, vdh @ gmwil.maam.............m.....m..
Recently, the empirical Neural Tangent Kernel (NTK) has become a prominent theoretical and practical tool for study- ing the behavior of neural networks both during training and at inference (Novak al., 2022; Hodgkinson 2025). The Gram matrix associated with the NTK has been used in lazy training (Chizat al., 2019) and shows promise tool to obtain uncertainty quantification esti-
al.,
Proceedings of the nd International Conference Machine Learning, Vancouver, Canada. PMLR 267, 2025. Copyright 2025 by the author(s).

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
mates (Immer et al., 2021; Wilson et al., 2025). Its log- determinant-hhihlyy sensitive small eigenvalues- also received recent attention, both as quantity of inter- est in model selection techniques that rely upon marginal likelihood approximation (Immer et al., 2023; Hodgkinson et al., 2023b) and as a way to quantify the complexity of a learning problem (Vakili et al., 2021). Similarly, both the NTK log-determinant and closely related quantities have re- cently appeared in quantification of generalization error via PAC-Bayes bounds (Hodgkinson et al., 2023c; Kim et al., 2023). While several approximations have been proposed (Mohamadi et al., 2023), convergence has only been shown in spectral norm, so they cannot be expected to capture the overall behavior of the full NTK.
Table I. Comparison of stochastic Lanezos quadrature (SI.Q: with degree 4, s Monte Carlo samples, and full re-orthogonalizttonn, MEMDET (Algorithm D.2), and FLODANCE (Algorithm 1) on dense 500,000 x 500,000 NTK matrix for a ResNet50 model trained on CIFAR-10 10 with 50,000 datapoints. MEMDET com- putes the exact log-determinant and serves as the benchmark, with relative errors of other methods measured against it. Costs and wall time are based on NVIDIA H100 GPU (S2/hour) and an 8- core 3.6GHz CPU ($0.2/hour) using Amazon pricing and include NTK formation from pre-trained network.
Method
Rel.
Est.
Wall
Name
Settings
TFLOPs
Frror
Cost
Time
SLQ
1= 100, 101
5203
55%
$83
1.8 days
MEMDET
LDL, n, 32
41,667
0%
S601
13.8 days
In practice, computing the log-determinant of the empirical NTK's Gram matrix is formidable task: besides suffer- ing from the pathologies that ill-conditioned matrices are subject to at scale, the NTK appears to have its own pe- culiarities that make the problem particularly challenging. Indeed, the task is considered so impenetrable as to be uni- versally avoided, since the NTK corresponding to most rele- vant datasets cannot be stored in memory. For example, stor- ing the NTK corresponding the relatively small MNIST dataset requires 2.9 terabytes of memory. The equivalent object for ImageNet-1k 1k requires 13.1 exabytes, an order of magnitude larger than CERN's s current data storage capac- ity (Smith, 2023). While is appealing to consider only a small subset of the training dataset when making necessary approximations, naïve estimates for the -og-determinant can incur significant bias. Furthermore, we will find that leading approaches for dealing with log-determinants of empirical NTK Gram matrices using sketching and other Monte Carlo approximations tend to be highly inaccurate.
FLODANCE R. = 500, = 0
0.04
4%
$0.04
1 min 1.5 hr
FLODANCE n, 5000. 9 = 4
41.7
0.02%
$4
rithm for determinant computation. This facilitates exact calculation of log-determinants of NTK Gram matrices cor- responding to neural networks with several million param- eters. We then provide novel approximation technique, FLODANCE, based on the scaling behavior of wide class of kernel matrices, by appealing to neural scaling laws. In detail, our main contributions are:
we derive a hierarchical for large-scale computation of -og-determnnanss, which we name MEMDET;
algorithm
under mild assumptions, we derive scaling laws for the ratio of pseudo peeudoddeterminants of kernel matrices contain- ing different subsets of the same dataset, enabling both a corresponding law of large numbers and central limit theorem for normalized log-determinants;
Scaling laws have played a prominent role in machine learn- ing theory and practice, providing insight into the asymp- totic behavior of generalization error (Li et al., 2023; Vakili
leveraging these scaling laws, we propose FLODANCE, a novel algorithm for accurate extrapolation of the log- determinant from small fraction of the full dataset;
et al., 2021), as well as guidance for compute-optimal re- source allocation when deploying deep learning models at scale (Kaplan et al., 2020; Hoffmann et al., 2022). It turns
we demonstrate the practical utility of our method by ap- proximating the NTK corresponding to common deep learning models; and
out that for large kernel matrices, the ratio of successive determinants can be cast in terms of the error of an associ- ated Gaussian process. This enables known scaling laws to be deployed in estimating the log-determinant, with surpris- ing accuracy. However in order to measure this accuracy against reliable baseline, a memory memory-constrane -conttrained method for exact out-of-core log- degermeermnnnnt computation is also required, and is of independent interest.
we provide a brigh-performance Python package detkit, it, which implements the presented algorithms and can be used to reproduce the results of this paper.
Crucially, we demonstrate that our approximation technique is able to obtain estimates with lower error than incurred by reducing the numerical precision in explicit computation. Our approximation techniques render an impractical task virtually routine, as shown in Table 1. To the best of our knowledge, this is the first time that the full empirical NTK corresponding to a dataset of this scale has been computed.
Contributions.
The central contributions of this work ad- dress the issues faced when computing log-determinants of large matrices. We are primarily interested in matri- ces that cannot fit into memory. To compute an accurate baseline, we derive MEMDET, memory-constrainee algo-
2

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
The remainder of this document is structured follows. In Section 2, we discuss the computational issues faced when computing determinants at scale, and we derive our MEMDET algorithm based on block LU computation of the log-determinant. Section 3 contains the appropriate scaling laws for pseudo-determinants of interest, the corresponding LLN and CLT for their logarithms, and the FLODANCE algorithm for their approximation. Numerical experiments presented in Section 4, and we conclude in Section 5.
do not generally incur significant error. However, in our setring, the matrix of interest is realized the product JJ where J is the Jacobian. The formation of the quadratic is well-known source of approximation error, and should be avoided if possible: paraphrasing Higham (2022), if 8 is the round-off used in choice of floating point arithmetic, then for 0 < 8 < vō we can consider the simple case of J
1
have
1
0
A summary of the notation used throughout the paper is provided in Appendix A. An overview of related work in linear algebra is given in Appendix B. Computational chal- lenges when dealing with NTK matrices are then discussed in Appendix C. Implementation and performance analysis of MEMDET are provided in Appendices D and E. Re- quired background for neural scaling laws, proofs of our theoretical results, and further analysis of FLODANCE ap- pear in Appendices F and G. Comparison of various log- determinant methods is given in Appendix H. Finally, an implementation guide for detkit appears as Appendix I.
1 1 1 11.
1 A(JJ))
1 1 1 1
JJ
where JJ rounds to the singular fI(JJ).. In general. sig- nificant precision loss should occur if coId(J.J") > 8-) where z 10-3 10-6, and 10-15 for 16-bit, 32-bit, and 64-bit precisions respectively. In the Gram matrices we con- sider, (JJd > 10 10 Due to the scale of the problems of interest here, the sheer size of the Jacobian (requiring at least hundreds of terabytes of space) makes directly oper- ating on J impractical. This necessitates breaking conven- tional wisdom and explicitly forming the Gram matrix. We will see later the incurred cost to accuracy when working in low precision.
2. Computing Determinants at Scale
The computation of determinants of large matrices, particu- larly those expected to be highly ill-conditioned, is widely considered to be a computationally "ugly" problem, which should be avoided wherever possible (Axler, 1995). How- ever, this is not always an option, e.g., when the determinant represents volume form, is required in a determinantal point process, or plays a role in training GPs via empirical Bayes. Although certain limiting behaviors under expecta- tion are well characterized (Hodgkinson et al., 2023a), they provide only a coarse approximation in finite-sample set- tings. The quadratic memory cost and cubic computational complexity of na've implementations make exact compu- tation prohibitive, necessitating the use of power-precision or randomized methods at scale. While approximate meth- ods can be useful in their own right, exact computation re- mains essential, at the very least to establish meaningful baselines. In our experiments, due to the pathological spec- tral behavior of the matrices we consider, we will see that our proposed exact method is comparable in speed to the state-of-the-arr Monte Carlo approximation techniques.
2.2. MEMDET: A Memory-Constrained Algorithm for Log- Dog-Determinant Computation
Conventional methods for computing the determinant of matrices include LU decomposition (for generic matri- ces), LDL decomposition (for symmetric matrices), and Cholesky decomposition (for symmetric positive-definite matrices). These methods also simultaneously provide the determinants of all leading principal submatrices [ of an m x m matrix (for k = 1...., m)) possibly after permuting M depending on the decomposition.
The LU decomposition has a computational complexity of approximately 3 while both LDL and Cholesky decom- positions have a complexity of approximately fm² While computational complexity is a concern for large-scale de- terminant computation, memory limitation poses an even greater challenge. These methods require substantial mem- ory allocation, either the size of the array if written in in-place, or twice the size if the input array is preserved.
To address this, we present MEMDET, a memory- constrained algorithm for large-scale determinant compu- tation. Below is a sketch of the algorithm, with details of its implementation, memory and computational aspects pro- vided in Appendix D, and description of the software im- plementing this algorithm provided in Appendix I. We fo- cus on the algorithm for generic matrices using IU decom- position, though the LDL and Cholesky decompositions fol- low similarly.
2.1. Low-Precision Arithmetic
Among the most common techniques for circumventing memory and computational bottlenecks when dealing with large-scale calculations in numerical linear algebra is to cast numerical values into lower precision, usually 32-bit, 16- bit, even 8-bit floating point values, instead of 64-bit (double precision) values that would otherwise be used.
Computations of the log-determinant in mixed precision
Consider the 2 x block LU decomposition of M (see for
3

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
k j
71b
1
1
k
B
k
A
C
B
i
S
i
S
91b
72
Figure I. Schematic diagram of MEMDET illustrating the efficient block processing order for I.U decomposition (Algorithm D.1, left panel) and LDL/Cholesky decompositions (Algorithm D.2, Algorithm D.3, right panel). The detailed ordering strategy is described in Appendix D.1.
instance, Dongarra et al. (1998b, Chapter 5.4)) with the first block, M11. of size b x b, < m, as
scratchpad space. In the latter case, a cache table tracks whether a block My should be loaded from the original matrix or from the scratchpad in future calls.
M1
M2
Lu L21
U 0
12 S
M
(1)
The computational cost of this procedure remains the same, independent of the number of blocks, n, (see Appendix E. E). However, increasing the number of blocks reduces memory usage at the expense of higher data transfer between disk and memory. Efficient implementation minimizes this by processing blocks in an order that reduces the loading of B and C. The left panel of Figure shows one such order, illustrating the procedure at the k-th iteration. Processing the blocks M,j starts from the last row of the matrix and moves upward. During the horizontal and vertical traverses on the path shown in the figure, the memory blocks C and B can remain in memory, avoiding unnecessary reloading. Once the procedure reaches the block (A= + 1, + 1). the memory block can be directly read into A instead of be- ing stored on disk as M initiating the block for the next iteration. Consequently, only the blocks shown in blue need be stored in the scratchpad space, while those in dark blue are already stored in the current state of the algorithm. This algorithm can be further modified to elimi- nate the need for A and use the memory space of S instead, but this would require storing the additional gray blocks on the scratchpad space in addition to the blue blocks. pseudo code and further efficient implementation details of the presented method can be found as Algorithm D.I.
M
M22
where L1 is lower triangular, Un is upper-triangular, and I is the identity matrix. The blocks of the decomposition obtained by computing Mn = [1 Un solving lower- Li M2 and upper-triangular triangular system U 12 1 and forming the Schur comple- system L21 = M21 UII
ment S = M L21 112- This procedure is then repeated on the (m- b) x (m - b) matrix S, treated as the new M, leading to a new b x upper-triangular matrix Un and a smaller Schur complement S at each iteration. This hierar- chical procedure continues until the remaining S of size b or less, at which point its LU decomposition is computed. The log-determinant of the entire matrix is the sum of the hog-determinants of all U, blocks at each iteration.
To make this efficie,, the algorithm is modified to hold only small chunks of the matrix in memory, storing intermediate computations on disk. Suppose M consists of ny Xn, blocks Myy of size bx b where i.j 1,..., n6. We pre-allocate four b b matrices A. B, C. and S in memory. The k-th stage begins by loading A + MR from disk and performing an in-place LU decomposition A = LU, with L and U stored in A. We then compute the Schur complement for all inner blocks Mij. i = + 1,..., 74g by loading B < My, and C +- M,, from disk, solving B + LLL-BBB B and C + CU-L in-place, loading S + M.j. and computing S < S - CB. The updated is then stored back to disk M.j + S, either by overwriting a block of the original matrix, or avoids doing so by writing to a separate
The algorithms for LDL and Cholesky decompositions fol- low a similar procedure with necessary adjustments for sym- metric and symmetric positive-definite matrices, as detailed in Algorithm D.2 and Algorithm D.3. These modifications include processing only the lower (or upper) triangular part
4

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
50 CIFAR-I0 10
ResNet50 CIFAR-10
10"h
10
95.3
de)(KK,)) det((K)))))
del(K..) ( del/(K,...))
10t1
10a1
-5.3
1020
1020
10
10a
10a
10
10
x
10
R
1
x10'
Figure Demonstration a scaling law for the ratios successive determinants det(KK,)// det((K,,)) (Assumption fo the empirical neural tangent kernel Gram matrix of a trained ResNet50 network on CIFAR- 10 with 50,000 datapoints.
of the matrix and handling permutations and diagonal scal- ing in the LDL decomposition. The right panel of Figure shows one possible block ordering for the algorithm for LDL and Cholesky decomposition: the ordering is optimal in minimizing the number of reads and writes, but it is not unique. Other block processing orderings with the same amount data access also exist, and these will be further discussed in Appendix D.1. A detailed analysis of the com- putational complexity and memory/disk data transfers of the algorithm is provided in Appendix E. Implementation of MEMDET algorithm can be found in Listing I.1.
cal results that networks initialization induce GP kernel in the large-width limit, enabling Bayesian inference for these infinite-width networks (Neal, 1996; Lee et al., 2018). Similarly, linearization of the gradient-flor dynamics dur- ing the late stages of training leads to the derivation of the NTK, whose infinite-wide analogue can be shown to be constant during this training phase (Jacot et al., 2018; Yang, 2020; Yang & Littwin, 2021).
The NTK was first derived in the context of neural networks. However, the quantity is well defined for more general class of functions. For continuously differentiable func- tion So : X' + Rd where d is the dimensionality of the model's output (e.g., the number of labels in a classification task), we define its (empirical) NTK as
Although the NTK matrices we work with are symmet- ric and positive definite (SPD), we do not solely rely on Cholesky decomposition, despite being the most suitable method for SPD matrices. This is because NTK matrices can lose their positive-defimitensss with even the smallest perturbations causing small eigenvalues to become negative, such as when converting from 64-bit to 32-bit precision for efficient computation. As a result, the Cholesky decom- position becomes unstable and fails, necessitating the use of LU and LDL. decompositions, suitable for more general matrices. We note that the block computations of LU de- composition can become unstable the matrices deviate from symmetry and positive-deffmitensss (Demmel et al., 1995), requiring pivoting of the blocks. However, in our empirical study of NTK matrices, we found that LU decom- position works well without block pivoting, though we do consider pivoting within each block.
Ke( x') =
(2)
where (z)) E Rx is the Jacobian the function So with respect to the flattened vector its parameters R... evaluated the point x. The assumption of continuous differentiability is often relaxed in practice.
Note that (z(x, dxd so computing the NTK across n datapoints yields a 4th-order tensor of shape (n, n, d, d). For computational purposes, this is typically flattened into a two-dimensional block matrix of size nd x nd, where each (i. j) block corresponds the x d matrix x(,,,,;))
The scaling laws we derive apply to kernel families satis- fying decay conditions on the eigenvalues of an associated integral operator (see Appendix F.I for precise statement drawn from Li et al. (2023)). This class of kernels includes NTKs, as shown by Bietti & Mairal (2019); Bictti & Bach (2021); Lai al. (2023); see Figure 2.
3. Scaling Law for the Determinant
We now turn our attention to the estimation of the log- determinant of Gram matrices corresponding to covariance kernels. Our primary motivating example is the NTK, an architecture-specific kernel associated with deep neural net- work models. The connections between neural networks and GP kernels are well-known, particularly the now classi-
3.1. Neural Tangent Kernels and Scaling Laws
Let K X' x X' + dxd be positive- kernel. For a sequence of inputs (r,) 0o C X and for each n E N, let
5

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
K, (M(2,, (x(2),)))) 1y)) be the corresponding Gram matrix for the first inputs, which is matrix of size nd x nd. Our objective is to estimate logdet(KK,) for large n, using computations involving only smaller values of n. Recall that f is a Gaussian process with mean function Ja :A:A - RR and covariance kernel K. : X' x d' + Rdxd denoted f ~~ GP(u, 8), if for every -i....;2in E X and d 2 1,
Res Net50
CIFAR-10
10P
95.3
tht del
100
(f(x,))== ((((<))))= (K(2,, JyD))=1) resultss founded on the following lemma, which we prove in Appendix F.2.
Normal
Rolling monn (o(n)
p(n) # m(u
Empirical
10 9
0%
9%
Lemma 1. Let f be zero-mean vector-valued m-dimensional Gaussian process with covariance kernel K. For each n 2 2, let
x10
PDF
Figure 3. Demonstration of the stationarity and second moment be- havior (Assumption 2) of the corresponding logarithmic process.
(1((=,)) f(x;) 1... n- E(n) = Bid =
Assumption 2 (STATIONARITY). Assume there exists a constants C' and exponent v > 0 such that the process 8n satisfying
denote the mean-squared error of fitting the f to the zero function using 3 Fn--- Then
det((kn) 8, 8n-1 = log (de((K))))) loog(Cn
pdet K,,)
(
d T/2
(On
2,3.....
B(n)",
pdet (
is stationary, ergodic, has finite second moment (E(o)) <
where r the rank of Cov(f(xn)
f(x) f
+00), and
1....., - 1). In the case where rank(K,))
(ank(K,-1)> d, this reduces to
Under Assumption and Assumption 2, we derive an ex- pression for the asymptotic behavior of the normalized log- determinants. The proof is given in Appendix F.2.
pdet(K,,)
S F(()), for any >>> 1,
pdet(Kn.))
Proposition 1. For larger numbers of imputs, letring
with equality if =d=1.. 1.
In 'lloddet(KK))
Ln 1===1 L1 11-1
log(n!)
Lemma is particularly interesting since it highlights a con- nection between the determinants of Gram matrices and error curves for GPs. This bounds the effect of adding or removing datapoint the determinant in terms of the prior variance of corresponding GP. Previous theoreti- cal (see Appendix F.1) and empirical studies (Spigler et al., 2020; Bahri et al., 2024; Li et al., 2023; Barzilai & Shamir, 2024) have established a power law relationship of the form E(m) (n 5) as n - oo for some & of this, and invoking Lemma we propose the following scaling law for determinants, demonstrated in Figure 2.
12 Co
n
law of large numbers (LLN) and central limit theorem (CLT) hold for the log-determinants L.a
(LLN) Under Assumption 1, there exist constants Co: > O such that Co.
L, = L, (1)
(4)
In view
(CLT) Under Assumption 2, there exist constants Co. V. o > O such that as n - 00,
Assumption 1 (SCALING LAw). Assume there exists constant C' > 0 and exponent > 0 such that - cc.
Vn (L, L..) N(,,00
(5)
det((,)) = - + (1).
We remark that log(n!) I is O(nlog(n) so the normalized log-determinant 4,, is (log((n) via (4).
(3)
det(K,,
This assumption allows for accurate estimation of the log- determinants of interest, and appropriate law of large numbers. However, to construct corresponding confidence intervals, we will need to assume some properties of the error term appearing in Assumption We further impose the mild assumptions of stationarity and bounded variance, demonstrated in Figure 3.
3.2. FLODANCE
As (4) is linear in the unknown parameters co M, these pa- rameters can be estimated by linear regression on a pre- computed sequence (I,,...,,,,) L.,.. This sequence can be
'This is generally more stable regression problem than esti- mating C and directly from (3).
6

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
computed using single pass of MEMDET on a subsample of the Gram matrix. After performing linear regression, by extrapolating to larger n, the normalized log-determinants of larger NTK Gram matrices can be estimated. The or- dering of the data is arbitrary, but it will affect the out- put of the regression task. We refer to this method as the Factorial-based Log-Determinant Analysis and Numerical Curve Estimation procedure, or FLODANCE. To improve performance, allow for non-asymptotic correction to the exponent V v(n)= == Vo + Ea= Mr -i (in practice, we find y S 10 works well). In light of (5) and discussion in Appendix F.3, for large n, we expect that approximately
Algorithm 1: FLODANCE: Pactorial-based Log- Determinant Analysis and Numerical Curve Esti- mation
Input
Precomputed partial NTK Gram matrix K,, of size m X m where m, n.d,
Model's output dimension d. Total number of datapoints n,
Data subsample size < n2 S n.
Burn-in length 1 - ng < n2g-
Number of terms in the Laurent series o Output: Estimated normalized log-determinant Ln of K, of size m X m where m = nd
9+1
// KKn,fokkkk) is the ls x l principal sub-matrix of K., 1 Run Algorithm D.2 on K,, to obtain (4x k= where lg x logabsdet (,,,.ttt) (,l,,,)))))
Σ
12 + 6n cn ~ ~((,00))), (6)
yn = nn.nn
where Un = Vn-I (Ln L1), Fn.0 vn and
Normalize and record every dl-th entry.
Tn,i VVa-T for 1,... m + 1 are the covari- ates. The corresponding numerical procedure is presented in Algorithm I and the implementation can be found in Listing 1.2.
Ing(n!)
2 Obtain (I (Tn)) for I.g f- 1m, and my + jd.
// Define design matrix X E R ( and response vector y
R00
3 for j = 1 to n, 7o do
In practice. we also observe that a burn-in period may be re- quired to obtain accurate estimates of cq and the 2 that appear in Proposition 1. Better performance achieved in our experiments by discarding the early determi- nant samples. effectively replacing the L term appearing in Un with different constants Inn,- for a burn-in of length no - 1. This seems to be due to the sudden emergence of very small eigenvalues that shift the model fit, and consti- tutes a consistent phenomenon that warrants further inves- tigation. We found that when needed, the burn-in required was always less than 500 terms, verified by cross- pross-valddttion
4
m -- no + j
5
y, + ny (1² 11) (Lng - -Lod)
often
X,, + (n, - 11}²
for i = 1 to + do
X,,51 s (n, nn 2n5 logT (n,+1)
8
The Log- gamma function computes log(n,)) !).
// Estimate regression coefficients B in y X,3 +
3 (x-x)='XXx XTy
1/ Estimate L, at larger value of n
10 (Co, Vo, U1) t- B
4. Numerical Experiments
11 i. + Lng + c - n 2 - g
-1
4--1 vn-in-log[(ne "log[(++1)
We now evaluate the accuracy of the FLODANCE algo- rithm for estimating log-determinants on large-scale prob- lems of interest. The test problems that we consider are NTK matrices corresponding to common deep learning models: ResNe19, ResNet18, and ResNet50 (He et al., 2016) trained on the CI dataset (Krizhevsky, 2009), and MobileNet (Howard et al., 2017) trained on the MNIST dataset (LeCun et al., 1998). Our experiments are split into two sections: first, the dataset size is reduced in order to en- able the matrices fit into memory on consumer device, in Section 4.1; and then larger subsamples and full datasets are considered, in Section 4.2.² As a baseline, we employ MEMDET to compute the relevant quantities, where accu- racy is limited only by numerical precision. Detailed run- time and performance diagnostics for MEMDET are pro-
12 return i.n.
vided in Appendix E. Additional experiments evaluating the performance of FLODANCE appear in Appendix G, in- cluding its robustness (Appendix G.1) and its application to Matérn kernel Gram matrices (Appendix G.2).
4.1. Smaller Data Sets
In order to demonstrate the scaling laws for ratios of suc- cessive determinants of NTK Gram matrices, a ResNet50 model was trained on a subset of 1000 datapoints from the CIFAR-10 dataset with d = 10 classes. The resulting ma- trix K1oon is of size 10,000 x 10,000. In Figure 2, we plot exact ratios det (e.) (e,, with the blue line
Experiments in this section were conducted on desktop-class device with an AMD Ryzen@7 5800X processor, NVIDIA RTX 3080, and 64GB RAM.
7

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Table 2. Comparison of approximations of the log- determinant é, with the exact computation e,, obtained in 64-bit floating point prccision (first row). Values represent average percentage relative errors over five trained networks, with standard deviations in parentheses. Bold
values indicate the closest approximation, with the next-best underlined. For the corresponding compute times, Table H...
Quantity
Model
Configuration
ResNet9
ResNet9
18
MobileNet
Dataset
CIFAR-10 10
CIFAR-10
CIFAR-I0
MNIST
Subsample Size
1000
2500
n 1000
2500
fn
Direct Computation (64-bit)
(Reference)
76538 (203)
181377 (649)
65630 (842)
183962 (7869)
Direct Computation (16-bit)
12.41% (0.12)
17.05% (0.13)
14.00% (0.24)
66.97% (2.13)
Direct Computation (32-bit)
3.67% (0.06)
6.77% (0.08)
5.25% (0.09)
14.27% (0.95)
SLQ
81.51% (0.16)
80.89% (0.24)
101.03% (1.61)
84.52% (1.51)
Relative Ervor
Block Diagonal
76.49% (0.12)
75.15% (0.16)
92.76% (1.55)
112.55% (1.22)
Pseudo NTK
118.35% (0.10)
35% (0.27)
122.95% (0.25)
75.32% (1.04)
FLODANCE FLODANCE FLODANCE
7g 1, n2
7.73% (0.77)
11.27% (1.10)
12.19% (0.30)
36.41% (2.53)
Rg = 1, = 100
5.61% (0.32)
8.54% (0.63)
8.09% (0.68)
35.51% 1% (1.46)
fgg 300, n1 500
1.34% (0.11)
1.37% (0.14)
2.9% (0.81)
23.19% (1.76)
representing the line of best fit under the scaling law As- sumption
20
For baseline comparison existing techniques, we com- pared FLODANCE approximations of both the matrices themselves and their log-determinants. For the matrix ap- proximations, we consider a plock-diagonal approximation ignoring between-data correlations, as well as the pseudo- NTK matrix studied in (Mohamadi et al., 2023). As an ap- proximate log- oeg-determinant technique, we consider stochas- tic Lanezos quadrature, often regarded the current state- of-the-art for large-scale log- degeeeermmnant estimation (Gard- ner et al., 2018), We also compare exact methods across 16-, 32-, and 64-bit (treated as exact) floating-poin precision, to assess the accuracy of our extrapolation against memory- saving maixed-precisie calculations.
Absolute Error li
15-
10-
5-
500
1000
1500
2000
2500
Block Diag
FLODANCE (Ng 1, 50)
SLQ
FLODANCE (ho 1, 100)
Direct Comp. (16-181)
FLODANCE (ho 200. n. 500)
Direct Comp. (32 bit)
To this end, ResNet9 was trained on a subsets of 1000 and 2500 images from CIFAR- 10, and R 1000 data- points. MobileNet was also trained on 2500 image subset of the MNIST dataset with d = 10 classes. A comparison of the different methods is presented in Table 2. We see that all existing NTK and log-determinant approximation techniques perform poorly when compared with 16- and 32-bit mrixed-precisionn calculations. On the other hand, the FLODANCE estimates that contained a burn-in phase con- approxi-
Figure 4. Comparison log-determinant accuracy for NTKs of ResNet9 trained on CIFAR-IO, 100 measured by absolute error, across variety approximation techniques for matrices of different sizes. Means across five trained networks are displayed, with shaded regions depicting one standard deviation.
4.2. Larger Data Sets
In our next experiment, we evaluate NTK matrices at an un- precedented scale, where exact determinant computation has not been previously reported. Due to memory con- straints, these matrices cannot be stored explicitly, making MEMDET essential for obtaining ground truth values.
sistently either outperformed the mixed-preciss mations or were competitive. When no burn-in was used, FLODANCE still outperformed the approximation meth- ods. This suggests that at this scale, the error in the scaling law approximation is less than the mixed- precisson errors discussed in Section 2. Extrapolating determinants based on expected behavior can bypass numerical issues scale. Comparisons for ResNet9 trained on 2500 images from the CIFAR-10 10 dataset are visualized in Figure 4. FLODANCE estimates consistently outperform the competing methods.
We consider two large-scale NTK matrices: K50,000+ dense matrix of size 500,000 x 500,000, for ResNet50 trained on CIFAR-10 with = 50,000 and d = 10, and an identical sized NTK for ResNer9. At this scale, comput- ing the full matrix in double precision poses a formidable challenge. To our knowledge, this is the first time an NTK
8

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
1
[a] ResNet8
10
(b) ReaNec50
3.0
1.0%
1.0%%
Laft Ordinate
Left Ordinats
Empirioal
Empiricai
2.5
Asymppott
Aqymptore
0..5%
2.0
Extropoaatio
2.5
Extrapolation
Right Ordimaae
-..9%
Sight Ohdinate
+0.0%
Rolativn Earur
Relative
<5
155
00.1%
1.5
0..1%
1.0
1.0
0.2%
0.2%%
0.5
0.5
0.0
0.0%
0.0
0.0%
-10*
x10
Figure 5. Log-determinant 6, for = ........00 50,000, corresponding to NTK submatrices where = nd and d = 10, from 64-bit NTK matrices of ResNet9 (a) ResNer50 (b) trained on CIFAR-10 with 50,000 datapoints. Values are computed using MEMDET (Algorithm D.2) with LDL decomposition (black curves, overlaid by colored curves). FLODANCE (Algorithm 1) is fitted in small region (shaded gray) and extrapolated over a much larger interval. The yellow curve the interval (me n..) = (10² x 10²) represents the fit, while the red curve in (n.. (5 10 5 10") shows the extrapolation. The blue curves, corresponding the right axis in each panel shows the relative error, reaching impressive 0.2% (a) and 0.02% in (b).
matrix of this size-at at full precision and over the entire CIFAR-10 dataset-has been computed, with the matrix it- self requiring 2 TB of memory (see Appendix C.I). This computation was carried out on an NVIDIA Grace Hopper GH200 GPU over 244 hours for ResNet50.
plicit computation. typically relying on stochastic approx- imations. However, in many problems, the sheer size of the matrices and their ill-conditioned nature, make not only approximation a difficult task, but forming the matrix itself to provide baseline becomes computationally intractable. We have addressed this problem two fronts. On the one hand, we defined MEMDET, memory-constrained algorithm for log- log-determinant computation, with different versions for general, symmetric, and symmetric positive- definite matrices. On the other hand, we derived neural scaling laws for large kernel matrices, and we introduced FLODANCE, a procedure for accurate extrapolation of log- determinants from small subsets of the data. The high speed and accuracy of our methods open the door for routine com- putation of interpolating information criteria and related diagnostic tools to enable principled model selection within deep learning frameworks (Hodgkinson et al., 2023b).
Having established this benchmark, we now evaluate the accuracy of FLODANCE at this scale. As shown in Fig- ure 5, FLODANCE with 9 6, Ro = 100, and n, = 5000 achieves an absolute error of just 0.2% for $00,000 on ResNet9, reducing computation time by a factor of (n/n,)² 1000. Similarly, for ResNet50, FLODANCE with q 4, ng 100, and n, 5000 achieves an absolute error of 0.02% with the same speedup. In contrast, SLQ exhibited relative error of 55% (see Table 1, also List- ing I.3 for implementation). Given their poor performance on smaller datasets, pSeudo-NTK and block-diagonal ap- proximations are omitted.
The ability to accurately compute and estimate matrices of this size further provides fascinating insights into the be- havior of the NTKs that considered in experiments, which treated square matrices of the size up to 500,000. Fur- ther, the memory constrained algorithms we described can be applied to other classes of matrices (Nguyen & Vu, 2014; Cai et al., 2015), where they can be expected to unlock sim- ilar insights into their scaling behavior. In terms of further computational tools, we are developing techniques ex- trapolate more refined spectral information about large trices from small sub-blocks (Ameli et al., 2025). Methods for blockwise decompositions of large Jacobian matrices would also circumvent the need to explicitly calculate JJJ", enabling higher resolution understanding of their behavior.
This experiment represents the first exact computation of an NTK determinant this scale, establishing new bench- mark for large-scale log-des degeeteermnnnt estimation. Additional experiments are provided in Appendix G, with Figure G.1 demonstrating the accuracy of the global FLODANCE fit, Figure G.2 illustrating sensitivity to subsample choice, and Figure G.3 examining sensitivity to subsample size.
5. Conclusion
The calculation of large matrix log-determinants is com- monly encountered but often avoided problem in statistical and machine learning applications at scale. Several tech- niques have previously been proposed to circumvent ex-

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Acknowledgments
nett (Eds.), Advances in Neural Information Processing Systems, volume 32: Curran Associates, Inc.
LH is supported by the Australian Research Council through a Discovery Early Career Researcher Award (DE240100144), FR partially supported by the Australian Research Council through the Industrial Transformation Training Centre for Information Resilience (IC200100022). MWM acknowledges partial support from DARPA, DOE, NSF, and ONR.
Billingsley, P. (1961). The Lindeberg-Levy theorem for
martingales. Proceedings of the American Mathematical Society, 12(5). 788-792.
Cai, T. T., Liang, T., & Zhou, H. H. (2015). Law of log determinant of sample covariance matrix and optimal estimation of differential entropy for nigh-dimensional Gaussian distributions. Journal of Multivariate Analysis, 137, 161-172.
Impact Statement
This paper presents work whose goal advance the field of Machine Learning. There are many potential societal consequences of our work, none which we feel must be specifically highlighted here.
Chizat, L.., Oyallon, E., & Bach, F. (2019). On lazy train- ing in differentiable programming. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc E. Fox, & R. Garnett (Eds.), Advances in Neural Information Processing Systems, volume 32: Curran Associates, Inc.
References
Demmel, J. W., Heath, M. T., & van der Vorst, H. A. (1993). Parallel numerical linear algebra. Acta Numerica, 2, 111-197.
Abramowitz, M. & Stegun, I. A. (1964). Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. New York: Dover, ninth dover printing, tenth GPO printing edition.
Demmel, J. W., Higham, N. J., & Schreiber, R. S. (1995). Stability of block LU factorization. Numerical Linear Algebra with Applications, 2(2), 173-190.
Ameli, S. & Shadden, S. C. (2023). A singular Woodbury and pseudo-determinant matrix identities and application to Gaussian process regression. Applied Mathematics and Computation, 452, 128032.
Dong, K., Eriksson, D., Nickisch, H., Bindel, D., & Wilson, A. G. (2017). Scalable log determinants for Gaussian process kernel learning. I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, Vishwanathan, & R. Garnett (Eds.). Advances Neural Information Processing Systems, volume 30: Curran Associates, Inc.
Ameli, S., van der Heide, C., Hodgkinson, L., & Mahoney, M. W. (2025). Spectral estimation with free decompres- sion. arXiv preprint arXiv:2506./199
Axler, S. (1995). Down with determinants! The American Mathematical Monthly, 102.
Dongarra, J., Hammarling, S., & Walker, (1998a). Key concepts for parallel out-of-core LU factorization. Com- puters & Mathematics with Applications, 35(7), 13-31. Advanced Computing on Intel Architectures.
Bahri, Y., Dyer, E., Kaplan, J., Lee, J., & Sharma, U. (2024). Explaining neural scaling laws. Proceedings of the Na- tional Academy of Sciences, (27, e2311878121.
Dongarra, J. J., Duff, L. S., Sorensen, D. C., & Vorst, H. A. V. (1998b). Numerical Linear Algebra for High Per- formance Computers. USA: Society for Industrial and Applied Mathematics.
Barrachina, S., Castillo, M., Igual, F. D., Mayo, R., & Quintana-Orti E. S. (2008). Solving dense linear sys- tems on graphics processors. In Euro-Par 2008 Parallel Processing (pp. 739-748). Berlin, Heidelberg: Springer Berlin Heidelberg.
Dongarra, J., Moler, C. B., Bunch, R., & Stewart, G. W. (1979). LINPACK Users' Guide. Society for Industrial and Applied Mathematics.
Barzilai, D. Shamir, O. (2024). Generalization in kernel regression under realistic assumptions. In Proceedings of the s1 International Conference Machine Learning, ICML'24: JMLR.org.
Engel. A., Wang. Z., Sarwate, A. D., Choudhury, S., Chiang, T. (2022). TorchNTK: A Library for calculation Neural Tangent Kernels of PyTorch Models.
Bietti, A. & Bach, F. (2021). Deep equals shallow for ReLU networks in kernel regimes. In International Con- ference on Learning Representations.
Fischer, S. & Steinwart, I. (2020). Sobolev norm learning rates for regularized least-squares algorithms. Journal of Machine Learning Research, 2105), 1-38.
Bietti, A. & Mairal, J. (2019). On the inductive bias of neural tangent kernels. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, & R. Gar-
Fitzsimons, J., Cutajar, K., Osbore, M., Roberts, S., & Filippone, M. (2017). Bayesian Inference of Log De- terminants. In Proceedings of the 33rd Conference
10

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Uncertainty in Artificial Intelligence, Proceedings of Ma- chine Learning Research Sydney, Australia: PMLR.
Models of Heavy Tailed Mechaniestic Universality. In Proceedings of the 42nd International Conference Ma- chine Learning.
Galoppo, N., Govindaraju, N., Henson, M., & Manocha, D. (2005). LU-GPU: Efficient algorithms for solving dense linear systems on graphics hardware. In SC '05: Proceedings of the 2005 ACM/IEEE Conference on Su- percomputing (pp. 3-3).
Hoffmann, J., Borgeaud, S., Mensch, A., Buchatskaya, E., Cai, T., Rutherford, E., de Las Casas, D., Hendricks, L. A., Welbl, J., Clark, A., Hennigan, T., Noland, E., Millican, K., van den Driessche, G., Damoe, B., Guy, A., Osindero, S., Simonyan, K., Elsen, E., Vinyals, O., Rae, J., & Sifre, L. (2022). An empirical analysis of compute- optimal large language model training. In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, & A. Oh (Eds.), Advances in Neural Information Processing Sys- Curran Associates,
Gardner, J., Pleiss, G., Weinberger, K. Q., Bindel, D., & Wilson, A. G. (2018). GPyTorch: Blackbox matrix- matrix Gaussian process inference with GPU accelera- tion. Advances in neural information processing systems, 31.
tems, volume 35 (pp. Inc.
Gelfand, A. E., Diggle, P., Guttorp, P., & Fuentes, M. (2010). Handbook of Spatial Statistics. Chapman & Hal- VCRC Handbooks of Modern Statistical Methods. CRC Press.
Howard, A. G.. Zhu, M., Chen, B., Kalenichenko, D., Wang. W., Weyand, T., Andreetto, M., & Adam, H. (2017). Mobilenets: Efficient convolutional neural net- works for mobile vision applications. arXiv preprint arXhv:1704.04860
Gelman, A., Carlin, J. B., Stern, Dunson, D. B., Vehtari,
A., & Rubin, D. B. (2013). Bayesian Data Analysis. Chapman and Hal//CRC, 3nd ed. edition.
Immer, A., Korzepa, M., & Bauer, M. (2021). Improving predictions of bayesian neural nets via local linearization. In A. Banerjee & K. Fukumizu (Eds.), Proceedings of The 24th International Conference gence and Statistics, volume 1 130 of Proceedings of Ma- chine Learning Research (pp. 703-711).: PMLR.
Gholami, A., Yao, Z., Kim, S., Hooper, C., Mahoney, M. W., & Keutzer, K. (2024). AI and memory wall. IEEE Micro.
Artificial Intelli-
Golub, G. H. & Van Loan, F. (2013). Matrix Compu- tations- - 4th Edition. Philadelphia, PA: Johns Hopkins University Press.
Immer, A., Van Der Ouderaa, T. F. A., Van Der Wilk, M., Ratsch, G., & Schölkopf, B. (2023). Stochastic marginal likelihood gradients using neural tangent kernels. In A. Krause, E. Brunskill, K. Cho, B. Engelhardt, S. Sabato, & J. Scarlett (Eds.), Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research (pp. (4333-14352), PMLR.
He, K., Zhang, X., Ren, S., & Sun, J. (2016). Deep residual learning for image recognition. In Proceedings of the IEEE conference on computer vision and pattern recog- nition (pp. 770-778).
Higham, N. (2022). Seven Sins of Numerical Linear Algebra. https: /ttpss///niihham.. ctw/c02/mmm//ccc wwwerraaalleeeeaaaeeeeee Accessed: 2025-01-29.
seven- s
Jacot, A., Gabriel, F., & Hongler, C. (2018). Neural tan- gent kernel: Convergence and generalization in neural networks. In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Biancee & R. Garnett (Eds.), Advances in Neural Information Processing Systems, vol- ume 31: Curran Associates, Inc.
Hodgkinson, L.. van der Heide, C., Roosta, F., & Mahoney, M. W. (2023a). Monotonicity and double descent in uncertainty estimation with Gaussian processes. In Pro- ceedings of the 40th International Conference on Ma- chine Learning (pp. 13085-13117).
Kaplan, J., McCandlish, S., Henighan, T., Brown, T. B., Chess, B., Child, R., Gray, S., Radford, A., Wu, J., & Amodei, D. (2020). Scaling laws for neural language models. arXiv preprint arXiv:2001.00 108311.
Hodgkinson, L., van der Heide, C., Salomone, R., Roosta, F,, & Mahoney, M. W. (2023b). The interpolating in- formation criterion for overparameterized models. arXiv preprint arXiv:2007.07785vl
Karlin, S. & McGregor, J. (1959). Coincidence proba- bilities. Pucific Journal of Mathematics, 9(4), 1141 1164.
Hodgkinson, L., van der Heide, C., Salomone, R., Roosta, F., & Mahoney, M. W. (2023c). A PAC-Bayesian per- spective the interpolating information criterion. In Ad- vances in Neural Information Processing Systems, Math- ematics of Modern Machine Learning Workshop.
Kim, S., Park, S., Kim, K.-S., & Yang. E. (2023). Scale- invariant Bayesian neural networks with connectivity tan- gent kernel. In The Eleventh International Conference
Hodgkinson, L., Wang, Z., & Mahoney, M. W. (2025).

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
on Learning Representations.
Processes for Machine Learning. MIT Press Cambridge.
Krizhevsky, A. (2009). Learning multiple layers of features from tiny images. Technical Report TR-2009, University of Toronto.
Rue, H. & Held, L. (2005). Gaussian Markov Random Fields: Theory and Applications. Chapman and Hal- I/CRC.
Kulesza, A., Taskar, B., et al. (2012). Determinantal point processes for machine learning, Foundations and Trends@ in Machine Learning, 5(2-3), 123-286.
Smith, T. (2023). An exabyte of disk storage CERN. https://homts httpss//hone.ccem//ees/////// comput Computing//e somputinggeeabbtteeeeee htgp/eeabytedddi--sttrraaa Accessed: 2025-01-29.
Lai, J., Xu, M., Chen, R., & Lin, Q. (2023). Generalization ability of wide neural networks on R. arXiv preprint arXiv:2302.05953
Spigler, S., Geiger, M., & Wyart, M. (2020). Asymptotic learning curves of kernel methods: empirical data versus teacher-student paradigm. Journal of Statistical Mechan- ics: Theory and Experiment, 2020(12), 124001.
LeCun, Y., Bottou, L., Bengio, Y., & Haffner, P. (1998). Gradient-based learning applied to document recognition. Proceedings of the IEEE, 86(11), 2324..
Stark, S. & Beris, A. N. (1992). LU decomposition op- timized for parallel computer with a hierarchical dis- tributed memory. Parallel Computing, 18(9). 959-971.
Lee, J., Sohl-dickstein, J., Pennington, J. Novak, R., Schoenholz, S., & Bahri, Y. (2018). Deep neural net- works as Gaussian processes. In International Confer- ence on Learning Representations.
Stein, M. L.. (1999). Interpolation of Spatial Data:
Some Theory for Kriging. Springer Series in Statistics. Springer New York.
Li, Y. Zhang, H., & Lin, Q. (2023). On the asymptotic learning curves of kernel ridge regression under power- law decay. In A. Oh, T. Naumann, A. Globerson, K. Saenko, M. Hardt, & S. Levine (Eds.), Advances in Neural Information Processing Systems, volume 36 (pp. 49341-49364): Curran Associates, Inc.
Steinwart, I., Hush, D. R., & Scovel, C. (2009). Optimal rates for regularized least squares regression. In COLT.
Ubaru, S., Chen, J., & Saad, Y. (2017). Fast Estimation of (r() via Stochastic Lanezos Quadrature. SIAM Jour- nal on Matrix Analysis and Applications, 38(4). 1075- 1099.
Matérn, B. (1960). Spatial variation. In Meddelanden från Statens boogsforskningsinsitttt, volume 49, No. 5. Almänna Förlaget, Stockholm. Second edition (1986), Springer-Verlee Berlin.
Vakili, S., Bromberg, M., Garcia, J., shan Shiu, D., & Bernacchia, A. (2021). Uniform generalization bounds for overparameterized neural networks. arXiv preprint 0et:2109906099.
Mézard, M. & Montanari, A. (2009). Infommtion, Physics, and Computation. Oxford University Press.
Venetis, 1. E. & Gao, G. R. (2009). Mapping the LU decomposition on a many-core architecture: challenges and solutions. In Proceedings of the 6th ACM Conference on Computing Frontiers, CF '09 (pp. 71-80). New York, NY, USA: Association for Computing Machinery.
Mohamadi, M. A., Bac, W., & Sutherland, D. J. (2023). A fast, fell-founded approximation to the empirical neural tangent kernel. In A. Krause, E. Brunskill, K. Cho, B. Engelhardt, S. Sabato, & J. Scarlett (Eds.), Proceedings of the 40th International Conference Machine Learn- ing, volume 202 of Proceedings of Machine Learning Research (pp. 25061-25081): PMLR.
Wang, K., Pleiss, G., Gardner, J., Tyree, S., Weinberger, K. Q., & Wilson, A. G. (2019). Exact Gaussian processes on million data points. Advances in neural information processing systems, 32.
Neal, R. M. (1996). Bayesian Learning for Neural Net- works. Berlin, Heidelberg: Springer-Verlag.
Wilson, J., van der Heide, C., Hodgkinson, I.., & Roosta, F. (2025). Uncertainty quantification with the empirical neural tangent kernel. arXiv preprini 250..2250
Nguyen, H. H. & Vu, V. (2014), Random matrices: Law of the determinant. The Annals of Probability, 42(1), 146-. 167.
Yang, G. (2020). Tensor programs II: Neural tangent kernel for any architecture. ArXiv, abs/2006. 145488
Novak, R., Sohl-Dickstein, J., & Schoenholz, S. S. (2022).
Fast finite width neural tangent kernel. In International Conference on Machine Learning (pp. 17018-17044). PMLR.
Yang, G. & Littwin, E. (2021). Tensor programs IIb: Ar- chitectural universality of neural tangent kernel training dynamics. In M. Meila T. Zhang (Eds.), Proceedings of the 38th International Conference on Machine Learn-
Rasmussen, C. E. & Williams, C. K. (2006). Gaussian
12

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
ing, volume 139 of Proceedings of Machine Learning
Research (pp. 11762-11772) PMLR.
13

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Appendices
Contents
Nomenclature
14
B Related Works in Numerical Linear Algebra
15
Memory and Computation Challenges of NTK Matrices
16
C.I Storage and Compute Requirements
16
C.2 Data Precisions in Our Computational Pipeline
17
Implementation MEMDET Algorithm
17
D.I Optimal Sequence Processing of Blocks
17
Complexity and Performance Analysis of MEMDET
21
E.1 Computational Complexity
21
E.2 Data Transfer and Memory Considerations
22
E.3 Empirical Performance Evaluation
23
E,4 Concluding Remarks of Performance Analysis
24
B.5
Computational Accuracy of MEMDET
25
F
FLODANCE and Scaling Laws
26
F1
Background Material for Neural Scaling Laws
26
F2 Proofs of Lemma 1 and Proposition
27
F3 Derivation of FLODANCE Parameterization
28
Further Empirical Analysis of FLODANCE
29
G.1 Sensitivity and Robustness on NTK Matrices
29
G.2 Extension to Matérn Kernel Gaussian Process
30
H Comparison of Log-Determinant Methods: Complexity and Runtime
32
I Implementation and Reproducibility Guide
33
A. Nomenclature
We use boldface lowercase letters for vectors, boldface upper case letters for matrices, and normal face letters for scalars, including the components of vectors and matrices. Table A.1 summarizes the main symbols and notations used throughout the paper, organized by context.
14

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Table A.I. Common notations used throughout the manuscript.
Context
Symbol
Description
Example Value
Dataset
π
Number of data points (c.g,, images)
50,000 for CIFAR-10
d
Number of model outputs (e-g., classes, labels)
10 for CIFAR-10
Size of NTK matrix is m m where m nd
500,000 for CIFAR-IO
MEMDET
Tib
Number of row/column blocks (i.e., n blocks in total) Block size; typically b R m/ng, for x b submatrices Available memory capacity (in bytes)
e-g., 16 e.g., 31,250
b
C
B
Precision (in bytes) a floating-pov point number
FLODANCE
n
Number of sampled data points from n Size of sampled NTK matrix m, n,d
e.g,, 2000
m,
e.g., 20,000
no
Start of fitting interval [no. n.
e.g., 100
q
Truncation order of Laurent series
e.g., 3
co Pos Vq
Regression coefficients
SLQ
Lanezos iterations (Krylov subspace size)
e.g., 100
8
Number of Monte Carlo samples
e.g., 100
Variables
K,
NTK matrix with n data points (matrix size m)
K
Sampled NTK matrix with n, data points (matrix of size m2)
bn
Log-determinant of NTK with n data points (matrix of size m)
E,
Estimated log-determinant
Ln
Normalized log- -eg-rmeerminant L, = n L Vectorization of neural network parameters Dimension of 0
0
p
Functions
pdet logabsdet
gecudoddttrrminan (product of nonzero eigenvalues) Natural logarithm of the absolute value of determinant
B. Related Works in Numerical Linear Algebra
The study of block decomposition methods in numerical linear algebra has a long history. Classical texts such as Golub & Van Loan (2013) and Dongarra et al. (1998b) provide foundational discussions on block LU, block Cholesky, and LDL. decompositions, detailing their computational advantages and numerical properties. These methods have been widely used to improve computational efficiency, particularly in nigh-performance computing (HPC) settings, where recursive block LU (Golub & Van Loan, 2013, Section 3.2. 11 1), parallel LU (Golub & Van Loan, 2013, Section 3.6), and block Cholesky (Golub Van Loan, 2013, Scction 4.2.9) play a central role in large-scale matrix computations.
Beyond theoretical foundations, numerous works have focused on efficient implementations block factorizations, particu- larly for parallel architectures. Stark & Beris (1992) optimized block LU decomposition for hierarchical distributed memory, aiming improve data locality while maintaining parallel efficiency. Dongarra et al. (1979) pioneered high- perffrmmnome implementations of block factorizations, laying the groundwork for modern HPC systems. More recent studies, such as Ga- loppo et al. (2005) and Barrachina et al. (2008), have extended block LU methods to GPU-based environments, leveraging parallelism but still assuming that intermediate submatrices in memory. While these approaches optimize performance in parallel settings, they do not address the challenge of computing factorizations when the full matrix size far exceeds available RAM. Traditional block methods typically assume that at least some intermediate submatrices can reside in memory, whereas our method (MEMDET) explicitly operates under constrained memory settings, using an out-of-core hierarchical block processing approach.
To address the issue of matrices exceeding main memory capacity, Dongarra et al. (1998a) introduced concepts for parallel out-of-core LU factorization, focusing on efficient data movement between disk and memory. While their work demonstrates how out-of-core computations can be applied to LU factorization, their approach does not extend to log-determinant computations or hierarchical block-wise processing. Similarly, studies on many-core architectures (Venetis & Gao, 2009) and hierarchical memory-aware LU factorizations (Demmel et al., 1993) have improved computational efficiency, but
15

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Table C.I. Memory requirements (for various floating-point precisions) empirical NTK matrices of common datasets. The memory is computed as (nd))² 3, where is the training set size (second column), d is the number of classes (third column), and 3 is the number bytes per floating point value.
Matrix Size
Dataset
Training Set
Classes
1
float32
float64
CIFAR-10
50,000
10
0.5 TB
1.0 TB
T.
MNIST
60,000
0.72 TB
1.5 TB
TB
SVHN
73,257
1.1 TB
2.2
4.2 TB
ImageNet-Ik
1,281,167
1000
3,282,778 TB
6,565,556 TB
13,131,111 TB
Table C.2. Estimated compute time (in hours using NVIDIA H100 GPU) for NTK matrix computation.
Compute Time (hrs)
Dataset
Model
float16
float32
float64
MNIST
MobileNet
25
50
CIFAR-10
ResNer9
6
24
70
ResNet18
14
63
65
ResNet50
37
177
297
ResNet101
107
442
1178
none have been designed specifically for computing log-determinants under extreme memory constraints, making our approach distinct.
C. Memory and Computation Challenges of NTK Matrices
C.1. Storage and Compute Requirements
The empirical NTK serves motivating example throughout this work, it encapsulates key computational challenges associated with large-scale matrix operations. Several software packages have been developed to compute NTK Gram matrices for various neural architectures using automatic differentiation frameworks (Novak et al., 2022; Engel et al., 2022). However, the full formation these matrices remains computationally prohibitive, even on common benchmark datasets.
Table C.1 presents the storage requirements for NTK matrices corresponding to various datasets, highlighting their enor- mous size. For instance, even CIFAR-10 requires terabytes of storage, while ImageNet- 1 exceeds exabyles, making full NTK computation infeasible for most practical applications. Despite its theoretical importance, the NTK Gram matrix is rarely used as practical tool, with approximations often employed to mitigate computational and memory constraints. Minibatching is one common strategy, and batch-wise NTK approximations have been explored for model selection (Im- mer et al., 2023). Yet, extending these estimates to full datasets remains an open challenge. Alternative approximation techniques (Mohamadi et al., 2023) have been proposed, but their convergence is only guaranteed in spectral norm, limiting their ability to capture the full spectrum the NTK. In contrast, the log-determinant-aa key quantity in this enckesss information from the entire eigenvalue distribution, making its computation particularly demanding.
Beyond storage limitations, the computation time for NTK matrices also presents major challenge. Table C.2 provides estimated compute times for NTK formation across various models and floating-point precisions on an NVIDIA H100 GPU. Even for relatively small datasets like CIFAR-10, NTK computation is expensive, with higher- calculations significantly increasing runtime. For large architectures such as ResNet101, pouble-precision NTK computation can require over a thousand hours, making exact evaluations impractical without algorithmic improvements like those introduced in this work.
16

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
C.2. Data Precisions in Our Computational Pipeline
Our computations follow a multi-stage pipeline, with each stage involving distinct data precision formats and practical constraints:
1. Model Training. All neural networks (e.g., ResNet9, ResNet50) were trained using 32-bit precision, which is the default and standard practice in most deep learning frameworks such PyTorch.
2. NTK Matrix Computation. The NTK matrix computed from the trained model and stored in various precisions (e.g., 16-bit, 32-bit, and 64-bit, from the same pre-trained model). The "precision" of the NTK matrix, as referred to throughout the paper, reflects the compute and storage format at this stage. Due to the high cost of forming these matrices, it is often tempting or necessary to compute and store them in lower precisions. Our low-precision experiments highlight the pitfalls of mixed-pree -xecissoe in these cases as per Section 2.1, regardless of the downstream use case.
3. Jog-Determinant Computation. Regardless of how the NTK matrix was computed and stored (16-bit, 32-bit, or 64- bit), all log-determinant computations were performed in 64-bit precision across all methods (c.g., MEMDET, SLQ, FLODANCE). This represents a "best-case" mixed- preeci-oeeiiionn setup.
Since MEMDET entirely eliminates memory requirement barriers, it became practical to perform high-precision computa- tions (e.g.. 64-bit in stage 3) even on large matrices- thus mitigating common concerns about the overhead associated with higher-preed pecber-prcciionn formats.
D. Implementation of MEMDET Algorithm
The pseudo-code of the MEMDET algorithm is given in Algorithm D.1 (generic matrices), Algorithm D.2 (symmetric), and Algorithm D.3 (symmetric pusitive-definit)) each computing the log-determinant of matrix M. The log-determinants of leading principal submatrices of a matrix M (possibly after permutation, depending on the decomposition) can also be readily computed. For example, for symmetric dositive-definite matrices, the Cholesky decomposition M = LLL, with lower-triangular L, gives logdet (ogdet(M(a..a))) ((...))) = Σ (L,),, where [[..:a] is the ls x k leading principal submatrix.
The memory requirements of MEMDET are determined by a user-defined parameter, allowing it to run on any system regardless of available memory. For an m X m matrix, the algorithm partitions the data into an my x n, grid of blocks, each of size b b, where b 1 ++ L(m - 1)/ngl. I. The computation requires either 3 or 4 concurrent blocks in memory: for n, = 2, only 3 blocks are needed, requiring 3b², B bytes, while for ni 2, 4 blocks are required, increasing the memory usage to 4b 3 bytes, where 3 is the number of bytes per floating point.
Given maximum memory limit (in bytes), the optimal number of blocks n, is determined by the parameter r = m m/3/c. with the following selection criteria:
I S 1, the entire matrix fits in memory, ni 1.
If v3* three blocks fit in memory, ny = 2.
Otherwise, m4 = =[2r].
D.1. Optimal Sequence of Processing of Blocks
is important to select an ordering of the blocks to minimize data transfer between disk and memory. The order which the block M,, is processed should be chosen to minimize the reading of the blocks B (corresponding the index j) and C (corresponding to the index i). Ideally, from processing one block to the next, one should update only one of the matrices B or C, but not both, to reuse one of the blocks already loaded in memory. We formulate this problem of finding the optimal sequence of blocks as follows for the case of L.I /CCholesky decomposition at the k-th stage of the algorithm. The case for LU decomposition can be formulated similarly.
Let G(V,E) denote a complete undirected graph with vertices V == (k+1,.... no}. where is the set of all possible edges e (v, w) between the vertices u, E V. Each vertex in V represents the event of loading one of the blocks B <- M or
17

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Algorithm D.1: MEMDET: Constrained- Memoory Comp. of Log-Det (Case Generic Matrix)
Input :Matrix M of size m x m,
// stored on disk, may not be loaded on memory
Maximum memory e in bytes
Output: ::: logarithm of the absolute value of the determinant (logabsdet) of M.
a: sign of the determinant M
1 r < m/B/c
If 3: number of bytes per floating plotting-point II nat number of row/column blocks, making na grid of blocks.
2 r S then 1-
3 else r - 2/v3 then 4-
4 else nh [r
b 1 [(m 1)/mb
If Size of each block is most b b Il Accumulates log-abs- abs determinant diagonal blocks Keeps track of the parity of matrix
+ 7 7a+1
// Allocate memory for block matrices
8 Allocate memory for b b matrix A
if na then Allocate memory for x matrices B, C x matrix S
n6 then Allocate memory
// Create scratchpad space on disk, large enough store (ng - 1) blocks
n, then Allocate empty file of the size (m(m - b) b²)33 bytes
// Recursive iterations over diagonal blocks for = 1 to ny do
13
if h- 1 then A - M PLU
Il Load from input array on disk I1 In-placee LU decomposition with pivoting (written A.
f L + logabsdet(U)
a + PP) (U)
17
if < no then // Iterate over row of blocks from bottom upward
18
for i 74, to le + step -1 do
C- My
Load from disk (from input array if or from scratchpad if >>
20
C-U
// Solve upper triangular system in-place
// herate over column of blocks in alternating directions per row if even then :Jstart:Jnnd 4- (I: 1,76)
22
else Jena) +- (Tg, +
23
for j = Jaart to Jend step (-1)² i- do
// Load B from disk (input array if = and = ni, otherwise from scratchpad)
24
if ; Tl,, or j * Jate then B <- Maj
25
if n,, then
26
B L L-LL-pp 'PBB
If Solve lower triangular system in-place
27
if, k> or * jend then Maj <- B
I/ Write to disk on scratchpud
28
fi
1 and = then
M,
II Load from disk (input array if 1 or scratchpad if >> 1)
from
30
A A C
// Compute Schur complement
32
S M,
I1 Load from disk (input array if = 1 or scratchpad if >>>>> 1)
33
S S C
Compute Schur complement
34
M,j < S
// Write to disk on scratchpad
35 return f, a
18

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Algorithm D.2: MEMDET: Constrained-Memnne Comp. of Log-Dee Det (Case II: Symmetric Matrix)
Input :Symmetric matrix M of size m x 77,
I| stored on disk, may not be loaded on memory
Maximum memory in bytes
Output: : (Tg)g=1) permutation {1.....1)) m)
I/ permitations induced by LDI. decomposition
L = (Ey 4=1 1:::4 E = logabsdet(M (p,, with the index set I = (TI....,TTT)
I/ log det()) )
(gg)g=1'
(g
If Sign of detf.)
mBB/c
#I 8: number of bytes per floating point Il number of rowdcolann blocks, making ni grid of blocks.
2 i I then n1
3 else i $ 2/v3 then ni + else + [2r]
b 1+ [(m 1)/msJ Initialize arrays E Ih" and π (l..... m)
II Size of each block is at most x b Hold diagonals and permutations, respectively
// Allocate memory for block matrices
Allocate memory for b x matrix A
7 if ng then Allocate memory for matrices B, C if n > then Allocate memory m, > 1 then Define pointers B, C.
b matrix S
/1 Used for swapping memory: (B,, C,), will refer to (B, C) or (C,B)
// Create scratekpad space on disk, large enough store (y(m, 1)/2 blocks
10 if n, then Allocate empty file of the size (m(m b)/2 4b²),3 bytes
/ Recursive iterations over diagonal blocks n2
for
12
if then A 4 MgA
Load from inpia array on disk
13
A PL
#1 In-pplace LDL! decomposition with pivoting (written AJ
14
d[L+ (1)6b6))) + diag(D)
If Accumulate diagonals of to Il Accumulate permutation indices
15
T (a--((---))))) 116648] (k- 1)b permutation(P)
16
mi then
Iterate over column of blocks backward (right to left)
17
for = n4 to k + step -1 do
18
if ng even then (B,, C.) (B, C) else (B,, C.) (C,B)
I1 swopp B and memories
19
j Th then B. My
20
Load from disk (inpua arnay if = or from scratchpad if Solve lower-triangular system in-place
21
B. + L LPL.
22
C B.
Deep copy of the memory pointed by memory pointed C,
23
B, + DD
// Processing order of rows: first process now then from now le downwarr to -1
24
R (j. k
2,j
25
for i = R(1) to R(j k) do
26
i i then C, + M,;
I1 Load disk (input array if 1, otherwise seratchpad)
27
j = n4 then C. e 'P'C
Solve lower triangular system in-place I/ Write to disk on scratchpad
29
ifm, > 2 and < then M., < C,
30
if k and k+ then
31
A + M,;
I1 Load from disk (inpat array if 1 scratchpad if
32
A - C.B.
If Comppte Schur complement
33
clse
Load from disk (input arroy if scratchpad if
35
S S C.B,
I/ Compute Schur complement
36
M, S
# Write to disk on scratckpad
37 (4g. .p) +- 10,1)
38 for y to do 0g) ) (fg-1 log(ld,D., sgn(ddd))))
I1 dhy the q-th element of
39 return f,
19

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Algorithm D.3: MEMDET: Constrained Memory Comp. of Log- Det (Case III: Symmetric Positive-Definit Definite Matrix)
Input :Symmetric positive-definte matrix M of size X
I/ stored on disk, not on memory
Maximum memory c bytes
Output: :l == (6 4 logdet M (snoo)
r < B/c
#1 3: number of bytes per floating-poits If Ray number of row/column blocks, making grid of blocks.
2 r S 1 then f-
else ifr 2/v then ni + 2 else na |2rl
b + 1(m 1)/mm] Initialize array d E I
Size of each block most b #1 Holds diagonals
// Allocate memory for block matrices Allocate memory for
matrix A
7 if > then Allocate memory for x b matrices B, C ru then Allocate memory x matrix S if nu > then Define pointers B,, C,
11 Used for swapping memory; C.) will refer to (B, C) or (C,
// Create scratchpad space disk, large enough to store nn,n, 1)/2 4 blocks
10 m, > then Allocate empty file of the size (m((m b)/2 1b²),3 )33 bytes
/ Recursive iterations over diagonal blocks
11 for 1 11, do
12
i: then A +- Ma
I1 Load from input array on disk #1 In-place Cholesky decomposition (written A. II Accumulate diagonals of L to d
13
<- LLT
14
dp (p-
diag(L)
if k < ng then
// Ieeaate over colann of blocks backward (right to left)
16
for ny to + step -1 do
17
if n is even then B. else B. C
II Alternate pointer B, to switch between and
if j = ng then B. - B
19
Load from disk (input array if = 1 or from scratchpad if
20
B. B.
Solve lower-triungular system in-place
// Processing order rows: first process row then from row & 4-1 downward
21
R (j, 1, 2, j.... 2, for = R(1) to R(j - kk
22
23
j- then C. = B, else
II Shallow copy of pointer C, pointing to B.
24
25
I1 Alternate pointer C, between C and if M, - j is even then C, = C else C Load C. from disk (imput array kf = 1, Whp otherwise scratchpad)
26
C, Ms
27
j n, then
28
C, L----C L-' C,
# Solve lower triangular system in-place Write on scratchpad
n, and <
then M C,
30
ifi = + and
then
31
A M2
I1 Load from disk (input array if scratchpad if
32
A A C.B.
// Compute Schur complement
33
else
34
S - Mij S S - C.B,
If Load from disk (input array if scrutchpad if
35
// Compute Schur complement Write to disk on scratchpad
My +
37 c
38 for
1 to m do E, <- f--1 + 2 2log(d,)
11 dry is the 9-th element of d
39 return e
20

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
C6
V.
e2
A
1/2
03
0/
44
C6
e
e9
k+1
C2
2+2
04
Cg
f/3
C10
C10
C8
3
4+4
C10
Figure D.I. Left: Example the processing order of blocks for a symmetric matrix at the k-th hierarchical step. In this step, to process S - My. the memory blocks B and C are selected from the set V = {vi, 12, Vs v. Middle: The corresponding complete graph G(V,)) Right: The corresponding line graph L(G), with one possible Hamiltonian path highlighted in red, starting from the node 610 and ending at the node e1-
C+MMa M,j k l.... ng. Each edge in E corresponds to the event of processing the block M At the k-th stage of the algorithm, eventually, all blocks i, le + 1, ,,...,,, will be processed, so E consists of all edges of a complete graph, including self-loops, with [E] [ To illustrate this concept, consider the example in Figure D.1. The left panel 2
ng
depicts the k-tth iteration of the algorithm for symmetric matrix, where the matrices B and C have four blocks to choose from the set V {u1, 12, U3 ). The corresponding graph G is shown in the middle panel of the figure.
The goal is to select an ordered sequence (ep). I,..., E| of edges such that each two consecutive edges ep and ep+1 in the sequence share a common vertex. This ensures that from processing one block to the next, only one of B or C is updated, while at least one block reused from the previous step.
To find such a sequence of edges, we define L(G). the line graph of G (also called the edge-to-vertex dual), where each vertex of T(G)) represents edge of G. Two vertices in (G)) are adjacent if and only if their corresponding edges in G share vertex. Thus, any Hamiltonian path in L(G) yields ordered edge sequence fulfilling our requirement.
As illustrated in Figure D.I, the right panel depicts the line graph of the given graph shown in the middle panel. with a possible Hamiltonian path highlighted in red. This path directly translates to the processing order blocks shown in the left panel. Notably, all valid Hamiltonian paths must terminate at the node é,, representing the block M This specific end point crucial as it allows for seamless transition to the next iteration (i.e., the k + iteration) without the need explicitly load the matrix A, would already be available in memory from the last processing block of the k-th iteration when S <- M was processed.
Given that G is complete and therefore Hamiltonian, it follows that its line graph T.(G) is also Hamiltonian. This implies the existence of at least one (but possibly many) Hamiltonian paths. Crucially, all Hamiltonian paths in L(G') have the same length. Consequently, any sequence of blocks derived from a Hamiltonian path constitutes an optimal solution to our problem. Thus, the block sequence presented Figure equivalent optimality to any other sequence obtainable from a Hamiltonian path.
The same problem can be formulated for the block LU decomposition, with the modification that G is complete and balanced bipartite graph G(V, V. E); however, the same logic and conclusion follow.
E. Complexity and Performance Analysis of MEMDET
E.1. Computational Complexity
Table E.1 provides a detailed breakdown of the computational complexity of MEMDET for generic matrices (second column, using LU decomposition) and symmetric matrices (third column, using LDI. or Cholesky decomposition). The
21

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Table F. Breakdown of computational complexity for MEMDET. The table presents the number of operations performed and the FL.OP
count per operation for generic matrices (LU decomposition) and symmetric matrices (LDL or Cholesky decomposition). The last row shows the total complexity, which remains independent the number blocks Rdy.
Generie Matrix
Symmetric (Positive-Deffnite) Matrix
Operation
Num. Operations
FLOOPs per Operation
Num. Operations
FLOOPs per Operation
Matrix Decomposition
12b
1/8²
1
+
b
71
102 + 12
Solve Lower Triangular System
2
202¹
102²
2 2n
m
1000
102
Solve Upper Triangular System
2
2
1
Full Matrix Multiplication
1
ng
+ !
b²
2 716
1
+
6²
Gramian Matrix Multiplication
2 216
Total Complexity
3
m
2
3
4 m
2
1
m
operations are categorized into matrix decomposition. solving triangular systems, and matrix multiplications used to form Schur complements. Each operation's 's complexity is given in terms of the number of times it is performed and the FLOP count per operation. In this analysis, one FLOP refers to fused multiply-add (FMA) operation one multiplication and one addition-as counted by modern GPU benchmarks. The table lists a unified complexity column for symmetric matrices, encompassing both LDL and Cholesky decompositions. While LDL includes additional operations such row permutations and diagonal scaling via D, these are excluded from the FLOP count due to their negligible cost relative to the leading terms.
The complexity of each operation is given by the number of times it is performed (a function of nn)) multiplied by the FL.OP count per operation (a function of the block size 6). Substituting 6 = m/n, into these expressions, the total complexity, obtained by summing across all operations, simplifies such that n, cancels out, as shown in the last row of Table E. Thus, the total computational complexity of MEMDET is independent of the number of blocks ng. and is identical that of conventional factorization algorithms where ni = 1.
Although the total computational cost remains the same, the contribution of individual operations shifts as n, increases. When n;- = 1, the entire computation consists solely of a matrix decomposition. As n, increases, the decomposition cost decreases while additional operations, such as solving triangular systems and matrix multiplications, account for larger fraction of the total complexity. At the extreme case of n4, = m, the algorithm consists primarily of matrix multiplications. This transition is illustrated in Figure E.1 (left panel). where the contributions of matrix decomposition, triangular system solving, and matrix multiplication plotted functions of ng. The total computational complexity, shown the black curve, remains constant, while the distribution of work among different operations shifts as ni, increases.
E.2. Data Transfer and Memory Considerations
While the total FLOP count is independent of ni. the number of data transfers between memory and disk increases with the number of blocks. Table E.2 summarizes the number of blocks read from disk to memory and written back to disk, as a function of ny. The actual volume of transferred data is obtained by multiplying the number of blocks by the block size, b² 3, where 3 represents the number of bytes per floating point. The right panel of Figure E. illustrates the total data transfer volume relative to the original matrix size.
For ni S 2, the entire computation performed in memory, avoiding disk I/O and thus requiring no scratchpad space. However, for n2 >> 2, the computation utilizes scratchpad space, and data transfer overhead increases approximately as On))). meaning that choosing an excessively high n, introduces unnecessary I/O costs. Despite this, the hierarchical design of MEMDET efficiently schedules block transfers, mitigating excessive data movement.
Table E.3 provides an analysis of the required memory and scratchpad space. The number of blocks stored in memory and on disk is determined by ny, with the actual space usage obtained by multiplying the number of blocks by the block size. By adjusting ng, MEMDET can be configured to run within any given memory constraint, making it adaptable to systems with limited memory.
22

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
[a] Complexity per Algehraic Operation
[b) Data Transfer
Symmeteic
Roal
Generie Symmetric
Matrix Decompustion
Sizee (bytes) / (m (m)
Solvo Thianguliz
ELOP m
Matris Multipiistion Total
2
1/2
0
173
R
Figure E. Theoretical computational complexity of MEMDET as function of the number of blocks na- The left panel shows the contributions of matrix decomposition, solving triangular systems, and matrix multiplication to the total complexity. The black curve represents the total computational cost, which remains constant, while the colored curves illustrate how the workload shifts across operations as ni, increases. The right panel displays the total data transfer volume (normalized by the original matrix size) for different n,, highlighting the increasing cost of disk 1/O the number of blocks grows.
Table E.2. Number of blocks transferred between disk and memory during MEMDET execution. The total data transfer volume is obtained by multiplying the number transferred blacks by the block size, b² bytes, where b = m/ng. Read operations occur all cases, while write operations to the scratchpad are only required for n, 2.
Operation
Generic Matrix
Symmetric Matrix
Read
] nn
nh
+ n0
n
2
+ n
My <
0,
7h
Write
5 1, Yb, > 3nd
1 +
2
Tu
7h
Table E.3. Number of concurrent blocks that must be allocated in memory (first row) and the total number of blocks allocated disk (second row) during MEMDET execution. The total required memory and disk space are obtained by multiplying number allocated blocks by the block size, 6² bytes. While the number of concurrent memory-residen blocks remains fixed, the total number of blocks allocated on disk increases with n, 2.
Hardware
Generic Matrix
Symmetric Matrix
Memory
or
3
0
nh S
0, 2 rig
Th S 2
Scratchpad
T4p
ns
nb
4
70 > 2
E.3. Empirical Performance Evaluation
To validate the theoretical complexity and memory analysis, we conducted empirical evaluations on SPD matrices of various sizes, ranging from m 20 to 2 The largest matrix tested was chosen to match the memory capacity of 64 GB system, allowing for direct comparison between MEMDET and conventional algorithms. For each matrix size, the algorithm was executed with different numbers of blocks, n& = 1,,, 8. where n,, corresponds standard full-matrix decomposition with the entire matrix loaded into memory. Fach experiment was repeated 10 times, and the mean and standard deviation of the profiling measures are reported.
Figure E.2 presents the experimental results. The left panel shows peak memory allocation. measured using a memory profiling tool, which precisely matches the theoretical predictions. As expected, when n, 1, the required memory
23

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
(a) Required Mernory
(b)
Performance Over Number of Blocks
(c)
Performance Over Matrix Size
10"1
10
10
Allocated Mennay (bytes)
1010
109
10a
Proness Time (sec)
Process Time (soc)
84-8
10"
10"
10
m
225
2
10"
mn 24
Ry-b
21
10'
21
10
212
109
210
210
10a
10P
212
255
Number of blocks da
Number blocks We
Matrix size
Figure E,2. Peak memory allocation (a) and CPU processing time (b, c) for MEMDET on symmetric positive-definite matrices of size m 200 ....,, 266, using Algorithm D.3. The matrices were processed using an n4 x ni, grid matrix blocks, where ni 1, 2, .,
2
2
2
2
Compute
0.3
Write
Prones Time ne (1s)
Proces Time me? (ts)
Process Time me (18)
Process Time m (us)
Reed
1.5
0.4
0.2
0.3
0.2
0.5
3
176
116
176
Figure E.3. Breakdown of MEMDET runtime (using Algorithm D.3) into computation and data transfer times, normalized by m for m 2 ...... 2 16 The total process time consists of reading from disk to memory (maroon), writing from memory to disk (light tan), and CPU computation (ochre).
equals the original matrix size, while increasing n, reduces the memory footprint. The middle and right panels display the measured process time as functions of ni, and respectively. At large m. the difference in process time across varying nn diminishes, indicating that the increase in data transfer cost does not significantly impact overall runtime.
To further analyze this effect, Figure E.3 decomposes the process time into computation time and data transfer time. At small m, the runtime dominated by disk I/O, but as m increases, computation time becomes the dominant factor. This confirms that for sufficiently large matrices, the performance of MEMDET approaches that of conventional in-memory methods.
E.4. Concluding Remarks of Performance Analysis
MEMDET maintains the same computational complexity as conventional factorization methods while distributing compu- tations across blocks. The total FLOP count remains unchanged, but increasing n, shifts the workload between operations (i.e., from matrix decompositions matrix multiplications). However, increasing n, also increases data transfer overhead, requiring balance between reducing memory usage and minimizing disk I/O.
The scheduling design of block operations optimizes memory usage while limiting unnecessary data movement, ensuring that MEMDET remains efficient under constrained memory conditions. By allowing users to specify a memory limit, MEMDET enables the processing of arbitrarily large matrices on systems with any limited memory size.
For large-scale applications, where conventional methods exceed memory capacity, MEMDET provides a practical alterna-
24

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Table E.4. Relative of MEMDET (using ni & blocks) and Numpy's a gg with respect NumPy's 8 sloget (used baseline)
across different models, dataset sizes, and precision formats.
Quantity
Model
ResNet9
ResNeto
ResNet18
MobileNet
Dataset
CIFAR-10
CIFAR-I0
CIFAR-10
MNIST
Subsample Size
n 1000
2500
1000
2500
MEMDET (16-bit)
2.2 x 10
5.4 10 10
5.4 10
8.7 x 10 16
Rel. Error
MEMDET (32-bit)
4.2 x 10
4.3 10
1.7 10 10
8.4 10 12
MEMDET (64-bit)
1.4 10-8
5.4 10
10-9
1.7 10-6
eigh (16-bit)
1.6 10
1.7 10
2.8 10
1.8 10
eigh (32-bit)
5.0 x 10
1.1 X 10
9.2 10
4.5 x 10
eigh (64-bit)
6.1 x 10
7.5 X 10
1.4 x 10
1.0 x 10
(a) Log-Determimancs Using MEMDET, eigh, and slogdet
Error of MEMDET and eigh Compared slogdet
10
10-6
32 bit (MEMDET)
32-bit (eigh)
52-bit (slogdet)
6d-bit (MEMDET)
10
6d-bit (eigh)
Relative Error
1
bit (slogdet)
logdet((K
10-10
10-12
bet (MEMDET)
32-bit
(eigh)
6d-bil (MEMDET)
6d-bit (eigh)
10a
10-14
29
20
241 212 28
29
210
11
Figure E.4. Log-determinants of growing NTK submatrices from ResNet9, computed in both 32-bit (red) and 64-bit (black) formats. Each submatrix has size m nd with = 10. (a) Comparison between MEMDET, s olddet, and e eigh for each input matrix. (b) Relative of MEMDET and gi with respect s logdett Regardless of the matrix precision, all log-determinant computations are performed in 64-bit precision, and errors remain well below 10"
tive. The experiments confirm that while data transfer overhead exists, it does not significantly impact runtime at large m. making MEMDET viable solution for large-matrix computations standard hardware.
E.5. Computational Accuracy of MEMDET
We validate the numerical accuracy of MEMDET by comparing its log-determinant output across various precision formats against two standard in-memory methods: numpy. winall.......... wwiiingg.slsgsee...... slogdet and numpy 1 linalg.a... eigh (from which the oog-determinant is computed as the sum of the logarithms of the eigenvalues). Table E.4 shows that MEMDET matches both methods with relative errors between 10" and 10 well within the margin of numerical agreement between the
baselines themselves.
To assess behavior at larger scales, we compute log-determinants growing NTK submatrices derived from ResNet9, with the matrices formed in both 32-bit and 64-bit srating-point formats. Figure E.4 confirms that MEMDET remains in tight agreement with bloddet and eigh across all scales. Notably, discrepancies between the 32-bit and 64-bit curves attributable solely to differences in the data precision of the underlying input NTK matrices.
Finally, we evaluate the effect of MEMDET's block size parameter nb, which determines the number of memory partitions used during computation. Figure E.5 shows that even for NTK matrices of size up to 100,000, increasing the number of blocks has no measurable impact on accuracy: all results remain within 10" to 10 relative error compared to the full in-memory LDI. decomposition.
16
25

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
x 100 (a) Using Various Number of Blocks fig
(b) Comparison with the Case of ny
10-6
2
n1 2
Tiy
2
Relative Error
10-10
T1y 24
logdet (K.)
2
fb 2
My 2
10 10-12
2
10
0
10-16 16
10
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
x10
x10
n
Figure E.5. Effect of the number of blocks used in MEMDET. (a) Log-dcterminants computed for NTK submatrices from ResNet50, with different values of my. (h) Relative with respect to the conventional LDL dccomposition (n, 1). All results match within 10 to 10 12 accuracy, indicating high numerical stability.
F. FLODANCE and Scaling Laws
In this section, collect supporting theoretical content related the FLODANCE algorithm. Appendix F.1 reviews background material neural scaling laws. Appendix F22 provides the proofs of our main lemma and proposition. Finally, Appendix F.3 presents the derivation of the FLODANCE parameterization.
F.1. Background Material for Neural Scaling Laws
FLODANCE builds recent theoretical developments in neural scaling laws that characterize generalization error in terms of kernel eigenvalue decay and source smoothness. Here, we review the key assumptions and results from Li et al. (2023) that underpin our bias analysis. These assumptions concern the spectral properties of the kernel, the embedding behavior of the corresponding RKHS, and the regularity the target function. We restate them below for completeness.
The following result from Li al. (2023), stated under a noiseless setting consistent with our framework, relies on the following assumptions:
Assumption F.1 (Eigenvalue Decay). There exists 3 > 1 and constants Cg. C3 such that
Czt X; S Cgt
(F.1)
where the X; are eigenvalues of the kernel k X x X IR under the decomposition guaranteed by Mercer's Theorem:
<( x' (=)(()))))1))
(F22
We need to define an embedding index associated to certain interpolation spaces that arise as the range of fractional powers of integral operators. In order to do so, we define the integral operator T I,² - I,2 that acts as the natural embedding of a RKHS H associated with our kernel K, precomposed with its adjoint. That is, T is the integral operator given by
(T)) / n(2, 6(2, (ffe) dd((x')
where f is the marginal distribution of on X, where is the source distribution on (X' x Y) underlying the dataset. The operator T' can be decomposed by the spectral theorem compact adjo-ntjoin operators via
==
Σ A e.
(F.3
26

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
For 2 O, this lets us define the fractional powers T LL + L2 of the operator T to satisfy
80
T() - Σ e
(F.4)
The interpolation space [7(]* associated to T can then be defined as
Ede < x Jo i=1
[H] range(T (3
aid ex
(F.5)
We now say that H has an embedding property of order a E (0, 1] if [H] can be continuously embedded into L Define then the operator norm, which has the form (see (Fischer Steinwart, 2020))
||[ + L
ess sup Σ {x)(2)) DEA i=1
We now have the following assumption on the embedding index, which is known to be satisfied if the eigenfunctions are uniformly bounded (Steinwart et al., 2009).
Assumption F2 (Embedding index). The embedding index a 1/3, where is the eigenvalue decay in (F.1), and og is defined as
ap inf {a [|[*]* + L = Ma < oo)
Finally, we have the following assumption on the smoothness of the source function fö = Eppy 1 æc), which is more precise characterization than requiring it to belong to some interpolation space.
Assumption F.3 (Source condition). There exists
0 and sequence (a,)<21 for which
00 S
a 2 ez:
1=1
and 0 < c S |a,1 < C' for some constants C.
These assumptions are required for the following theorem, taken from (Li et al., 2023).
Theorem F.1. Under Assumptions F.I to F.3, fix > 1 and suppose that x n n" "frr 0 2 B. Then
B = En)) = Opoly (n min(s.2)8).
F.2. Proofs of Lemma and Proposition
Proof of Lemma 1. Fix 0 and write K = Kn + almm where Imm is the x identity matrix. In block form, K contains K according to
K x(x xn)) | |=( T =("n;"))) ald
K
1,,,,,,2, ,.,)) = (x(2, = E (ndddxx where x( determinants differ by the Schur determinant:
l)dxd Consequently, since K; and K are both positive-definite, their
(ct() det det((K det (((F,,),) + aI. "T "
Combining the Sylvester rank inequality with Corollary 20 from (Ameli Shadden, 2023), we take a 0 and observe that
pdet(K,,)) < pdet(K,..) pdet(Cov(ffr,)) f(x;)) for i I,..., 1)).
(F.6)
27

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
This lets apply the AM-GM inequality and then bound the Frobenius norm in terms of the nuclear obtain
(d,)((K,)) pdet
< (nde((oo(((e)))]((t)) for - 1,.... 1)) -
Ay r
Σ X3
pdet(K,..) K 1,,,))
j=1
trace (Cov(f(zz, f(z;) for i 1,..., 11))
1/2
F(n)",
where the A are the non-zero eigenvalues of the covariance matrix, and r is its rank.
Proof of Proposition 1, From equation (3),
logdet(K, (I.)) - logdet(K,,11 = logg C' - // 2/log n + log/1 1 + (1))
= llog C - v log n + 0p()),
(F.7)
and
logdet(KK,) logdet(K,) (n 1) log C v log(n!) + op()).
Letting Co = log C' logdet(K,) and dividing by implies the first result. For the second result, we replace the (I())
term in (F.7) with - Consequently,
n
In In
- cn log(n!) (n 1) 1/2
(1-2)
8,
vn
and from (Billingsley, 1961), (n - 1) 11-1/2 Ei1 8; converges weakly to normal random variable with zero mean and variance o² =E[
F.3. Derivation of FLODANCE Parameterization
Here derive the equation (6), which leads to the numerical procedure in Algorithm To show that the coefficients co and Vo...., Vy
be obtained using standard linear regression procedures, is necessary determine the form of the covariates n. Letting Vn = vo + Ei=1 V:n" Proposition implies that, asymptotically in n,
1 I, = 1 + 1 (1- n2
9
1/2 + Vn
log(n!)
iid o(0,03).
co
1/
Cu.
Cn
~
n
n
for some Co: Vo...., y and Rearranging, there is
n (Ln L1)
n
9
1/
log(n!)
CD
1/0
+
+ En
Vn
i=1
Equivalently, the above relation can be recast linear regression problem:
++1
yn nO +
v VI + En-
where co. Uo. M, are the regression coefficients for the covariates F. defined as
-
= 0, = 1,...g. g.
E Fnd
log
11
28

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
10
(a) ResNet9
x10
(b) ResNetã0
110
3.0
3.5
Left Onlivate
Laft Ordinate
2.5
Empirica:
Asymptore
Enpicioal
Asymggrate
3.0
Fit
2.5
2.5
2.5
2.0
Right Ondinate
Right Ordraate
Abeolute Ernor
Absolute Ferorr
<5
1.55
1.5 +
1.5
1.5 K5
1.0
1.0
0.5
0.5
0.0
0.0
0.0
0.0
x10+
x10*
Figure G.I. Log-determinant f, for n2 1,..., 50,000, corresponding to x NTK submatrices where nd and d 10. from 64-bit NTK matrices ResNet9 (a) and ResNet50 (6) trained on CIFAR-IO with 50,000 datapoints. Values are computed using MEMDET (Algorithm D.2) with LDL decomposition (black curves, overlaid by colored curves). The orange curves represent theoretical based on the parametrization derived in Algorithm 1. Fitting performed globally over the entire interval (ng. n) (1,5 10"), demonstrating the accuracy of the thcoretical model. The blue curve, corresponding the right axis (scaled to the-thousandth of the left axis), shows the absolute crror.
The target variable for regression is given by
n
yn
Vn- (Ln - L1).
We note that in numerical implementations, the log(n!) should be evaluated using the log-gamma function: log I'(n
1) log(n!).
G. Further Empirical Analysis of FLODANCE
This section investigates the accuracy, robustness, and generality of FLODANCE through two complementary studies. Ap- pendix G.1 revisits NTK matrices, quantifying the method's ' sensitivity to fitting-intn interval length and to random subsampling. Appendix G.2 then applies FLODANCE to a multi-output Gaussian process with a Matérn kernel, demonstrating that the same scaling-law machinery applies well beyond neural tangent kernels.
G.1. Sensitivity and Robustness on NTK Matrices
Global Fit of FLODANCE on NTK Matrices. To validate the theoretical parameterization given Appendix F.3, we compute the fog-determinants NTK submatrices of size I,...,7 n from ResNet9 and ResNet50 trained CIFAR-10, with 50,000. The degermmiinsssss presents the empirical results: black curves (largely overlaid by orange) show the computed degemmeermnnnts, while the orange curves show the theoretical fits derived from the parameterization in Algorithm 1. The fitting is performed globally over the entire interval, demonstrating the accuracy of the theoretical model in capturing the log-determinant behavior. This global fit complements the extrapolation-but small subset and extrapolated to the full range. The near near-perfect overlap between the curves highlights the quality of the fit, with errors remaining below 0.05% for most of the interval.
obtained using MEMDET (Algorithm D.2) with LDL decomposition. Figure G.1
application shown earlier in Figure 5, where FLODANCE is trained on a
Uncertainty Under Subsampling. To assess the robustness of FLODANCE predictions, we evaluate how subsampling affects degermmermnnnt estimates. In Figure G.2, we generate 15 independent subsamples of NTK matrices of size nd, with d = 10 and 10,000, from ResNet50 trained CIFAR-10. Panel (a) shows the exact log-determinants f, computed using MEMDET for each subsample. The black curve denotes the ensemble mean, while the gray shading around the mean (barely visible) indicates the standard deviation across subsamples. The red curve (right axis of panel) shows the normalized standard deviation, which remains below 0.1% throughout, indicating remarkable stability of the -egedeetrssnants under subsampling.
29

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
X10º
(a) Variation of Empiriaal Ensemble
x10
(b) Variation of Prediction
1%
2.0%
Laft Onlinate
Left Orlinate
Mean STD
Empirioal (Masa))
Asyuptote
(Mran)
1.5%%
Right Ordrrats
Vaar (al) [
Extrapplation (Maan;
Normollaed
Right Ordinate
-5
0.15%
Red. Riror (Mean)
1.00%
Rel. Esvor (STD)
2
0.5%%
0.01%
0.0%
0.0
0.8 1.0
0.8
1.0
X10*
x10*
Figure G.2. Sensitivity of log- determmennn estimates to subsampling variation. (a) Exact log-determinants f, of 15 randomly subsampled NTK matrices size m X m, m nd, with number of classes d' and data points 10,000, computed using MEMDET on ResNe150 trained on CIFAR-10. The black curve denotes the mean, and the shaded gray (barely visible) shows the standard deviation
across subsamples. The right ordinate shows normalized standard deviation (red), which remains below 0.1%. (b) Predicted log- determinants using FLODANCE fitted over small interval (No, n.) (1, 10") (yellow), and extrapolated to (n,, n) (10", 10* (red).
The left axis shows the predicted mean; the right axis shows the relative error (bluc), with mean error under 1% and variation across ensemble (shaded blue).
Panel (b) evaluates the predictive performance of FLODANCE the same ensembles, where the model is trained on the interval (No: n..) = (1, 10²) and extrapolated to (n,, n) (10", 10*). The predicted mean log-determinant is shown in red (extrapolated) and yellow (fitted), closely tracking the ensemble mean of the exact log-determinants (black, largely obscured). The right axis displays the relative error: the blue curve denotes the mean, which stays below 0.5%, and the shaded region shows the standard deviation across ensembles, which remains comparably tight. These results confirm that FLODANCE exhibits both accuracy and robustness under random subsampling.
Effect of Fitting Interval Size. We investigate the sensitivity of FLODANCE to the choice of the fitting interval size ng which governs the trade-off between computational cost and extrapolation accuracy. Recall that FLODANCE fits model to log-determinants computed over the interval [no n,,] and extrapolates to the full range [M,, nJ. where 7o and 50,000 in this experiment. Figure G.3 evaluates this trade-off on NTK matrices from ResNet50 trained on CIFAR- 10. As n, increases, more data is used for fitting (increasing the cost), but less extrapolation is required (increasing the accuracy).
Panel (a) shows the root-mean-sque square (RMS) error of the fit the training interval [mo, n,). which increases with n, due the growing number of points match. Panel (b) displays the RMS error in the extrapolation region [n., ni]. which 50,000,
decreases as the extrapolation range shrinks. Finally, (c) plots the relative error of the prediction the endpoint which drops sharply from over 50% to under 0.1% n, increases. These results illustrate the central design principle of FLODANCE: by choosing moderate n, one can achieve accurate predictions while keeping the cost of computing exact hog-determinants limited to smaller submatrices.
G.2. Extension to Matérn Kernel Gaussian Process
Finally, we demonstrate the generality of FLODANCE beyond NTKs by applying it to multi-output Gaussian process with a Matérn kernel, widely used in spatial statistics due to its tunable smoothness. The isotropic Matérn correlation function of Matérn (1960) (see also Stein (1999, p. 31)) between two spatial points x, x' E IR. is given by
2--V p(x, x' aa = r(v) (V22 |
(vzyller
x a''le) K, a
x x 112
22
where I(-) is the Gamma function and K.()) is the modified Bessel function of the second kind of order v (Abramowitz & Stegun, 1964, Section 9.6). The hyperparameter v modulates the smoothness of the underlying random process, and the hyperparameter a > is the correlation scale of the kernel. We construct the multi-output covariance using the
30

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Fitting Interval [1,
(b) Extrapolation Interval (R,,
104)
[e] Prediction Ervor
10a
Relative Eevor
RMS Error
RMS Erre
10%
101
10a
10-1%
10
14 5.x 10"
10
10P
R4 5xx 101
10
10
R2 5.6 10"
10
Figure G.3. Sensitivity of FLODANCE to the choice of fitting interval size n., Based on ResNet50 trained on the full CIFAR-10 10 dataset with n = x 10 data points and 10 classes, resulting in NTK matrices of size m = rad = 10 FLODANCE extrapolates log- determinants by fitting a model on submatrices of size no = to n,,, and extending this fit to larger sizes up to n. (a) Root-mean-squal square error (RMSE) of the [n., n.. which decreases with n,, (c) Relative error of predicting the log-determinant
interval [1. n,). showing increasing fitting error as ni, grows. (b) RMSE of extrapolation in the interval
5 10", again decreasing with n.,
(a) Scale Law
x 111
(b) Prediction Based Scale Law
10a
1.0%
4.22
Left Ordinate
0.0
Emsiiicaa
Asymplale
0.8%
(et( /eett(K..)))
0.2
Fit
Extrapolation
0.4
Righi Ohilinate
0.6%
10-1
Relativ
0.6
0.4%
0.2%
1.0
10-9
-1.2
0.0%
10
10
10
10a 10 0.0
0.2
0.4 0.6
n 0.8
1.0
72
x10*
Figure G.4. Application of FLODANCE to multi-output Gaussian process with a Matérn kernel. generate n 10,000 spatial locations in I and assume d 10 10-dimensional output per location, resulting in covariance matrix of size m nd 100,000. The covariance structure follows Matérn kernel with smoothness v = 1.5 and seale parameter a = 0.04, combined with a lincar model of coregionalization (LMC) for output covariances. (a) Scale law illustrated by the ratio of successive determinants
increasing submatrix sizes. (b) Log- iegddttrrminnt prediction using FLODANCE. The black curve (left axis, largely obscured) is the exact log- ogg-determinant f, computed by MEMDET. FLODANCE is fitted on [1 n, 10 (yellow) and extrapolated to [n., n 10") (red). The blue curve (right axis) shows the relative error of prediction, which remains below 0.4%.
(non-homogeneous) linear model of coregionalization (LMC) (Gelfand et al., 2010, Section 28.7) given by
<(æ, x' aa'la,, (a,v) (((x)+((()) p(æ x' v),
where R dddd is the local covariance of the model's 's vector output of size d. In this model, we use the matrix square root, 2 to project the scalar Matérn correlation into the x d coregionalization space while preserving positive- definiteness.
In our experiment, we generated n 10,000 random spatial points dimension p for Gaussian process with output dimension d 10, resulting in covariance matrix of size nd 100,000. We set the Matérn correlation scale to 0.04 and the smoothness parameter to 1.5. The local covariance fields a (<x) were instantiated as random symmetric positive-definite matrices drawn from a Wishart distribution.
Figure G.4 illustrates the effectiveness FLODANCE on this Matérn-basee covariance structure. Panel (a) shows the empirical scale law via the ratio successive log degermmrmmnstss, which exhibits a smooth trend consistent with the theoretical behavior observed for NTKs. Panel (b) evaluates the extrapolation accuracy of FLODANCE: the black curve denotes the
31

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Table H.J. Computational complexity of log-determinant estimation methods. The first row corresponds to the exact method (MEMDET), while all other methods are approximations. One FLOP is counted as a fused multiply-ade add (FMA) operation.
Method
Approach
Complexity
Description
MEMDET
Direct factorization
m
1m1² + m2
m Full matrix size, = nd
FLODANCE
Submatrix extrapolation
ng Number of data samples from n fm; mm + mmm + (a 3)+
mg Sampled matrix size, m, = n,d
Laurent series Iruncation order
SLQ
Stochastic trace estimation
(m²t ms)
Krylov subspace size
Number of Monte Carlo samples
Pseudo NTK
Cross-class block reduction
(
22
+
+ m
m Full matrix size
d Number of model outputs
Block Diagonal
Class-wise block approx.
ma 1md + m
m : Full matrix size
d Number of model outputs
Table H.2. Wall-clock runtimes (in seconds) for various log-determinant approximations é,, using the same models and confi
as in Table 2,
Quantity
Model
Configuration
ResNet9
ResNe19
ResNet18
MobileNet
Dataset
CIFAR-10
CIFAR-10
CIFAR-10
MNIST
Subsample Size
n = 1000
n = 2500
n = 1000
n = 2500
Direct Computation (16-bit)
6.69
49.35
5.50
70.90
Direct Computation (32-bit)
7.08
49.39
5.87
54.05
Runtime (seconds)
Direct Computation (64-bit)
7.09
51.57
5.97
51.22
SLQ
8.904
67.70
22.18
148.8
Block Diagonal
0.003
0.009
0.003
0.008
Pseudo NTK
0.015
0.077
0.014
0.082
FLODANCE FLODANCE FLODANCE
ng 1, ng 50
0.008
0.008
0.008
0.009
1, n. 100
0.017
0.018
0.016
0.017
no = 300, n2 500
0.377
0.350
0.354
0.376
exact log-determinant é, computed via MEMDET, while the colored curves show the (yellow) over the small interval [1, ng 1001 and the extrapolation (red) to in,, 10*). The right axis displays the relative error (blue), which remains below 0.4% throughout. This example demonstrates FLODDNCE''s flexibility in handling more typical structured kernels.
H. Comparison of LogeDeterminant Methods: Complexity and Runtime
We summarize the computational characteristics of the log-determinant approximation methods evaluated this work. Table H.1 compares their theoretical computational complexity, highlighting how each method accesses or approximates the full matrix whether through exact computation with full matrix access (MEMDET), or through approximation strategies such as subsampling (FLODANCE), matrix-vector product oracles (SLQ), or blockwise approximations (Pseudo NTK and Block Diagonal).
Table H.2 reports wall-clock runtimes (in seconds) for several NTK datasets, using the same models and configurations as in Table 2. It compares direct 64-bit computations performed on matrices stored in 16-, 32-, and 64-bit precisions to various approximations, including FLODANCE, SLQ, and others. All measurements were conducted on the same hardware under comparable conditions. FLODANCE consistently achieves sub-second runtimes even on subsamples of size n = 2500, outperforming other approximations while avoiding the instability issues associated with SLQQ.
32

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
I. Implementation and Reproducibility Guide
We developed Python package detkit' that implements the MEMDET algorithm and can be used to reproduce the numerical results of this paper. A minimalistic usage of the detkit. memdet function is shown in Listing I.1, where the user can specify various parameters: the maximum memory limit meax_ the structure of the matrix via the assume argument-set to gen for generic matrices (Algorithm D.1), sym for symmetric matrices (Algorithm D.2), and spd for symmetric positive-defimtt matrices (Algorithm D the data provided in full in its lower/upper triangular form (tri angle), the arithmetic precision used during computation (mi imixed_ xeeieeeeellleeeeeeeeeeeeeeeeeeee s sion on), the location of scratchpad space on disk (scratch) ch_dir),, and enabling parallel data transfer between memory and disk (parallel_io) The function in this example returns the log and sign of the determinant of the full-size matrix (1d, s gn), along with the diagonal entries of matrix D (diag) and the array of permutation indices (pezm) for the permutation matrix P from the LDL decomposition PTMP- = LDLT.
Listing I.I. A minimalistic usage detkit package. The function memdet: computes logabsdet MM) using the MEMDET algorithm.
Instail detkit with "pip install detkit"
from detkit import memdet
import zazr
NTK matrix M disk M oarr, filerame... maara', mode='r" F()
Compute lM) and (g)( with Algorithn D.2
Assume M is symmetric and only its upper triangle part is referenced.
1d, sign, diag, perm, info memdet (M, max max_mem=' 32GB', assume=' sym', triangle='u', overwrite-False, mixed__rrcciioon flcatt44, acratch_dir- /tap', parallel "eessoratttess.rreeeee... verbose-Irue, ieturn_info-True, flops-True)
The next example, shown in Listing I.2, demonstrates the use of the FLODANCE method via the F Fit class in detkit. This method fits a scaling law to a small subset of log-determinants computed from submatrices and extrapolates to larger submatrix sizes using Algorithm 1. Specifically, we use the diag array from Listing I.I to compute the log- determinants l, = ΣKa log Dil. k = nop ................................ y, for submatrices of size kd x kd, where d is the output dimension of the model gc.g., number of classes in CIFAR- -1). These submatrices correspond to a permuted ordering of the original matrix during LDL decomposition, i.c., M == PMP.. While the sampling of submatrices could, in principle, be performed in any order, the LDL decomposition conveniently provides the deg-determsnants of successive principal submatrices formed by selecting the first kd rows and columns-at no additional cost; advantage exploit in this approach. The resulting sequence & is then fit over the interval & E [,,, and the fitted FLODANCE model is used to predict log-determinants in the extrapolation range n..
We have also concurrently developed separate high-pertn perforrrformacce Python package, imate,* which implements stochastic Lanczos quadrature (SLQ), randomized method for approximating the log-determinant at scale. This package is imple- mented with a C++/UUDA backend and supports execution on both CPU and multiple GPUs. Listing I.3 demonstrates a minimalistic usage of the imate. 1 logdet function.
'detkit is available for installation from
http:///pppp.ooogppprroeettttttt the documentation can he found
https: / /mme gi thub. io/det t, and the source code is available at https://gtthh com/ ame / t. "imate is available for installation from PyPI (https: htttps::/ypp... http://pppppppppp../// wttpg//rrojcc//mateee the documentation can be found https://amee i.gi thub. o immtt, and the source code is available at attttps://gtthii ccm/ amel 1 imate.
33

---

Determinant Estimation under Memory Constraints and Neural Scaling Laws
Listing L2. The class tlogeer. fits and extrapolates log-determinants using FLODANCE in Algorithm
import nampy as
from detkit import FitLogdet
Range of datapoints (n) and nomber of labels (d) for CIFAR-10 d = np. ra.re (500C0), 10
Compute & logabsdet(Mpuas)) for sub-matrices of the size = 1,..., Here, diag is an array of length m obtained from Listing 1.1 ell np, cumsum inp. log abs (diag)))
Keep every d-th element ell ell [1][(6e [d-1)
Choose a fit interval, such as (mg: n.) = (10",5 5 10")
15 no, ns = 1e2, 5e3
fit_mask = (n > n0) & (n < ns)
Fit using Algorithm with order truncated Laurent series flodet: FitLogdet (q-1)
20 fiodet.fit (n [fitts__aask]
Extrapolate in a larger interval, such as in (n,, = (5 10 xx 109) n_eval - geomspaae (ns, 2)
ell_eval flodet. [n_eval]
Listing I.3. A minimalistic usage of imate package. The function logdet computes logdet(M) using the stochastic Lanezos quadra- algorithm.
Install imate with "pip install imate"
from imate import logdet import numpy
Number of data (n) and labels (d) d 50000, 10
NTK matrix M on disk M numpy, memmap ( filename. anpen'meepy', mode="'', dtype= float32', shape- (nxd, n=d))
Compute (Mg using stochastic Lanczos quadrature (SLQ) method
Assume M symmetric.
1d, info = logdet (M, method=' slg', min_num_sample==100,, nax mak_nun_cappes==20,
lanozos_dest [inozos_degree=100,, error_rtol=0.01, 01, confiidence_level-0. 95,
wtt.gggiiiianceeeeeee....... 001, orthegonalize--1, sur_threads-O, num num_gpu_devices-0 gou-True, verbose-True, return_info-Trul
plot-True)
34


Reviews:
[{"review_id": "K8jSq5mc9j", "reviewer": "Reviewer_EHZs", "summary": "This work proposes a method to scale the calculation of matrix determinant to extremely large matrices, especially in context of ill-conditioned matrices such as the empirical NTK matrix. They first develop a memory-constrained algorithm which computes determinant exactly and can serve as a baseline for their later experiments. The core idea behind their main method is to derive an scaling law for the determinant and then use the obtained fit to extrapolate to determinants of even larger matrices. The approach is quite interesting and, besides helping with the particular problem, might also inspire techniques that utilize extrapolation for computation of spectral quantities pertaining to neural networks.", "questions": "", "limitations": "", "rating": 3, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We would like to thank the reviewer for taking the time to review our paper, and for their constructive feedback.\nRegarding general applicability of our methods: we stress that MEMDET computes the exact (64-bit) determinant, and is a general method that works for arbitrary matrices. \nOn the other hand, our neural scaling law approach, FLODANCE, is only suggested to work for Gram matrices satisfying Assumption 1. \nIn our writeup and experiments, we focus on NTK matrices, since they are known to be a particularly difficult class of kernel matrices. Other algorithms (e.g. SLQ) already perform very well on better behaved, more traditional classes of matrices. \nFor modern neural networks, the corresponding NTK kernels are non-stationary, and the Gram matrices tend to be highly ill-conditioned and very large, often not being feasible to compute explicitly: hence the need for our algorithm. \nWe believe that this makes NTK matrices particularly useful examples to highlight the utility and efficacy of our method. \nHowever, we do recognize that other Gram matrices are commonly used in ML applications. In particular, we computed the Gram matrix associated to the linear model of coregionalization [1], with 10 outputs, and Matern kernel over 10,000 data points. Plots demonstrating the performance of FLODANCE on this dataset, along with experimental details, can be found at http://anonymous.4open.science/r/memdet-E8C1/notebooks/matern.pdf .\n\nAs the reviewer has pointed out, the $\\log$-determinant often appears as one term in more elaborate quantities of interest. \nExamples of this are the quadratic-form term in the $\\log$-marginal likelihood of a Gaussian process (along with their gradients for the purposes of training), and the curvature term that appears in the IIC. \nEstimating these quantities is a clear goal. \nHowever, the other terms appearing in these quantities are highly challenging to estimate as well, requiring separate novel techniques that we believe will each be of independent interest.\nThis work is part of a broader research program under active development, and the current submission already contains two novel techniques for computing and estimating a crucial term, with corresponding software. \nWe believe this strikes an appropriate balance of utility, impact and digestibility for an outlet such as ICML. \n\nTo clarify the sensitivity of FLODANCE to subsample size and seed subset, we have conducted the following experiments.\nFirst, http://anonymous.4open.science/r/memdet-E8C1/notebooks/param_selection.pdf shows the effect of the number of training samples on prediction accuracy, and describes the experimental setup. \nAs is to be expected, the general trend is that using more samples results in better prediction. \nThere is a natural tradeoff between computation time and accuracy that is up to the practitioner to determine.\n\nSecond, we tested the sensitivity of the estimate to the training sample. To this end, we sampled an ensemble of 15 submatrices of the NTK for ResNet-50 trained on CIFAR10. Each submatrix contians 10,000 datapoints with 10 classes, forming a $100,000\\times100,000$ NTK matrix. Experimental details and plots can be found at http://anonymous.4open.science/r/memdet-E8C1/notebooks/scale_law_pred_uq.pdf\n\nThird, we upgraded our computations on both architectures to use the full CIFAR-10 dataset with $(n = 50,000)$, resulting in NTK matrices of size $(m = 500,000)$. For ResNet9, we also recomputed log-determinants using 64-bit precision via MEMDET. The updated results are provided in http://anonymous.4open.science/r/memdet-E8C1/notebooks/scale_law_fit.pdf (was Figure H.1) and http://anonymous.4open.science/r/memdet-E8C1/notebooks/scale_law_pred.pdf (was Figure 4). The 64-bit computation for ResNet50 is underway and will be completed in a few days, with updated results included in the final version.\n\nWe note that while the new plots are derived from the same NTK kernel and methodology, their appearance differs from earlier versions due to a technical adjustment. To improve quantization in 32-bit precision, we had previously scaled the NTK matrices, giving an additive term proportional to $n$ in the $\\log$-determinant. Hence, the curves differed from the current ones by a linear function of $n$. This shift does not affect the fitting behavior or predictive quality, but it alters the shape of the curve. In the updated 64-bit computations, this scaling was no longer necessary, resulting in cleaner and more interpretable visualizations without any transformation applied.\n\nAs is evident by these additional studies, FLODANCE is quite robust to these factors. \nWe thank the reviewer for suggesting these experiments, as we believe their inclusion in the final version of the document will improve the quality of our work. Please let us know if anything else needs verification or clarification, as we are eager to engage further if necessary to improve our work.\n\n[1] https://doi.org/10.1007/BF02066732", "claims_and_evidence": "- The experimentation is all specific to NTK. But the method is pitched as a method which is generally applicable. While that's true, I don't think there is strong empirical evidence to support this claim. Especially the scaling law approach may not work in other contexts. \n\nAlso, see methods and evaluation criteria section", "methods_and_evaluation_criteria": "- It would have been more interesting to see the use of determinant in some specific context, and then observe that a better estimation leads to better improvements. But, you know, often it's not the value of determinant in isolation that is needed, and rather some other quantities too. And then, in deep learning, things can always turn out a bit weird, so would have been nice to ground it somehow.\n\n- Another thing I am not fully sure is if you just want to calculate determinant, how many samples of the dataset must be used to do so. This is very unclear. It could be that all estimates using 1000 samples are terrible, but you really need like 50K samples to get to the most accurate estimates. Therefore the gold standard is a bit hard to gauge. Perhaps the authors can do some testing in an appropriately downscaled setting, where say the determinant is exactly computed, in 64 bit precision, on all datapoints in the training set. Right now it is unclear how much to read in their FLODANCE estimates closer to the memdet baseline. \n\n- Finally, it is also unclear how specific are the determinant estimates to the 'seed'-subset used for the scaling laws. See the questions section below for elaboration.", "theoretical_claims": "--", "experimental_designs_or_analyses": "In general, the experimental design seems fine but I have a few questions about them:\n\n- what is the floating point precision used for the baselines such as block diagonal, SLQ, Pseudo NTK? \n\n- do all the methods under comparison use the same sized subset? \n\n- what is p in the table 2? also, i suppose memdet baseline is the L_n one? if so, please label it more explicitly, you don't remember all possible symbols in your first few reads :)", "relation_to_broader_scientific_literature": "I think the use of scaling laws could be an interesting direction for numerical estimation, which I believe is relatively unexplored.", "essential_references_not_discussed": "--", "other_strengths_and_weaknesses": "Overall the paper is well written, and the idea of using scaling laws is clever. But how accurate and robust these determinant estimates are and how much the downstream applications gains from it are far from clear.", "other_comments_or_suggestions": "It's not clear to me how much the details of MEMDET are relevant to a ML audience. Maybe put some of that in appendix, and describe the scaling law and its theory a bit better.", "questions_for_authors": "- How specific are the determinant estimates to the subset used for deriving the scaling laws? Can you present some experimental evidence which would support that your estimate is robust to this change? On a similar note, how big must the subset be to allow for this robustness?\n\n - I believe the scaling laws has been shown in the context of LLMs trained in an online manner. Is there evidence that similar scaling laws also hold in the classical vision settings as used here? Otherwise that should be stated as an assumption as well.", "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": "i8Th00dpQ8", "reviewer": "Reviewer_zrLg", "summary": "This work addresses the memory and computation bottlenecks in estimating log-determinants of large matrices such as the Neural Tangent Kernel (NTK) for large models and datasets. The paper introduces MEMDET, a memory-constrained algorithm for exact log-determinant calculations for matrices too large to fit in memory by using block LU decompositions with an efficient block ordering to minimize data transfer. The paper then proposes FLODANCE, an accurate algorithm for extrapolating log-determinants computed on a small subset of the dataset to the much larger full dataset by exploiting the scaling behavior of the log-determinants for a class of kernels.", "questions": "", "limitations": "", "rating": 5, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their constructive feedback and positive review of our work. This work is indeed part of a broader active research program aimed at providing computational tools for estimating linear algebraic quantities at scale. We will make this more clear in the final version of our paper, and cite the related work  for hyperparameter selection that you suggested.\n\nRegarding the assumptions underlying the scaling law behavior, we acknowledge that these technical assumptions have been identified but not verified.\nTo do so for the empirical neural tangent kernel matrix may not be realistic, however, we have come across [1], which examines the spectral decomposition of the analytic NTK (an approximation of the empirical variant we consider) on the sphere. In [1, Theorem 4.8], it appears that our Assumption F.1 is verified over a bounded eigenbasis (Assumption F.2), and Assumption F.3 is automatically satisfied as we consider $f_p^\\ast = 0$. If desired, we are happy to include this discussion in the final manuscript. \nMore generally, however, the theory is instead used to establish plausibility of the scaling law, since there is a well-characterized class of kernels for which these assumptions hold. We then use the scaling laws as an empirical theory in much the same way that others do in the ML literature. \n\nIn terms of this empirical behavior, a plot of the residuals of the successive $\\log$-determinant increments after removal of the scaling law can be found at http://anonymous.4open.science/r/memdet-E8C1/notebooks/stationarity.pdf. \nThis figure demonstrates that these values are approximately normally distributed, empirically validating the scaling law in Assumption 1 (as well as Assumption 2).\nThis figure is included in the updated version of our work, and clarification has been added regarding the distinction between this empirical verification and the assumptions stated in Appendix F.\n\nRegarding downstream applications, we note that the $\\log$-determinant is typically only one of multiple terms appearing in quantities of interest. For this reason, we believe it is premature to include studies computing these partial quantities as approximations, until we have developed the toolkit to treat the other terms.\n\nThe non-asymptotic terms that appear in the exponent model sublinear behavior that occurs in the pre-limit. \nTo be clear, they model the function $\\nu = \\nu(n)$ by allowing for subsequent terms beyond the constant in its Laurent series. \nThese terms were included for flexibility as we found that they do improve accuracy for matrices at the scale of our experiments. \nHowever, they are not strictly necessary: when the exponent is treated as a constant in $n$, we found that our method still outperforms SLQ (see http://anonymous.4open.science/r/memdet-E8C1/notebooks/relative.pdf for a version of Table 2, updated to display relative error for reviewer ikyx). \nErrors and cost also now appear in an updated version of Table 1 found at http://anonymous.4open.science/r/memdet-E8C1/notebooks/comparison.pdf . Algorithm 1 has also been updated to explicitly state how the terms are computed.\n\nPlease let us know if anything else needs verification or clarification, as we are eager to engage further if necessary to improve our work.\n\n[1] Murray, M., Jin, H., Bowman, B., \\& Montufar, G. (2023, February). Characterizing the Spectrum of the NTK via a Power Series Expansion. In International Conference on Learning Representations.", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "Yes", "theoretical_claims": "I checked the proofs in Appendix F and G.", "experimental_designs_or_analyses": "All experiments in Section 4 appear sound.", "relation_to_broader_scientific_literature": "In contrast to conventional numerical methods for estimating log-determinants such as SLQ, I find FLODANCE's approach of exploiting the scaling behavior of the log-determinants of a specific class of kernels for efficient estimation quite creative and novel. Analyzing the scaling behaviors under various limits has shown to be a powerful tool in machine learning broadly, such as predicting optimal hyperparameters at scale [1]. This paper demonstrates that this idea can inspire similar algorithmic progress in numerical linear algebra. \n\n[1] Tensor programs v: Tuning large neural networks via zero-shot hyperparameter transfer. Yang et al. 2022", "essential_references_not_discussed": "None that I'm aware of.", "other_strengths_and_weaknesses": "Strength: Overall, the paper is very well presented. The conceptual insights and experimental results are both strong.\n\nWeakness: Three nontrivial assumptions are required to prove the scaling law (Theorem F.1), but they were not verified in the experiments. In this case, it's unclear how the theory presented in Section 3.1 can actually explain the observed scaling behavior of the log-determinants.", "other_comments_or_suggestions": "While FLODANCE performed strongly in the experiments relative to alternatives, it's unclear how useful the absolute level of accuracy it achieves is for practical applications (e.g. model selection). Can the authors compare FLODANCE with alternatives in some downstream applications that require estimating the log-determinants?", "questions_for_authors": "1. Can you motivate the form of the non-asymptotic correction to the exponent $\\nu$? How much does it affect the accuracy of the estimate?", "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": "GSo9pqSwW2", "reviewer": "Reviewer_ikyx", "summary": "The paper designs two algorithms for computing log-determinant of large PSD matrices under memory constraints. The first algorithm is named MEMDET, designed based on block LDL decomposition. The other one is named FLODANCE, designed based on neural scaling law assumptions. Empirical results show that the proposed methods achieve significant speedup while maintaining high accuracy.", "questions": "", "limitations": "", "rating": 4, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their constructive feedback and positive review of our work. The reviewer has raised a few points that deserve addressing, and we believe that doing so will improve our work.\n\nThe first of these is our dependence on the scaling law assumption. In the submitted version of the document, we provided criteria on the kernel which guarantees the scaling-laws hold in Appendix F. However, these are technical conditions, and their satisfaction is non-trivial to check (although we have discussed some possible verification in the response to Reviewer zrLg). This means that like most of the literature exploring scaling laws, we are relying on empirical behavior to check that the scaling laws are satisfied. To this end, http://anonymous.4open.science/r/memdet-E8C1/notebooks/stationarity.pdf contains a figure that compares the ratio of successive determinants to the fitted scaling law. Note that the log-scale on the $y$-axis turns this into a residual plot. These residuals are approximately normally distributed, which is evidence that Assumpsions 1 and 2 are satisfied. This plot has been included in the updated version of our document, along with discussion.\n\nRegarding the focus of our numerics on the NTK: this was simply done because these matrices are known to be a particularly difficult class of kernel matrices to work with. For modern neural networks, the corresponding kernels are non-stationary, and the matrices tend to\nbe highly ill-conditioned and very large, often not being feasible to compute explicitly: hence the\nneed for our algorithm. We believe that this makes NTK matrices particularly useful examples\nto highlight the utility and efficacy of our method. \nHowever, we do recognize that other Gram matrices are commonly used in machine learning applications. In particular, we computed the Gram matrix associated to the linear model of coregionalization [1], with 10 outputs and a Matern kernel over 10,000 data points. Plots demonstrating the performance of FLODANCE on this dataset, with experimental details, can be found at http://anonymous.4open.science/r/memdet-E8C1/notebooks/matern.pdf .\nWe have included this experiment in the appendix of an updated version of our work.\n\nOther reviewers have also drawn attention to the lack of the distinction between our methods and their competitors. To this end, we have included a brief section in the  the appendix outlining the stochastic Lanczos methods, as well as the Pseudo-NTK. We believe this will help our exposition, in addition to a brief higher-level summary in the main document. \n\nThank you for your suggestions regarding Table 2. An updated version of this table containing relative error can be found at http://anonymous.4open.science/r/memdet-E8C1/notebooks/relative.pdf . We remark that MEMDET is an exact method that allows for determinant computation of large matrices that do not fit into memory. It achieves this through block-decompositions, relying on an underlying algorithm (LDL$^\\top$, LU, or Cholesky) that is amenable to this type of decomposition. As such, the use of MEMDET in our FLODANCE approximations reduces to LDL$^\\top$, since the training sample is small enough to fit in memory. We have emphasised this in the final version of the document.\n\nPlease let us know if anything else needs verification or clarification, as we are eager to engage further if necessary to improve our work.\n\n[1] https://doi.org/10.1007/BF02066732", "claims_and_evidence": "The authors claim that the proposed methods achieves fast and accurate log-determinant estimation under memory constraints. The claims are well-supported through descriptions on algorithm design, complexity analysis and empirical studies.", "methods_and_evaluation_criteria": "The proposed methods are evaluated under the NTK estimation scenario for deep learning. the evaluation criteria focus on accuracy, running speed and memory usage, aligning well with the studied problem.", "theoretical_claims": "I checked the proofs of the main theoretical results and they are rigorous. However, the proof is largely relying on the scaling-law assumption, which may require further justification.", "experimental_designs_or_analyses": "The proposed methods are evaluated across various neural networks and sample sizes. They are practically sound but may benefit from more experiments beyond the NTK scenario.", "relation_to_broader_scientific_literature": "The proposed methods can facilitate relevant researches on the practical or theoretical behavior of neural tangent kernels. They may also inspire new matrix-based algorithms for, e.g. kernel approximation.", "essential_references_not_discussed": "The paper may benefit from additional discussions about the Lanczos family of matrix trace estimation techniques.", "other_strengths_and_weaknesses": "Strength:\nThe proposed algorithm designs are promising, and the empirical performance gain of the proposed methods is significant.\n\nWeakness:\nThe neural scaling law assumption may need further justification, e.g. empirical results on some ill-conditioned / noisy scenarios.", "other_comments_or_suggestions": "No comments.", "questions_for_authors": "* From Table 2, the estimation seems very unstable since the results vary a lot across different methods. How is the numerical stability of the proposed MEMDET method? If the ground truths in Table 2 are acquired by MEMDET, it seems not very fair since FLODANCE is using MEMDET as the backbone. How are these results compared to the native method? Also, I believe the absolute value of the log-determinant does not matter a lot, instead we are more interested in its relative value across different settings. From this point, showing the relative behavior of each method (e.g. verifying the neural scaling law) may give a more solid comparison.", "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": "pY7PzJyfZs", "reviewer": "Reviewer_QjXZ", "summary": "This paper proposes a scalable way to compute Neural Tangent Kernel log-determinants of dense matrices which may arise in training deep neural networks. The empirical Neural Tangent Kernel has been shown as the effective tool to study the behavior of neural networks during both training and inference. In particular, the Gram matrix (i.e. the kernel matrix) can be used in lazy-training and in obtaining uncertainty quantification estimates. Computing log-determinants of such NTK kernels for a large $n$ is extremely hard and sensitive to small eigenvalues. The challenges include: 1) dealing with storing the large kernel matrix; 2) a scalable way to compute functions of the kernel matrix (in this case, the log determinant). The authors use an existing block $LU/LDL^T$ decompositions to efficiently load blocks of pre-computed Gram matrix from disk. However, realizing that this still poses a problem for a large value of $n$, the authors propose regressing over precomputed log determinants up to $p \\times p$. For a large $n$, the authors propose an asymptotic approximation. The authors test the proposed approximation on the NTK matrices corresponding to widely used models such as ResNet9, ResNet18, and ResNet50.\n\n## update after rebuttal\nWe thank the authors for their rebuttals. For the most part I am satisfied with the replies. I will maintain weak accept for the submission. In order for the submission to get a higher score (for this conference and for other follow-up work based on this submission), the authors are encouraged to:\n- demonstrate some empirical comparison in computation time; the submission focuses on approximation quality mainly.\n- tabulate symbols used in the writeup, as at least one other reviewer also mentioned difficulties in keeping track of the symbols.", "questions": "", "limitations": "", "rating": 3, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their careful review of our paper and positive feedback on our work. We acknowledge that some points need clarification and aim to do so below.\n\nRegarding the implicit value of $m=10$: this is the number of classes used in all the classification datasets we used in experiments. We will state this explicitly in Section 4, as we understand the confusion.  We further thank the reviewer for their suggestion regarding plot and variable labeling, as well as algorithm naming. We have included these improvements and corrections in the updated version of the paper.\n\nRegarding Section D: this is included as algorithmic analysis of MEMDET. Many of these were conducted for our own benefit to choose appropriate block sizes, and to predict how long the computations would take. MEMDET is an exact algorithm (with no approximation occurring) and was required as a baseline for our experiments, since we could not fit the full NTKs into memory on the machines we were using. We believe that exact methods will always have a place, and so hope the information in Appendix D will be useful to others.\n\nTo clarify the difference between algorithms: Pseudo-NTK summarizes the data by summing over rows and columns $m\\times m$ block, then dividing by $m$. \nThis gives a good approximation with respect to the operator norm. \nHowever, it should not be expected to give refined determinant estimates. \nBoth SLQ and FLODANCE build their approximations from cheaper quantities. \nIn the case of SLQ, matrix vector products are used in a Monte Carlo approximation. \nThis gives a low bias estiamte with high variance, particularly when the condition number of the matrix is large. \nOn the other hand, in FLODANCE we use extrapolation from matrix minors. \nThis gives a low variance estimate in our experience (see http://anonymous.4open.science/r/memdet-E8C1/notebooks/scale_law_pred_uq.pdf), while the magnitude of the bias is unclear, and subject to the stationarity of the process in the training region (see http://anonymous.4open.science/r/memdet-E8C1/notebooks/stationarity.pdf). \nWe agree that including more details of these methods will improve our work, and have added a short description in the main body, and a more in-depth section in the appendix of the updated version of the document.\n\nPlease let us know if anything else needs verification or clarification, as we are eager to engage further if necessary to improve our work.", "claims_and_evidence": "The claims made in the paper seem to be sound and based on well-known linear algebraic techniques for computing decompositions which may be used to compute log determinants. A clever way to scale this up to $n$ by an asympotitic approximation is interesting.", "methods_and_evaluation_criteria": "The proposed benchmark datasets are well chosen and suitable for evalating the problem. The comparison against competing methods (such as SLQ and Pseudo NTK) are also included.", "theoretical_claims": "I am familiar with the correctness of MEMDET algorithm since these are based on already-existing block decompositions of large matrices. I did check the appendix for proof of Proposition 1 and it seems to be correct.", "experimental_designs_or_analyses": "Experimental designs are more or less sound, although for a non-expert, it is not too clear why the implied value of $m = 10$ is used for all of the experiments (I see this from the fact that $K_{1000}$ results in $10,000 \\times 10,000$ matrix; see Section 4.1).", "relation_to_broader_scientific_literature": "The authors provide broader connection to other machine learning applications such as model selection, quantification of generalization bounds, etc. The considered problem of log-determinant computation and the proposed methodology is quite relevant to many of the interesting problems in machine learning.", "essential_references_not_discussed": "I am not aware of other essential references which should have been included. However, it may be enlightening to include how the proposed method is different algorithmically from SLQ and Pseudo NTK approaches.", "other_strengths_and_weaknesses": "Strengths:\n1. Considers an important problem of interest in machine learning. Log determinant computation is inherently hard.\n2. Scaling law consideration for approximation is interesting.\n\nWeaknesses:\n1. Some of the plots in the paper are not easy to understand (see below).\n2. The difference between the proposed method and the competing methods such as SLQ should be highlighted. I could not find any information in the writeup itself.", "other_comments_or_suggestions": "1. Table 1: SLQ  = Stochastic Lanczos Quadrature, but this is nowhere mentioned in the writeup until this point.\n2. Section 3.2: shouldn't the precomputed sequence be: $\\{ L_1, ... L_p \\}$, not up to $L_n$?\n3. In Algorithm 1: I was reading the index for all variables used in the regression equation for $y_n$. It was strange to see the index shifted by $+1$, but the formula for $y_n$ in Section 3.2 was shifted by $-1$. This needs to be consistent.", "questions_for_authors": "1. Why is $m = 10$ assumed in for all datasets in the experiment?\n2. Are there any way you can use other variable names to help the readers, e.g. remembering the definition of $m$ across the entire writeup is challenging. The definition of $p$ is okay enough though.", "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:

Meta-review: This paper proposes a method to calculate the determinants of large matrices which arise as Neural Tangent Kernel for large neural networks. The reviewers agree that the method is sound and most praised the originality and creativity of the work. After the discussion period, all reviewers now give positive scores. Although reviewer QjXZ initially complained of the lack of run time comparison to some other methods, the authors have provided a very convincing rebuttal with additional experiments, which should be included in the camera ready version. Reviewer EHZs also complains of the application being limited to NTK type matrices by the application of the scaling laws, this is an issue that can be solved by better writing and contextualization in the revision. 
Since two reviewers complained of the difficulty following the notation, I recommend adding a table of notation. 

Overall, this paper makes a very solid contribution and should be **accepted**.
Decision: Accept (poster)