OPML feed of all feeds.
Subscribe to the Atom feed, RSS feed to stay up to date.
Thank you to arXiv for use of its open access interoperability.
Note: the date of arXiv entries announced right after publication holidays might incorrectly show up as the date of the publication holiday itself. This is due to our ad hoc method of inferring announcement dates, which are not returned by the arXiv API.
Powered by Pluto.
Source on GitHub.
Maintained by Nima Anari, Arnab Bhattacharyya, Gautam Kamath.
from ECCC Papers
from CCI: jobs
The Algorithms and Complexity group of IRIF is seeking excellent candidates for postdoctoral positions in classical and quantum computing, with starting date in October 2027 (negotiable).
Applications with a CV including list of publications, a summary of research, and names and email addresses of at least 3 references should to be sent to algocomp-apply@irif.fr and received by Oct. 24, 2026.
Website: https://www.irif.fr/en/equipes/algocomp/index
Email: adiro@irif.fr
from ECCC Papers
from Francis Bach
A common task in machine learning is to estimate or optimize “log-sum-exp” functions with (potentially continuously) many terms such as $$ \log \Big( \int_{\mathcal{X}} e^{v(x)} dq(x) \Big),$$ where \(v: \mathcal{X} \to \mathbb{R}\) is some potential function, and \(q\) is a probability distribution on the set \(\mathcal{X}\). This has many applications throughout data science, often through the normalization of probabilistic models, but also as a smooth approximation to the maximum, in transformers through its derivatives, or in reinforcement learning when using entropy regularization [19]. Sometimes the set \(\mathcal{X}\) is finite (potentially big) and the integral can be done by explicit summing, but often an exact computation is infeasible, and sampling from the probability distribution \(q\) is used instead.
The key difficulty comes from the variance of such estimates, in particular when \(v\) takes large values. In the simplest example, for \(z_1,\dots,z_n \in \mathbb{R}\) independent and normally distributed with mean \(\mu\) and variance \(\sigma^2\), the relative squared error for estimating \(\mathbb{E}[e^z]\) is $$\frac{ {\rm var}\big( \frac{1}{n} \sum_{i=1}^n e^{z_i} \big) }{( \mathbb{E}[ e^{z} ])^2} = \frac{1}{n} \frac{ {\rm var}(e^z) }{( \mathbb{E}[ e^{z} ])^2} = \frac{ e^{\sigma^2}-1}{n}.$$ It converges to zero when \(n\) grows (as can be expected from the law of large numbers), but explodes exponentially when \(\sigma\) grows. Even taking the logarithm does not change the exploding variance, that is, \({\rm var}\big( \log \big( \frac{1}{n} \sum_{i=1}^n e^{z_i} \big)\big)\) also can be shown to grow asymptotically similarly in \(\frac{ e^{\sigma^2}-1}{n}\) (when \(n\) is large, as can be obtained from the delta method).
While difficult to estimate, the log-sum-exp function comes with many nice properties (and that’s why people love it); I particularly like the fact that (1) it is a smooth approximation to the maximum (see, e.g., this earlier post), and (2) it is a way to normalize probabilistic models that is adapted to maximum likelihood estimation, in particular in hierarchical probabilistic models, where (conditional) independence assumptions lead to separability of associated loss functions (as thoroughly used in probabilistic graphical models).
The main question I try to answer in this post is:
Can we keep the advantages of optimizing log-sum-exp functions while being less exposed to their computational / statistical disadvantages?
The magic of least-squaresAt the other end of the spectrum sits least-squares regression, with essentially the exact opposite features:
Can we reconcile them? In other words, is least-squares really all I need? (my colleagues sometimes mock me for my love of least-squares).
Note that there is another (classic) attempt at seeing the world through least-squares: doing it in series through Newton’s method, leading in this context to iteratively reweighted least-squares, but this is for computations only, with no statistical improvement. What we are aiming at is stronger: can we get least-squares-based closed-form estimators for maximum-likelihood problems that typically require optimization of a convex function (such as logistic or softmax regression)?
Interestingly, my new attempt can be summarized in one integral equation $$ t \log t\, – t + 1 = \int_0^1 \!\! \frac{ (t-1)^2}{\rho t + 1-\rho} (1-\rho) d\rho,$$ which can be checked by usual integration tricks. Let’s see why and how!
Relative density estimation as a testbedIn this post, I look at a simple fundamental problem where we can study and compare various estimation frameworks, noting that it can be extended in several ways (in particular, through mutual information, see below).
We consider two probability distributions \(p\) and \(q\) on \(\mathcal{X}\); our goal is to estimate the logarithm of the relative density \(\log \big(\frac{dp}{dq}(x)\big)\). This turns out to be equivalent to estimating the Kullback-Leibler (KL) divergence because of the variational formulation [1] $${\rm KL}(p\|q) = \int_{\mathcal{X}} \log \big(\frac{dp}{dq}(x)\big) dp(x) = \sup_{v: \mathcal{X} \to \mathbb{R}} \int_{\mathcal{X}} v(x) dp(x) + 1\, – \int_{\mathcal{X}} e^{v(x)} dq(x). \tag{1} $$
This is one particularly important instance of an \(f\)-divergence (see, e.g., [2]), with the following definition and variational formulation based on the Fenchel conjugate \(f^\ast\) of \(f\): $$D(p\|q) = \int_{\mathcal{X}} f \big( \frac{dp}{dq}(x) \big) dq(x)= \sup_{v: \mathcal{X} \to \mathbb{R}} \int_{\mathcal{X}} v(x) dp(x) \, – \int_{\mathcal{X}} f^\ast(v(x)) dq(x),$$ the representation being a consequence of \(f(t) = \sup_{ u \in \mathbb{R}} ut-f^\ast(u)\) applied to each \(t = \frac{dp}{dq}(x)\). The KL divergence corresponds to \(f(t) = t \log t \, – t + 1\) and \(f^\ast(u) = e^u \, – 1\).
Note that for the particular case of the KL divergence, when optimizing with respect to a constant on top of \(v\), we obtain the Donsker-Varadhan representation [3] $${\rm KL}(p\|q) = \sup_{v: \mathcal{X} \to \mathbb{R}} \int_{\mathcal{X}} v(x) dp(x)\, – \log \Big( \int_{\mathcal{X}} e^{v(x)} dq(x) \Big). \tag{2}$$
We see the log-sum-exp function appearing explicitly. To estimate the potential \(v\) from i.i.d. samples from \(p\) and \(q\), the traditional variational approach corresponds to replacing integrals with empirical averages. For \(q\), this leads to a potentially unstable empirical average when only samples are available. The goal of this post is to explore another way (I present the main principles behind this new framework; see [7] for more details).
The framing through \(f\)-divergences is really key to the new approach, as other divergences will be instrumental in the definition. Before we move on, we state another (equivalent) variational formulation with two potentials \(v\) and \(w\), which we will need later: \(D(p\|q)\) is equal to $$\sup_{v,w: \mathcal{X} \to \mathbb{R}} \int_{\mathcal{X}} v(x) dp(x) + \int_{\mathcal{X}} w(x) dq(x) \mbox{ such that } \forall x \in \mathcal{X}, w(x) \leqslant -f^\ast(v(x)). \quad \tag{3}$$ At optimum, we get \(w(x) = -f^\ast(v(x))\), and we recover Eq. (1) in the KL case. The constraint is convex, but for \(f(t) = t \log t – t + 1\), it is far from what traditional convex optimization methods typically allow. This formulation appears in [25, Theorem 4.4] and has the nice property of preserving the symmetry of the problem (that is, if \(p\) and \(q\) are swapped, this is equivalent to replacing \(f\) by \(t \mapsto t f(1/t)\), and this corresponds to swapping \(v\) and \(w\).) In what follows, we will obtain candidates for functions \(v\) and \(w\) that satisfy the constraint \(\forall x \in \mathcal{X}, \ w(x) \leqslant -f^\ast(v(x))\), typically without equality.
Weighted chi-square divergencesAnother relevant function \(f\) for \(f\)-divergences is, for a parameter \(\rho \in [0,1]\), $$ f(t) = \frac{1}{2} \frac{ (t-1)^2}{ \rho t + 1-\rho}. $$ It leads to a weighted chi-square divergence $$ D(p\|q) = \frac{1}{2} \int_{\mathcal{X}} \frac{ \big(\frac{dp}{dq}(x)-1 \big)^2}{ \rho \frac{dp}{dq}(x) + 1-\rho} dq(x).$$ It has been used in various areas of applied mathematics [4, 5], and comes under several names for special cases, such as Pearson chi-square divergence for \(\rho =0\), or Neyman chi-square (or reverse Pearson) for \(\rho=1\), or Le Cam divergence for \(\rho=1/2\).
The function \(f\) above, which is of the form “quadratic over affine” has the variational representation through “quadratic plus affine” functions: $$ \frac{1}{2} \frac{ (t-1)^2}{ \rho t + 1-\rho} = \sup_{u \in \mathbb{R}} \ (t-1) u \, – \frac{1}{2} ( \rho t + 1 – \rho) u^2, $$ with the optimal \(u = \frac{t-1}{\rho t + 1 \, – \rho}\) (this, by the way, is not the Fenchel representation).
Thus, applying this for each \(x \in \mathcal{X}\) to \(t = \frac{dp}{dq}(x)\), for this function \(f\), we have $$ D(p\|q) = \!\! \sup_{u(\rho,\cdot):\mathcal{X} \to \mathbb{R}} \int_{\mathcal{X}} \Big\{ \big(\frac{dp}{dq}(x) -1\big) u(\rho,x) \, – \frac{1}{2} \big( \rho \frac{dp}{dq}(x) + 1 – \rho\big) u(\rho,x)^2 \Big\} dq(x). $$ This is exactly a quadratic variational problem since the function \(u(\rho,\cdot): \mathcal{X} \to \mathbb{R}\) only appears quadratically. The optimal variational function is then \(\displaystyle u(\rho,x) = \frac{ \frac{dp}{dq}(x) \, – 1}{ \rho \frac{dp}{dq}(x) + 1-\rho}.\)
Note that the quadratic cost function that we just defined is explicitly a least-squares prediction problem, to predict \(y\) given \(x\), where with probability \(\rho\), \(y\) takes the value \(1/\rho\) and \(x\) is sampled from \(p\), and with probability \(1-\rho\), \(y\) takes the value \(-1/(1-\rho)\) and \(x\) is sampled from \(q\) (this is thus reminiscent of noise contrastive estimation [26]). We thus see a potential instability at \(\rho=0\) or \(\rho=1\) (we will see later that this is not the case).
This variational formulation of \(D(p\|q)\) can be rewritten by assembling expectations with respect to \(p\) or \(q\) together, as $$ \sup_{u(\rho,\cdot):\mathcal{X} \to \mathbb{R}} \int_{\mathcal{X}} \big( u(\rho,x) \, – \frac{\rho}{2} u(\rho,x)^2 \big) dp(x) + \int_{\mathcal{X}} \big( -u(\rho,x)\, – \frac{1-\rho}{2} u(\rho,x)^2 \big) dq(x) .$$ We then get exactly a two-potential formulation with $$v(x) = u(\rho,x) \, – \frac{\rho}{2} u(\rho,x)^2 \ \ \mbox{ and } \ \ w(x) = -\, u(\rho,x) \, – \frac{1-\rho}{2} u(\rho,x)^2 $$ as the two potentials (the rigorous reader can check that \(\forall x \in \mathcal{X}, w(x) \leqslant -f^\ast(v(x))\), which is not obvious, so that they are the optimizers for Eq. (3)).
To summarize, for the function \(f: t \mapsto \frac{1}{2} \frac{ (t-1)^2}{ \rho t + 1-\rho}\), we have exactly what we want, that is, a variational formulation through least-squares. How can it be extended to more general \(f\)-divergences?
Extension by integrationIf we can write \(\displaystyle f(t) = \frac{1}{2} \int_0^1 \frac{ (t-1)^2}{\rho t + 1-\rho} d\nu(\rho)\) for some non-negative measure \(\nu\) on the interval \([0,1]\), then we can directly use the developments above to get a representation of \(D(p\|q) \) as $$\!\!\sup_{u:[0,1] \times \mathcal{X} \to \mathbb{R}} \int_{0}^1\!\!\!\! \int_{\mathcal{X}} \!\! \Big[\! \, (\frac{dp}{dq}(x) -1) u(\rho)(x) \, – \frac{1}{2} ( \rho \frac{dp}{dq}(x) + 1 \, – \rho) u(\rho)(x)^2 dq(x) \!\Big] d\nu(\rho), \tag{4}$$ with now a function \(u\) from \([0,1] \times \mathcal{X}\) to \(\mathbb{R}\).
We then get exactly a two-potential formulation, as presented in Eq. (3), with $$v(x) = \int_0^1 \big[ u(\rho,x) \, – \frac{\rho}{2} u(\rho,x)^2 \big] d\nu(\rho) \tag{5} $$ and $$w(x) = \int_0^1 \big[ -u(\rho,x) \, – \frac{1-\rho}{2} u(\rho,x)^2 \big] d\nu(\rho), \tag{6}$$ (which satisfy the constraint \(w(x) + f^\ast(v(x)) \leqslant 0\), which is, again, not straightforward) where \(u(\rho,\cdot)\) is a maximizer of Eq. (4). This corresponds to performing a continuum of least-squares problems in parallel.
These developments are valid for all \(f\)-divergences with an integral representation, and in particular the KL divergence, since we have $$t \log t\, – t + 1 = \int_0^1 \!\! \frac{ (t-1)^2}{\rho t + 1-\rho} (1-\rho) d\rho,$$ that is, \(d\nu(\rho) = 2 (1-\rho) d\rho\). This integral representation does not come out of nowhere; in fact, it comes from the theory of operator convex and operator monotone functions that we explored in an earlier post. It includes KL, obviously all weighted chi-square divergences (with \(\nu\) a Dirac measure), and all the \(\alpha\)-divergences [14], but unfortunately not the total variation.
Now that we have a generic framework to learn potentials \(v\) and \(w\) as integrals of least-squares estimates \(u(\rho,\cdot)\) for each \(\rho \in [0,1]\) and their squares, we can start to use function spaces to parameterize them, starting from linear models (see below for more general models).
Linear models with closed-form spectral estimationIf each function \(u(\rho,\cdot)\) is modeled as linear in some feature vector \(\varphi: \mathcal{X} \to \mathbb{R}^m\), that is, \(u(\rho,x) = \theta(\rho)^\top \varphi(x)\) for some \(\theta(\rho) \in \mathbb{R}^m\) (a family of parameters indexed by \(\rho\)), the optimization problem in Eq. (4) leads to $$ \sup_{\theta(\rho) \in \mathbb{R}^m} \ (\mu_p – \mu_q)^\top \theta(\rho)\, -\, \frac{1}{2} \theta(\rho)^\top ( \rho \Sigma_p + (1-\rho) \Sigma_q) \theta(\rho), $$ with the moments of \(\varphi\) with respect to \(p\) and \(q\): \(\mu_p = \mathbb{E}_p[\varphi(x)]\), \(\mu_q = \mathbb{E}_q[\varphi(x)]\), \(\Sigma_p = \mathbb{E}_p[\varphi(x)\varphi(x)^\top]\), and \(\Sigma_q = \mathbb{E}_q[\varphi(x)\varphi(x)^\top]\).
The optimal \(\theta(\rho)\) is then obtained by solving a linear system: $$\theta(\rho) = ( \rho \Sigma_p + (1-\rho) \Sigma_q)^{-1} ( \mu_p \,- \mu_q), $$ with an optimal value $$ \frac{1}{2} ( \mu_p \, – \mu_q)( \rho \Sigma_p + (1-\rho) \Sigma_q)^{-1} ( \mu_p\, – \mu_q).$$
Doing this for all \(\rho \in [0,1]\) and integrating, this leads to a new divergence that depends on the distributions \(p, q\) and on the feature map \(\varphi\): $$F(p\|q,\varphi) = \frac{1}{2} \int_0^1 ( \mu_p \, – \mu_q)^\top( \rho \Sigma_p + (1-\rho) \Sigma_q)^{-1} ( \mu_p\, – \mu_q) d\nu(\rho), \tag{7}$$ and a new candidate for \(v(x)\) to estimate \(f'(dp/dq(x))\) using Eq. (5). By construction, we obtain a lower bound on the \(f\)-divergence \(D(p\|q).\)
This is closed-form but requires integration with respect to \(\rho\), which is not practical, in particular since one should expect to need many quadrature points if using quadrature to estimate integrals, since the least-squares problems may diverge when \(\rho\) tends to 0 or to 1.
Eigenvalues to the rescue! Given that we need to invert matrices \(\rho \Sigma_p + (1-\rho) \Sigma_q\) for all \(\rho \in [0,1]\), using an eigenvalue decomposition, as done for ridge regression when solving for multiple values of the regularization parameter [6], seems natural. In our situation, we need a generalized eigenvalue decomposition of the pair \((\Sigma_p,\Sigma_q)\), that is, a basis \((v_1,\dots,v_m)\) of \(\mathbb{R}^m\) such that $$ \forall i,j \in \{1,\dots,m\}, \ v_i^\top \Sigma_q v_j = 1_{i=j} \ \mbox{ and } \ \Sigma_p v_i = \lambda_i \Sigma_q v_j.$$ A few lines of algebra (see [7]) then lead to $$F(p\|q,\varphi) = \sum_{i=1}^m \frac{1}{2} \int_0^1 \frac{ \big( ( \mu_p \, – \mu_q)^\top v_i \big)^2}{ \rho \lambda_i+ 1-\rho } d\nu(\rho) = \sum_{i=1}^m \big( ( \mu_p \, – \mu_q)^\top v_i \big)^2 \frac{f(\lambda_i)}{(\lambda_i-1)^2}. $$ There is a similar “simple” formula that is summing over eigenvalues for \(\theta(\rho)\), and the potentials \(v(x)\) and \(w(x)\) as quadratic-linear forms $$v(x) = \varphi(x)^\top M \varphi(x) + 2c^\top \varphi(x)\ \mbox{ and } \ w(x) = \varphi(x)^\top N \varphi(x) \, – 2c^\top \varphi(x),$$ with detailed formulas for \(M,N,c\). See [7] and, for convex analysis aficionados, the section below the references for more details, in particular a nice link with sum-of-squares optimization.
Although individual least-squares problems may have instabilities around \(\rho=0\) or \(\rho = 1\), after integration, estimation remains stable for all functions \(f\) such that \(t \mapsto f(t)/(t-1)^2\) remains bounded (for the KL, it is decreasing with \(f(0)=1\)).
Note that throughout this blog post, the concept of “simple closed-form formula” is quite subjective: by it, I mean stable routines from numerical linear algebra with explicit guarantees: this includes inverting linear systems and (generalized) eigenvalue decomposition [8].
Computational complexity. Using classical numerical linear algebra routines, the running-time complexity is \(O(m^2n +m^3)\) to compute the divergence and find estimates of \(v\) and \(w\), which is problematic when \(n\) or \(m\) is large. Feature learning as explained briefly below and more thoroughly in [7] allows one to learn an \(r\)-dimensional linear representation that is shared across all \(\rho\)’s, with iterative algorithms that have iterations of complexity \(O(r^3 + rmn)\), which is efficient when \(r\) remains small.
Benefits of spectral estimationNow that we can solve all these \(\rho\)-dependent least-squares problems in one shot, does it hold its promise in reducing variance? The main competitor here is the direct variational approach that maximizes Eq. (1) or Eq. (2).
Given data, for the spectral method, we simply (and classically) replace expectations with empirical averages (which corresponds to using empirical moments) and potentially add regularization, i.e., replace \(( \rho \Sigma_p + (1-\rho) \Sigma_q)^{-1} \) with \(( \rho \Sigma_p + (1-\rho) \Sigma_q + \lambda I)^{-1} \). This corresponds to performing ridge regression for all \(\rho\)-dependent least-squares problems.
The benefits can be measured either with theoretical arguments or by simulations. We provide both below.
High-dimensional evaluation on a Gaussian model. The simplest possible set-up is the Gaussian case with common covariance matrices (which I thoroughly explore in [9]). Since our unregularized estimator is invariant under affine transformations, we can consider \(p\) Gaussian with mean \(\Delta \in \mathbb{R}^m\) and covariance identity and \(q\) Gaussian with mean \(0\) and covariance identity. In the high-dimensional limit where the dimension \(m\) and the number of samples \(n_p\) and \(n_q\) grow to infinity with fixed ratios, the performance of the variational and spectral estimators only depends on \(s = \| \Delta\|^2\) and the “aspect ratios” \(\alpha_p = \frac{m}{n_p}\) and \(\alpha_q = \frac{m}{n_q}\). We consider linear features.
This setup is favorable to the variational estimator because the true log-density ratio is affine in \(x\), while the new spectral estimator incurs a bias (which can be explicitly characterized, see [9]). Is the increased bias compensated by the reduced variance? This can be precisely analyzed in the high-dimensional regime where \(n_p, n_q, m\) tend to infinity with fixed ratios \(\alpha_p = \frac{m}{n_p}\) and \(\alpha_q = \frac{m}{n_q}\), using random matrix theory [11] or the convex Gaussian min max theorem (CGMT) [12]. This allows us to compute asymptotic performance for the (unregularized) variational and the spectral approach, with the following performance for fixed \(s=1\) below, for all values of \(\alpha_p\) and \(\alpha_q\) (see [9] for all details).
Differences in performance between the spectral and variational estimators. Negative: spectral wins, positive: variational wins.
As expected, for large numbers of observations (small \(\alpha_p,\alpha_q\)), the variational method leads to better performance due to a reduced bias, but for smaller numbers of observations, its increasing variance makes the spectral method preferable.
Simulations. We consider a simple situation with data in two dimensions (\(d=2\)) and a non-linear log-density which is learned using random features based on ReLUs, that is, \(\varphi(x)_i = (w_i^\top x + b_i)_+\) for \(i \in \{1,\dots,m\}\) for randomly chosen \((w_i,b_i) \in \mathbb{R}^{d+1}\), with an increasing number of observations \(n = n_p = n_q\). The problem of estimating KL divergence between generic distributions is a non-parametric problem with convergence rates that exhibit the curse of dimensionality, and unless \(n\) is very large, or special sparsity assumptions are made, we can only estimate accurately in small dimension (see below for higher dimension when feature learning is used).
We see that KL is better than using the Pearson divergence (which is the traditional way [22, 23, 24] to use least-squares for density estimation, but suffers from the improper geometry in particular in the way it deals with positivity of densities), and better than variational, except for a large number of observations, where variational and KL spectral are the same.
Comparison of estimators of \(v\) using the criterion \(D(p\|q)\), with data in \([0, 1]^2\) with \(q\) uniform and \(p\) with independent components such that \(\log(dp/dq)\) is a sum of a few cosines. We consider ReLU random features, with \(m = 512\).
Within information theory, the KL divergence is often used as a measure of independence between random variables, leading to the mutual information: given a product space \(\mathcal{X}_1 \times \mathcal{X}_2\) and a joint distribution \(p(x_1,x_2)\), we can consider \(q(x_1,x_2)\) as the distribution with independent components that have the same marginals as \(p\), which we write \(q(x_1,x_2) = p(x_1) p(x_2)\), following the usual graphical model convention. Then the optimal log-relative-density is $$ \log \frac{ p(x_1,x_2)}{p(x_2) p(x_1)} = \log \frac{ p(x_2|x_1)}{p(x_2)}.$$
Hence, our framework for closed-form estimation allows us to perform conditional density estimation \(\log p(x_2|x_1)\) (with the additional need for \(\log p(x_2)\)). This can be done in general for any \(\mathcal{X}_2\), but when \(\mathcal{X}_2\) is finite, our closed-form estimation is exactly a way to perform softmax regression.
Note, however, that this new closed-form estimate is not equivalent to least-squares estimation for classification directly on one-hot encodings, which is essentially what the potential obtained from the Pearson divergence would do (and then leads to classical canonical correlation analysis). In the figure below, we see that for simple Gaussian data in two dimensions, the new spectral estimator for the KL divergence behaves “similarly” (but not identically) to softmax regression.
Comparison of estimators of conditional densities, by plotting the surfaces where one class dominates, learned from classification data (based on data with the same colors). Left: softmax regression, right: new closed-form spectral estimator.
It is also interesting to consider adding quadratic features (because the log-density here is quadratic since the class-conditional covariance matrices are not equal), and also compare to the classical square loss (which corresponds to using Pearson divergence, that is, \(\nu\) is a Dirac at \(\rho=0\)). With more features, all methods tend to have more similar classification regions (if the set of features is big enough to model all real-valued functions, they are identical). Note that on the top right, we see the masking problem of least-squares where some classes totally disappear (this is solved by adding features in the bottom right plot, but indicates an unnatural cost function).
Comparison of estimators of conditional densities, by plotting the surfaces where one class dominates, learned from classification data (based on data with the same colors). Top: using linear features, bottom: using quadratic features. From left to right: softmax regression, spectral estimation for KL divergence, spectral estimation for Pearson divergence.
Computational complexity. To obtain estimates with machine precision for \(k\) classes with feature vectors in dimension \(d\), then softmax regression using Newton’s method would take \(O( nd^2 k^2 + d^3 k^3)\) per Newton iteration, while the closed-form estimator takes only \(O(d^2 n + kd^3)\) for one eigenvalue decomposition, which is a significant gain. When gradient-based algorithms are used with feature learning, both computation times can be reduced.
Feature learningLinear models are great, but if one lesson has been learned since deep learning took over, it is that we need to learn features in a more end-to-end way, and large sets of predefined features are not enough for various reasons, in particular adaptivity to unknown (linear or non-linear) latent variables. In our variational framework where we have a lower bound on the KL divergence, this is simply maximizing \(F(p\|q, \varphi)\) in Eq. (7) with respect to \(\varphi\), and leveraging the fact that \(F(p\|q, \varphi)\) is convex in moments of \(\varphi\).
Indeed, as a convex function of the moments \(\mu_p-\mu_q,\Sigma_p,\Sigma_q\), it is lower-bounded by a constant plus $$ {\rm tr} \big( M \mathbb{E}_p [ \varphi \varphi^\top ] \big) + {\rm tr} \big( N \mathbb{E}_q [ \varphi \varphi^\top ] \big) + 2c^\top \big( \mathbb{E}_p [\varphi] – \mathbb{E}_q[\varphi]),$$ for matrices \(M,N\) and a vector \(c\) which can be computed from any given \(\bar\varphi\) (this is exactly what was needed to compute the potentials \(v\) and \(w\), and chosen so that the lower bound is tight for \(\varphi = \bar\varphi\) (see the section below the references). This allows for a minorization-maximization algorithm [15, 16] similar to the expectation-maximization (EM) algorithm [17]. The link with EM allows us to reuse many of the computational tricks developed there, such as online EM [18], to make the algorithm scalable to large networks to parameterize \(\varphi\) and large numbers of observations. More on this in a next post.
ConclusionIn this blog post, I introduced a new framework for relative density estimation that circumvents the exploding variance of means of exponentials. This was obtained by a continuum of stable least-squares problems, and made computationally feasible through a single generalized eigenvalue decomposition. Beyond tackling the exploding variance problem for models that need to be normalized, the unintended consequence for normalized models (where the sum/integral can be computed) was to obtain a closed-form estimator for softmax regression.
The recent paper [7] explores other consequences, in particular in terms of rates of estimation for the KL divergence, sometimes with minimax rates and partial adaptivity to linear latent variables. Overall, there is a long way to go, but I see this new framework as a potential replacement for the last layer of neural networks, where cross-entropy loss and log-sum-exp dominate. More on this in the next post.
Acknowledgements and tool usage disclosure. I would like to thank Frederik Kunstner and Nicolas Flammarion for helpful clarifying suggestions. Frontier LLM models were used to produce figures and correct typos.
References[1] XuanLong Nguyen, Martin J. Wainwright, and Michael I. Jordan. Estimating divergence functionals and the likelihood ratio by convex risk minimization. IEEE Transactions on Information Theory, 56(11):5847–5861, 2010.
[2] Yury Polyanskiy and Yihong Wu. Information Theory: From Coding to Learning. Cambridge University Press, 2025.
[3] Monroe D. Donsker and S. R. Srinivasa Varadhan. Asymptotic evaluation of certain Markov process expectations for large time-III. Communications on Pure and Applied Mathematics, 29(4):389–461, 1976.
[4] László Györfi and Igor Vajda. A class of modified Pearson and Neyman statistics. Statistics & Risk Modeling 19(3): 239-252, 2001.
[5] Lucien Le Cam. Asymptotic Methods in Statistical Decision Theory. Springer Science & Business Media, 2012.
[6] Gene H. Golub, Michael Heath, and Grace Wahba. Generalized cross-validation as a method for choosing a good ridge parameter. Technometrics 21(2):215-223, 1979.
[7] Francis Bach. A Spectral Framework for Closed-Form Relative Density Estimation. Technical report, arXiv:2605.10668, 2026, to appear in Advances of Neural Processing Systems (NeurIPS).
[8] Gene H. Golub and Charles F. Van Loan. Matrix Computations. Johns Hopkins University Press, 1996.
[9] Francis Bach. Regularized Variational and Spectral Log-Density-Ratio Estimation in the Gaussian Location Model. Technical report, arXiv:2607.01895, 2026.
[10] Didier Henrion, Milan Korda, and Jean Bernard Lasserre. The Moment-SOS Hierarchy: Lectures in Probability, Statistics, Computational Geometry, Control and Nonlinear PDEs. World Scientific, 2020.
[11] Zhidong Bai and Jack W. Silverstein. Spectral Analysis of Large Dimensional Random Matrices. Springer, 2nd edition, 2010.
[12] Christos Thrampoulidis, Samet Oymak, and Babak Hassibi. The Gaussian min-max theorem in the presence of convexity. Technical report, arXiv:1408.4837, 2014.
[13] Trevor Hastie, Robert Tibshirani, and Jerome Friedman. The Elements of Statistical Learning. Springer, 2009.
[14] Andrzej Cichocki and Shun-ichi Amari. Families of alpha-beta- and gamma-divergences: Flexible and robust measures of similarities. Entropy, 12.6:1532-1568, 2010.
[15] David R. Hunter and Kenneth Lange. A tutorial on MM algorithms. The American Statistician, 58(1):30–37, 2004.
[16] Julien Mairal. Stochastic majorization-minimization algorithms for large-scale optimization. In Advances in Neural Information Processing Systems, 2013.
[17] Arthur P. Dempster, Nan M. Laird, and Donald B. Rubin. Maximum likelihood from incomplete data via the EM algorithm. Journal of the Royal Statistical Society: series B (methodological) 39(1):1-22, 1977.
[18] Olivier Cappé and Eric Moulines. On-line expectation–maximization algorithm for latent data models. Journal of the Royal Statistical Society Series B: Statistical Methodology, 71(3):593–613, 2009.
[19] Brian D. Ziebart, Andrew Maas, J. Andrew Bagnell, and Anind K. Dey. Maximum entropy inverse reinforcement learning. AAAI Conference on Artificial Intelligence, 2008.
[20] Francis Bach. Sum-of-squares relaxations for information theory and variational inference. Foundations of Computational Mathematics, 25.3: 865–903, 2025.
[21] Keiji Matsumoto. A new quantum version of f-divergence. In Nagoya Winter Workshop: Reality and Measurement in Algebraic Quantum Theory, pages 229–273. Springer, 2015.
[22] Masashi Sugiyama, Taiji Suzuki, and Takafumi Kanamori. Density Ratio Estimation in Machine Learning. Cambridge University Press, 2012.
[23] Zaid Harchaoui, Francis Bach, and Eric Moulines. Testing for homogeneity with kernel Fisher discriminant analysis. Technical Report 0804.1026, arXiv, 2008.
[24] Mónica Ribero, Antonin Schrab, and Arthur Gretton. Regularized \(f\)-divergence kernel tests. Technical Report 2601.19755, arXiv, 2026.
[25] Michel Broniatowski and Amor Keziou, Minimization \(\varphi\)-divergences on sets of signed measures, Studia Scientiarum Mathematicarum Hungarica, 43(4):403–442, 2006.
[26] Michael U. Gutmann and Aapo Hyvärinen. Noise-Contrastive Estimation of Unnormalized Statistical Models, with Applications to Natural Image Statistics. Journal of Machine Learning Research, 13(11):307−361, 2012.
We defined our new divergence in Eq. (7) as $$F(p\|q,\varphi) = \frac{1}{2} \int_0^1 ( \mu_p \, – \mu_q)^\top( \rho \Sigma_p + (1-\rho) \Sigma_q)^{-1} ( \mu_p\, – \mu_q) d\nu(\rho). $$ It is convex in \(\Sigma_p\), \(\Sigma_q\), and \(\mu_p – \mu_q\), and there is a nice duality theory here, showing that \(M,N,2c\) defined above are in fact derivatives with respect to the parameters above. The main result here, shown in [7], is: $$F(p\|q,\varphi) = \ \sup_{M,N,c} {\rm tr}(\Sigma_p M) + {\rm tr}(\Sigma_q N) + 2c^\top (\mu_p – \mu_q) \qquad \qquad \qquad \qquad$$ $$\qquad \qquad \qquad \qquad \mbox{ such that } \ \forall \lambda \geqslant 0, \left( \begin{array}{cc} \lambda M + N & (\lambda-1)c \\ (\lambda-1)c^\top & -f(\lambda) \end{array} \right) \preccurlyeq 0.$$ The semi-definite constraint above exactly implies that the potentials \(v(x)\) and \(w(x)\) are such that \(\forall x \in \mathcal{X}, \ w(x) + f^\ast(v(x)) \leqslant 0\), in exactly the same way as semi-definite programming can be used for optimization through sums-of-squares (see, e.g., [10]). Indeed, for any \(\lambda \geqslant 0\), we have $$\lambda v(x) + w(x) -f(\lambda) = { \varphi(x) \choose 1}^\top \left( \begin{array}{cc} \lambda M + N & (\lambda-1)c \\ (\lambda-1)c^\top & -f(\lambda) \end{array} \right) { \varphi(x) \choose 1} \leqslant 0,$$ which is equivalent to \(w(x) + f^\ast(v(x)) \leqslant 0\).
In earlier work [20], I had developed a similar framework for feature maps that were normalized so that \(\| \varphi(x)\|=1\), leading to the lower bound on \(D(p\|q)\) equal to $$\sup_{M,N} {\rm tr}(\Sigma_p M) + {\rm tr}(\Sigma_q N) \ \mbox{ such that } \ \forall \lambda \geqslant 0, \ \lambda M + N + f(\lambda) I \preccurlyeq 0, $$ which is then equal to an expression common in quantum information theory [21], that is, $$ {\rm tr} \big[ \Sigma_q f( \Sigma_q^{-1/2} \Sigma_p \Sigma_q^{-1/2})\big].$$ The key novelty is that we no longer need normalized features, which makes feature learning significantly easier.
from Gil Kalai
Tobias Boege and Geva Yashfe have posted a remarkable paper, Recognition of algebraic matroids is undecidable. It brings together matroid theory, algebraic geometry, model theory, and undecidability.
An algebraic matroid abstracts algebraic independence in a field extension
. If the elements of
are represented by elements of
, the rank
of a subset
is the transcendence degree over
of the field generated by its representatives. In characteristic zero, Ingleton proved in 1971 that every algebraic matroid is linear over an appropriate extension field. In positive characteristic, however, algebraic matroids include examples that are not linear over any field. In 1975, Ingleton and Main showed that non-algebraic matroids exist: their example was the Vámos matroid
, which violates an extension property of algebraic matroids. The class of algebraic matroids is closed under matroid union and truncation; Lindström later constructed infinitely many excluded minors.
Boege and Yashfe prove that:
The contrast with characteristic zero is striking: recognizing algebraic matroids in characteristic zero is decidable, because there algebraic and linear matroids coincide.
Recognizing whether a matroid has a particular sort of realization is a central question in matroid theory. One precursor is Mnëv’s universality theorem: realization spaces of oriented matroids can model arbitrary primary semialgebraic sets, and the associated realizability problem is complete for the existential theory of the reals. Another is work of Lukas Kühne and Geva Yashfe showing that multilinear representability is undecidable: there is no algorithm to determine whether a matroid can be represented by a -arrangement of vector subspaces for some positive integer
.
Here is a very rough glimpse of the new proof. Classical von Staudt constructions turn incidences of points and lines in a projective plane into equations for their coordinates. The group configuration theorem of Hrushovski and Zilber, together with work of Evans and Hrushovski, allows Boege and Yashfe to recover suitable projective planes from patterns of algebraic dependence. A central difficulty is to identify the Frobenius map using only information recorded by a matroid. They do this by connecting the additive and multiplicative algebraic groups through the affine group. In the resulting structure, the operations commuting with Frobenius give a copy of the rational function field
, with
distinguished. The authors then translate solvability of equations over this field into algebraic realizability of finite matroids. The necessary undecidability theorem for equations over
is due to Pheidas for odd
and Videla for
.
The Vámos matroid has rank four on eight elements. Here is a concrete description. Divide its elements into four pairs
. Every set of at most three elements is independent. Among the 70 four-element sets, exactly five are dependent:
The conspicuous exception is , which is independent. Every other four-element set is independent as well, and therefore a basis. The five dependent four-element sets are called circuit-hyperplanes.
Why does this example matter here? Ingleton and Main showed that this pattern cannot arise from algebraic dependence among elements of a field extension: under the stated rank conditions, dependence of the five listed unions would force to be dependent too. Thus
is not algebraic over any field.
Conspicuous means easy to notice or standing out clearly. So “a conspicuous exception” is an exception that draws attention.
Here are some beautiful figures from the paper.
Authors: Samruddhi Pednekar, Supartha Podder
A longstanding open question in Boolean function complexity asks whether approximate degree composes multiplicatively under block composition. Although a general multiplicative upper bound is known, matching lower bounds have previously been established only for restricted classes of functions. We resolve this question for all total Boolean functions by proving the matching lower bound. Together with Sherstov's upper bound, our result shows that, for every pair of total Boolean functions $f:\{0,1\}^n\to\{0,1\}$ and $g:\{0,1\}^m\to\{0,1\}$, \[ \widetilde{deg}(f\circ g) = Θ\!\left( \widetilde{deg}(f)\,\widetilde{deg}(g) \right), \] where $\widetilde{deg}$ denotes constant-error approximate degree.Authors: Aleksandr Storozhenko
For an arbitrary isotropic log-concave distribution $P$ on $\mathbb{R}^d$, we prove that the polynomial $(Cm)^m\|v\|_2^m - \mathbb{E}_{X\sim P}\langle X,v\rangle^m$ is a sum of squares for every even $m\ge2$, where $C>0$ is a universal constant. This removes the dependence on the Poincaré constant in the theorem of Kothari and Steinhardt (arXiv:1711.07465), recovering the optimal moment bounds for log-concave distributions. As an immediate corollary, we obtain computationally efficient algorithms with dimension-free error guarantees for a wide range of high-dimensional statistical estimation problems. Our proof uses stochastic localization to decompose $P$ as an average of random strongly log-concave measures, whose centered moments admit the subgaussian certificates of Diakonikolas, Hopkins, Pensia, and Tiegel (STOC 2025; arXiv:2410.21194). With a covariance-adapted choice of localization, we show that a fourth-moment certificate derived from Letwin's variance inequality for quadratic forms (arXiv:2607.24164) suffices to control this averaging at every even degree.Authors: Aaron Potechin, Jeff Xu
It is well known that the \Lovasz-Theta function of a random graph $G(n,\tfrac{1}{2})$ is $Θ(\sqrt{n})$. More precisely, it is tightly concentrated in the interval \( [\sqrt{n},\, 2\sqrt{n}], \) where the upper bound follows from an explicit dual witness for the associated semidefinite program. Numerical evidence and heuristic arguments suggest that the true value is $(1+o(1))\sqrt{n}$. However, closing this gap has remained a longstanding challenge, resisting existing techniques even in light of recent progress on sharp algorithmic thresholds and non-asymptotic free probability. In this work, we resolve this question by proving that the \Lovasz-Theta function of $G(n,\tfrac{1}{2})$ is $(1+o_n(1))\sqrt{n}$ with high probability, determining its asymptotic value up to vanishing relative error.Authors: Taotao He, Runfa Hu
We study the computational complexity of peak-oriented rationalizability, a survey based revealed-preference test introduced by Seror (2026). We provide a polynomial-time algorithm for testing rationalizability and recovering a utility function, establishing that peak-oriented preference elicitation is computationally tractable. In contrast, we show that computing the peak-oriented Houtman-Maks index is NP-hard. These results delineate the precise computational boundaries of peak-oriented revealed-preference analysis.Authors: Martin Bichler, Abheek Ghosh
We study equilibrium computation in discrete multiplayer Colonel Blotto games with player-specific battlefield values. In the two-player model with common battlefield values, equilibria can be computed in polynomial time. We show that this tractability breaks down in the multiplayer model with player-specific values under the standard uniform tie-breaking rule. In particular, computing a $(c/n)$-approximate Nash equilibrium is PPAD-hard for some constant $c>0$, even when every player has three resources, where $n$ is the number of players. The main technical step is PPAD-hardness for computing a constant-approximate well-supported Nash equilibrium. In contrast, under uniform tie-breaking, a pure Nash equilibrium can be computed in polynomial time when every player has one resource. We also prove PPAD membership for computing $\varepsilon$-approximate Nash equilibria for inverse-exponentially small $\varepsilon$. Finally, for non-uniform monotone tie-breaking, we show PPAD-hardness even when every player has one resource and all players have identical battlefield values.Authors: Sebastian Ben Daniel
Valiant and Vazirani isolate a satisfying assignment of a circuit with probability $Ω(1/n)$. Dell, Kabanets, van Melkebeek, and Watanabe showed that success above $2/3$ implies $\mathrm{NP}\subseteq\mathrm{P/poly}$ and asked about the range in between. We show that every positive constant already implies the collapse: if a randomized nonuniform polynomial-size pruning procedure succeeds with probability $ε$ on affine circuit inputs with at most $2^{\lfloor 2/ε\rfloor}$ satisfying assignments, then $\mathrm{NP}\subseteq\mathrm{P/poly}$. Success $10/\log L$ on affine inputs with at most $L^{1/3}$ satisfying assignments suffices, where $L$ is the description length, and on inputs with one or two satisfying assignments the threshold $2/3$ drops to $3/5$. No cryptographic assumption is used, and the procedure may read the entire circuit. The proof compiles a pool of circuits into one circuit whose satisfying assignments are indexed by tags in $\mathbb{F}_2^d$. Each member is assigned an affine region of tag space, and if one member is unsatisfiable, the satisfying set shrinks to that member's region. Because regions may overlap and have different dimensions, the collapse reduces to a combinatorial bound: no set of tags meets more than a $2/d$ fraction of an equally weighted family of affine subspaces of all dimensions below $d$ in exactly one point. This regional counting cannot go below order $1/\log L$. The range between $Θ(1/n)$, achieved by affine hashing, and $O(1/\log n)$ remains open.Authors: Kainan Zhou, Chuhong Xu, Gangzhen Qian, Zhaoyi Li
A generative search answer can cite a supported passage yet omit a source relationship that changes its interpretation. We specify a claim-gated audit of the query-source-answer tuple. An omission is resolved only when relationship evidence, answer adoption, materiality, and disclosure are all observed; incomplete evidence remains unresolved rather than being treated as independence. The specification separates this endpoint from citation support and review priority, and binds decisions to versioned evidence spans. A reference checker makes the record contract executable. On an exhaustive synthetic suite, it reproduces all 81 three-state predicate combinations and rejects 192 deliberately malformed records. Common-guard baselines and predicate ablations isolate endpoint logic from missing-evidence handling, while controlled transitions check support separation and evidence removal. These are finite contract-conformance results, not detector accuracy or evidence of improved user outcomes. We define the independent annotation, held-out evaluation, and paired utility tests still required to establish semantic validity and deployment benefit.Authors: Daqing Wan, Jun Zhang
For an $[n,K]$ Reed--Solomon code, the covering radius is $n-K$. Gandikota, Ghazi, and Grigorescu proved deterministic NP-hardness of bounded-distance decoding when the decoding radius is $d$ below the covering radius for every $1\le d\le c\log n/\log\log n$, where $c>0$ is an absolute constant. We prove that, for every fixed rational $0<α<1/2$, bounded-distance decoding is NP-complete under deterministic polynomial-time many-one reductions over explicitly represented finite extension fields for the additive gap $d=\lfloor n^α\rfloor$ below the covering radius. The hard codes have odd block length~$n$, dimension $K=(n+1)/2-d$, decoding radius $(n-1)/2$, and rate tending to $1/2$. The alphabet size is subexponential in the evaluation set size: for a fixed $0<η<1$ depending only on $α$, it is $2^{Θ(n^η\log n)}=2^{o(n)}$. The proof passes through moments subset sum on $n-1$ nonzero field elements, with required subset size $(n-1)/2$ and $d$ prescribed moments. The arithmetic ingredient is a uniform positive-completion theorem over prime fields $\mathbb{F}_q$ with $q\ge d^{2+ρ}$, for any fixed $ρ>0$. A sharper form follows from a higher-dimensional point-count estimate based on Deligne's theorem; the weaker form used in our reduction is proved more elementarily using additive-character orthogonality, the one-variable Weil bound, a moment identity of order $2d$, and Newton identities. A universal completion pool, an extension-field quotient construction, and a deterministic linear-size simultaneous power condenser complete the reduction.Authors: Julien Baste, Simon Dreyer, Clara Marcille, Mikaël Rabie, Ronan Toullec-Streicher
The Zero Forcing (or corruption) of a graph is the problem of finding a minimum-size ``corrupting'' set. It corresponds to a subset of its vertices that can corrupt the whole graph by iterating the following rule: if a corrupted vertex has exactly one neighbor that is not yet corrupted, the neighbor gets corrupted. The iteration of this process comes from the fact that the corruption of a vertex might enable new corruptions (from itself or some of its neighbors). For this reason, one can consider a step of corruption, where all the possible instances of the corruption rule are applied at once. This paper investigates Zero Forcing on temporal graphs, where the topology of the graph evolves throughout the experiment. At each time step (or snapshot) of the graph, a step of corruption is resolved wherever possible. We study the problem of finding a minimum-size corrupting set such that the whole (temporal) graph is corrupted at the end of the experiment. We present a panorama of results, including NP-hardness in some not-so-restrictive scenarios, polynomial algorithms, and a solution to an open question when the whole graph must be corrupted in a single step.Authors: Eshan Chattopadhyay, Pooya Hatami, Chin Ho Lee, Shachar Lovett, Avishay Tal, Emanuele Viola
We prove that the XOR of $k$ majorities on disjoint blocks of \(\ell\) bits has correlation at most \((2d/\sqrt{\ell})^k\) with every degree-\(d\) polynomial over \(\mathbb F_2\). By known techniques, this implies pseudorandom generators with polylogarithmic seed length for low-degree polynomials over $\mathbb F_2$ and for alternating circuits with parity gates.Authors: Zhi-Long Chen, Nicholas G. Hall
We study nonpreemptive scheduling on a single machine with release dates, due dates, positive job weights, and a common processing time. The objective is to minimize total weighted tardiness. Although closely related equal-processing-time problems admit polynomial-time algorithms, the complexity of this problem has remained open in the literature since 2010. We prove that its decision version is strongly NP-complete, even when every job can meet its due date if processed immediately upon release. The reduction is from unweighted MAX-CUT and uses a quadratic number of jobs with polynomially bounded numerical data. Its main ingredient is a constructive normalization theorem that converts every sufficiently inexpensive feasible schedule into a binary choice for each graph vertex; after normalization, total weighted tardiness equals a constant minus a scaled cut value. We also give a deterministic polynomial-time phase-grid assignment algorithm for the shifted objective $Φ=F+p\sum_jw_j$, where $F$ is total weighted tardiness. The algorithm enumerates at most $N$ release-date residues modulo $p$, solves one minimum-cost assignment problem for each residue, and returns the best phase-grid schedule. It runs in $O(N^5)$ arithmetic operations and achieves the tight ratio $3/2-1/(2N)$ for this algorithm. Because the added term $p\sum_jw_j$ is independent of how the jobs are scheduled, the shifted and original objectives have exactly the same optimal schedules. However, the approximation guarantee applies to the shifted objective; for the original objective, the analysis provides an additive bound. Thus, the paper both resolves the long-standing complexity question and provides a complementary worst-case guarantee for the phase-grid assignment algorithm.Authors: Simon Mackenzie
A communication problem can have far more possible inputs than its communication cost would suggest. Must its difficulty already be present on a much smaller set of inputs? We prove that every finite total Boolean matrix of deterministic communication complexity $c\ge4$ has a submatrix on $2^k$ of its original rows and $2^k$ of its original columns, with $k=Θ_\varepsilon(c)$ and complexity at least $(1-\varepsilon)(k+1)$, for every fixed $0<\varepsilon<1$. Since $k+1$ is the maximum possible cost on such a square, the retained problem can be arbitrarily close to maximally hard. This answers affirmatively the lossless condensation question of Hamed Hatami; Göös, Newman, Riazanov, and Sokolov (STOC 2024), who recorded it as Open Problem 2, conjectured a negative answer. Hrubeš previously guaranteed input length $Ω(\sqrt c)$. The same argument gives an original $2^{c-2}$-by-$2^{c-2}$ square retaining at least $c/3-O(\log c)$ bits of communication complexity. The proof builds on Hrubeš's counting and covering argument. We count submatrices equipped with short communication protocols: a player names a covering submatrix, then the players run its protocol. This avoids the loss from converting rectangle partitions into protocols. A recursion on rectangles makes the argument constructive. For fixed rational $\varepsilon$ and any target depth $d\ge4$, a deterministic algorithm returns either a protocol of depth below $d$, or a square of original inputs at input length $Θ_\varepsilon(d)$ with the same near-maximal guarantee. Its running time is $2^{O(2^d)}$ times a polynomial in the table size. If the original complexity is at least $d$, the algorithm necessarily returns the square. An extension gives constant-factor condensation for any fixed number of number-in-hand players, with bounds independent of the finite output alphabet.Authors: Xiao Fan, Jingyuan Li, Yubo Han, Hongbin Guo, Guanya Li, Yang Hu, Wenchao Zhang, Weibin Ji, Yi Zhang
Resting-state functional connectivity (FC) is widely used to classify brain phenotypes and disorders. Most pipelines use the full connectome and seek gains through model design. We instead examine how FC geometry constrains classification and cross-site transfer. Across-subject FC variation concentrates in a small effective subspace, suggesting substantial redundancy in nominal dimensions. Across cohorts, these subspaces may differ in orientation even when their effective ranks are comparable, potentially limiting transfer. Across 2,330 subjects from HCP, ABIDE, and ADHD-200, effective-rank analysis reveals strong spectral concentration. Projection onto leading components at the effective-rank scale recovers most of the full-FC classification performance. In ABIDE, site-specific effective subspaces are weakly aligned, and their principal-angle overlap predicts pairwise transfer after covariate adjustment despite comparable per-site effective ranks. Controlled rotations that alter subspace orientation while preserving the mean and covariance spectrum drive transfer toward chance, whereas displacement-matched label-orthogonal rotations do not. These results identify subspace orientation as a key factor in transfer degradation under controlled perturbations. This study offers a geometric diagnostic of FC generalization and suggests evaluating cross-site harmonization by its ability to align effective subspaces alongside classification accuracy.Authors: Evanthia Papadopoulou, Martin Suderland, Zeyu Wang
We study a \emph{collapse process} to construct the farthest Voronoi diagram of lines in three dimensions, given a spherical map of the diagram's unbounded features. The collapse process sweeps through the diagram in order of decreasing distance from the farthest lines. It follows the \emph{shrinking map}, a cell complex on a topological sphere that encodes the locus of points with a fixed farthest distance. We show that, in three dimensions, the collapse process has exactly four non-terminal local event types that change the structure of the shrinking map: \emph{deletion, swap, local minimum}, and \emph{local maximum} events; plus one terminal event. This list is complete. We give intrinsic three-dimensional geometric descriptions of the four non-terminal events. First, we classify the two vertex-related events, deletion and swap events, by the spherical convex hull of the four tangent points from a vertex to its four defining lines. Then, we analyze the events related to the local extrema of the distance function along the trisector of three lines. We show that the distance function along a trisector has at most $4$ local maxima and $8$ local minima, and that both bounds are tight. The extrema can be found via a polynomial of degree $12$. As a byproduct, this gives a direct method for finding the smallest sphere tangent to three given lines. At each local extremum, the tangent sphere touches the three lines at points lying on a great circle. The collapse process and the completeness of the event list also apply, under similar general position assumptions, to the farthest Voronoi diagram of convex sites under strictly convex distance functions.Authors: Nandana Ghosh, Ankush Acharyya, Rakesh Gupta, Supantha Pandit
We study exact one-center optimization for axis-parallel segments using axis-parallel squares under endpoint-based coverage, where a segment is \emph{$1$-covered} if the square contains at least one of its endpoints. In the monochromatic problem, we seek a minimum-side-length square that $1$-covers all $n$ input segments. We obtain $O(n\log n)$-time algorithms for both unrestricted and segment-constrained centers. The unrestricted bound matches the known bound implied by the two-representative color-spanning-square problem, whereas the segment-constrained result is new. We also prove matching $Ω(n\log n)$ lower bounds for both center models in the fixed-order algebraic decision-tree model. In the bichromatic problem, an admissible square must fully contain all $m$ blue segments, minimize the number of red segments with an endpoint in its interior, and, subject to this minimum, maximize its side length within a prescribed bounding box. We give deterministic $O(m+n\log^2 n)$-time and $O(m+mn\log n)$-time algorithms for unrestricted and blue-segment-constrained centers, respectively.Authors: Ignat Romanov, Andreas Hadjipieris, Neofytos Dimitriou
Smart maritime infrastructures provide continuous access to heterogeneous sensing streams, enabling repeated experimentation, digital-twin development, and AI-based maritime services. However, sensing hardware alone is not sufficient for scene-specific model development: historical video streams must also be spatially indexed, contextualized, and reduced to informative subsets for annotation. This paper presents a frame-to-panorama localization and context-aware sampling pipeline for ship detection in historical PTZ maritime video lacking reliable pan, tilt, and zoom metadata. The main contribution is an end-to-end data-curation approach that recovers camera-view information from historical PTZ video and combines it with environmental context and visual diversity to construct compact, scene-specific training sets. Specifically, frames are localized on a reference panorama using SuperPoint and LightGlue, enriched with weather and solar-state metadata, and selected through diversity sampling to preserve variation across camera view and environmental conditions. A second context-aware stage targets under-represented distant-vessel cases near the horizon using tile-level visual embeddings and Gaussian Mixture Model clustering. Applied within the CMMI MDigi-I Smart Marina testbed, the proposed pipeline reduces 40,718 candidate frames to 220 images for annotation, corresponding to a 99.5% reduction. A YOLO26-m detector fine-tuned on this subset achieves a mean AP50 of 94.78% $\pm$ 0.51% and a mean AP50-95 of 75.10% $\pm$ 1.73% under sequence-grouped five-fold cross-validation. These results demonstrate that highly redundant infrastructure video streams can be transformed into compact, spatially and contextually diverse training sets for scene-specific detector adaptation while substantially reducing annotation effort.Authors: Jerrin Bright, John Zelek
3D foundation models recover video cameras and geometry in one forward pass, but some of the strongest are up to scale. Joint people-scene reconstruction then requires two missing outputs: metric scale and persistent person identity. We ask whether one up-to-scale foundation representation can support both through lightweight adaptation. Exact metric labels are scarce, but unlabeled in-the-wild video is abundant. We use people in curated web video to initialise the solution: a posed metric body and 2D keypoints give an approximate, closed-form scale pseudo-label. These pseudo-labels pretrain a Scale Readout, which is then fine-tuned together with a lightweight adapter using exact metric supervision from standard real-video training splits. At inference the head predicts metric scale from foundation-model tokens, without the ruler or its teachers. For person identity, we probe the pretrained foundation model alone and find evidence that its intermediate query-key features encode person correspondence across frames. In most evaluated moving-person clips, a mid-layer token prefers that person over the vacated location and other people. A tiny projection reads this correspondence; together with metric pelvis motion and proposal confidence, it drives dustbin-aware Sinkhorn association of per-frame bodies. WildHSR combines both readouts to reconstruct metric cameras, scene and people from monocular video. Each window is predicted feed-forward; analytic association and Sim(3) composition connect windows. On EMDB-2, WildHSR is the first feed-forward method in the published comparison to beat the best optimization-based WA-MPJPE and RTE while leading feed-forward methods on all three world-frame metrics. On RICH, it leads feed-forward people-and-scene methods on WA-MPJPE and W-MPJPE. The complete pipeline runs at 10.1 fps on one GPU.Authors: Ruoxi Shi, Xinyue Wei, Fanbo Xiang, Zexiang Xu, Hao Su
Approximate convex decomposition (ACD) converts triangle meshes into small sets of convex parts and is a standard preprocessing step for physics simulation, collision detection, and large-scale robot learning. The majority of modern ACD methods produce high-quality decompositions through an expensive search over candidate cutting planes, with per-mesh runtimes of tens of seconds that force game pipelines into overnight bakes and keep articulated-object datasets on CPU clusters for days. Prior work has accelerated isolated stages, most recently VisACD's GPU-based visibility metric, yet the dominant costs -- search, mesh cutting, and convex hull construction -- have remained on the CPU because their natural decomposition into many small homogeneous phases trails off in a fading last wave at every kernel boundary, and the variable-sized output of each phase forces a host round trip simply to allocate the next launch's input. We address these obstacles by adopting the warp, rather than the thread or thread block, as the unit of algorithm design, an idea introduced in the graph-processing community for a different pathology and which we adapt here to fuse the many heterogeneous phases of a computational-geometry pipeline into single warp-resident kernels, paired with a device-side heap allocator that lets the buffers between fused phases be sized and allocated on the device. Building on this template, we present CuACD (CUDA ACD), the first fully GPU-resident ACD system, together with a suite of reusable GPU components, released as open-source standalone CUDA modules that drop into any search-based ACD pipeline. On the V-HACD benchmark, PartNet-Mobility, and an Objaverse subset, CuACD achieves more than an order of magnitude of speedup over CoACD at matched or better quality.Authors: Ainesh Bakshi, Apoorv Vikram Singh, Xinyu Tan
We prove that the spectrum of an unknown $d$-dimensional quantum state can be estimated to error $\varepsilon$ in total variation distance using \[ O\!\left(d^2\min\left\{ \frac{1}{(\varepsilon\log d)^4},\; \frac{1}{(\varepsilon\log d)^2} \right\}\right) \] copies. This matches the recent lower bound of Wang. When restricted to unentangled measurements, we give an algorithm with an additional factor of $d$ in copy complexity, which we conjecture to be optimal. We develop a framework for recovering the small eigenvalues of a quantum state by matching Chebyshev moments. We bound the variance of each Chebyshev moment estimate in terms of scalar derivatives of the corresponding polynomial, using classical and quantum Efron--Stein decompositions. Different rescalings of the Chebyshev polynomials balance approximation error and variance, yielding two regimes in our copy complexity bound.Authors: Max Göttlicher, Lennard Hofmann, Christoph Niederbudde
Zero Forcing (ZF) and Power Dominating Set (PDS) mark vertices in a graph based on a common forcing process starting from a problem-specific set of initially marked vertices. ZF initially marks the selected vertices while PDS additionally marks their neighbors. In the forcing process, a marked vertex with only one unmarked neighbor may force that neighbor which then becomes marked, too. A solution marks the entire graph by exhaustive application of this rule. One variant generalizes the forcing threshold; vertices may force when they have a fixed number of $k$ unmarked neighbors. Another variant limits propagation to a fixed number of rounds. We classify the parameterized complexity of the problem variants obtained by combining these choices of initialization, round limit and forcing threshold. We show that with appropriate choices, these variants range in parameterized complexity from fixed-parameter tractable to complete for every even layer $W[2\ell]$ of the $W$-hierarchy, and up to $W[P]$-complete. Our results demonstrate that small changes in any one of these three dimensions can lead to a sharp change in problem complexity.Authors: Bo Liu
Choosing Shellsort gaps is a well-known open problem. For over sixty years, successful sequences have relied on human-designed formulas, numerical searches, or number-theoretic constructions. Although stronger general bounds exist for dense or mainly theoretical families, the worst-case upper bound for a short, sparse, and practically competitive construction has not advanced beyond $N^{4/3}$ for decades. We ask whether the sequence itself can instead be learned from execution. We present an RL-driven, self-supervised system that searches over executable gap generators. Every proposal is valid by construction, and executed candidates return exact comparison and move counts; no classical sequence is used as a target. Across five independent searches, the system discovers a common rational-geometric family. A second self-supervised stage tunes only a finite prefix, producing the practical sequence $1,3,8,20,47,116,300,585,1416,3303,\ldots$. Once frozen, it obtains the lowest equal-task average operation count among seven classical baselines on 25 large tasks with $10^7Authors: Jesse Beisegel, Ekkehard Köhler, Robert Scheffler, Martin Strehler
We study a variation of the classical Maximum (Minimum) Leaf Spanning Tree problem. In many applications, Depth-First Search (DFS) is used to compute a spanning tree of a graph. Such a search tree is constructed by connecting each vertex $v$ with the last vertex the search has visited before $v$ and we call this a last-in tree. By restricting the Maximum (Minimum) Leaf Spanning Tree problem to last-in trees of a graph search, we ask for a search ordering that leads to the largest (smallest) number of leaves in its search tree. Recently, Bergougnoux et al. [Journal of Computer and System Sciences 154 (2025)] have studied the parameterized complexity of these problems for DFS. They showed that the minimization problem is para-$\mathsf{NP}$-hard and the maximization problem is $\mathsf{W}[1]$-hard when parameterized by the number of leaves. When parameterized by the number of internal vertices, both problems have polynomial kernels. Here, we examine whether these results also hold for the variant Lexicographic DFS (LDFS). We show that the hardness results of DFS can be transferred to LDFS. We also present exponential kernels for the number of internal vertices as the parameter. We complement this by showing that polynomial kernels do not exist, unless $\mathsf{NP} \subseteq \mathsf{coNP} / \mathsf{poly}$. We also consider last-in trees of searches that do not follow the DFS scheme. In contrast to (L)DFS, minimizing the number of internal vertices is para-$\mathsf{NP}$-hard for several searches including Breadth-First Search.Authors: DongYun Byun, Akira Matsubayashi
The $k$ vertex-disjoint paths problem asks whether, given a graph $G$ and $k$ pairs of vertices $(s_1,t_1)$, \ldots, $(s_k,t_k)$, $G$ has $k$ pairwise vertex-disjoint paths connecting $s_i$ and $t_i$ for all $1\leq i\leq k$. If $G$ is undirected, then this problem is NP-complete, but there exist FPT algorithms parameterized by $k$.Since these algorithms involve an extremely large function on $k$, algorithms for restricted graphs have also been investigated. In particular, a $2^{2tw\log tw+O(tw)}\cdot n$ time algorithm for undirected graphs with $n$ vertices and treewidth $tw$ is proposed by Scheffler (Technical Report 396, TU Berlin, '94), and it is proved by Lokshtanov, Marx, and Saurabh (SIAM J. Comput. '18) that, under the ETH, there exists no $2^{o(pw\log pw)}\cdot n^{O(1)}$ time algorithm for either directed or undirected graphs with pathwidth $pw$ and for $k=Ω(pw^4)$. It has not been known whether the lower bound also holds for a smaller $k$. In this paper, we prove that, for both the directed and undirected cases, there is an algorithm faster than Lokshtanov et al.'s lower bound for $k=tw^{o(1)}$ by proposing a $2^{O((tw+k)\log k)}\cdot n$ time algorithm. Besides, we prove a lower bound that, under the SETH, there exists no $(2-ε)^{pw\log pw}\cdot n^{O(1)}$ time algorithm for directed graphs and for a general $k$. This lower bound is tight because, with slight modifications, Scheffler's algorithm runs in $2^{pw\log pw+O(pw)}\cdot n$ time also for directed graphs.Authors: Kyungjin Cho, Eunjin Oh, Sebastian Wiederrecht
We study the \textsc{$k$-Disjoint Paths} problem on a graph embedded on a surface with bounded Euler genus. Given a graph $G$ with $n$ vertices and $k$ vertex pairs embedded on a surface of Euler genus $g$, we present a $2^{O(k^2+g^2)}n$-time algorithm that computes $k$ pairwise vertex-disjoint paths connecting the given vertex pairs if such paths exist. Our approach relies on the decomposition of $G$ into $O(k+g)$ planar subgraphs while bounding the complexity of the boundaries between these subgraphs. This approach enables the use of techniques for compressing linkages in planar graphs. Moreover, our techniques yield two kernels of size polynomial in $k$, $g$, and the treewidth of the graph, and of size $2^{O(k+g)}$. These results extend recent advances on \textsc{$k$-Disjoint Paths} on planar graphs [Cho et al. SODA 2023] and [Włodarczyk and Zehavi FOCS 2023] to surface-embedded graphs.Authors: Santosh S. Vempala
We prove nearly quadratic lower bounds for randomized algorithms for linear optimization and uniform sampling over convex bodies in the membership oracle model. For linear optimization, this matches the known nearly quadratic upper bound up to a polylog factor in the dimension. For uniform sampling, this improves on the previous linear lower bound. Our construction also implies the same lower bound for volume estimation.Authors: Jonas Ellert, Lukas Nalbach
The Lempel-Ziv (LZ77) factorization decomposes a text into the least possible number $z$ of phrases that each refer to an earlier occurrence. It is this phrase count, rather than the encoded size, that governs the size of LZ-based compressed indexes, and computing a factorization with few phrases is a time and space bottleneck in their construction. In practice, computing LZ77 quickly has so far required building a suffix array. Ellert [SPIRE 2023] gave algorithms that compute the exact LZ77 factorization, and a 3-approximation of it, in sublinear working space. They have remained unimplemented, because two of their components resist a direct implementation: a lookup table that degenerates to patterns of length at most two for realistic inputs, and an orthogonal range reporting data structure that is impractical. We replace both, fine-tune every remaining stage, and obtain the first practical implementation, which runs in space close to the text rather than to the suffix array. On one thread, our 3-approximation factorizes 12-19x faster than the classical LPF algorithm while using 14x less memory; on 32 threads, even our exact algorithm is 1.4--2.9x faster than parallel LPF, at 9x less memory. In practice the approximation ratio stays far below 3. As a side result, passing only its perfect phrases to a downstream compressor yields a precompressor that is on par with the state of the art [Dinklage, SEA 2026] in compression ratio, and better in memory consumption and parallel throughput.Authors: Johannes Fischer, Lukas Nalbach
Approximate pattern matching (APM) on highly repetitive texts is a central task in bioinformatics. Bi-directional r-indexes support left- and right-extension of a pattern and thereby accelerate APM algorithms based on search schemes, but existing variants -- br-index and b-move -- suffer from two bottlenecks: $O(σ)$ character-predecessor/-successor queries on the run-length-encoded BWT per extension, and predecessor queries on sparse bit vectors to maintain a value in the suffix array interval and to access the PLCP array while locating. We present Move-rb, a bi-directional r-index built on the optimized r-index Move-r. Although Move-rb is 2x larger than br-index, it is up to 24% smaller than b-move, answers APM queries 1-4 orders of magnitude faster than br-index, 1.9-10x faster than b-move and up to 5.5x faster than the state-of-the-art bi-directional (uncompressed) FM-index columba, which is 36-42x larger. Memory usage (including index size) during APM locate queries is reduced by 1.5x (up to 6.5x) for Hamming distance and 2.5x (up to 7.4x) for edit distance. Move-rb can be constructed 3-14x faster while using 19-141x less memory than br-index, b-move and columba. A variant using a relative Lempel-Ziv-encoded suffix array locates up to 10x faster while being 1.2-2.4x larger. We achieve these speedups by optimizing index operations and search scheme APM algorithms: Without a direction switch, Move-rb computes an all-$k$-character extension in output-optimal $O(k)$ time and a single-character extension in the same time. Augmenting Move-rb with $O((r+\overleftarrow{r})\log(n/r_{\min})\logσ)$ bits reduces a single-character extension to $O(\logσ)$ time, where $r$ and $\overleftarrow{r}$ are the numbers of runs in the BWT of the text and its reverse, and $r_{\min}=\min(r,\overleftarrow{r})$. A direction switch incurs only $O(\log\log_ω(n/r_{\min}))$ additional time.Authors: Pravesh K. Kothari, Andrew D. Lin, Peter Manohar
In this paper, we consider the problem of refuting $\mathbb{F}_2$-linear equations with random right-hand sides. Formally, we give a sub-exponential $2^{O(n/\log n)}$-time randomized algorithm that takes as input an arbitrary $m \times n$ matrix $A$ and a uniformly random vector $b \in \mathbb{F}_2^m$, and outputs a witness showing that no assignment satisfies more than a $\frac{1}{2}+ε$ fraction of the equations provided that $m \geq 2^{O(n/\log n)}$. The setting above is the semirandom refutation variant of the famous work [BKW03] that gives a $2^{O(n/\log n)}$-time search algorithm for the learning parity with noise (LPN) problem with $m \geq 2^{O(n/\log n)}$ equations. Building on the search algorithm of [Lyub05], we also give a $2^{O(n/\log \log n)}$-time refutation algorithm that succeeds with only $m \geq n^{1 + γ}$ equations, for a small constant $γ$. Finally, we prove that our algorithm is not captured by the sum-of-squares hierarchy by proving a degree-$Ω(n)$ sum-of-squares lower bound, showing that ''[BKW03]-style'' algorithms achieve better runtime than can be done under sum-of-squares. We thus obtain a natural example of a noise-tolerant signal recovery problem that exhibits a nontrivial gap between the performance of efficient algorithms and that of those based on the sum-of-squares hierarchy.Authors: Xiaoyu Li
The Komlós conjecture, now a theorem, asserts that whenever the columns of a matrix $A\in\mathbb{R}^{m\times n}$ have Euclidean norm at most one, some signs $\varepsilon\in\{-1,1\}^n$ make every coordinate of $A\varepsilon$ bounded by an absolute constant. Guo, Fang, and Lu gave the first polynomial-time algorithm for finding such signs, a deterministic spectral signing procedure with discrepancy $8272$ and running time $O((mn^9+n^{10})\log(m+n))$. We give a deterministic algorithm that finds signs with $\|A\varepsilon\|_\infty<99$ using $O(mn+n^{ω+2}\log^3 n)$ arithmetic operations, where $ω>2$ is any fixed attainable matrix-multiplication exponent; with the current bounds on $ω$ this is $\widetilde O(mn+n^{4.372})$. Our algorithm uses the same framework: it rounds a single fractional coloring and watches all rows through the top eigenvalue of a Gram matrix of energy-corrected barriers. Steps follow flat directions, rescaled so that no barrier near its threshold moves faster than a constant, and the regularizer grows as coordinates freeze; together these bound the number of updates by $O(n^2\log n)$. A weak quadratic charge on the tracked row sums leaves $O(n\log^2 n)$ rows to evaluate at any time, and a motion clock bounds when any other row could approach its barrier. Each update is a short sequence of matrix products. Its direction is read off by conditional expectations from a polynomial soft projector, its small constraint residual is repaired in affine row representations, and one identity accounts for every change of representation. For rational input the algorithm has polynomial bit complexity.Authors: Marius Bächler, Markus Chimani, Henning Jasper
In this paper, the parameterized complexity of the multiplicative $α$-spanner problem with independent weights and lengths on undirected graphs is considered for the first time. All prior FPT results (except one on DAGs) assume basic instances (i.e., with unit weights and lengths) and are parameterized in the stretch factor $α$ and the (in practice typically non-constant) number of removed edges. We show that several parameterizations do not allow FPT algorithms. However, our exclusion approach generalizes an existing algorithm for basic instances to arbitrary weights and lengths. It is parameterized by the total removed weight and a new tightness parameter. The latter is more precise than $α$ and allows us to also improve the best known result for basic instances. Our second algorithm, called inclusion approach, uses the natural parameterization in the spanner's total weight. We prove that this sole parameter leaves a W[2]-hard problem, but also show FPT algorithms exist when augmented with secondary parameters.Authors: Tomohiro Koana, Soh Kumabe
We give a deterministic polynomial kernel for Odd Cycle Transversal, derandomizing the randomized kernel of Kratsch and Wahlström (TALG 2014). Our algorithm uses a deterministic polynomial-time construction of almost multilinear representations of gammoids. Such a representation assigns a block of columns to each element so that, for every subset of elements, the normalized matrix rank approximates its matroid rank to within a prescribed additive error $δ$. The construction builds on recent breakthroughs in NC algorithms for matching. Our kernelization algorithm then computes the required representative families from these representations.Authors: Yihan Zhang
We initiate the study of approximating the top eigenvalue and eigenvector of a random symmetric matrix $ A \in \mathbb{R}^{n\times n} $ using $ q(A)b $ where $q$ is a degree-$d$ polynomial and $b$ is a standard Gaussian vector independent of $A$. For spiked GOE $ Y = λvv^\top + X $, we identify $ d_\star = \frac{\log(n)}{2\log(λ)} $ to be the critical degree threshold above which accurate approximation of the top eigenvalue and eigenvector is possible. This sharpens the common belief that spectral methods can be implemented by $ O(\log(n)) $-step power iterations and offers a precise connection between spectral methods and low-degree polynomial algorithms, a popular proxy for all polynomial-time algorithms. For GOE $X$, we identify $ d_\star = n^{1/3+o(1)} $ to be the critical degree threshold for top eigenvector approximation, whereas constant degree suffices for top eigenvalue approximation. Moreover, in the limit where $ d/n^{1/3} $ converges to a positive finite constant, we compute the exact asymptotic eigenvector approximation accuracy in terms of the expected squared overlap. These results significantly improve upon predictions made in randomized numerical linear algebra for deterministic data matrices that the iteration count of power methods is governed by the inverse spectral gap. Technically, our analyses leverage extremal properties of Chebyshev polynomials and draw upon the rich literature of random matrix theory.Authors: Ilan Doron-Arad, Hadas Shachnai, Gilad Shmerler
We study the following class of matroid optimization problems with a linear constraint (P-MOL). Given a matroid M=(E,I), two weight functions $v,w:E\to R_{\ge 0}$, and a threshold $L\in R_{\ge 0}$, find $opt v(S)$ where S is either an independent set or a base of M satisfying a budget-type constraint: $w(S)\le L$ or $w(S)\ge L$, and $opt\in\{min,max\}$. P-MOL provides a unified representation for a broad family of NP-hard optimization problems, including budgeted matroid independent set, constrained minimum-basis, and knapsack-cover variants with a matroid constraint. Also, it naturally extends to multiple matroid constraints. In particular, we consider the matroid intersection cover (MIC) problem, where feasibility is defined by the common independent sets of two matroids and one seeks minimum $v(S)$ subject to $w(S)\ge L$. Our main result is a unified EPTAS for all nontrivial P-MOL variants, obtained by generalizing a technique of Hassin and Levin (SIAM J. Comput., 2004) for solving the constrained minimum spanning tree problem. Specifically, for any fixed $ε>0$, we present an algorithm running in time $|E|^{O(1)} (1/{ε^2})^{O(1/ε)}$ that outputs a feasible solution S whose value is at most $(1+ε)OPT$ for minimization variants and at least $(1-ε)OPT$ for maximization variants. This resolves the complexity status of all members of P-MOL, as none of these problems admits an FPTAS (Doron-Arad, Kulik and Shachnai, ICALP'24). Finally, we separate the P-MOL family from its extension to matroid intersection. We show that an EPTAS is unlikely to exist for the matroid intersection variant of P-MOL under a covering constraint, whereas an EPTAS is known to exist under a budget constraint. This highlights a qualitative difference between these two types of linear constraints that does not arise in the single-matroid setting.from TCS+ Seminar Series
The next TCS+ talk will take place this coming Wednesday, September 30th at 1:00 PM Eastern Time (10:00 AM Pacific Time, 19:00 Central European Time, 17:00 UTC). Sepehr Assadi from University of Waterloo will speak about “Greedy is Optimal for the Semi-Streaming Matching Problem” (abstract below).
You can reserve a spot as an individual or a group to join us live by signing up on the online form. Registration is not required to attend the interactive talk, and the link will be posted on the website the day prior to the talk; however, by registering in the form, you will receive a reminder, along with the link. (The recorded talk will also be posted on our website afterwards) As usual, for more information about the TCS+ online seminar series and the upcoming talks, or to suggest a possible topic or speaker, please see the website.
Abstract: We prove that no single-pass semi-streaming algorithm (deterministic or randomized) can achieve a better-than-half approximation to the maximum matching problem. This implies the optimality of the naive greedy algorithm, answering a longstanding open question in graph streaming literature since the introduction of the model. Our proof consists of two main parts:
1. Blueprint framework: reducing the problem of proving lower bounds for semi-streaming matching to constructing certain combinatorial objects which we call blueprints; and,
2. Blueprint construction: an optimal construction of such blueprints usable within this framework.Putting these two parts together implies our semi-streaming matching lower bound.
Based on joint work with Max Jiang and Mars Xiang in https://arxiv.org/abs/2607.14644 (STOC 2026) and https://arxiv.org/abs/2607.14656 (arXiv; July 2026)
Authors: Tanmay Inamdar, Pallavi Jain, Pranjal Pandey
Envy-freeness is a fundamental notion of fairness in the allocation of indivisible goods. In this paper, we study envy-free allocation under Leontief preferences, which model perfect complements. Although Leontief preferences have been extensively studied in the context of allocating divisible goods and market equilibria, they have received comparatively little attention for the allocation of indivisible goods. We show that, unlike additive valuations in cardinal preferences, an envy-free allocation always exists for Leontief preferences when there are at least two goods. In contrast, envy-free allocations may fail to exist when there is a single good, however it can be decided in polynomial time. We next study the problem of computing a welfare-maximizing envy-free allocation. We prove that this problem is NP-hard in general, whereas it is polynomial-time solvable when there is only a single good or agents have identical demands. Finally, we investigate the parameterized complexity of this problem.Authors: Alex Bortolotti, Monaldo Mastrolilli
The Ideal Membership Problem (IMP) asks whether a polynomial f belongs to an ideal of Q[x_1, ..., x_n]. Polynomial Calculus (PC) certifies membership by deriving f from the generators, and a degree-d derivation needs at most n^O(d) steps. We write PC-IMPd for the problem of producing a degree-bounded PC certificate, and call it solvable when one is guaranteed to exist and can be found in time n^O(d). Over Q, unlike over finite fields, a derivation may need exponentially many bits. We study PC-IMPd on instances arising from constraint satisfaction problems, and ask for which constraint languages L it is solvable. Our main contribution is a reduction framework for PC-IMPd, based on pp-definitions, pp-interpretations, and pp-encodings, that mirrors the algebraic approach to CSP complexity. Solvability is preserved by these constructions and, in the language of algebras, by passing to subalgebras, finite direct powers, and homomorphic images. We obtain new tractable classes over ternary and larger domains: every language closed under the median operation on a finite chain has solvable PC-IMPd, by reduction to the Boolean majority algebra, and in particular so does every language over {0, 1, 2} closed under a fixed-value majority. This also places IMPd(L) in P for such languages, advancing the classification of IMPd over ternary domains. In the process, we settle the last open case of the Boolean dichotomy for IMPd(L) and complete the Boolean classification of PC-IMPd(L) with an unconditional lower bound for an instance of PC-IMP1. A recent PC-to-SoS simulation reduces degree-automatability of Sum-of-Squares (the open problem of finding a degree-d SoS proof in time n^O(d) when one exists) to solvability of PC-IMPd. Each new tractable class therefore yields a family of constraint systems on which SoS proofs are degree-automatable.Authors: Yeyuan Chen, Miryam Mi-Ying Huang, Yinchen Liu, Er-Cheng Tang
For every fixed integers $r \ge 4$ and $2 \le k \le r-2$, we construct $r$-dimensional cubical sheaf complexes whose degree-$k$ CSS codes have positive constant rate, linear distance, and constant soundness, with bounded row and column weights. Taking $r=4$ and $k=2$ gives a family of asymptotically good binary qLTCs. At the core of our construction is a uniform product-expansion theorem for explicit Reed-Solomon codes on norm-one evaluation sets. The key point is that the expansion constant stays bounded away from zero as the local code lengths grow. We place these codes on arithmetic cubical complexes, obtaining constant local expansion for both the resulting sheaf and its dual. Together with the local-to-global framework of Dinur, Lin, and Vidick (FOCS 2024) and sheaf duality, this gives linear distance and constant soundness, while an asymmetric choice of local code dimensions gives positive rate. The resulting codes are explicit and polynomial-time computable.Authors: Tianrun Zhao
We prove quantum soundness of the total-degree diagonal line-vs-point test using the individual-degree soundness theorem of Ji, Natarajan, Vidick, Wright, and Yuen. A random change of coordinates yields projective polynomial decoders of total degree at most $md$. The uniform-line slice of the test bounds the weight of outcomes of degree greater than $d$, which are removed by a common relabeling. This reduction does not yield a dimension-independent soundness bound: the $\operatorname{poly}(m)$ dependence of the individual-degree theorem persists, as discussed in Section 1.2 of arXiv:2009.12982.Authors: Andrea Serani, Kevin Maki
This paper presents a representation-aware framework for geometry-based screening of hull-form producibility at early design stages. The proposed signature combines dimensionless total and signed developability deviation with curvature-class area fractions, distributed fields, metric-specific validity, and representation provenance. These descriptors characterize surface features relevant to plate forming and developability, but are not calibrated predictors of fabrication cost, forming effort, or process feasibility. Native IGES/STEP boundary representations (BReps) are evaluated through direct differential geometry and trimmed-domain quadrature, whereas triangulated surfaces use discrete curvature recovery and area-weighted aggregation. Analytical and semi-analytical controls verify the formulation, while matched-face BRep-to-mesh tests assess discrete curvature recovery. Application to DTMB 5415, KCS, JBC, and KVLCC2M shows that curvature intensity and areal extent provide complementary information and that derivative-based outcomes can be representation sensitive. KCS, for example, exhibits approximately 24% greater developability deviation than DTMB 5415, while double-curved regions occupy 72.9% of its valid surface versus nearly the entire DTMB valid surface. The resulting quantities provide an early geometric screening layer for subsequent use as objectives, constraints, surrogate responses, or design-space features. HullProd, the companion open-source software, implements the signature, distributed fields, validity, and provenance.Authors: Sergey Pupyrev
A graph is 1-planar if it has a drawing in which every edge is crossed at most once. We show that the smallest cubic non-1-planar graphs have $30$ vertices. Two such graphs are the Tutte-Coxeter graph of girth eight and a graph of girth seven that we call the Byte graph. Every subcubic graph with fewer than $30$ vertices is 1-planar. Our proof is computer-assisted, but directly testing all relevant graphs is impractical. To establish non-1-planarity of the two graphs, we extend a SAT-based solver with a custom clause propagator based on separating cycles and a case split based on graph automorphisms, allowing independent cases to be solved in parallel. To show that all smaller subcubic graphs are 1-planar, we introduce the concept of $k$-flexibility: every set of at most $k$ prescribed edges can remain uncrossed in some 1-planar drawing. We use this property to reconstruct 1-planar drawings of larger graphs from drawings of smaller $k$-flexible graphs. This replaces exhaustive testing of more than forty billion cubic graphs with computations on far fewer graphs of smaller order.Authors: Shubhada Aute, Fahad Panolan, Geevarghese Philip
An edge weighting of a graph induces a coloring of its vertices in which the color of a vertex is the total weight of the edges incident with it. Such an edge weighting is proper if adjacent vertices always receive distinct colors. Deciding whether a graph admits a proper weighting is known to be NP-complete for the weight set $\{0,1\}$, and also for $\{1,2\}$. In recent work (arXiv:2604.12363) we showed that both problems are FPT parameterized by the vertex cover number $k$, but it was open -- to the best of our knowledge -- whether either parameterized problem had a polynomial kernel. In this work, we show that both problems have polynomial kernels when parameterized by $k$. We also show that both problems are W[1]-hard parameterized by treedepth, answering another question from our earlier work. We then study the pre-weighted versions of the two problems, in which the weights of some edges are fixed in advance, and the task is to extend the assignment to a proper weighting of the whole graph. We show that both pre-weighted problems are FPT parameterized by the vertex cover number $k$. For the $\{1,2\}$ version the running time is $2^{O(k \log k)} \cdot n$; for the $\{0,1\}$ version we obtain the same running time when every pre-weight is $1$, and a slower FPT algorithm in the general case. We also show that both pre-weighted problems are W[1]-hard parameterized by either of (i) the feedback vertex set number or (ii) the treedepth of the input graph. Since a graph with no pre-assigned weights is a special case, our algorithms for the pre-weighted versions solve the two original problems as well, in time $2^{O(k \log k)} \cdot n$, significantly improving on the bound of $2^{O(k^4)} \cdot n^{O(1)}$ from our earlier work.Authors: Siddharth Pritam
Existing link-based combinatorial preprocessing methods speed up the computation of persistent homology by removing a vertex or edge only when its link remains a cone. We replace this condition with a quantitative homological certificate. The reduced homology of the filtered link of a generator (a vertex or edge) defines a local obstruction module whose future part describes the effect of deleting that generator. Its barcode certifies either exact deletion or an explicit bound on the bottleneck error, and a conflict colouring extends this guarantee to families of generators. Our implementation, HomTrim, removes an additional 13% to 41% of the input edges beyond domination-only preprocessing and reduces backend persistence time by factors ranging from 1.55 to 4.61 on weighted flag filtrations. The same module also yields regularity diagrams that measure how far a generator can move before becoming visible to homology, together with a multiscale stability result.Authors: Gregor Diatzko, Félix Lasseux, Sabine Storandt
When designing algorithms for geometric graphs, exploiting structural parameters can lead to significantly improved bounds. Two prominent parameters in this context are $c$-packedness and $λ$-low density, both of which locally restrict graph complexity. Parameterized algorithms based on these parameters have been developed for computing well-separated pair decompositions, balanced separators, as well as distance oracles. Nevertheless the practical applicability of algorithms parameterized by $c$ or $λ$ remains unclear. While $c$-packed and $λ$-low-density graphs have been proposed as realistic models for road networks, the actual parameter values of large real-world instances have so far remained unknown, and existing theoretical guarantees are partially too loose for practical usage. In this paper we first devise scalable implementations for the approximate computation of $c$ and the exact computation of $λ$. Our experiments on road networks with millions of edges reveals a significant gap between the two parameters. On the theoretical side we prove that $c\in O(λ\sqrt n)$ which complements the known result that $λ\in O(c)$. Furthermore we present improved parameterized algorithms for balanced separator computation that reduce the separator size in theory and practice. We also show how to compute a tree decomposition with a width linear in the respective parameterized balanced separator size in polynomial time. This structural result yields a variety of new algorithmic consequences. Among them is an exact distance oracle with query time $O(c)$ for $c$-packed graphs after polynomial-time preprocessing, which improves upon the previous $O(c\log n)$ bound. Our experiments show that the proposed techniques efficiently produce small balanced separators and enable the construction of concise exact distance oracles on large road networks.Authors: Alexandr Andoni, Manuel Paez, Krish Singal, Tian Zhang
We study the problem of constructing metric spanners in general metric spaces in subquadratic time when given blackbox access to a fast algorithm for batch approximate nearest neighbor search. In particular, we show the following results for any metric space $\mathsf{M} = ([n], \mathsf{d})$ with aspect ratio $Δ$ admitting a $c$-approximate batch nearest neighbor search algorithm with runtime $τ_{\mathsf{M}}(n)$, (1) There exists an algorithm that, for any $k \in \mathbb{N}$, constructs an $O(c k)$-distortion spanner with $\tilde{O}(kn^{1+1/2k} \log Δ)$ edges and runs in time $\tilde{O}(τ_{\mathsf{M}} \cdot k n^{1/k} \log Δ)$. (2) Any algorithm that learns at most $o(n^{1+1/k}/k)$ pairwise distances by querying a distance oracle and a blackbox batch nearest neighbor search oracle necessarily incurs $Ω(c k)$ distortion. Our results entail that (truly) sub-quadratic time algorithms for spanner construction is equivalent to subquadratic time BANN (up to constant-factor losses). As a further application, we use our fast spanner constructions to obtain a fast algorithm for approximating the Wasserstein distance $\mathsf{W}_q$, for all $q > 1$, over any metric space admitting an efficient batch approximate nearest neighbor search algorithm. Together with recent new efficient algorithms for approximate nearest neighbor search in $\ell_p$ spaces, for $p > 2$, our results entail the first subquadratic time algorithms for spanner construction (with the stated size-distortion tradeoff) and $\mathsf{W}_q$ distance approximation over these metric spaces.Authors: Tyler Chen, Diana Halikias, Christopher Musco, David Persson
The Girard--Hutchinson estimator provides an extremely simple randomized estimate of the Frobenius norm of a matrix $A$ that can only be accessed implicitly via matrix-vector products. In particular, if $Ω$ is a random Gaussian matrix with $r = O(1/\varepsilon^2)$ columns, than $\frac{1}{r}\|AΩ\|_F^2$ provides a $(1\pm \varepsilon)$ multiplicative approximation to $\|A\|_F^2$ with high probability. In this work, we introduce a closely related estimator, given by \begin{align*} {\frac{1}{r}\|AΩ\|_F^2 + \frac{1}{r}\|Ψ^T A\|_F^2 - \frac{1}{r^2}\|Ψ^T AΩ\|_F^2}, \end{align*} where $Ψ$ is a second, independent random Gaussian matrix with $r$ columns. We prove that this estimator yields a $(1\pm\varepsilon)$ multiplicative approximation to $\|A\|_F^2$ when $r = O(1/\varepsilon)$, a quadratic improvement over Girard--Hutchinson. This dependence on $\varepsilon$ is optimal. Our method, which we call Hutch# (pronounced ``Hutch sharp''), matches the complexity of the Hutch++ algorithm [Meyer, Musco, Musco, Woodruff, 2021]. However, unlike Hutch++, Hutch# uses only \textit{non-adaptive} matrix-vector products with $A$ and $A^T$ and requires no orthogonalization or other adaptive linear algebra steps. Thus, Hutch# combines the simplicity of the Girard--Hutchinson estimator and the optimal query complexity of Hutch++.