Log-concave Sampling from a Convex Body with a Barrier: a Robust and Unified Dikin Walk
Authors: Yuzhou Gu
Nikki Lijing Kuang
Yian Ma
Zhao Song
Lichen Zhang
Venue: NeurIPS 2024 poster (2024)
Keywords: Log-concave sampling
Dikin walk
Primary area: probabilistic_methods
TL;DR: We design a Dikin walk for log-concave sampling over polytopes and spectrahedra with the fastest mixing time and efficient per iteration cost.

Abstract: We consider the problem of sampling from a $d$-dimensional log-concave distribution $\pi(\theta) \propto \exp(-f(\theta))$ for $L$-Lipschitz $f$, constrained to a convex body (described by $n$ hyperplanes) equipped with a barrier function, contained in a ball of radius $R$ with a $w$-warm start. 

We propose a \emph{robust} sampling framework that computes spectral approximations to the Hessian of the barrier functions in each iteration. We prove that for the polytope constraints, sampling with the Lee-Sidford barrier function mixes within $\widetilde O((d^2+dL^2R^2)\log(w/\delta))$ steps with a per step cost of $\widetilde O(nd^{\omega-1})$, where $\omega\approx 2.37$ is the fast matrix multiplication exponent. Compared to the prior work of Mangoubi and Vishnoi, our approach gives faster mixing time as we are able to design a generalized soft-threshold Dikin walk beyond log-barrier.

We further extend our result to show how to sample from a $d$-dimensional spectrahedron, the constrained set of a semidefinite program, specified by the set $\{x\in \mathbb{R}^d: \sum_{i=1}^d x_i A_i \succeq C \}$ where $A_1,\ldots,A_d, C$ are $n\times n$ real symmetric matrices. We design a walk that mixes in $\widetilde O((nd+dL^2R^2)\log(w/\delta))$ steps with a per iteration cost of $\widetilde O(n^\omega+n^2d^{3\omega-5})$. We improve the mixing time bound of prior best Dikin walk due to Narayanan and Rakhlin that mixes in $\widetilde O((n^2d^3+n^2dL^2R^2)\log(w/\delta))$ steps.

Log-concave Sampling from a Convex Body with a Barrier: a Robust and Unified Dikin Walk
Yuzhou Gu New York University yuzhouguunyu. edu
Nikki Lijing Kuang
University of California, San Diego 1 11kuan@uussd. edu
Yi-An Ma
University of California, San Diego yianma@ucsd, edu
Zhao Song
Lichen Zhang MIT CSAIL
Simons Institute for the Theory of Computing, UC Berkeley magic. Limuxkde@gmail.com com
1 Lichenzocsail mim edu
Abstract
We consider the problem of sampling from a d-dimensional log-concave distri- bution x(f) x exp(-f(0)) for L-Lipschitz f. constrained to convex body with efficiently computable self-concordant barrier function, contained in ball of radius R with a w-warm start.
We propose robust sampling framework that computes spectral approximations to the Hessian of the barrier functions each iteration. We prove that for polytopes that are described by hyperplanes, sampling with the Lee-Sidford barrier function mixes within (( dT.2 R2) log(w/ō) steps with a per step cost of Ond--1),, where w Rs 2.37 the fast matrix multiplication exponent. Compared to the prior work of Mangoubi and Vishnoi, approach gives faster mixing time able to design generalized soft-threshold Dikin walk beyond log-barrier.
We further extend our result to show how to sample from a d-dimensional spec- trahedron, the constrained set of a semidefinite program, specified by the set {o Rdd Di= x;A; z с} where AI....,Au: C are n X n real symmet- ric matrices. We design a walk that mixes in O(Ad + d..22 log(w steps with per iteration cost of O(n (w- We improve the mixing time bound of prior best Dikin walk due to Narayanan and Rakhlin that mixes in Ö(n228 + n² dL R2) log(w/)) steps.
1 Introduction
Given convex body, generate samples from the body according to structured densitics a fundamen- tal problem in computer science and machine learning. It has extensive applications in constrained convex optimization [Lovász and Vempala, 2006, Narayanan, 2016), differentially private learn- ing [Wang et al., 2015, Lin et al., 2024) and online optimization [Narayanan and Rakhlin, 2017]. A central theme in the theory of sampling is to leverage stochasticity and reduce per iteration costs without having proportionally increase the number of iterations. This theme has played out in continuous optimization for both first and second order methods. For the first order methods, focus is and
on reducing the variance of the gradient estimators [Johnson and Zhang, 2013, Zhang, 2013, Defazio et al., 2014]. For the second order methods, matrix sketching is often used to
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

reduce computation and storage related to the Hessian matrix [Lee et al., 2019, Jiang et al., 2020b, 2021, Song and Yu, 2021, Qin al., 2023].
In Markov chain Monte Carlo (MCMC), much of the progress is on the first order methods. Non- asymptotic analyses performed for the stochastic gradient Langevin algorithms and their variance reduced extensions [Dubey et al., 2016, Raginsky et al., 2017, Brosse et al., 2018, Chatterji et al., 2018. Zou et al., 2018, Dalalyan and Karagulyan, 2019, Li et al., 2019, Ding and Li, 2021J. For second order methods, there has been a paucity when it comes to applying the sketching techniques.
In this paper, we focus on sampling from a log-coneave distribution constrained convex body IC that is described by n constraints, where second order information is proven essential for capturing the geometry of the convex body and consequently for achieving fast convergence rates [Narayanan, 2016, Narayanan and Rakhlin, 2017, Chen et al., 2018, Laddha et al., 2020, Mangoubi and Vishnoi, 2023, 2024]. In particular, we associate the convex body with self- concordant barrier function [Nesterov and Nemirovskii, 1994, Renegar, 1988, Vaidya, 1989] and utilize the Hessian matrix H (z) of the barrier function in the sampling algorithm. We then use the Hessian matrix to propose samples and compute the probability to accept or reject the proposed samples. We note that this problem involves subtleties beyond the scope of constrained optimization. In continuous optimization, which focuses on finding the descent directions, unbiased estimators with reasonable variance oftentimes suffice [Vaidya, 1989, Lee et al., 2015, Huang et al., 2022]. In MCMC, on the other hand, the target probability distribution needs to be maintained along the entire trajectory. This poses significant challenges to speed up the sampling algorithms. In the scenario of uniform sampling over polytope, Laddha et al. [2020] shows that for simple logarithmic barrier', an unbiased estimator in fact suffices. However, more complicated barriers such as the Lee-Sidford barrier [Lee and Sidford, 2014, 2019] can only be approximately computed an unbiased estimator is extremely hard to be obtained. Moreover, for sampling from more sophisticated log-concave distributions and convex bodies, more general approach is needed.
d-dimensional
We therefore propose to obtain high precision estimator the desired acceptance rate with an improved per step running time and without sacrificing the rapid mixing time. In particular, we ask the following question:
Can we significantly reduce per iteration cost, while preserving the convergence rale of the log-concan sampling algorithms over various convex bodies?
We answer the above question in the affirmative. To this end, we present slew of results regarding log-concave sampling. For polytopes, we provide a walk that mixes in Ö + dLL² R²) steps with per iteration cost 0(nd"~1). )) Prior statc-of-the thett-rrrttttt results are due to Mangoubi and Vishnoi [2023, 2024], for which their walks mix in Ö(nd + d R²) with per iteration cost* ((nz((A) + d²). Our walk mixes faster whenever n 2 Our result partially resolves an open problem Mangoubi and Vishnoi [2023] as they asked whether it's possible to design a Dikin walk whose mixing time is only depends on and independent of I., R. We remove the and R dependence on the dominating term d², hence for the case of sampling from uniform distribution (T. = 0), we recover the state-of-art result of Laddha al. [2020] which mixes in (dd)) steps with per iteration cost O(dd"1)).
We also consider sampling from convex bodies beyond polytopes. In semidefinite programming (SDP), one often focuses on the dual program of an SDP, where the constraint set is defined as spectrahedron K- fa E Rd Σd A; 2 C)S for n symmetric matrices A, and C' [Jiang et al., 2020a, Huang et al., 2022]. We propose a walk that samples from a log-concave distribution over a spectrahedron IC' using the Hessian information of the barrier.... The walk mixes in ((d + dL R²) steps. For log- -arrierrier of an SDP, explicitly forming the Hessian matrix takes a prohibitively large
'Logarithmic barriers. or log- barrierss have been extensively studied in the context of sampling and convex op- timization. For dddimensiioal polytope with constraints, it typically leads to convergence rate polynomially depends on n.
LLc--Sifoodd barriers are the first polynomial-time computable barriers with complexity parameter depends only logarithmic on n and in extension, the convergence rate.
We use (-) to hide polylogarithmic dependence on n, d and 5 where 8 is the TV distance between the target distribution and our Markov chain Jr.
A Wee use nnz(A) to denote the number of nonzero entries matrix A.
We use A 2 0 to denote A is positive semidetinite matrix, and A > B to denote A- B > 0.

---

O( d (nddddd-)))) time. We utilize our robust sampling framework to approximately compute this Hessian via tensor-based sketching techniques, achieving runtime of (n n dddw-o). As long as n
(which is usual setting for SDP, where d < n), this provides a speedup.
1.1 Related works
In this work, we focus on Dikin walk [Kannan and Narayanan, 2012], refined variant of ball walk [Lovász and Simonovits, 1993]. Roughly speaking, ball walk progresses by moving to a random point in the ball centered at the current point with the obvious downside that when the convex body is flat, ball walk progresses slowly. Dikin walk overcomes this problem by instead moving to a random point in good ellipse centered at the current point, in an effort to better utilize the local geometry of the convex body.
For sampling over polytopes, a number of Dikin walks use the ellipse induced by the oog-barrier functions [Kannan and Narayanan, 2012]. The work of Laddha et al. [2020] shows that uniform sampling over polytope can mix in ((vd) steps, where v is the self self-concordant parameter of the barrier function Going beyond uniform distributions, Narayanan and Rakhlin [2017] proposes a walk that samples from a log-concave distribution on a polytope in ((²a)) steps. This bound is later improved by Mangoubi and Vishnoi [2023], in which they show that for logarithmic barrier, mixing time of (nd dI,2 R22) ) is attainable, where L is the Lipschitz constant of J and R is the radius of the ball that contains K., A more recent work by Kook and Vempala [2024] has obtained better bound when the function is a-relative strongly-convex and the density π is 3-Lischhitz They manage to obtain a mixing bound of (vv33 log(w /6)). However, their algorithm requires stronger assumptions on f and hence is incomparable to our result. In most previous works, the focus has been on improving the mixing time, rather than on per iteration costs.
For any convex body, it is well-known that a universal barrier with = d exists [Nesterov and Nemirovskii, 1989, Lee and Yue, 2021]. however it is computationally hard to construct the Hessian of the universal barrier. In short, the universal barrier requires to compute the volume of the polar set of a high dimensional body, which is even hard to approximate deterninistically [Furedi and Barany, 1986. Nesterov and Nemirovskii, 1994]. The seminal work of Lee and Sidford [Lee and Sidford, 2014] presents nearly-universal barrier with v = (dd log" r2) and the Hessian can be approximated in Ond==1) time. The Lee-Sidford barrier was originally designed to solve linear programs in O(V)) iterations, and it has been leveraged lately for Dikin walks with rapid mixing time. For uniform sampling, Chen et al. [2018] gives an analysis with a walk that mixes in (dd255 steps. Subsequently, Laddha et al. [2020] improves the mixing to (dd) steps. In pursuit to better leverage the local geometry, Lee and Vempala [2017] proposes a walk relying on Riemannian metric that mixes in (ndd/4)) steps. The mixing rate is later improved to (nn22/3 via Riemannian Hamiltonian Monte Carlo (RHMC) with log-barrier [Lee and Vempala, 2018] 1018 and 0nn /1/3 with mixed Lee-Sidford and log-barrier [Gatmiry et al., 2024]. In this work, we show that log-coneave sampling can also leverage Lee-Sidford barrier to obtain a mixing time of Ö(d² + R²). In comparison, the hit-and-run algorithm [Lovász and Vempala, 2007] mixes in Odd R /r2) steps where is the radius of the inscribed ball inside IC. For the regime where I. = 0(1) and r = 0(1), our walk mixes strictly faster than that of hit-and-run.
We note that for sampling over more general convex bodies, a recent work [Chen and Eldan, 2022] proves that for uniform sampling over an isotropic convex body, the mixing time bound Ö /v3) where tha is the KLS constant [Kannan et al., 1995]. However, it is unclear how to generalize their result to non-isotropic convex bodies or log-concave densities.
1.2 Our results
Our results concern log-concave sampling from polytopes and spectrahedra. For polytopes, we state the result in its full generality.
Technically, what they have shown is that if the Hessian is p-symmetric self- concorann barrier function, then the walk mixes in (udd) steps, and they subsequently prove for various barriers of interest,
3

---

Theorem 1.1 (Robust sampling for log- concave distribution over polytopes), Let ó E (0, and R 2 1. Given a constraint matrix A E Read with a vector b E TR", let k: {x E TR Ax S bj be the corresponding polytope. Suppose K. is enclosed ball of radius R with non-empty interior. Let f: K. - IR be an L-Lipschitz, convex function with an evaluation oracle. Finally, let π be the distribution such that π x x.
Suppose we are given an initial point from K. that is sampled from a w-warm start distribution" with respect to π for some
O. then there is an algorithm (Algorithm 1) that outputs point from distribution iI where TV S
Let g - lR be u-self- u-selffconcorlont barrier function such that in time Cg. spectral approximation Hg the Hessian of g denoted by Hg can be computed satisfying
(1-c) Hg 3 H (1+
for = (L/dd) Then, Algorithm 1 takes at most
Ō(vvd dL R²) log(w/6))
Markov chain steps. It uses function evaluations and an extra C2 + d time at each step.
Let us pause and make some remarks regarding Theorem 1.1. As long as the Hessian matrix can be approximately generated with an error inversely depends on d, then our algorithm is guaranteed to converge. Moreover, if the approximate Hessian can be generated more efficiently, then this directly implies an improvement of our algorithm. On the mixing time side, Theorem 1.1 is nearly-optimal up to polylogarithmic factors and the dependence on L and R. In particular, the vd mixing time bound is also achieved by Laddha et al. [2020] in the case of uniform sampling. We instantiate the meta result of Theorem 1.1 into the following theorem.
Theorem 1.2 (Robust sampling with nearly-univense barrier). Under the conditions of Theorem 1.1, let g be the Lee- Lee-Sidford barrier with v Od n). Then, we have
C Ö(nd
The algorithum takes at most
Ō((d + dLL log"n))oog(((0))
Markov chain steps.
The Lee-Sidford barrier [Lee and Sidford, 2014, 2019] is the first polynomial-time computable barrier with a nearly-optimal self-concordance parameter. Several prior works [Chen et al., 2018, Laddha et al., 2020] utilize this barrier for uniform sampling. For log-concave sampling, the walk of Kook and Vempala [2024] requires strong-conveeitt-lik assumption on f in order to attain an (dd) mixing time. Our work is the first to obtain an (dd) mixing time for log-concave sampling over polytopes, when L and R are small.
References
Mixing time
Per iteration cost
Lovász and Vempala [2006]
d' /T²
nd ndd=1
Narayanan and Rakhlin [2017]
d' + RR.
Mangoubi and Vishnoi [2023]
nd R.
ndw-I
Mangoubi and Vishnoi [2024]
nd - dR.
nuz(A + d²
This work
d² + dL²R²
nd
Table I: Comparison of algorithms for sampling from a log-concave density e -J over a d-dimensional polytope with n constraints, where ff is L-Lipschitz. We omit Ö(-) to mixing time and per iteration cost. We assume the evaluation of can be done in unit time. The first row is the hit-and-run walk. Our algorithm has the fastest mixing time among all Dikin walks, and for = 0(1), r = O(1), our walk mixes faster than hit-and-run.
We say distribution p is w-warm start of with respect to distribution if sup
pla) S

---

Our next result is closely related to semidefinite programming. Given the the set K = fr E R Ei TiA; = C} for symmetric n. matrices A1.... Ad, C, we consider sampling from log-concave distribution over K. We use the popular log-barrier for K. studied in Nesterov and Nemirovskii [1994], Jiang et al. [2020a], Huang et al. [2022]. To further utilize our robust sampling framework, we develop a novel approach for generating a spectral approximation of the Hessian. Compared to Jiang et al. [2020a]. Huang al. [2022] in which the Hessian is maintained in very careful way in conjunction with the interior point method, our spectral approximation does not require control the changes over iterations and is highly efficient when n >> d. which is popular regime for semidefinite programs,
Theorem 1.3 (Robust sampling for log-concave distribution over spectrahedra). Let å E (0, 1) and R 2 Given a collection of symmetric matrices A1...., Ad E R and a target symmetric matrix C' Raxn, let == {x E IR4 "i=1", A, - C} be the corresponding spectrahedron. Suppose K. is enclosed in ball of radius R with non-empty interior. Let f : - IR be L-Lipschitz, convex function with an evaluation oracle. Finally, let a be the distribution such that π x e
Suppose we are given an initial point from K. that is sampled from a w-warm start distribution with respect to x for some w >> 0, then there an algorithm (Algorithm 1) that outputs a point from distribution y where TV,,, T) S 8.
Let g KC R be the logarithmic barrier function with n. There an algorithm that uses
Ö((nd + dLL. R²) log(w/8))
Markov chain steps. In each step uses function evaluations and an extra cost of
Ön + +nddw-55
Prior work due to Narayanan and Rakhlin [2017] provides walk that mixes in Ö(n + n= steps and each step could be implemented in O'dn" n22d=1) time, our algorithm improves upon both mixing and per iteration cost. On the front of per iteration cost, computing the Hessian of the log-barrier of K exactly takes time O(dn + n d-ddd1), 1), for the regime where n > the dominating term is dn Note that whenever >> d d
algorithm is more efficient than that of exact computation. This a common regime in semidefinite program in which the number of constraints small compared to the size of the primal solution.
References
Mixing time d=
Per iteration cost n1 + ndd
Lovász and Vempala [2006]
Narayanan and Rakhlin [2017]
n'd d'dL n² d d n² n² dow-
This work
nd RL
Table 2: Comparison of algorithms for sampling from a log-concave density eJ over a d-dimensional spectrahedron with mixing time and per iteration cost. We assume evaluation of f can be done in unit time. The first row is the hit-and-run walk. Our algorithm mixes faster than Narayanan and Rakhlin [2017] in all parameter regimes, and for L = 0(1), 0(1), da R² 2 n, our walk mixes faster than hit-and-run.
symmetric matrix constraints, where f is T.-Lipschitz. We omit Ō(-) to
2 Technical Overview
Before diving into our techniques, first illustrate an informal version of our algorithm. The algorithm itself is a variant of the Dikin walk [Kannan and Narayanan, 2012] called the soft-threshold Dikin walk [Mangoubi and Vishnoi, 2023]. Starting with w-warm start point, the algorithm approximately computes the Hessian of the barrier function, then proposes a new point z from Gaussian with mean at the current point and variance ((x)-1 where T is a multiple of the approximate Hessian plus a multiple of the identity. The proposal is then accepted via a Metropolis filter. Note that compared to the standard Dikin walk framework, our algorithm crucially allows the Hessian to be approximated with a good precision. This means any improvements on generating spectral approximation leads to runtime improvement of our algorithm in a black-box fashion.
5

---

Algorithm Sampler for log-concave distribution exp(-I) over convex body with a barrier function g. is the description of the convex set, function / an I.-Lipschitz concave function and g is -seeo-connnrrantt barrier function, 2o E K is w-warm start point, and ó is the total variance distance between our Markov chain and the target distribution. Let S I enclosed by a ball of radius/ R.
procedure SAMPLINGCOWVEEBBOYYAA..g....
2: Initialize x to be point in K
3: T + ( R²² log(w/6))
4:
for t 1 - do
5:
Sample a point (N, Ia)
6:
(z) is the exact Hessian matrix at x, we never explicitly form
7:
Compute (x) where (1 - 1/d)
(1 (+)
8:
((x) d. (x) + dT2... Ia
9:
- x -0(x))1122
10:
if E Int(K)) then
11:
11 (z) is the exact Hessian matrix
12:
Compute (2) such that (1 - (/d)H(() = (z) (14)))))22)
13:
(2) = d (2) 2 Ia
14:
exp(
15:
Our walk is lazy, i.e., it only moves to new point with half probability accept + with probability 0.5 m0in(,
16:
17:
else
18:
reject
19:
end
20:
end for
21:
end procedure
2.1
Soft-threshold Dikin walk for log-concave sampling via barrier-based argument
For uniform sampling from structured convex bodies such as polytopes, Kannan and Narayanan [2012] introduce refined ball walk that utilizes self-concorlarl barrier function associated with the convex body. To extend the Dikin walk framework for log-concave sampling, Narayanan and Rakhlin [2017] adds an additional coefficient exp(-f(x))
exp(-/(#]] to the acceptance probability of Metropolis filter. This achieves a sub-optimal mixing time bound of Ö(v2)). By introducing soft-threshold regularizer to the Hessian, Mangoubi and Vishnoi [2023] proves that for log log-barrier function, is possible to approach the optimality with an Ö(nd + dL ((2) mixing time. Unfortunately, the method used by Mangoubi and Vishnoi [2023] is specialized to log-barrier, making it hard to generalize to other barriers. They show that under proper scaling, the true Hessian H (x) under a sequence of polytopes converges to the regularized Hessian under the given polytope in the limit. They define the j-th polytope as
A
b
La
j
la
A
by
IJ
1a
the number of copies of Ia and j 1, is function of j. While this polytope- based method works well for log-barrier, it functions poorly for more intricate barriers, such as volumetric barrier and Lee-Sidford barrier. Take volumetric barrier as example, recall that H (x) A (1) A))(2) A where 2(S(x))-1 -1 A) is the statistical leverage score [Drineas et al., 2006, Spielman and Srivastava, 2011]. Leverage score is a numerical quantity that measures how important row of matrix is, compared to other rows. This means that, if one duplicates row infinitely many times, the row will be assigned score of 0. More concretely, suppose we are given an n x matrix with n identical rows, then the leverage score of each row 4. Taking n -> oo, it's easy to see that all rows have score 0, and the Hessian of volumetric barrier will zero out the rows contributed to the

---

infinitely many copies of identities. Thus, for volumetric barrier, extra constraints of the polytope sequence vamish. On the other hand, the Mangoubi and Vishnoi [2023] construction heavily relies on the infinitely many occurrences of La for it to converge to the regularized Hessian.
To circumvent this issue and provide simpler argument, we develop a proof strategy that is barrier- We example a class of popular barrier functions for polytopes that can be classified as weighted log-barrier functions [Vaidya, 1989]. Roughly speaking, these functions take in the form of g(z) = - - wi (og(S(r))) with w = 1, being the log-barrier, being the volumetric barrier [Vaidya, 1989] and if we choose w to be a proper power of the Lewis weights, then the barrier is precisely the Lee-Sidford barrier [Lee and Sidford, 2019, I Laddha et al., 2020]. These weights are stable, meaning that if two points are close in certain local norm induced by proper ellipses given by these weights, then their weights are also close [Lee and Sidford, 2019]. This property proves to be particularly useful when proving the mixing rate of our Dikin walk. In addition, we carve out three sufficient conditions for a weighed log-barrier barrier function to work for log-concave distribution, and has (vv)) mixing rate:
based, instead of polytope-basee
w = (S(x)))
v-symmetry. The Hessian (xx) is v-symmetry, that is, for any 2 E E. (,(1)
(2x KC) S (VD) where (r) = {y E Ted (y - 2) Hx)((y x) S r²) We show that D = V for log-barrier, volumetric barrier and Lee-Sidford barrier for polytopes, and log-barrier for spectrahedra.
Bounded local norm.
The variance term [z V log (et((())))] S doly log n.
We show this condition holds for all barriers of interest in this paper.
. Convexity of regularized barrier. The function F(x) log (ee((H())) + Ia) is convex at x for any x E K. We note that for log-barrier, volumetric barrier and Lee-Sidford barrier, has been shown that log (x))) is indeed convex at x. We further prove that this still holds when a copy of identity is blended in.
We want to highlight that our approach is more robust and generic than that of Mangoubi and Vishnoi [2023], as it solely depends on the structure of the Hessian matrix. Our argument should be treated as a generalization of Sachdeva and Vishnoi [2016], in which they utilize the bounded local norm and convexity of log det H (x)) for log-barrier, together with the Gaussianity to conclude that the difference between log det are small. We go beyond log-barrier for polytopes and uniform sampling.
2.2 Approximation preserves the mixing time
The core premise of our algorithm is we allow the Hessian to be approximated, so that each step can be implemented efficiently. This poses significant challenges proving mixing, as even we have a good approximation of one-step in the Markov chain, the approximate chain does not necessarily mixes as fast as the original chain. For example, if every step is approximated within e-TV distance to the original Markov chain, then in T steps, the TV distance between the resulting distribution under approximate walk and the distribution under original walk could be as large as Te. This means we have take e very small (O(Tm))) to guarantee same convergence property as original chain, which is unacceptable. A related issue is that the stationary distribution under approximate walk may not be the same as the target distribution. Therefore we need way to control the mixing properties of the approximate walk. We resolve this problem by exactly computing the acceptance probability under approximate walk. That is, after proposing the next step from x, ((z).
we sample (=), and accept with probability min x())G_(x)
(x) (where (z) is the probability of proposing 2 starting from ø under the approximate walk). Recall that in the exact walk, we accept xx))G_(x) } (where (z) is the probability of proposing 2 stating from x with probability min 11, (=))..
under the exact walk). If we use this acceptability for the approximate walk, then we only have one-step approximation guarantee, and will suffer from the problem mentioned previously.
By using this modified acceptance probability, we can prove reversibility of the approximate walk, and that stationary distribution of the approximate walk is indeed our target distríbution. For the mixing rate, we bound the conductance of the approximate walk, which requires us to prove the following properties:
(1) The proposed step is accepted with decent probability.
7

---

(2) Starting from two points moderately close to each other (i.e., 1|- - 2)(0(a) I c). the approximate steps taken Pz and Ps are close in TV distance (i.e.. TV(P,, Ps) S 1 - e') for some absolute constants e >> 0.
We are able to prove these properties by comparison with the original chain. That is, assuming the original chain satisfies properties (1)(2), then the approximate chain also satisfies the same properties, as long the approximation is good enough. This enables prove these key properties in hassle-free way. To prove these comparison results, we need prove, for example, Gz (z).). Because both Gz(z) (z) and (=) have nice factorizations
1
(2) x dct( exp(
~11 - = 212(,))
(x) x det( (o())) 12 exp(
2 || = (
we only need to prove det der()(x)) = det 4 (z) and exp(- 2لمم - -=((aa)))
separately. Both of these can be handled via (ar) = ((x). using our approximation procedure.
2.3 Approximation of Lewis weights for rapid mixing
Given our robust and generic Dikin walk framework, we realize with the Lee-Sidford barrier whose complexity v = d'log' n. This ensures the walk mixes in Ō + d2.2 steps, improving upon the log-barrier-based walk of Mangoubi and Vishnoi [2023, 2024] that mixes in O(nd + dL" steps, as the log-barrier has complexity = n. To implement the Lee-Sidford barrier, it is imperative to compute the Ep (for p > 0) Lewis weights, which is defined as the following convex program:
n
min 0
log det M, subject to Σ (a Ma.)p/2 <d.
and the Ep Lewis weights is vector up Rdd where (w,,): = a M,a, and M, is the optimum of Program (1). Solving the program exactly is not efficient, and long line of works [Cohen and Peng, 2015, Lee and Sidford, 2014, 2019, Cohen et al., 2019, Fazel al., 2022, Jambulapati et al., 2022] provide fast algorithms that approximate all weights up to (1 + ef-ftct The Hessian of Lee-Sidford barrier has the form (up to scaling) (x) = A
where W, (,( A) is the diagonal matrix of fp Lewis weights with respect to 5(z)-1A. A. When clear from context, we use W, p as shorthand for (p(S(x))=AA) A). As W, could be approximated in Ö(nd~=1) time, we could then perform subsequent operations exactly to form the approximate Hessian. We instead provide a spectral approximation procedure that runs in (mnz(A) + Tmaa (d, d,, d)) given the approximate Lewis weights where Tmaa b, c) is the time of multiplying an X b matrix with a x e matrix and Tmar(d, d²,d) ~ d=.22 This is important as even though it doesn't 't lead to direct runtime improvement, if approximate Lewis weights can be computed in ((nzz(A) + (d)) time, then we obtain a per iteration cost upgrade a black-box fashion. Our approach is based on spectral sparsification via random sampling: given matrix B E =x one can construct matrix B E exd consists of re-scaled rows of B, such that (1 - e) B B BT B 3 (1 +) BTA B. If we sample according to the leverage scores [Drineas et al., 2006, Spielman and Srivastava, 201 IJ. one
can set O 2d d). In our case, set B However, computing leverage score exactly is as expensive as forming the Hessian exactly. which requires 01n time. To speed up the crucial step, we use randomized sketching technique to reduce the row count from n to poly(d), then approximately estimate all leverage scores with this small matrix. One either to approximate these scores in time O (nnz( (A) log + Taaat e 2d, d)), or use the subsampled randomized Hadamard transform [Lu et al., 2013] in time O(nd log n + Tiatt (d, E 2d, d)). We remark that approximating leverage scores for subsampling has many applications in numerical linear algebra, such as low rank approximation [Boutsidis et al., 2016, Song al., 2017, 2019]. Utilizing this framework, we provide an algorithm that approximates the Hessian of Lee-Sidford barrier in time (() + Ó(nnz(A) + Tmaa (d, d', dd)). ). as advertised.
'A, and perform the sampling.
the sparse embedding matrix [Nelson and Nguyên, 2013]
8

---

2.4 Approximation of log-barrier Hessian over a spectrahedron
For a spectrahedron induced by the dual semidefinite program, is also natural to define its log-barrier with its corresponding Hessian being
Hing (+) = A(S(a)-* & ()a)-
where A E RBdxn with the i-th row being vec(A.) for y. x n symmetric matrix A,, and S(x) = n symmetric matrix. Nesterov and Nemirovskii [1994] shows
Zim1 T; - C' where C' is also that by smartly arranging the organization of (1)) it can be computed exactly in time (dn"" d==In?). n Beyond exact computation, Jiang al. [2020a], Huang et al. [2022] present algorithms that approximately maintain the Hessian matrix under low-rank updates. These maintenance approaches are crafted towards solving semidefinite programs in which the trajectory can he carefully controlled. For our thft-threshold Dikin walk, though the proposal generates point that is relatively close to the starting point, due to its Gaussianity nature, much less structure can be exploited and thus maintained. Instead, we propose maintenance-frre approach that uses randomized sketching to generate a spectral approximation of Hog (1).
d
To illustrate the algorithm, let's 's' first define a matrix
B = A(S(z)) 1/2 @ S(2)-112)
is not hard to see that Hogg (2) = BB Due to the Kronecker product structure of (5(x)-1/2 @ 5(x)-1/2 it is natural to consider sketches for Kronecker product of matrices [Diao et al., 2018, 2019, Ahle et al., 2020, Song al., 2021]. Following the standard procedure for using sketch to generate spectral approximation, we can choose a TensorSRHT matrix [Ahle et al., 2020, Song et al., 2021] T and compute == (S(x) 11/2 1/2 8 -1/yyT then form AR. This unfortunately, not efficient enough, as multiplying A with R might be too slow. To further optimize the runtime efficiency, we use a more intrinsic approach based on the matrix T and B. It is well-known that the i-th row of B can be first computed as (xx -1/2 A,S(x)) then vectorize. If we can manage to apply the sketch in a row-by-row fashion, then fast matrix multiplication can be utilized for even faster algorithms. We recall that T = P- (H 1 @ D2) where H is the Hadamard matrix and P is a row sampling matrix. To compute a row, we can first apply HD, and D2 individually to S(x) in Ōn?)) time since H is a Hadamard matrix. We can then form two matrices X and Y with rows being the corresponding sampled rows from HD1 and H Finally, we compute XA, to form one row. Using fast matrix multiplication, we can form each row in Tma n, n). = O(n2 time. Applying this procedure to d rows leads to an overall 2(n288w-5) time, which beats the (dn")) time as long as n >> d. Our runtime is somewhat slow when n is not that larger than d, this is due to we have to set = O(1dd). For more popular regimes where it suffices to set E ((), our algorithm outperforms exact computation as long as n 2 d.
3 Applications
There are many applications of our Dikin walk, such as differentially private learning, simulated annealing and regret minimization. For a more comprchensive overview of the applications, we refer readers to Narayanan and Rakhlin [2017].
Differentially private learning. Let X [aI..... Fm. denote an m-point dataset, associate a convex loss function l; to point Ti, the learning problem attempts to find an optimal parameter 0, IC that minimizes {()) = m (0). We could further enforce privacy constraints learning: we say a randomized al algorithm M : X + R is cDP if for any datasets X, X' E 2mm differ by single point and for any S C IR. (x[ S S] < e Pr[M(X') C S]. The differentially private learning seeks to solve the learning problem under DP guarantees, and it has been shown that if one allows for approximate DP, then it can be achieved via the exponential mechansim and sampling from the distribution (0) x e ((()) [Wang et al., 2015]. If K. is a polytope or spectrahedron, and in addition l is I.-Lipschitz, one can implement our walk to obtain an Ō + dI 2 R2) mixing for polytopes or Ö(nd + dlL,2 R²) mixing for spectrahedra. Via a standard reduction from TV distance bound to infinity distance bound [Mangoubi and Vishnoi, 2022], our walk can also be adapted for a pure DP guarantee [Lin et al., 2024].

---

Simulated annealing for convex optimization. Let f: R + R be an L-Lipschitz and convex function, we consider the problem of minimizing f over K. when we can only access the function value of f. This is also referred to as zeroth-order convex optimization. One could solve the problem via the simulated annealing framework, in which one needs to sample from the distribution (a) x e [a))TT where T is the temperature parameter [Kalai and Vempala, 2006]. Our walk can be adapted to the framework for polytopes, it mixes in Ö + d' L2R)) steps and for spectrahedra, it mixes in O 15 + d R²) steps, improving upon prior best algorithms that mix in minfc + d3.5 .dl + d L L² steps for polytopes and d2 +
R2R² steps
for spectrahedra [Narayanan and Rakhlin, 2017, Mangoubi and Vishnoi, 2023].
Online convex optimization. Consider the following online optimization problem: let K. be convex set that we could choose our actions from and let f,...., fy be a sequence of unknown convex cost functions with & K + R. The goal is to design good strategy that chooses from IC so that the total cost is small, compared to the offlin optimal cost, in which one could choose the strategy after seeing the entire sequence of 4''s. The algorithm works as follows: round
choose a distribution 1t- supported on IC and play the action Y, ~ pt-1- Our goal is to minimize the expected regret, defined as Regy (U) = E 4(())-- E' C(U)) with respect to all randomized strategies defined by distribution pt-- It turns out if one sets (x))) 1) E'_ (x(x) where 1/ > is step size parameter and sets Pr e e then this update rule reflects the multiplicative weights update algorithm [Arora al., 2012]. Moreover, one could prove nearly-optimal regret bound of O(VT) the KL divergence between (40 and all possible pe 's are bounded. The algorithmic problem of sampling from f4 is equivalent to the log-concave sampling problem study in this paper, and we could efficiently generate strategies with good regrest, similar to Narayanan and Rakhlin [2017]. Let 4's be I-Lipschitz linear functions, Narayanan and Rakhlin [2017] uses their varying Dikin walk to achieve regret of O VT), while we could significantly improve this bound to (dvTT).
4 Conclusion
In this paper, we design a class of error-robust Dikin walks for sampling from a log-concave and log-Lipschitz density over convex body. The key features of our walks are that their mixing time depends linearly on the complexity of self-concordant barrier function of the convex body, and they allow computationally intensive quantities to be approximated rather than computed exactly. For polytopes, our walk mixes in O + dLL2R3 R²) steps with a per iteration cost of Ö(nd" 1). and for spectrahedra, our walk mixes in Ö(nd + 'dL' R²) steps with a per iteration cost of Ö(n + ndw For polytopes, our walk is the first successful adaptation of the -Sdffod barrier for log-coneave sampling under the minimal assumption that f is Lipschitz, improving upon the mixing of prior works [Narayaman and Rakhlin, 2017, Mangoubi and Vishnoi, 2023, 2024]. For spectrahedra, improve the mixing of Narayanan and Rakhlin [2017] from nzdd3 + n²d R22 to nd + d nē, for the term that depends on only and d, we obtain improvement of ndf and for the term that depends on T. and R, we shave off the quadratic dependence on n. Moreover, we adapt our error-robust framework and present a sketching algorithm for approximating the Hessian matrix in Ö(n" + 2 5) time. Our results have deep implications, as it could be leveraged for differentially private learning, convex optimization and regret minimization.
While our framework offers the fastest mixing Dikin walk for log-concave sampling, its mixing time has a rather unsatisfactory dependence on the radius of the bounding box, R. It would be interesting if one could design walk that does not depend on R. We also note that we require the density of be log-Lipschitz rather than Lipschitz. Kook and Vempala [2024] shows that if f in addition satisfies the relative strongly-convex property, then there exists a walk whose mixing does not depend on R and one only needs the density of be Lipschitz. Another important direction is to further improve the mixing rate of sampling over a spectrahedron using barrier functions such as volumetric barrier and hybrid barrier [Anstreicher, 2000], while maintaining an small per iteration cost. Finally, extending RHMC to log-concave sampling would be essentially, as has the potential to give the fastest mixing walk by utilizing Riemannian metrics instead of Euclidean. As our work is theoretical in nature, we don't foresee any potential negative societal impact. It has the potential positive impact to reduce energy consumption and carbon emission when deployed in practice.
10

---

Acknowledgement
We would like to thank Jonathan Kelner, Yin Tat Lee, Santosh Vempala and Yunbum Kook for valuable discussions and anonymous reviewers for helpful comments. The research is partially supported by the NSF awards: SCALE MoDL-2134209, CCF-2112665 (TILOS), is also supported by the U.S. Department of Energy, the Office of Science, the Facebook Research Award, well as DOC-RRA--T--23 20-2009--999-0 from the CDC's Center for Forecasting and Outbreak Analytics. Lichen Zhang supported by NSF awards 1C5-2195217 and DMS-202242 For more information related to the paper and adjacent topics, see https / //www youtube... com/ /0zhaosong2031 and https: //space bilibili. com/ /55455877766566611
References
Thomas D Ahle, Michael Kapralov, Jakob BT Knudsen, Rasmus Pagh, Ameya Velingker, David Woodruff, and Amir Zandieh. Oblivious sketching of high-degree polynomial kernels. In Pro- ceedings of the Fourteenth Anmual ACM-SIAM Symposium on Discrete Algorithms (SODA), pages 141-160, 2020.
Josh Alman and Virginia Vassilevska Williams, refined laser method and faster matrix multiplica- tion. In Proceedings of the 2021 ACM-SIAM Symposium Discrete Algorithms (SODA), pages 522-539. SIAM, 2021.
Kurt M Anstreicher. The volumetric barrier for semidefinite programming. Mathematics of Operations Research, 2563):66--380 2000.
Sanjeev Arora, Elad Hazan, and Satyen Kale. The multiplicative weights update method: meta- algorithm and applications. Theory of Computing, 8(1):121-664 2012.
Christos Boutsidis, David P Woodruff, and Peilin Zhong. Optimal principal component analysis in distributed and streaming models. In Proceedings of the forty-eighth annual ACM symposium on Theory of Computing (STOC), pages 236-249, 2016.
Nicolas Brosse, Alain Durmus, and Eric Moulines. The promises and pitfalls of stochastic gradient langevin dynamics. Advances in Neural Information Processing Systems, 31, 2018.
Niladri Chatterji, Nicolas Flammarion, Yian Ma, Peter Bartlett, and Michael Jordan. On the theory of variance reduction for stochastic gradient monte carlo. In International Conference on Machine Learning, pages 764-773. PMLR, 2018.
Yuansi Chen and Ronen Eldan. Hit-and-run mixing via localization schemes, 2022.
Yuansi Chen, Raaz Dwivedi, Martin Wainwright, and Bin Yu. Fast meme sampling algorithms on polytopes. J. Mach. Learn. Res., 2018.
Michael B. Cohen and Richard Peng. Lp row sampling by lewis weights. In Proceedings of the Forty-Seventh Annual ACM Symposium on Theory of Computing. STOC '15 2015.
Michael B Cohen, Ben Cousins, Yin Tat Lee, and Xin Yang. A near-optimal algorithm for ap- proximating the john ellipsoid. In Conference on Learning Theory, pages 849-873. PMLR. 2019.
Amnak Dalalyan and Avetik Karagulyan. User-friendly guarantees for the langevin monte carlo with inaccurate gradient. Stochastic Processes and their Applications, 129(12):5278-53 10
Aaron Defazio, Francis Bach, and Simon Lacoste-Julien. Saga: A fast incremental gradient method with support for non-strongly convex composite objectives. Advances in neural information processing systems, 27, 2014.
Huaian Diao, Zhao Song, Wen Sun, and David Woodruff. Sketching for kronecker product regression and p-splines. International Conference on Artificial Intelligence and Statistics, pages 1299- 1308. PMLR, 2018
11

---

Huaian Diao, Rajesh Jayaram. Zhao Song, Wen Sun, and David Woodruff. Optimal sketching for kronecker product regression and low rank approximation. Advances in neural information processing systems, 32:4737-4748, 2019.
Zhiyan Ding and Qin Li. Langevin monte carlo: random coordinate descent and variance reduction. J. Mach. Learn. Res., 22:205-1,
Petros Drineas, Michael W. Mahoney, and S. Muthukrishnan. Sampling algorithms for 12 regres- sion and applications. In Pruceedings of the Seventeenth Annual ACM-SIAM Symposium on Discrete Algorithm, SODA '06, page 1127-1136, USA, 2006. Society for Industrial and Applied Mathematics.
Ran Duan, Hongxun Wu, and Renfei Zhou. Faster matrix multiplication via asymmetric hashing, In FOCS, 2023.
Kumar Avinava Dubey, Sashank Reddi, Sinead A Williamson, Barnabas Poczos, Alexander Smola, and Eric P Xing. Variance reduction in stochastic gradient langevin dynamics. Advances in neural information processing systems, 29, 2016.
Maryam Fazel, Yin Tat Lee, Swati Padmanabhan, and Aaron Sidford. Computing lewis weights to high precision. In Proceedings of the 2022 Annual ACM-SIAM Symposium on Discrete Algorithmss (SODA), pages 2723-2742, SIAM, 2022.
Z Furedi and Barany. Computing the volume is difficult. In Proceedings of the Eighteenth Annual ACM Symposium on Theory of Computing, STOC '86, page 442 447, New York, NY, USA, 1986. Association for Computing Machinery.
François Le Gall and Florent Urrutia. Improved rectangular matrix multiplication using powers of the gyppersmith-winngreee tensor. In Proceedings of the Twenty-Nnth Annual ACM-SIAM Symposium on Discrete Algorithms, SODA '18 page 1029-10000 10061046,220 2018.
Francois Le Gall. Faster rectangular matrix multiplication by combination loss analysis. In Pro- ceedings of the Thirty-Fifth Annual ACM-SIAM Symposium on Discrete Algorithms, 20 2024.
Khashayar Gatmiry, Jonathan Kelner, and Santosh S. Vempala. Sampling polytopes with riemannian hme: Faster mixing via the lewis weights barrier. In Shipra Agrawal and Aaron Roth, editors, Proceedings of Thirty Seventh Conference on Learning Theory, volume 247 of Proceedings of Machine Learning Research, pages 1796-1881. 100.-RR,..MMR, 30 Jun-03 Jul 2024.
Roger A. Horn and Charles R. Johnson. Matrix Analysis. Cambridge University Press, 1990.
Baihe Huang, Shunhua Jiang, Zhao Song, Runzhou Tao, and Ruizhe Zhang. Solving sdp faster: robust ipm framework and efficient implementation. In FOCS, 2022.
Arun Jambulapati, Yang P. Liu, and Aaron Sidford. Improved iteration complexities for overcon- strained p-norm regression. In Proceedings of the 54th Annual ACM SIGACT Symposium on Theory of Computing, STOC 2022, 2022.
Haotian Jiang, Tarun Kathuria, Yin Tat Lee, Swati Padmanabhan, and Zhao Song. A faster interior point method for semidefinite programming. In 2020 IEEE 6lst annual symposium on foundations of computer science (FOCS), pages 910-918. IEBE, 2020a.
Haotian Jiang, Yin Tat Lee, Zhao Song, and Sam Chiu-wai Wong. An improved cutting plane method for convex optimization, convex-concave games, and its applications. In Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing, pages
2020b.
Shunhua Jiang, Zhao Song, Omri Weinstein, and Hengjie Zhang. A faster algorithm for solving general lps. In Proceedings of the 53rd Annual ACM SIGACT Symposium on Theory of Computing, pages 823 832, 2021.
Rie Johnson and Tong Zhang. Accelerating stochastic gradient descent using predictive variance reduction. Advances in neural information processing systems, 26, 2013.
12

---

William B Johnson and Joram Lindenstrauss. Extensions of lipschitz mappings into a hilbert space. Contemporary mathematics, 26(189-206):1 1984.
Adam Tauman Kalai and Santosh Vempala. Simulated annealing for convex optimization. Mathemat- ics of Operations Research, 31 2006.
R. Kannan, L. Lovász, and M. Simonovits. Isoperimetric problems for convex bodies and a localiza- tion lemma. Discrete Comput. Geom., 1995.
Ravindran Kannan and Hariharan Narayanan. Random walks on polytopes and an affine interior point method for linear programming. Mathematics of Operations Research, 37(1):1-20, 2012.
Yunbum Kook and Santosh Vempala. Gaussian cooling and Dikin walks: The interior point method for logconcave sampling. In Shipra Agrawal and Aaron Roth, editors, Proceedings of Thirty Seventh Conference on Learming Theory, volume 247 of Proceedings of Machine Learning Research, pages 3137-3240. PMLR, 30 Jun-03 Jul 2024.
Aditi Laddha, Yin Tat Lee, and Santosh Vempala. Strong ceff-concordance and sampling. In Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing (STOC), pages 1222-1222, 2020.
Beatrice Laurent and Pascal Massart. Adaptive estimation of a quadratic functional by model selection. Annals of Statistics, pages 1302-1338, 2000.
Yin Tat Lee and Aaron Sidford. Path finding methods for linear programming: Solving linear programs in ((sqr((rank) iterations and faster algorithms for maximum flow. In 55th IEEE Annual Symposium on Foundations of Computer Science, FOCS 2014, Philadelphia, PA, USA, October 18-21, 2014, pages 424-433, 2014.
Yin Tat Lee and Aaron Sidford. Solving linear programs with sqrt (rank) linear system solves. arXiv preprint arXiv:1910.000.0 2019.
Yin Tat Lee and Santosh S. Vempala. Geodesic walks in polytopes. In Proceedings of the 49th Annual ACM SIGACT Symposium Theory of Computing, STOC 2017, 2017.
Yin Tat Lee and Santosh S Vempala. Convergence rate of riemannian hamiltonian monte carlo and faster polytope volume computation. In Proceedings of the 50th Annual ACM SIGACT Symposium on Theory of Computing, pages 1115-1121, 2018
Yin Tat Lee and Man-Chung Yue. Universal barrier is n-self-comcordlan. Mathematics of Operations Research, 46(3): 129- 2021.
Yin Tat Lee, Aaron Sidford, and Sam Chiu-wai Wong. A faster cutting plane method and its implications for combinatorial and convex optimization. In Foundations of Computer Science (FOCS), 2015 IEEE 56th Annual Symposium on, pages 1049-1065. IEEE, 2015.
Yin Tat Lee, Zhao Song, and Qiuyi Zhang. Solving empirical risk minimization in the current matrix multiplication time. In Conference on Learning Theory, pages 2150-2157. PMLR, 2019.
Zhize Li, Tianyi Zhang, Shuyu Cheng, Jun Zhu, and Jian Li. Stochastic gradient hamiltonian monte carlo with variance reduction for bayesian inference. Machine Learning, (08(8))10111172)) 2019.
Yingyu Lin, Yian Ma. Yu-Xiang Wang, Rachel Emily Redberg, and Zhiqi Bu. Tractable MCMC for private learning with pure and gaussian differential privacy. In The Twelfth International Conference on Learning Representations, 2024.
László Lovász and Miklós Simonovits. Random walks in a convex body and an improved volume algorithm. Random structures & algorithms, 4(4):359-412, 1993.
László Lovász and Santosh Vempala. Hit-and-run is fast and fun. preprint, Microsoft Research, 2003.
László Lovász and Santosh Vempala. Simulated annealing in convex bodies and an O* O*(n²) volume algorithm. Journal of Computer and System Sciences, 72(2):392-417, 2006.
13

---

László Lovász and Santosh Vempala. The geometry of logconcave functions and sampling algorithms. Random Structures & Algorithms, 30(3):307-358, 2007.
Yichao Lu, Paramveer Dhillon, Dean P Foster, and Lyle Ungar. Faster ridge regression via the subsampled randomized hadamard transform. In Advances in neural information processing systems, pages 369-377, 2013.
Oren Mangoubi and Nisheeth Vishnoi. Sampling from log-concave distributions with infinity- guarantees. S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, Advances in Neural Information Processing Systems, volume 35, pages 12633-12644 Curran Associates, Inc., 2022.
Oren Mangoubi and Nisheeth K Vishnoi. Sampling from structured log-concave distributions via soft-threshold dikin walk, In Advances in Neural Information Processing Systems, NeurlPS'23, 2023.
Oren Mangoubi and Nisheeth K. Vishnoi. Faster sampling from log- concccc densities over polytopes via efficient linear solvers. In The Twelfth International Conference on Learning Representations, 2024.
Hariharan Narayanan. Randomized interior point methods for sampling and optimization. Ann. Appl. Probab., 2016.
Hariharan Narayanan and Alexander Rakhlin. Efficient sampling from time-varying log-concave distributions. The Journal of Machine Learning Research, 18(1):4017-40455 2017.
Jelani Nelson and Huy L Nguyên. OSNAP: Faster numerical linear algebra algorithms via sparser subspace embeddings. In 54th Annual IEEE Symposium on Foundations of Computer Science (FOCS), pages 117-126. IEEE, 2013.
Yurii Nesterov and Arkadii Nemirovskii. Self-concordant functions and polynomial-tmme methods in convex programming. USSR Academy of Sciences, 1989.
Yurii Nesterov and Arkadii Nemirovskii. Interior-point polynomial algorithus in convex program- ming. SIAM, 1994.
Lianke Qin, Zhao Song, Lichen Zhang, and Danyang Zhuo An online and unified algorithm for projection matrix vector multiplication with application to empirical risk minimization. In AISTATS, 2023.
Maxim Raginsky, Alexander Rakhlin, and Matus Telgarsky. Non-convex learning via stochastic gradient langevin dynamies: nonasymptotic analysis. In Conference on Learning Theory, pages 1703. PMLR, 2017.
James Renegar. A polynomial-time algorithm, based newton's method, for linear programming. Math. Program., 1988.
Sushant Sachdeva and Nisheeth K Vishnoi. The mixing time of the dikin walk in polytope simple proof. Operations Research Letters, 44(5):630-634, 2016.
Shai Shalev-Shwartz and Tong Zhang. Stochastic dual coordinate ascent methods for regularized loss minimization. Journal of Machine Learning Research, 14(1), 2013.
Zhao Song and Zheng Yu. Oblivious sketching-based central path method for solving linear program- ming problems. In 38th International Conference on Machine Learning (ICML) 2021.
Zhao Song, David P Woodruff, and Peilin Zhong. Low rank approximation with entrywise 11-norm error. In Proceedings of the 49th Annual ACM SIGACT Symposium Theory of Computing, pages 688-701, 2017.
Zhao Song, David P Woodruff, and Peilin Zhong. Relative error tensor low rank approximation. In Proceedings of the Thirtieth Annual ACM-SIAM Symposium on Discrete Algorithms, pages 2772-2789. SIAM, 2019.
14

---

Zhao Song, David Woodruff, Zheng Yu, and Lichen Zhang. Fast sketching of polynomial kernels of polynomial degree. In International Conference on Machine Learning, pages 9812-9823, PMLR, 2021.
Daniel A Spielman and Nikhil Srivastava. Graph sparsification by effective resistances. SIAM Journal on Computing. 10(1):1111))110199999 1926, 2011.
Joel A Tropp. Improved analysis of the subsampled randomized hadamard transform. Advances in Adaptive Data Analysis, 3001n02):115-126 2011.
Pravin M Vaidya. new algorithm for minimizing convex functions over convex sets, In 30th Annual Symposium on Foundations of Computer Science (FOCS), pages 338-343, IEEE Computer Society, 1989.
Santosh Vempala. Geometric random walks: survey. Combinatorial and computational geometry, 52(573-612):2, 2005.
Yu-Xisn XuaXin Wang, Stephen Fienberg. and Alex Smola. Privacy for free: Posterior sampling and stochastic gradient monte carlo. In Francis Bach and David Blei, editors, Proceedings of the 32nd International Conference Machine Learning, volume 37 of Proceedings of Machine Learning Research, pages 2493-2502, Lille, France, 07-09 Jul 2015. PMLR.
Virginia Vassilevska Williams. Multiplying matrices faster than coppersmith-wwnnoraad. In Proceed- ings of the forty-fourth annual ACM symposium Theory of computing (STOC), pages 887-898. ACM, 2012.
Virginia Vassilevska Williams, Yinzhan Xu, Zixuan Xu, and Renfei Zhou. New bounds for matrix multiplication: from alpha omega. In Proceedings of the Thirty-Fijh Annual ACM-SIAM Symposium Discrete Algorithms, SODA'24, 2024.
David P. Woodruff. Sketching as a tool for numerical linear algebra. Foundations and Trends in Theoretical Computer Science, 10(1-2):1-177, 2014.
Difan Zou, Pan Xu, and Quanquan Gu. Stochastic variance-reduced hamilton monte carlo methods.
In International Conference on Machine Learning, pages PMLR, 2018.
15

---

Appendix
Roadmap. Since there are many technical details in the appendix, we provide a roadmap. The appendix can be divided into 5 parts: the first part conveys some preliminary information and states the sufficient conditions we require barrier functions to have, in Section A and B. The second part proves the mixing time of the Dikin walk when barrier functions satisfy the conditions, and the proofs are divided into Section C and D. The next part focuses on the runtime complexity of sampling from polytopes, including how to generate a spectral approximation of the Hessian, approximate leverage scores and Lewis weights to high precision, in Section E. and how to incorporate these al algorithmic prototypes to implement an efficient sampling algorithm the log-barrier for sampling from a spectrahedron, as it is relatively less explored before. Finally, we prove the convexity of the function log det (x) + Ia) for log-barrier, volumetric barrier and in extension Lee-Sidford barrier in Section H and I.
Section F. We dedicate Section to study
A Preliminaries
For a vector x E R", we use ||=]|2 to denote its lg norm, i.e., ||æ||2 == (Ei=1 2? m/?11/2 We use ||ق||| to denote its & norm, ||ك||| := Di=1 |xs]. We use || to denote its & norm, i.., || == maX |aa]] For a random variable X, we use E[X] to denote expectation. We use Pr[.] to denote probability.
For any positive integer, we use [xx] to denote the set {1, 21,..... , .n..
We use O₂ to denote a length-d vector where every entry is O. We use 1a to denote length-d vector where every entry is 1. We use to denote identity matrix which has size d x d simply / when dimension is clear from context.
For matrix A, use A to denote the transpose of matrix A. For square and non-singular matrix A, we use A to denote the inverse matrix A. For a real square matrix A, we say it is positive definite (PD, i.e., A > 0) if for all vectors E IR" (except for 0,,), we have .Az > 0. For a real square matrix A, we say it is positive semi-definite(PSII,, i.e., = 0) if for all vectors x E TR" we det(A) to denote the determinant of matrix A. For a matrix A, we use |||] to denote its spectral norm, use ||A|| F to denote its Frobenius norm, i.e., ||A||e == CE j=1 A2 1//2 and use unz( A)) to denote the number of nonzero entries in A.
have Ax 2 0. For square matrix,
Given function f K. + R. we say it's convex if for any x, E K., f(x) > (wv + fyy) We say it's L-Lipschitz if (xx) - [(y) S L- || t - ulla for a fixed parameter L > 0.
T (a v).
For two distributions P1 and P2, we use TV(P,, P2) to denote the total variation (TV) distance between P1 and Pz.
For convex body K. C Rdd we use (K.) to denote the interior of K.
We use to denote Kronecker product. We A @g B to denote A @ I + I B. Given matrix A, we use vco(A) to denote its vectorization, since we will always apply vec(-) to symmetric matrix, whether the vetorization is row-major or column-major doesn't I't matter. We to denote Hadamard or element-wise product.
For two vectors x, y E Rdd we use (x. y) TTTy y to denote the standard Euclidean inner product over Rdd and for two symmetric matrices A, B E Rdxdd we use (A, B) tr(A" B] to denote the trace inner product.
In Section A.I, we present several definitions and notations related to fast matrix multiplication. In Section A.2, we provide some backgrounds about convex geometry. In Section A.3, we state several basic probability tools. In Section A.4, we state several basic facts related to trace, and Kronecker product. In Section A.5, present the standard notion about self-concordance barrier.
A.1 Fast matrix multiplication
Definition A.1 (Fast matrix multiplication). Given three positive integers a, b, c, we use Tmatt (a, b, c) denote the time of multiplying x matrix with another b x matrix.
For convenience, we also define the w(.,..) function as follows:
16

---

Definition A.2 (Fast matrix multipicaation, an alternative notation). Given z, V. we use d=((..0.2) to denote the time of multiplying d' du matrix with another d d' matrix.
Lemma A.3 (Gall and Urrutia [2018]). We have the following bounds of :):::));
11) w.
31) 4.199712 (see Table 3 in Gall and Urrutia (2018]).
Here, denotes the exponent of matrix multiplication, currently = 2.373 [Williams, 2012, Gall and Urrutia, 2018, Alman and Williams, 2021, Duan et al., 2023, Williams al., 2024, Gall, 2024].
A.2
Convex geometry
We define B(x, R) C Rd as
Blx, R) = tu e R 11 2112 s RJ.
Define polytope IC as
: ExER Ar sb.
Define spectrahedron K. as
K: to Rd
Σ A C},
where Ai...., Ad. C E R"x" are symmetric matrices.
We will often use the notion of Dikin ellipsoid:
Definition A.4 (Dikin ellipsoid). We define the Dikin ellipsoid D9 C. Id follows
Do = {w E Edd HT (0) 1)
where H(0) E Bdxd is the Hessian of self-concordant barrier 0/
We define the standard cross-ratio distance (see Definition E.1 in Mangoubi and Vishnoi [2023] for example).
Definition A.5 (Cross-ratio distance). Let u, K. If & v, p. q be two endpoints of the chord in K. which passes through and such that the four points lie in the order of p. u, v, q, let
a(u, w) = 1 u112
1p
4/12
pp u 0
We additionally set (u, v) = 0 for = v. For two subsets S1 S2 K., define
S2) = min v) : u E S1, Sz).
We state the standard isoperimetric inequality for cross-ratio distance.
Lemma A.6 (Lovász and Vempala [2003]). Let Rdd RR. be log-concave density, with support on a convex body K. Then for any partition of Rd into measurable sets S1, S2, Ss. the induced measure x² satisfies
**s) (2) x² (Si)
n (S2)
A.3 Probability tools
We state some useful concentration incqualities.
17

---

Lemma A.7 (Lemma on page 1325 of Laurent and Massart [2000]), Let X ~ Xx be a chi-squared distributed random variable with k degrees of freedom. Each random variable has zero mean and o² variance. Then
Pr[X - k 2 (2v + 2t)o² S exp(-t),
Príko X 2 2v % S exp(-t))
Lemma A.8 (Matrix Chernoff bound [Tropp, 2011]). Let XI,...,X. X, be independent copies of a symmetric random matrix X Raxd with [[[] 0, ||||| S y almost surely and | X]] S o². Let W = Lieel X₁. For any e E (0,1). 1).
SE²
Pr[l[W| 2 c) < 2d exp
( + ne/33
A.4 Basic algebra facts
Lemma A.9. Let A and B be m matrices. Then have
tr[AB] tr]
If A is symmetric, then tr tr(AB] +t[AAB].
If A and B are PSD, then (A, B) 2 O, and (A, B) 0 if and only if AB =0mxm mm
If = 0 and B 2 C. then (A, > (A.C).
Lemma A.10. Let A, B, C, D be conforming malrices. Then
(A (C @ D) AC 8 BD.
(A0s (C @s D) = (AC @g BD + AD 85 BC).
(A@ B) A B
If A and B are non-singular, then A B is non-singular, and () A-1 B-1
vec(ABC) = (CT @ A) vec(B).
A.5
Self-concordant barrier
We provide standard definition about eeff-concordance barrier,
Definition A.11 (Self-concordant barrier). A real-valued function F Int(K)) + R. is a regular self-concordant barrier if satisfies the conditions stated below, For convenience, if x & Int((C), we define F(x) co.
1. (Convex, Smooth) F is a convex thrice continuously differentiable function on Int(XC).
2. (Barrier) For every sequence of points {x:} E Int((C) converging to point x o Int((CC), lim;-ro f(x))= 6c.
3. (Differential Inequalities) For all Ted and all x E Int((C),, the following inequalities hold.
(a) F(z) is 2-Lipschitz continuous with respect the local norm, which is equiv- alent to
F'x)),,, h, k] < 2222FF(F)ll hl)
(b) F(x) is v-iischhiitz continuous with respect to the local norm defined by F,
[V ||
We call the smallest positive integer v for which this holds, the self-concordance parameter of the barrier.
18

---

B Sufficient Conditions for Log-concave Sampling
To prove that our algorithm works for log-concave distribution, we state three key assumptions (see Section B.1) which are the sufficient conditions to prove the mixing rate of log-concave distribution, for different barrier functions.
B.1 Conditions: v-symmetry, convexity and bounded local norm
We state three key conditions in order to lower bound the conductance for log-concave sampling. Assumption B.1. Given a convex body K, let H K Rdxd be a self-concordant matrix function, we assume that
(i) v-symmetry. For any x K, we have (1) S Kn (2x K) S Ey(DD) where (x(r) fy Rd xy- Hz) z) S r3.
(ii) Convexity, Let F K. - IR be defined as F(x) Log(det(H(()) + Ia)), then F is convex in &.
(iii) Bounded local norm. Ler V log(det(H(z))) denote the gradient of (og((ct(H(()))) in x. For any K, we have
((H(2))-1/2 {og(det(H(<)))
For the 7-symmetry assumption, we will prove that for barriers of the concern (log-barrier, Lee- Sidford barrier for polytopes and log-barrier for spectrahedra), 7 = This characteristic is proved in Laddha et [2020] for polytopes, and we prove that for spectrahedra Section G. For convexity, due to the large amount of calculations, we defer to Section H and I. For bounded local norm, the results for polytopes are similarly shown in Laddha al. [2020], and we will prove for log-barrier over spectrahedra.
B.2 Local PSD approximation for any self-concordant matrix function
In this section, we prove a generalization of the matrix function self-concordance regularized matrix function. When H/ is Hessian of the log barrier, this fact was proved in [Mangoubi and Vishnoi, 2023, Lemma E.3]. We generalize this fact general barriers.
Lemma B.2. Let (0, 1). Let H KC + I. be self-concordant matrix function and d(u) == a a- + +1 - For any u, v E K. such that ||u - u v(l I αu/z we have
(1 1/2 ||u
w (m))))
(()) (u) (1 + 1/2 ||u 22 b(v).
Proof. By [Laddha et al., 2020, Lemma 1.1], we have
(1 ||u ell =l H(v) (u)) N (1+ ||u vlIH(w)) HH()).
Because Huu) ((u), we get
(1 1/2 ||u ulla(w))) H(v) Hu)) 1 (1+ 1/2 I/u ullo(w))) H()).
Using (u)) = a -1 H (u) + 17 La again we get
(1 - a 1/2 /2 ||u - w/ (x.)) (v) = (2) = (1 a12 ||u ollo(w)) D(v).
This completes the proof.
B.3 Bounded local norm for regularized Hessian function
We note that the bounded local norm condition generally holds for standard barrier functions, but in our core argument, we will instead rely the bounded local norm of the regularized Hessian function. We prove that the bounded local norm condition implies the bounded local norm on the regularized Hessian function.
19

---

Lemma B.3. Let H K. -t I be matrix function. define F(z) = log (x) + 12)
L R. then we have the following inequality:
(x) L21a) { log [C
(eg
Proof. Let (x) UAUT be the eigendecomposition of the matrix ( (H)),, then
(22 = [| + [|(H(=)+LII)) +
= -|U + 1,2,,)=//
II(A LLdaa
the given condition implies
|H(x))1/2v log det (z)) ||((x)) 3/2 3/2 (
=IUA 3 V=
=|AA -|AA~
As 1 2 O, we can expand the regularized squared & norm as
II(A 1/12-3122111
Σ
L2 v3
(A:
i=1
d 1
- Σ
i=1
|| 1113//11/
This completes the proof.
C Key Tools for Robust Sampling
In this section, we list several tools which can be shared in the core proofs of all the main theorems. This section is organized as follows. In Section C.I, we provide tool for Gaussian concentration. In Section C.2, we show that PSD approximation implies the determinant approximation (by losing factor of d). One of the major idea in this work is instead of using exact accept probability, we will only need to use approximate accept probability. In Section C.3, we prove several useful robust properties for approximate accept probability. In Section C.4, we present general lemma for bounding the TV distance between exact process and approximate process. In Section C.5, we provide some specific choices for parameters and then bound the TV distance. Section C.6, we provide a definition which will be used in proof of lower bound on conductance. In Section C.7, we show how lower bound the conductance by cross ratio.
C.1 Gaussian concentration lemma
Lemma C.I. Let & (V, Id), thenfor any t /2d, we have
Pr[li[]22 2 1 S exp(-(12 d)/8).
Proof. We consider the squared lg norm of &: ||<| Σd_ S?, which is x² distribution with degree of freedom d. By Lemma A.7, we know that
PrIlee| 2
2vdk 2k² d] s exp(-1?),
set 2.2 = (t2 -d/88, obtain the desired bound.
20

---

C.2 Spectral approximation to determinant approximation
Lemma C.2. Ler cx E (0, 1). Given four psd matrices satify the following conditions
(1 cx (L))
(1+
(1 - eΦ ( (1+
Then we have
11-10-Ep d) det (xx))())))
det (xx))())) det(
I (1 10.60 cx d) det(F(x))
(e))) 2
(et()(z))
Proof. Since ((x)) E (1 + c (xx), have the following spectral bound:
(1 61- 4a =
= (1 +
we can take the determinant of the middle term:
det( (z) = dct( 1/2) detp))))))(())))))))) det((()))
det(
(x)(
e det((1 (1 co
E (1 (o)
E + 10 10. d.
Thus, we complete the proof.
C.3 Approximate accept probability
We start with some definitions:
Definition C.3. We define
Gz(z) == (et(@(u)) exp(-0.55 x =1(3
Gu(u) det (x))) 1/2 ex((-05 =((2(z))
(=(u) == det (x)))) 1/2 exp( 0.55| x =
Definition C.4. We define
(u) = minn (x exp(-s(u))
(.)(a) (()
G. (u) exp(- ((xx)
min Gz(2) exp( f(()) }
paccept (u) = E
(x)(x))
{(a)
G. exp( exp(-f(x))
We will use Pz. Ps shorthand for the above notations, and we call Pa the exact process and P the approximate process.
Lemma C.5. Fix x and (x). Suppose ((x) is (CH, 6H) good approximation (xx). Then
=
E
2,0(a) (u) paccept (u)l S 0.001.
(x)(x)
Proof. Note that
Gz(u)
(det 1/x)) exp( -1]l1=11--- u (1)(m))
G5 (u)
(det 1/2) exp(- - u (p)
21

---

Because ((x) is a good approximation to ((x), by Lemma C.2 we have
(det (z)))
(1 cx)d//2 -
I (1 + (x) 4/2
(det (x)))//2
and
| 2 [(()))))
(1 cH)² v
I (1 + cH)²
| 2 (xx)
By choosing step size sufficiently small, we have I|u
x x=li((0)
< d with probability 0.9999. Therefore
Pr ||| uI ul(c0)) - || - u1 [
S 3rHd 2 0.9999.
Then
0.5|-0
u| (1))
Pr
3cdd
I 1 + 3cHd
Σ 0.9999.
exp(-0.5)
((a))
For edd small enough, we have
paccept (u) 2.)(())
Pr
0.9999
< 1.0001 2 0.9998,
paccept (u) ((x))
and therefore
Pr[ Pe (u)
paccept (u) S 0.0001] 2 0.9998. xx)(2)
(.)(
So
E (x (xx)
Pacceps 2,4(2)
-
paccept
2, (w)
[u)l
- 0.9998 0.0001 + 0.0002
S 0.001.
This completes the proof.
C.4 TV distance between exact process and approximate process
Lemma C.6 (TV distance between exact process and approximate process). Let denote the (1 approximation of 4. Let dh denote the failure probability. Let Cp be defined as Lemma C.5. Let Pz denote the exact process. Let Pg denote the approximate process. If E Ic, then
TV(P,. Pa) S 88 + cp + dden-
Proof. We say ((x) is good if it is an (CH, 84) approximatin to ((x), and ((x) is bad otherwise. We can upper bound TV (P,,,, P.) as follows:
TV(P,, Pa) E[TV(P (0)), Pz. P,,
E(TV(P2, Ps. P2,ū(z)) ((x)) is bad]
E(TV( (o))), Pr P2))(a)) ((x)) is good]
< Pr(0(x) is bad + EETV (P Pa (()) ((x) good]
s Sn + ETT P. (()) is good].
where the third step follows from TV.:.) :) < 1, the last step follows from Lemma E.3.
It remains to compute the second term in the above equation. By Pinsker's 's inequality, when ((x) is good, have
TV(N(z, - (z)), -- (x))))
22

---

0.5 (z, (2) (x))
det ((z) det (xx)
-.25
1/2 - log
< 0.25 ((1 + E - - )
S 0.5 dcg-
where the first step follows from Pinsker inequality.
Now, we can upper bound TV P2) in the following sense,
(Vs Pa)
S 8n
TT
Pa a. b(x): P Pu,5(a)) ((x) is good]
S 8h +
E
[ (u) - paocept (u)l 1fu n/) 1fu n'Hl
(m,
(())
n. (,)(()
(a,, (xr)
any coupling
5H
E a~N(z, (a))
ppcceept(1)
2.7(2) (u) paccept
TV (z,) (-,
(x)(())
= =Ō + 6p + TV(N(a, -- (x)), t (x)))
I Öh + Ep + deH
The forth step follows from definition of cp. the fifth step follows from Eq. (2).
C.5 TV distance between exact and approximate process: instantiation
Lemma C.7. df 0.001 and бн 0.001. Then, have
(V,,,, P.) 0.01
Proof. We can upper bound it as follows
TV(P,, P.) < 811 Ep + de
< 0.001 + + 0.001
< 0.01
where the first step follows from Lemma C.6, the second step follows from choice En, 8n, and the last step follows from Lemma C.5.
n
C.6 Lower bound on conductance, definition
Definition C.8. Let 3 E (0,0.1) denote some fixed constant. Let S C and 82 = K3YS1 and T(i)) 1/2. We define St and Sy as follows
$ = (a S): <
S = (a E S2 ((() S B).
C.7 Lower bound on conductance, lemma
Lemma C.9 (Lower bound the conductance by cross ratio), Let S and St, be defined as Definition C.8. The conductance & satisfies
ス 3 (S(
16
Proof. The proof follows the general format for conductance proofs for geometric Markov chains (see e.g. Section 5 of Vempala [2005]).
Let S1 C and let S2 KA, Let Si and 5% be defined Definition C.8. We define S
5 = (KC)5)))
23

---

By Lemma 12, we have that
f #ES1
(x
P
(3)
DES,
Thus, by Lemma A.6, we have
(5a)) 2 ( (2)=" (5%)
(4)
Case 1. First, assume that both x *5 ( and T' (52) 2 nT (S2). In this we have
J S1 ((()))(())de
/ (z))) + 1 o Ps ()))))
2 2.s.
s.
S2
(S,))
8
(8), S2) (S')='(()))
2 08) S2) min( ' ( x" (S2))
3
> oSS 8.) min(m"in (Si)), π =² (S2)).
(S)
16
where the first step follows from Eq. (3), the second step follows from Definition C.8, the third step follows from Eq. (4), the fifth step follows from one of x (S)) and x* (S%) at least 1/2.
Case 2. Now suppose that instead either x TT(s) < (Si) or
x (St) < (S2).
Case 2a. If nr (81) < (Si) then we have
1 / (a(s))
/ (a())
f (r()))=())dz
2 S.s
S,
1
2 Si\S
1 3
2 (TT
3
8 min(( (Si), (S2)).
(6)
where the first step follows from Eq. (3), and the third step follows from Definition C.8.
Case 2b. Similarly, n* (8%) 1 (S2) we have
[ 5
((()))x(z)dz + + / (()))
/
(2))) = 2 Ss
((()))((z)dz
SSAS
3
IV IV
4 ²T (S2)
( 3 min( (Si), x (S2)).
where the first step follows from Eq. (3), and the third step follows from Definition C.8.
Therefore, Eq. (5), (6), and (7) together imply that
min( n (S2)) / (2()))m(z)dr' >
8 ())
(8)
16
for every partition S, U US = K.
Hence, Eq. (8) implies that
inf SCK:K (5)5 T (S) / (=(()))))()dd > OS S))
1
3
16
24

---

D Correctness for General Barrier Functions with Regularization
In Section 1. show how to bound the density ratios. In Section D.2, show how to bound the determinant. In Section D.3, we show how to bound the local norm. In Section D.4, we explain how bound the TV distance between the TV distance between two approximate processes. In Section D.6, we show how to bound the TV distance between two Gaussians. In Section D.7, we prove reversibility and stationary distribution. In Section D.8, we prove a lower bound on cross ratio distance. In Section D.9, we relate the cross ratio distance to local norm by utilizing 7-symmetry.
exact distributions. In Section D.5, we show how to bound
Algorithm Our algorithm for sampling from polytope with log-bau bog-barrier (formal version of Algo- rithm 1). As the only differences between log-barrier for polytopes and others the choice of and how to generate the approximate the Hessian, we only present this algorithm.
procedure MAIN(A E Raxd b R" 2,8 (0, 10 F. 0. 0)
Theorem 1.1
2:
K. (xr E Rd Ar S b)
(11/d) 3: Let be the log-barrier for K. with v = n
4:
5:
- ((1/(dL2))
6:
T (va
11²)) 1og(w/8)
7:
I + Za
b We are given to E Int(K)
8:
e @(1)
9:
+ O(e/d)
for -
11:
Sample point ~ W(, Ta)
Let H denote the Hessian of barrier function >> Algorithm
13:
((x) + SUBSAMPLE(A, b, x, n, eH)
14:
H (z) can be written Eies o,a,a, where |S| = E d dlog d 2 > The above step takes nnz(A) log + Tmai (d, CH d. d) time
15:
((x) - - (x:) +
Ia
17:
+ x + ((x) -1/2g
18:
if E Int(K)) then
19:
(1)) < SUBSAMPLE(A, z, n, n,
B Algorithm
20:
(2) + a + a
accept x with probability
exp(- [(=)) (det (de))) (exp)))eexxxx(003. exp(-0.5| =
min {
exp(-f (x)) (det (x))) exp(-0.5||x =1
22:
else
reject 2 end if
24:
end for
26:
return x
end procedure
D.1
Bounding the density ratio
Lemma D.1 (Bounding the density ratio). Let E Int(K).. Suppose f is L-Lipschitz and S 1/(10dT2), then
Pr exp(-f(z)) 2~
(x)) Lexp(= f(x)) 2 0.999.
in addition, f is B-smooth and 1j 1/(10)3), then we further have
exp((f(xx) 2 2 0.499.
Pr
x~ (()) (xxp(-f(2))
25

---

Proof. Proof of Part 1. Since 2 ~ M(x, 4(x) (() we have that
== +x) -1/2g
z + (a- H (z) + (a) 1a)
for some { ~~ (~)
Since a (1(2) + 11 la - la. and (2: and La are both positive definite, we have that
n la (a- (1) + (a)
Thus, upper bound
|| 2||2 follows:
|| - x l| = II(a 'H(x) + 17 (a) 1a)-1a1)11)1
= Hx() + + -aa
VVeTnaat
v llEll2
where the third step follows from Eq. (9),
Recall that & is sampled from (0, Ia), using Lemma C.I, show
Prilllla > t| < exp(-(12 - -d)8), Vt > 22. Combining the above two equations, we have
19
Pr[l - x112 > 20/20d/] S exp(
d) < 0.001.
(10)
8
Using 7 S 1/(80dL2), we have
Pr[l - æ212 > 1/(2L.) < 0.001.
(11)
Since f is L-Lipschitz, then we have
exp(-f(x)) = cxp(-(((e)-f(x)) 2 exp(-Li - x112)). exp(-f(x))
Therefore,
Pr exp(-[(=)) 1/2 2
Pr
exp(-L]1: =112) 2 1/2)
2~ exp( (ff(x))
z~ N=, (())
Pr
|| x|12 S log(2)/L
N (xx))
Pr
|2- æ|12 S 1/(2L)
N(z ((x))
0.999
where the last inequality holds by Eq. (11).
Proof of Part 2. Moreover, in the setting where f is differentiable and 3-smooth, we have that, since x is multivariate Gaussian random variable,
VV(===>)V\(>><>00 5 2
If (z (z-xx) (z( S O, we have that
f(x)
< (<x)
B| x x01
<1811=
26

---

Therefore,
(2)
Pr
2- (=
(x)) [(x)
2
(T2)
1
Pr
and (z-
=0
Pr
[2= - 2) Vf((x) > 0]
(~ (x))
T(x)
2
z~N(z,bb1(s))
2
Pr
[e (e
-
0.5
2~ Ma, (x))
2
v 12
Pr
112--211 v
0.5
vB
2~ (())
2 0.999 0.5
= 0.499,
where the last Inequality holds by Eq. (10) since 17 S 2043
Lemma D.2. Let || 2 2=))(=)) < 0.001. Ler n - 13. Then we have
If(x) f()) 0.1
Proof. We have
15(x) f(=) < L. | - =||2
< L v || - 2||
< 0.1
where the second step follows from ((x) 2 17 La. the last step follows from 7j S 0.01/L².
D.2 Bounding the determinant
Lemma D.3 (Bounding the determinant). Consider any e Int(K)) and & N(O, Ia). Ler
x (#(x)) 1/2gg. Then
Pr [og(det(((z))) log(det 2 -.1 2 0.999.
Proof. First, note that
(gg(det(((2)) - log(det(((x))) = log( det (()))
det dct(((z))
log( det(d H(2) + dI? Ia)
del(d. H(x)) + dI2 Ia)
del (+(()) + L.2 Ia) det(H()) + I Ia)
log(
thus, it suffices to consider the function F(x) = (og(de)(H()))))))))))))))))))))) + L Ia)) and bound F(z) F(x).
By Assumption (ii), we know that F(x) is convex.
Thus,
F(2) - F(x) 2 (2 - x) VF(x).
We know that
== (+(x)))
where ~ N((,
Thus,
F(z) (23-F(a) F(x) 2 E (T{)))) (0(x)) 1/2 1/2
27

---

By property of Gaussian distribution, we know that
< (())))1127((())
is a Gaussian with mean 0 and variance [[(((=))
By definition of & and Assumption (iii). have
[1(H())++ + L² (2) O(d)
and therefore
Thus, rescaling the function and applying the standard concentration inequality, we complete the proof.
Li
Lemma D.4 (Bounding the determinant, shifted). Let a < 0.001/d. For fixed x, 2 E K. such that ||2 (a) 0.001. Then, we have
log(det(((2))) log(det((()))) 2 -0.01.
Proof.
Let F(x) = log del (H) + Ia). Then
log(det( (((det(((=))) - Log(det = F(z) F()) 2 (z (2=2)
VF(x)
where the second step is because F(x:) is convex (by Assumption (ii)).
We have
- xz)-x) VF(x)) = 2) ( (0))))((=)//(((()())))) ((t)))3)
S [[(= - x) T)))
- 0.001 [0==))
$ 0.001 0.1.
where the second step is by Cauchy-Schwand the third step is by assumption ||> ==| 0.001, the fourth step is by 0.001/d /d and Assumption (iii).
D.3 Bounding the local norm
Lemma D.5 (Bounding the local norm). Consider any a (n), and g ~ NNO, Ia). Let
7 $ 0.001/d. Let = r + (#(x))-114 11/2.. Then
Pr [||= 100==1100 =110()) -I|= - x112(2)) -11=-x11((3) S 50.01] 2 0.999
Proof Sketch. Our proof is generalization of [Sachdeva and Vishnoi, 2016, Proposition 7J, where they prove the result for non-regularized log-barrier. Note that their proof also works for regularized barriers because two points are close the (xx) norm indicates that they are close in the H (æ) norm. The only difference is that when computing the Gaussian polynomial as in Sachdeva and Vishnoi [2016], we need to handle non- uniform weights instead the uniform weights of log-barricr. This could also be handled by observing that if two points are close in (x)) norm, then their weights close in the sense that [|w,(() (z)))((p(()) - (x)))) S Co where Cp < 1 some small constant that depends p (see [Lee and Sidford, 2019, Lemma 34]). One could then check that the argument of [Sachdeva and Vishnoi, 2016, Proposition 7] is robust under small perturbation to the weights, and conclude similar conclusions.
Lemma D.6 (Bounding the local norm, shifted). Consider any x In),K), and & NN(0, 1a). Let " 0.001/d. Let ||z x x=]((z)) < 0.001. Let (+xx)) -1/2g. Then
Pr |||u
=1
||u - x==11 < 0.1 > 0.999
28

---

Proof. Consider the following inequalities:
= [(())) ||ue x ==(2(w))
S u x =1(2(()) + + (2x - v, x - 2) 2))
S 21/u x (())))) + 41(u 2 2) 2)
2 n x (12(z) + 4/wl
where the second step follows from Lemma B.2, in the last step ((x)))/ 2)112).
First, using 1) is small, we can show ||u x ((x)) 0.01 with probability 0.9999.
Second, using || xz))(()) 0.001 and Gaussian concentration, we can show that
Pr[[w] mol 0.01] 2 0.9999.
Put them together and union bound, we obtain the desired result.
Lemma D.7 (Bounding the local norm, shifted). Consider any Int(K),, and & (0 Ia). Let " 0.001/d. Lel || x=)((x)) < 0.001. Let
(0()))11/25 We have
11/2g.
Pr[l]|u 2 l()tz) - x =1( 0.1] 0.998
Proof.
Let w' ~ + (0))((=2/) 11/2g' where e' ~ ~N, I2). We have
|||u (|0(z) - ||u x (())
< |||u = [la(=)) - ||ul - -2((a(=) + I||u - -2((b(>) 1|1l - ==(l(ur)
+ I1|u - 21(2(ur) - =|u ---2pp ftu=2llacc -1) + ||/u - = =lli(w) - ||u - (
So
Pr[lllu =1ba(=) I|u x x| 0.1]
Pr[lllu' =1lb(2> ||u' = 21 (w)))
> 0.05]
Lemma D.5
Pr[ll|u =1(b(w) ||u 2 [() 0.05]
Pr(u + u'
Lemma D.6
Lemma D.10
0.001 + 0.001 + 0.001 0.003.
where first step is by union bound, second step is by Lemma D.5, Lemma D.6, and by TV(u, u') S 0.001 when | - =l))())) 0.001 (Lemma D.1
D.4 Bounding TV distance between exact processes
Recall that define PS(S) Pr[æ SJ. We prove a bound between exact distributions of two
points & and z.
Lemma D.8 (TV distance between the exact distribution). For any x, 2 E such that || - 2|l S 0.001, we have that
TV(P,, 5 5 0.99.
Proof. We define
Gn, = (dot( hxtpp///_2ppp_311100000
(20 0.55lu =))b()))
Then
TV(Pa. Ps.)(e))
Pmin:
where
[min(1, Gz(u) exp(-f(u))G,() (x)()f(w)G.(() G.(u)
Pmin
(2, 1(x)) (u cxp(=f(=))G. (xp(-f()))G. G.(4)
29

---

For convenience, we define
(u)
01 (²)
02 cxp(-f(u))G..(t)
exp( (xp(-/(u)GG.(t)
G2(u)
0g
exp( (f(()(())G G.(u)
Let 0 E (0, 1). We have the following:
TV(P,, P.)
Pmin
E Emmin {1, a1 02 03)l
min(1, 0} Pr(min (1, 02 az) 2 0 -minfa1, 42, 42, fr,,infa1,,2.1 as} < <
min 0}
(atpLL,0)PP[mmm((((e Pr[minf(1, 02 az} 2 0]
S 0 Pr[min{a1, 02 as) 2 8]
0 0g
0.01 0.95
0.99
where the sixth step follows from Eq. (12), and the last step follows from choice of 0 and 0p.
Let a; =e- bb, It suffices to prove that for some 0 E (0, 1), we have
Pr[min(a1 a2,03) az} 2 0] 2 0g-
(12)
Suppose for some 7 E [1,5] we can show that
Pr[hi < r| 2 01
Pr[b < r| 2 02
Pr[by < r| 2 0g
Then by union bound, we can know that
Primax{(i, b2, b3} r]
Pr]min(a1; 02 03) > 6 e
Primin(a1, 02, 0z) > 0]
2 (1 01) (1-B2) (1 fz)
2 0.95
where the second step follows from - = 0 2 0.01, and the last step follows from 0 02 = 03 0.99. In the following, we will establish the three bounds of interest.
Part 1. By definition, we have
(u)
a1 (u)
(d)) exp(-0.5]|lu =112(a) + 0.5| =1);(a)) (det( (de)()
exp(p0000 (oop--05||lu uu 21(,(.) + 0.5||u x =1) 0.5 log det (x))) 0.5 0.5 (ogg(dd(((((()))
Using Lemma D.7 (on the term 0.5|/u 211ا(=) 0.5| u x =(م) and Lemma D.5 (on the term 0.5/Iu x (w) 0.5| x=122 we have
Pr (r.50/ 2 2l 0.5|lu 2 (1a(w)) + 0.5 x (B) 0.5|/u x 21 S 0.17] 2 0.999.
Lemma D).
Lomma D.5
30

---

From Lemma D.4, we know
0.5 (oe(((((((())))) - 0.5 (L)))) & 0.17.
Thus,
Pr[61 <] 0.99
Part By definition, for 02 we have
exp( (()
02 cxp( (x)
exp(-f(u)) (d))()
exp(- ox)() (de)( exp(-0.5||== -ll( + 0.5|/u - x x1(2()))
exp(-f(u) + f(x) 0.5 log
0.5 (og)((de((((())) - -0.5llx | -u12 0.5|/lu - ==1((())
By Lemma D.I, we have
Pr[f(u)- - f(x) < 0.17] 2 0.999
By Lemma D.3, we have
Pr(0.5log(dee(((2))) 0.5 (.5 < 0.1r] 2 0.999
By Lemma D.5, we have
Pr|0.5||z ull3 ul(3(u)) 0.5||u-=/le - =1) S 0.1r] 2 0.999
Then, combining the above three equations, can get the following result:
Pr(b2 S r] 2 0.99
Part
exp( (() G.(w)
ag exp( (xp(- (=(u)
exp(-f(u)) (det
exp(-f(=)) (d())) exp(-0.5] -uli 0.5|lu =112(z))
exx(-f(uu) + f()) 0.5 (og(de((((u))) - 0.5 log(det(@(z))) - 0.5|| ull + +..5 x =l)(()))
By Lemma D.1 (on the term f(u) - fxx) and Lemma D.2 (on the term f(x) f(=)), we have
Pr(f(u) - f(x) f(x) f(z) < 0.17] 2 0.999.
Lemma
By Lemma D.3, we have
Pr[0.5 (og((e((0(z)))) - 0.5 (ollode(((((())) < 0.1r] 2 0.999.
By Lemma D.5 (on the term 0.5| u =1la
0.5|/u x x112(,.) and Lemma D.6 (on the term
0.5l|= u| 2l 0.5 Iu x
1122 we have
Pr (r.555 - u u 0.5 |u 21 ( + 0.5||u
x: 2
0.5 x <112(0) S 0.2rr] 2 0.999
Lemma
Lemma D.5
Thus, combining the above three equations, we obtain the following:
Pr(bs < r] > 0.99.
31

---

D.5 Bounding TV distance between approximate processes
Lemma D.9 (Robust version of Lemma E. 10 in Mangoubi and Vishnoi [2023]). If 8n < 0.001, da 0.001, de < 0.001 and for any x;, = K, || =1) 0.001, then have
TV (P... P.) - 0.99.
Proof. We can upper bound TV(Pa, Pa) as follows:
TV(P,, P₂) - TV(P,, P.) +V)
+ + Ps) < TV + P.y
< 0.01 + 0.9
- 0.99
First step is by triangle inequality. Second step is by Lemma D.8. The third step is by Lemma C.7.
D.6 Bounding TV distance between Gaussians
Lemma D.10 (Lemma E.6 in Mangoubi and Vishnoi [2023]). For any 2, E such that || x1|0(2) I 0.001, we have
TV(N(x, x- ------())))) ((2, & /3da + 1/2. || =1(
Further, ad < 0.001, then we have
(v(N(x, + - (z)) (x))) < 0.001.
Lemma D.11 (Robust version of Lemma D.10). For any x, K. such that ||> - 2l(0(a) S 40
with probability at least 1 1/1000, we have
Tv(N(x, T6 ((,))),(2, /2dcH + / 3da 1/2. || 2 21(
Proof. By triangle inequality,
TV(N(>, (V (x)), :6 (z)))
& TV(N(c, $- (x), (x))
+ TV(N(z, (z), (z))
TV(N(2, ((, -z)())
Second term is bounded using Lemma D.10. First and third term are bounded using similar ways. We have
TV(W(<, क- (x). $ g (x)))))
< D(V(,, -- क (x)))
(1(((()) 1/100(=)1//000000000 1/2 d + log det ((x) det (xx)
< ((1 d (og
< den-
Thus, complete the proof.
32

---

D.7 Reversibility and stationary distribution
For any I E K, we define the random variable Zy to be the step taken by the Markov chain in Algorithm from the point z, that is, set = x ((x))_//2 where M/, Ia). If K. set Ze with probability
min exp(-I(=))
2(0)) 13
exp(-f(x)) (det (det( =16/
Else, set
We provide modified version of Proposition E.12 in Mangoubi and Vishnoi [2023] where our Markov chain differs from theirs.
Lemma D.12 (Reversibility and stationary distribution). For any Si, S2 S K. we have that
(())
yES, PSi))
ZES,
Proof. Let a, b be two i.id. random variables (one can view and as random coins to generate the corresponding sparsifers) such that aa.. is function of and Lb.y is function of b,
Let Pa.=() be pdf of (0,, T and (-)y()) be the pdf of NO,, T. Because and are i.i.d, they interchangeable, i.e., (a, b) has the same distribution (b, a).
JE P2
[e T(Φ) S q(a)) f Pa,(()) J. q(b)) minf T())pa,(())
d dy da dz
SESS,
JoES.
J. g(a) / ((b) { I min min(T(y)Pr., (z) r Jy dr db da IZES,
DES
g(a) ((b)
ZES, Joes. min(T())Pa
da db da
DES9 I g((a) [
/ (e)
dx da dy
Joes
T(V) P,
JEES.
The first and fifth steps are by definition of P.. The second and fourth steps are by changing order of integration. The third step is by interchangeability of a and b.
This proves reversibility and stationary distribution.
D.8 Lower bound cross ratio distance
Definition D.13. define
F(n, ar R) := 0.1/ h.1//na-1 1 122
For simplicity, use F to denote F(n, 3. RR. Lemma D.14. Let S and S be defined Desinition C.8. Le F be defined as Definition D.13. We have
(S, > 1000F.
Proof.
Using Lemma D.17, we have that for any u, v E K.,
alu,v) F. | vl]
33

---

We will try to prove (u, 2) > 1000F in the following. To do that, we will prove it by making contradiction.
Suppose (u, v) S 1000F, then we have to require ||u - u v(lo(w) S 0.001.
Once we had that ||u (w) I 0.001, using Lemma D.9, for any u, v E K. we have that
TV(P,, P.) - 0.99
(13)
On the other hand, Definition C.8 implies that, for any u E S1, S/2 have that TV(Pu. P.) (S)) (S 2 2 0.99.
(14)
where the last step follows from 1 - 2/3 > 0.99. Thus, obtain contradiction, which means that S)) 2 au, v) > 1000F.
D.9 Bounding cross ratio distance through y-symmetry
In Lemma E.2 of Mangoubi and Vishnoi [2023], they only prove the cross ratio distance bound for log-barrier function. Here, we generalize it to arbitrary barrier function. The key concept our proof relies on is the notion of E-symmetry (Assumption (i)).
Lemma D.I.. Let H E Rddd be a functtion,, and u, v E K. Consider the chord that passes through u, v and intersects with K. Let p, be the two endpoints of the chord, with the order p, u, V. q. Then, we have
1p - ull =|l I VD.
Proof. Suppose can show that p E K' (2u K.), then it naturally follows that p E E, E.(vV) and lpp -ull S Vb.
p E K. trivially follows by the construction of the chord. To see p 2u K. is enough show that for some K. we have p = 2u. y. We claim 2u p is such a choice. To see 2u K., we need to show that 2Au Ap S We partition the constraints into two sets.
Case 1. Suppose for as we have a. u S a. p. then 20, S 2a S a, p + b, therefore, 2 u a. p $ b;. Case 2. Suppose otherwise, a. u > a. p. then consider
20, a, p b. (a, a p) + (a, u b,)
(a u b,) (a] b,) (b. (4-aTu))
41
d2
We shall show that d1 S d2 via a geometric argument. Consider the projection both p and u onto the hyperplane a,r bj. and the two chords that pass through p and respectively, orthogonal to the hyperplane. Let us denote them with ln Lu respectively. Observe that lp and l, are parallel, and ||, - Lu112 (a, u - b.) (a, b.) = d, It should then be obvious that di S d2: let proj(u) denote the orthogonal projection of onto the hyperplane, then
|/- proj(u))l2 = 1|4. - loll + +11 proj(p)ll2
since proj()()) 2 0, we have that di S d2, and consequently,
20 u a. parpe by
This proves that 2u E IC, and we can conclude that E 2u - K. which yields our desired result.
Remark D.16. The above argument be extended to spectrahedra and more general convex sets.
34

---

Let u, v E K. be two arbitrary points, and consider the chord that passes through u, with two endpoints being p. q. Note that p. ŌK.
Lemma D.17. For any u.u K. let H be any matrix function satisfying Assumption (i), (ii) and (iii). We assume that K. is contained in ball of radius R and has nonempty interior.
For any parameter Then, we have
(0, 1), (0, we define matrix function ((u) == a-1, H(()) I
1 2vv + 17 R
v), 2
|u vl)
Proof. Without loss of generality, assume that
lp - u212 s llu 4112-
(15)
We can lower bound 0² du, v) as follows:
1
0² (u,
||
0
IV
max
V 9
Iu v
IV IV
max
| q
7
||u 0
max
|/u q
14
||u r
m
x1
|| 9
1 ||u 11
22 v 2
2 ||p 1
2
p q
1 ||u R 21i
|/ 1
2 || u
R²
1 uu 1
2
| 2112
+
2 V
2
R²
(u v) T
1
X xa
xa-lE
H
1
IT)
2va
2R
1
v112(u)))
2va 11 + 21
R2
where the second step is by 1p g/12 2 max [lp u112: ||0 - g/12]. the third step is by ||u glla <
||u gll2, the fourth step is by ||p 4 p u/|2, the fifth step follows from ||p u||2 S ||u-glla
(see Eq. (15)), sixth step is by max [A, B} 2 0.5 (A + B), where the seventh step follows from [pp 2/12 R and Eq. (16). The eighth step is due to Lemma D.15.
remains to show Eq. (16). Due to the fact that p, u, are on the same chord, vector u v and p on the same direction. This means we write
(p u) for some E IR, and
Iu v H1l
||c (p 2) ||
Ip w H
||p - u||# H
1u 2 2
(16)
1p u
This completes the proof.
35

---

E
Computing Approximate Lewis Weights and Leverage Scores
In Section E.1, we present an algorithm to approximate the Hessian. In Section E.2, show how to compute leverage score to low precision. In Section E.3, we show how to compute Lewis weights to high precision.
E.1
Subsampling to approximate the Hessian
Algorithm 3 Subsampling Algorithm
procedure SUBSAMPLE(A e R b IR" ,,,,,,,,,d cH (0,
2: Eo < 0.001
3: m 4 ezj'd d
4: Form B E nRa where i-th row of B is a, ((a;, x) b,)-1
5: Compute the O + appproxximaii to the leverage score of B
6: Sample matrix B Rmxd such that (1 EH) <H) B 3 B B 3 (1 B
7: return B
8:
end procedure
In this section, we provide a sparsification result for (w) using leverage score sampling.
To better monitor the whole process, it is useful to write H (w) as A (ww) 2A, where A E Faxd is the constraint matrix and S(w) is diagonal matrix with S(w); a by. The sparsification process is then sample the rows from the matrix S(w)- A.
Definition E.1. Let B E Raxd be a full rank matrix. We define the leverage score of the i-th row of B
oB(B) = by (B by
where by is the i-th row of B.
Definition F.2 (Sampling process). For any E K, let H(w) = AT S(w)-22 Let 2 B-ar(S-' Aw) suppose we sample with replacement independently for s rows of matrix ((w) A, with probability pt of sampling row for some 1, Let ilj) denote the index of the sampled in the j-th trial. Define the generated sampling matrix as
1 [(w) = PIO) (0 b:()))²
Lemma E.3 (Sample using Matrix Chernoff), Let e, 8 E (0, 1) be precision and failure probability parameters, respectively. Suppose (w) is generated as in Definition E.2, then with probability at least 1 8, we have
(1-e) Hww) - (ww) (1+ H(w).
Moreover, the number of rows 0(3 -
2d log log(dd)))))
Proof. The proof will be designing a family of random matrices X. Let Yt = (A ( Note that
a; be the i-th sampled row and set Y, p. V.v. Let X, = Y, Ia-
n
yo Σ (AT
a.a (A -(w)-2
(ATS(ww -2A) 1/2
S(w) u.0 (A 5(w) -2 AA)
1=1
(AT)(w)) AA (A2 S(w) AA) S(w) 2A) -1/2
Ia-
(17)
36

---

Also, the norm of ys connects directly to the leverage score: 11 = S(w),,, a. (AT
a
(18)
We use i(j) denote the index of that has been sampled during j-th trial.
Unbiased Estimator. Note that
E[ = E[y] Ia
11
Σ
pi
ucy -
La
Pr
0.
Bound on | To bound
|X|, we provide a bound for any ||x;|| follows:
| |x|| = || - -aal
S 1 ||&ill
-1+ + 11
Pi
11:112
S1+ 8 o. (S (w)A)
d
(W)SS((W)A)) (S(w)A)
+1.
1+
8
Bound on B[XT X]|I- We compute the spectral norm of the covariance matrix:
Mi( Mit))
E
(G) x = La + R
-
P
Pi
n
oe (Sw)) A)
Ia
yall
21a
Pi
d
8 VeVi
Ia
1=
1)Td:
the spectral norm is then
d B
| (1) X
1
Put things together. Set 7 = 1 $ and 0² = 1, we apply Matrix Chernoff Bound as in Lemma A.8:
SE
Pr[ 2 {] 2d exp
d/3 - 1 + d/3)c/3 2d exp(-se² (33dd)
where we choose s = 0(B - -2 Finally, we notice that
1
1
W-
wiyy - 1a)
Pi())
37

---

(AW)(w)) -2 AA) 1/22 1/2
Σ
(A ((w) -A)--AA--- 1/2 La
j=1 Pi() (1/g) (w())-bb(() b(())
-(w)) (1o Therefore, we can conclude the desired result via ||] 2 c.
E.2 Computing leverage score
In this section, we provide an algorithm that approximates leverage scores in near-input sparsity time. Our algorithm makes use of the sparse embedding matrix of Nelson and Nguyên [2013] that has (log(nn/e) nonzero entries per column.
Lemma E.4 (Approximate leverage scores). Let A E Baxd, there exists an algorithm that runs in time ( (((=((nn(A))) A) dd)))) thut oulputs a vector Wea (A) E IR", such that wz((A w2 A).
Proof. We follow an approach of Woodruff [2014], but instead we use a high precision sparse embedding matrix Nelson and Nguyên [2013] and prove two-sided bound on leverage score.
Let 8 be a sparse embedding matrix with 2d poly log(d/(<d))) rows and each column has ((c-"lo((nd) nonzero entries. We first compute S in time O nnz(A)), then compute the QR decomposition SA = QR in time 0(e-2(') Note that R E I dxd hence R-1 can be computed in (dd)) time.
Now, let G' E Rddxt matrix with t log(n/8)) each entry of G is i...d. (00, random variables. Set q: = ||e: ARR for all E [n]. We argue is good approximation to (w> (A));
First, with failure probability at most 6/n, we have that qi Re, ||e| AR- 113 via LindeonnLinddenstaauss lemma [Johnson and Lindenstrauss, 1984]. Now, it suffices to argue that l|e, AR 112 approximates lle U well, where U E Raxd is the left singular vectors of A. To this, first observe that for any Rd.
|AR 212 (1 e) ||SAR 2212
(1 +) ) |Qx1
(1 + c) 11=112
where the last step due to Q has orthonormal columns. This means that all singular values of AR- are in the range [ - e, 1+. Now, since U is an orthonormal basis for the column space of A, AR and U has the same column space (since R is full rank). This means that there exists a change of basis matrix T' E ddd with AR-1T U. Our goal is to provide a bound on all singular values of T. For the upper bound, we claim the largest singular value the contradiction that the largest singular is larger than + 2 and let be its corresponding (unit) singular vector. Since the smallest singular value AR is at least e, we have
most + 2e, to see this, suppose for
| 2 (1 (l-e)l[Ta TT
(1 (1
> 1,
however, recall AR-IT U, therefore ||AR 1 112 = ||00 = = a contradiction. One can similarly establish a lower bound of 1 2c. Hence, the singular values of T are in the range of [1 - 2e, 1 + 2e]. This means that
||e7 AR
1 112 [leT UT 1/
= (1 + 2c) le, |e,,,, U
(1 + 2c)(w2
as desired. Scaling to e/2 yields the approximation result.
Now, regarding the running time of computing 4,, note that we can first multiply R-1 with in time Ō(e-2 this gives a matrix of size d x t. Multiplying this matrix with A takes Ō mmn(())) time. Hence, the overall time for computing W (A) g R is Ō(e-1 mnz(A) + +e- -2d)). D
38

---

Remark E.5. It suffices to choose 1 for generating approximate leverage scores and subsam- pling the target matrix (x).
10
E.3
Computing Lewis weights
Before stating the main technical tool for this section, we want make several remarks regarding the Lewis weights.
Since its introduction by Cohen and Peng [Cohen and Peng, 2015] as algorithmic tool for lp row sampling. Lewis weights, as its natural formulation follows from a convex program, is not known to be solved exactly in polynomial time. All known algorithms [Cohen and Peng. 2015, Lee and Sidford, 2019, Fazel et al., 2022, Cohen et al., 2019, Jambulapati et al., 2022] can only compute c-approximate Lewis weights.
Definition E.6 (fp Lewis weights), Let A Ruxd The Lewis weights of is vector W wp(A) R" satisfying
(x(A) - (W,(())
where (A(A) IR is the diagonal matrix that puts wp(A) on the diagonal, a Ruxxd - R is the operation that oulputs all leverage scores of inpul matrix.
are now posed to state main tool from Lee and Sidford [2019].
Theorem E.7 (Theorem 46 in Lee and Sidford [2019]). Lei P = fx: Ax > b} denote the interior of
non-empty polytope for non-degenerate A E Raadd There an O(d n)-self concordant barrier db defined by lp Lewis weight with p ((logn) satisfying
A W = V (x)) = (q W,A,
where A₂ diag(Ax b) and We is the E, Lewis weight of the matrix Ag. Furthermore, we can compute or update the Wx. (ib(z) and (2)((x) as follows:
Initial Weight: For any I E Rdd, one can compute a vector Wa such that (1 - -e)W> S Wa S
(1 + e)w,, in O(nd~11/2 log" n log(n/e)) time.
Update Weight and Compute Gradient/Hessian: Given vector w. such that w. (1 + 100)u, for any y with || - yllaT W.A. S log² for some small constant 0, we can compute wy v and H such that wy = (1 + E)W
(=(())lv<) 1 S, and (1 (l=e)(2(()) 3 (1-
in Ond log n log(n/e)) time.
F Fast Approximate Sampling from Polytopes: Complexity
In this section, we provide the runtime analysis for sampling from polytopes via log-barrier and Lee-Sidford barrier. The result for log-barrier is in Section F.1 and for Lee-Sidford barrier is in Section F.2.
Log-barrier in nearly-linear time
Throughout this section, let E (0, ()) and n (0, 12(ddL)))
Lemma F.1. If g is log-barrier function with n, then each iteration of Algorithm 2 can be implemented in O(nnz(A) + Tmaa(d d3, dd) time plus O(1) calls to the oracle of f-
Proof. We go through each step of Algorithm 2 and add up the time and oracle calls for each step:
We first need to sample dd-dimensional Gaussian random vector ~ N(0 Ta), which can be performed in O(d) time.
39

---

At each iteration, by the definition of log- berrier function, we can write the Hessian as follows:
aza
((x)
a x b,)
We define define matrix B E Raxd as follows
: A
where 5; = ((a,, 2) b;), Vi E [n]. Then we can write H (z) = BTB E Rddd
For matrix B, we can compute a constant approximation of its leverage score, then we can sample according to leverage, and generate diagonal matrix such that
(1 (-=)BB' B BT DDB (1 + + This step takes (An + d=.3..))) time.
Computing (xx) can be done in d time and computing the vector
2 (+)) &
requires to invert the matrix (xx). taking the square root and multiplying with vector &. Inversion and square root can be performed in time (dd)) by computing its spectral decomposition, and the matrix-vector product can be done in time (dd)) Hence, the overall time of this step is (dd))
To determine whether E K, one can simplify verify A= in time ((nz((A))
We also need compute (z) and (z). This similar computing (x) and ((x).
We then need compute the determinant det (+()) and det (+(() which can be done in ((d)) time via spectral decomposition.
We then need (1) oracle queries to the function values of f.
Therefore, adding up the number arithmetic operations and oracle calls from all the different steps of
Algorithm 2, we get that each iteration of Algorithm can be computed in Ó(nnz( + Tmat (d, d,, d) arithmetic operations plus (1)) calls to the oracle of f-
Remark F.2. We note a recent work Mangoubi and Vishnoi [2024] has provided an algorithm that runs nnz(Al A) d² time per iteration, by adapting techniques from Laddha et al. [2020]. Their per iteration cost improvement only works for log barrier, and they could only manage to obtain an O(nd + d R²) mixing time. Our framework the other hand generalizes to other barriers such as Lee-Sidford barrier, which provides a nearly-optimal mixing rate. We also adapt the framework to sampling from spectrahedra.
F.2 Lee-Sidford barrier via approximation scheme
Lemma F.3. If g is the Lee-Sidford barrierfunction with d d'oog n, then, each iteration of Algorithm 2 can be implemented in Óndd=1) time plus (11) calls to the oracle of f.
Proof. For Lee-Sidford barrier, suffices to use (xx) A[ Ar as where Wz is the Ep Lewis weights for (log m). However, as Lewis weights cannot be computed exactly, we shall use the algorithm of Lee and Sidford [2019] compute an e-approximation.
To invoke Theorem E.7, we need to make sure that = and x satisfy || x x:IAIW,A.
log We
can achieve so by scaling down a by factor of log" n. This only blows up the convergence by Õ(1) factor, so it is acceptable.
The runtime analysis then similar to Lemma F.I. Computing the initial Lewis weights takes O(nd=11/2 log(1/e)) time owing to Theorem E.7. As the algorithm requires Õ(d)) iterations, this time can be amortized. For each iteration, it needs query the Lewis weights data structure that outputs an e-Lewis weights in 0nd 1og(1//lgg/// (og)) time. we will choose as O(1d)) this time is Ond We can then use this approximate Hessian to progress the algorithm. All subsequent operations take O(d)) time. Thus, the cost per iteration is ((dd=1))
40

---

G Approximate Sampling from a Spectrahedron
We present algorithm that efficiently and approximately samples from spectrahedron, which is popular convex body utilized by semidefinite programs.
G.1 Definitions
Definition G.I. Let ..,..., Aa E Ruxa be collection of symmetric matrices and C E R be symmetric. define corresponding spectrahedron as
K [=
CJ.
We define the log-barrier for spectrahedron and its corresponding Hessian. Definition G.2 (Nesterov and Nemirovskii [1994]). Let K. be a spectrahedron described by symmetric matrices {A Ad} " and C R= The log barrier &log (z) is defined
(log(1) log (S(xx))
and its corresponding Hessian defined as
H1og(²) = A(S(z) & S(z))1AA
where S(x) == Ei a,A, CE Ruxn and A E Bdxnn and the i-th row of A is the vectorization of A, E Rux=
To compute (1)g it is handy to define matrix B E Ødxx²
Definition G.3. We define a matrix B Rdxn² B A(S(x) & 5(x)-1/22) Consequently, Hoog (z) = BB
state a useful lemma for computing the matrix B.
Lemma G.4. Lel matrix B E RRxxn be defined as in Def. G.3. Then, the i-th row of can be computed vec(S(x))1 112 (,
Proof. The proof relies on a simple fact of Kronecker product and vectorization:
vec ((x))()))4) 1/2 A,S(z) 1/3) (S(x)-1/2 ((x) vec(A:), which is the definition of the i-th row of B.
G.2 n-symmetry of Hog for spectrahedra
In this section, we prove that Hoog is n-symmetry for spectrahedra. Lemma G.5. Hoog(r) r--symeery,, that is, for any K.,
- K.) (2
Ez(Vn).
Proof. Pick point (1). then we know that (y - x) TA(((=)) @ S(x) S 1. We note that the set IC n (2x KC) can be characterized as the set of all points E R such that
u,
d
Σ
(2r u);A; 2 C.
These two conditions imply
S(xx) 3
[w: x;)A, 3 S(x),
41

---

d [1 - 14= Im
-1p
We will prove y satisfies these conditions. Note that y E (1) implies that
y1 x4,) 1/2
1
yt - AA,{)() <
==1
we let M denote Ei= (w: x, (. then the above condition implies that M M 3 In meaning all the eigenvalues of M lie between [-1, 1], thus we have shown that (1) C K)).
For the other direction, we know that
d
In
(3) -= In.
1=1
and henceforth
(y x) (A((()) ( ((s))
(Ue (:))
$ 1 n,
where we use the fact that - < M < Jn hence all eigenvalues of M lie between [-1, 1]. and the squared eigenvalues of M have its magnitude most 1. We then the squared Frobenius norm is equivalent to squared & norm of singular values of M, as desired.
G.3 Bounded local norm of (z))(r)
We prove that H1og(2) has the bounded local norm property.
Lemma G.6. Let Hing (2) E IRdxd be defined in Def. G.2. Let F(x) log (Hing (x)). Then, F(x) is convex in x.
Proof. The proof is by observe that the function F(x) is actually the volumetric barrier function for the SDP spectrahedron. Thus, the function F is convex, for more details, see Nesterov and Nemirovskii [1994].
Lemma G.7. Let F(x) be defined as in Lemma G.6. Then, have
(=2((((())12AA,5t))) Øs (545)=122
((((((()11/AAA55()) 05 (144)
VF(x)
os (S(x) 1/2 S(x) ((m -1/2))
where P = Hx)
Proof.
For simplicity, we let H to denote Hogg For each coordinate (d), we have
0r, (x) - (z))
OF
Or,
-r( 0(S(z) S(z) AT] dz;
42

---

(r(H())~\A( -5(x)= (1050)
dr; - +5()))++(()))+ 0S(x)= )
dz,
- -(((())))A(((())) 08(x) + ((x) (S(x)) 00(())
dx, S(x) "))A
dz,
=
examine the following term:
= (TATTT1007(((50444)))11111111110
(166""""243336""7"""""""""00"0"3 - ((P((3100044))8555)))111111"1""" ((x))
(=P())(ee-""AAA5())118)))) 1,)]
where P E 18n2 xn²
is the projection of Bf ie. P P (2) ( Thus,
VF(x)=
8s (Sn))
This completes the proof.
Lemma G.8. The matrix function (x) A E ddd has the following bounded norm property: for any direction E Rd.
Proof. Let Zt = x + t. h for some fixed vector h E Rd, S. = E *1 C E Ruxn
pdx=" P. AT (AAAT)) A, g=²xn" and E, [e= & L.) & In) >n (4n & e P. & e;). Note that for any matrix M gaxn tr(Pi (M x
Ar A(S, 8 S. S-1/3
[,]] = tr(p (I, & M)] tp MJ. Also note that tr(E.] = tr(P,] = d.
Let H, = AS," & ST A E ddd We have
H H.) H, 1/
10
-t((H,' H.) H ( H.)
0t
Ot
OS 0 S-
USS @ ST Al
=tr(H,
A H.IIA AT
Øt
ōt
= trrAAH; AT OS. S-1
Ø(S-1 Ø S.
AH,
AT
0t
0t
0S
@ S:1
tr(PP(()) 1/2
S. 1/2
SS
1/2 S. 1/2
Øt
-P 1/2
8 S5 1/2 05,-
S
(S
1/2 S. 1/2 1
Ot
1/2 1/2 -1 05
8 S
05
8 S,
tt(( (S
1/2 S.
(S.
2 S 1
(S.
1/2
1/3)
dt
Øt
= tr( 1/2
h: A:A 1/2
8g 2
43

---

d
tr(P: ((S 1/2
, 1/20 8s In + (5, 1/2
S 1/2 @ S 1/2
h, 1/
==1
-2tr(s, 1/2
h: A 1/2 1/2. St -, 2
h4 1,,,))511
1=1
tr(Pt. (5
1/2 8 &SS, S, 14 A
1 1
d
d
<4tr(S,
Et.S S,, 1/2 1/2
h:A:)5, 1/3
d
I t4) 1/2
hgA,)S 1/3
<4tr|E]]r[(S tr( 1/2
S, 1/399)
=4dlh
the fifth step is by for PSD matrix M, we have MPM 3 M where P is projection matrix, the sixth step by the below derivation, the seventh step is by tr(P: (M² @s I,) = 2r (rrMM, the eighth step is by for PSD matrix M, 2M @s I, 2 M & M, the ninth step is by Hölder's 's inequality where (A, B) S || [[B]|,, for B s1 be its Schatten-1 norm. B is PSD, |B|] = Lr(B), the
tenth step is by for PSD matrix A, ||| & tr[A].
Note that
0
aS Ot
S,²
St
Ōt
= - 5,
A:)S,
i=1
Thus,
0s,
@ 8.
d
S
h A.); 8 S.1
01
1=1
S. S
and we have
(S-1/2 s S. 1/11
1/2
05, S..1
(S, 2 S-
05,
@ S,
(S- 83 S. 1/123-11
1/2
Ot
Øt
d
=SS-"/ 1/2
@ S. 1/2 - 2
MAA:) & + S S.
S ((, 1/2
@ //12/--11
(S- -S-1//2 1/2
SS, 1/2
1 hAA)),; & S. + S. 8 S
h, (S, 1/2
1/2 @ S. 1/1
8
h4 h4 -(S-1/2 0g 1n)
1/2
1=1
plug in the above derivation obtain the desired result. Finally, we send t + 0 and conclude the proof.
Lemma G.9. Let (x) be the matrix function A(S-1 1 S-1)AT, then we have
|1 Ō(d).
44

---

Proof. We have
|
max (M()))1/V<<)))
w:||a||a=1
max
==]le]a==1
max
max
S 2d ||u|l
2d.
where the fourth step is by tr [M] S vd. ||M]| for d x PSD matrix M, the fifth step follows from Lemma G.8, the last step follows from || H()) = 1.
G.4 Fast Hessian approximation
We need a particular type of sketch for Kronecker product of matrices. Definition G.10 (TensorSRHT [Song et al., 2021]). The TensorSAHT IR" x R - TRe defined as S == J (HD, & HD2,, where each row of P E 13, contains only one 1 at random coordinate and one can view P as sampling matrix. H is a x n Hadamard matrix, and D1, D2 are two n. X n independent diagonal matrices with diagonals that are each independently set to be Rademacher random variable (uniform in {-1,1})
Lemma G.11 (Lemma 2.12 in Song et al. [2021]). Let IT be a TensorSRHT matrix defined in Definition G.10. If O(e 2d log (nd/(<6))), then for any orthonormal basis U E I xd have that with probability at least 1-8, 8, the singular values of IIU lie in the range of [1 e, c).
The following result provides efficient embedding for B.
Lemma G.12. Let B E ıxxnn be defined as in Def. G.3, e (0, denote an accuracy parameter and 8 (0,1/10) denote failure probability.
Let Raxo? be TensorSRHT matrix with s O(E 2ofe-?dloo 2d (nd/(c8))) then we have
Pr[ T x|12 = (1 B)| æ1|2, Vo E Rd 2 1 1-6. 6.
Moreover, IIB can be computed in time
O(d Tmae (n, n, s)).
Proof. The correctness part follows directly from Lemma G. It remains to argue for the running time. We need to unravel the construction of both S and (see Definition G.2). Recall that P. (HD1 & HD2) and
((S(x))1/2 & S(x) voc(A,)
1 P.HHD & HD2) (S(x) 1/2 @ 5(x) -1/2
voc(A,)
V8
- (H @D (IID)a"")""""""DDaaa""""""""" )ecc]AA
1
V8
1 1
since P is a row sampling matrix, the product can be computed as follows:
First compute HDrS(z)-1)) -1/2 and (2)))) Since H is Hadamard matrix, this step can be carried out in On lo time.
45

---

Applying P to the vector can be interpreted as sampling s coordinates from the ma- Let (i., (d., sampled by P. We form two matrices X, Y R2x= where the k-th row of X It is easy to
trix ( (, A,5(x)-1/2D,11
denote the coordinates
is (HDyy(x))11/))))
and the k-th row of Y is (HD=5(=)-1/)))>> verify that the (e) entry of XAYT Rxxs is the corresponding entry of (,55)))))
This step therefore takes Tima (n, n, 8) time.
As we need apply the second step to d rows, the total runtime
O Tmat (n, n, s)).
Remark G.13. For comparison, compute B straightforwardly using the Kronecker product identity takes (dd)) time. As long as the row count is independent of n. Lemma G.12 approximately computes B in time O(n² poly(d)). For the regime where d & n, our algorithm is much more efficient.
H Convexity of Regularized Log-Barrier Function
In this section, we prove the convexity of log det Hog (x) + Ia) for both polytopes and spectrahedra.
H.1
Convexity log det (Hog(() ) Ta): polytopes
a.a We prove that for Hxx) = Ei= (a, x-b,)² the Hessian of the log barrier, we have the convexity.
Definition H.1 (Ridge leverage score). Let E Badd we define the X-ridge leverage score of as X (B) = b; B Aaa b,.
If X we abbrieviate as (B).
For convenient of the proof. we also define
( : by (B B+ bb.
We also define xx (B) as follows:
a (B) = 6i (B Xa.) 2b.
Fact H.2. Let X be Rdxd be invertible and symmetric. Then have
0X - 0x 0t
X-1
Ot
Lemma H.3. Let E Raxd be invertible and symmetric. Then, Ölog det X = -x.
Proof.
The proof is by chain rule:
i log det X
det x det X,
to compute det X, recall the adjugate of X, denoted by adj(X) where X-1 det x adj(X). and remains to show that 0 det X adj(X) this can derived using cofactor expansion of matrix, which states that fix row i, we have
det X = Σ Xi adj(X)
by product rule,
d X ax,s
0Xii ox,s
θ Dadj(X)a im
Σ
adj(X)
aXi,
46

---

adj(X)
where a = 0 for all & since the cofactor of i, k is the principal minor by removing and ax column k. Therefore,
det X = ad)(X)
and consequently,
i log del X
ddet X
det X
del X adj(X) =X =x
The goal of this section is to prove the d x d matrix F2F(2) > O. We start with computing Fx(x) E R
Lemma H.4. We define 82, = a, for each E [n]. Let S₂ denote a diagonal matrix such that
Sg = diag(8,). Let (x) = = S 2 E Rdxd, Let F(x) = log (e) + Ia), then
n
VF(x)
Σ
(S A)
Proof. First, we know that
0s 2
-201.j
(19)
0z, 52,,
For each j E [d], we can write Jz, as follows:
OF (r(( Ta) (zH(2) + Ia)
dz,
0zj
(r(((())) Ia) a1)
n
08 2
-r((([(x)) + 1a)
0:0
0æj
T 0411
= (r((H(()) + 1a)
20,a
1=1
as a (H (x) + La)
aq
= - 2
tr
i=1
aij S2,i
(S A)
where the forth step follows from Eq. (19), the fifth step follows from tr(AB) = tr(BA], and the last step follows from o.
We use chain rule:
OF
Dwa
DEA
0. A)
Sx,i
Fact H.5. We have the following partial derivative:
((
a)
z.
2
(S-
A) Ojd
dm
8
47

---

Proof. We analyze the partial derivative term
dm a (M(x) 1a) =a! (((()) + La) -1
①
)a
(x)))
=a, (-H(()) + Ia)
(H()) + Ia) -)a;
Ozy
(2(H(( 2 + 1a) 1a)-111 n
aga
aj
(H (z) Ia) 1a)=
j=1
( (H(r) + Ia) 1 0 )2 دnd
2
88
j=1
(5z A) O1
2
1.j
(20)
82.
j=1
where the second step follows from Fact H.2.
Fact H.6. We have
,, - 5, A) - Σ (S,
j=1
Proof. Note that
(H(x) + Ia) Ia) (H(=) (H) (z) + +1a
(1(()) + 1a) 1a)
+ (1()))
= (H(x) + 1.) 12)
22
j=l
therefore,
0, A) = a. (H(x) Ia) Ia)_'aa -'a;
a. (1) + 1a)
20.
(z) +1))+++0000 18) -1))=aaa) 2a,z
1=1
= 5. A)+ +
(S A)
(S-
Thus we complete the proof.
Lemma H.7. Let F(x) log det (H + Ia). then
V²F(x) = 0.
Thus, F(x) is convex.
Proof. Note that
0 F(x)= = Σ a T (H + -a)='a,
a;a; Om
i Oi.k
+
0i.k
0 (H())) + Ta) -a)=
dz, Ork
$ dm
i=1
3 a (H(x) 1a) a, d,
a,
a (H(x) + Ia) a.
(21)
$ dz
Plug in Eq. (20) into Eq. (21), we have
3 0S A)
2
Σ
~
SS A) Ci,kaj.
$2
S S S
i=
==
48

---

and
2
A)
0:0;
2
2 8-1
(5 A)
ajaj
2
82
i=1
1=1
j=1
3 i=1
(5 A) + Σ
(S
A))
0:0;
2
(ST A)
aaa
Su,j
j=1
40
40
a;a}
(S
A)
+ 3
(S
A)
2
S A)
2
8
i=1
f= j=1
t=1 j=1
Σ (SF+ A)
0:01
Σ
0. ( ST=
0:0; 8z
(ST
A).
aiay
a;a]
8
i=1
i11=1
== Bi B2 + B3.
where the first step follows from Fact H.6, the last step follows from
0:a 2 2
Bi = ( A)
10 B2
a;a; 8z1
(S A)
a;a} Sz,j S0,
B3 2
(S- A)
(
(=1 j=1
is obvious that B1 and B2 0. is not hard to see that
0g0
0ja Sng
0,a
n,
82,5
Thus we have Bg ž 0 and we complete the proof.
H.2 Convexity of log det (Hog (x) fa): spectrahedra
In this section, we prove that log det ( (Hng(2) + Ia) is convex function on æ. Combined with the bounded variance property, we provide algorithm that samples from the SDP spectrahedron.
Definition H.8. We define the
matrix g=² as
:= 15/1/2
1/2 53 1----115 1/2)ATTTTTA (H + Ia) Aa)- 1/2 S-1/2..
We define projection matrix P E R. IR as follows
: (S-1/2 @ (A-P))AA--ASS-11/2---- S-1/2).
Remark H.9. It is not hard to see that P is a PSD matrix. Note that the matrix H + Ia is PSD, since H. is PSD, therefore (H + -a) is also PSD. Thus, the ridge- cetje-proeection PSD, as its two "arms" " symmetric.
Lemma H.10. Let P. E xg" be the projection matrix and let P E xn? be the ridge-projection matrix defined in Def. H.8. Then, we have P.P.
Proof. Let H = T be its eigendecomposition, we need to show that
(H + -a) --1
Let us expand the LHS:
(UEUT -- (U(E I = UUE + ((E++)) -UT
3UN-'UT,
where the last step follows from for any non-negative number X, (1+ A) <X-1
49

---

The following lemma is generalization of Theorem 4.1 of Anstreicher [2000].
Lemma H.11. Let Hog (a) E Idddd be defined as Definition G.6. Let F F(z) = log det (oog (x) + Then, (x) is convex in x.
Proof. Using standard algebra and chain rule computations, we have
F²F(2) 2Q(x)) R(x) 2T(x).
For simplicity, we drop x. The R, E IR ddd can be defined as follows
Qi = A(H + La) A A [ A
R1g A(H + La) AA - - A5-
Ti == A(H + AAT (8-'A,,5 ' A(l + AA) (S- 05
Using Lemma H21 we know that
0 03)(2) 3 V2 3
Thus, F(x) is convex.
Lemma H.12. For any if (xx) > 0, then have
Part I. Q 2 0;
Part T > 0;
Part 3. T {Q + R);
Part 4. F² ž 0;
Part 5. R Q;
Part 6. F2F(x) = 3Q(())
Proof. Proof of Part and 2. Let E I2 and & & O.
Then we have
T = Liseie
is
A(H AAT -B-
= (B* 2 + Ia)
where B = B(E) = Σ A,A, and B = 5-1/2
Similarly, we have
RE A(H -S-
P. B B
We can rewrite g TTT as follows:
< <T AT AHH - 1-
= (B P
is obvious that Ia, and D are all PSD matrices.
Using Part of Lemma A.9 and Part 6 of Lemma A.10, we have
Q 2
T 0
50

---

Since & is arbitrary, then we know that
= 0
T = 0
Proof of Part Note that P is a projection matrix, and 3 P by Lemma H.10.
(Bel)) Ia) (B PB
= (
1 - + Ia) + (
= 2
where the last step follows from Part 2 of Lemma A.10. Applying Part 1 of Lemma A.9, we have
(P. ( e (BB & (P. @1)) (B
2
The above equation implies the following
E
< T 2 (Q R)
Note that, here & is arbitrary, thus we know that
(Q+
Proof of Part 4. We have
V2 =Q + - 2T
2 2Q + R 2 (Q+T)
Q
Proof Part Let u; orthonormal eigenvectors of with corresponding eigenvalues A
Then using [Horn and Johnson, 1990, Theorem 4.4.5), B + / has orthonormal eigenvectors i 8 Uj
with corresponding eigenvalues (A] )3), while (see Horn and Johnson (1990]) B 8 B has the same eigenvectors U; 69 Uj, with corresponding eigenvalues Az2z-
It then follows from (A, >2)² for each i,j that
B
Using Part of Lemma A.9, we know
PP 1a) > B)
Which is
E Q 2 & RE
Since it's arbitrary &, then we have
Q= R
Proof of Part We have
V²F(x) =2 + R - 2T
2Q + R
33
where the second step follows from T' > 0, the last step follows from Part 5.
51

---

H.3 Discussion
is well-known that for popular barrier functions, such as volumetric barrier and Lee-Sidford barrier [Lee and Sidford, 2014, 2019], the function log det (x) is convex in æ. However, given a PSD matrix function (x) for which log det H (x) is convex in a, is generally not true that log det (z) + Ta) is convex.
Fact H.13 (Folklore). Suppose that (x) positive semi-definite matrix for in the domain. Suppose log det H is convex. It is possible log det (H (z) + I) is not convex.
Proof. log (z) is convex in x, then for any fixed x E K, E R log (x) convex in t (for t sufficiently close to 0).
For sufficiently closely to 0, define ft) == (H + tu)). Then
log det (s (e) Σ log
iE(d]
and the condition that log det (z) is convex is equivalent to
0²
(t) - f!(1)²
log det H(x tv)
Σ
2 0.
ōt2
f:(t)²
iE|d]
Now
log (e tv) + Ia)
Σ log(fi(t) + 1). sE|d]
The condition that log (el((H()))) + Ia) is convex is equivalent to
8² dr3
(((A())) + 1) 51(1 (f((t) 1)2
log degdet((g [[ La)
2 0.
iEld
- Suppose d 22) = 1100 = - 4,0) = (0) = 3. Then
82
log det H(x + tv) = (1)
f1( + (2(
0.
0t
t=0
f2
0²
(0) log det (x + (w) + I) = f1(0) +1 f2 +1 0.
(0)
σt
This gives an example for which log det H(x) is convex in region but log det (x + Ia) is not.
I Convexity of Regularized Volumetric Barrier Function
In this section, we prove the convexity of log det of the regularized volumetric barrier function. This crucial, as the convexity proof for regularized Lee-Sidford barrier is identical up to replacing leverage score by Lewis weights.
I.1 Definitions
Definition I.1. Let A E Raxd, Let b E R", For each i E [n]], we define
81, (a, x by) S₂ = Ar b
For each i [ni], define
(A,)AA)) TT, Ax) -aa,i -
For each i for each e [n]
o.l(As) = T (A. Ax) - az.1 '
Let Σ (A,) In-
52

---

Definition 1.2. We define (()) E ddxd as follows
H(x) = A [(Ar) Ag
dxn nxn axd
Definition I.3. For each j [d], define (z(2) E R as follows
Pạ(œ) = aA{)
dr,
nxn
For each j E [d], for each le E [d], we define (z,((z) Raxx as follows:
(2)(As))
P3A(2)
nxn
In the previous sections, we mainly use notation xx (Az). For simplicity, from this section, we will use notation (1)) instead.
Definition 1.4. We define Q(x) E Ruxa as follows
Q(x) == Ox (Ax).
Fact 1.5, We have
(i,((Aa) (Ax); 1n)
Proof.
The proof is straightforward from definition.
We list number of standard calculus results.
Fact I.6. We define the following quantilies:
Let S. Ar R":
Let S diag(sz) Ruxa denote the diagonal matrix by putting su the diagonal;
Let Ax,j denote the j-th column of matrix A Rand
Let Ag = =S=1A:
Let aI, denote the i-th row of A. for each i E [n].
Then, we have for each j E [d]
Part 1.
Øsz
A.
Or
nx1
Part 2.
8
Azgj
0x
nxl
nxl
Paart 3.
08
-2 s A.
0x;
axl
ax1
53

---

Part 4.
-
2 2diag
Sx o An.))
- -S A)
dr,
axn
Part 5.
0A S
2 2 diag(-S; A
Oxj
nxn
Part 6.
0A A2 -2A__ diag( (oog(-Aspe Az dx,
Part 7.
A A,)
([TAA)) AT (a, (A) Az) Part 8. For each i [n]]
Dan dzz
Aus. azii scalar
Part 9, For each i [n]
dag -2 Azis ansaaa T
dz,
scalar
dxd
Part 10.
AArii
A2 A
dx,
scalar
scalar
Part 11.
dArks Ozk
Angi Aniip scalar
scalur
scalar
Part 12.
8A -A02
Part 13.
-Ag,s.s O Ag, d.s.a..
Oxk.
Part 14.
0A.
diag( An,as) Ag nxd
dx
nxn
nxd
54

---

Definition 1.7. We define Σ(Az) E R and (A2) E R as
Σ(Aa) = (A((H(()) + 1a)-1A}1 In
-., (Ax) = (x) + la) A
Definition I88. We define Fx) as follows:
F(x) = log(det H (z) + Ia))
Definition 1.9. We define
Raxd to be
dm dz
OF
0 OF
① OF
+
O Ox
Ox Ox
Ox dx
where
0 OF
(1)
Or. Or 1ij.k = +t[(((()) + 1a)
0r,
0 OF
(H(x)
(H(z)
Ör Or Tr 12. -t(((((r)) + 1a)
(H(r) + (a)
Ook
dx,
1.2
Gradient of otj
Fact 1.10 (First derivative of leverage score). We define the following quantities:
For each j [d], let Ax, I= denote j-th column of A;
Let Ag ST1150 A Raxd,
Let Ag,u. E R.a denote the j-th column of Ag Raax,
Let L(Ao) E paxx denote diagonal matrix where (i. entry is (A,((A))
Let (Az) = 0 E R"X";
Let Leee), B= denote column vector of on, (As) Raxn.
Then we have for each j E [d]
Part For each E [m]].
0a,,; 2(0.,(AA) Anns)) -00,;(A)) Arji. Or,
Part 2. For each E [n]) 1 E [n]
do,,((A)) = 2( (Ax) (l,, oc.(Ae) (Auss + Arss) dz;
Proof. We know
(o,,((x) dx,
da Ta A,) azzi -1ar.
Ōx,
aa J J, (AT ATTA))
Ox
a
0(A Ax)-
(A Aa)
(0 0x,
dx
For the first term in the above, we have
-2Am.i (a T (A] Ab) -2A -..((Ae)
55

---

For the second term, we have
(0 T 0(Ag Ax)
2(az,50 T (4, Az) -1 A Ag AA)))) A2)
dzz
-2 (Ag Az) -1 A diag( A (Az, Ar (A₂ Az) - ax,i
-2 (A An) )
Gs 1,,00,,A,, (A As) -'u.i
T1. (Ax) A 1=1
Az,,,).
Similarly, we can prove for 00.a 823
1.3 Gradient of o Lemma 1.11. Let a and E defined as Definition I.1. Then, we have
Part
(A,,(As)) 2 (Az) (on,s As,e.) -2 (0x,(Ar) Ag nx1
o
dx,
nXn
nx1
ax1
Part 2.
do (0..-(A)) dx,
2 20. (An) diag( Ar, (A.))) diag(Ar,eg) (A,((An) diag(An.)))
Part 3.
(E(Az) diag(o o0 (Ax) Ac,,.)) 2dagg(22 E(Ax) Ag,e.j
Ozg
nXn
nxl
2 diaag(a (A - E(Aa))
Proof. It follows from Fact I.10.
I.4 Gradient for (x:)
Lemma 1.12. Recall the definitions of the following quantities:
Let (xx) E Raxdd be defined us Definition
For each j [d], let P5(x) Raxn be defined as Definition L.3.
For each j [d].fo for each ls E [d], let (x)((r) Raxa be defined as Definition I.3.
Then, we have
Parr 1. For each j E [d]
0H())
2A diag ((aa(E(Az))))) (E(AA))
0a. + (z)
56

---

Part 2. For each j E for each &o E [d]
a ( + 6A. diag( { (2,,
0zk Ozj
-2A Pa diag( (AAA A
22A (a) diag (A,,,,))
+ (,,
Proof. Proof of Part 1. We have (H(x)
3(A ((A))A))
0zj
dx,
DAY -(OA)) 0zz
E(A,)A, + A (
AAz
AA + A E(A,)
dzj
0xz
-A, diag((z,, )) (A,,)
P(r)((
A,E(A,) (an,g(A,n,jj)A
= - 2AA diag E A2) A2,
+ Az)
Proof of Part 2. Then, have
0 1 (1(2))
0 diagAA + A (x))A))
Oak 0z,
00k
δ
(-2A E(A,) (A,g(A,,,))))) Ōak P)()A))
(22)
Dak
where the first step follows from Part 1, the second step follows from the sum nule.
For the first term of Eq. (22), we have
24, Σ diag( (aog((A,xx.)) AA.,.))AA))
0zk
④ (E(A)) (iag(Ar...))
- 2 (A [{,) diag( (n,, 2A
Oak
Osk
2.A -(A,)
(diag( Ag (ax, )) I(A, diag Ac,t,))
E(AA)
Ooo (Ax)
0rk
② (A (EA,)) diggAA))aa((AA (Ar,ag(A,,,)))))) 2A_
à (E(Az) diag) (Ar, + As
2 -2
Osk
Osk
E(Ax) diag( 0
0
dok (Au,.))) - ((A)) diag( (A,,,) dok
(Aa)
++ 2AA diag( (ong(AA,,,)))(((( diag(A (A,..))
-2 Pa (x) diag( A Ar
+ 2.4 Σ Az) diag( Az,,.) diag( A Az,a,k) Ar
24 E(A,) diag( Ag (,) diag( A b)
(23)
For the second term of Eq. (22), have
(A]
Ork
(AT)) + AT A (P)(())) +A]P((2) 8
(An)
Ozy
dzk
0r:
- (A,) -AAg)dia(AA...... Aia.(A,,,, + A. ())())). P,) (Au,,
57

---

(A diag(A..,.)PP + A P,a(x)A, A ding(A..A)P311)44
-((A)) ( P.
(24)
By combining Eq. (22), Eq. (23), and Eq. (24), we have
+ 24 diag( As, (A(a) diagí (Az,,,) Ax
Ozk dr,
24 Ps. (2x) diag( (Ax,-.)
+ 24 ((AA)) diag (Az,e.)) diag( A2,,,k) Az
+ 2 E(Ax) diag( (Az,, diag( (2,,, Az
2AA (3)
+ A A2)
= + ++6AJ diag(Az, ) diag(Az, (A,,,.)) Ar
2AA (x( diag A
2AA P3 (z) diag (A,,,
+ P.
Thus, we complete the proof.
1,5
Hessian of (x)
Lemma I.13. Recall the definitions of following quantities:
Let == A]] E Baxd be definition 1.2.
Q(x) J=,- (Ax) E R be defined as Definition 1.4.
Then, we have
Part 1. For each j E [d] and each k E [d]
82 H²)(()) = C. + C2 0rzzk
where define some local d size matrix variables
C1 +2 diag (
C2 = _
5),
C3 _.A__ diag(A,, diag @02 (x) A
C4 = _ diag(A,, diag( Q02 (x) A,, Ar
C. +84, ((Q(x) diag( Ag,, . )Q diag( A ( Qx()) A,
Proof. By Part 2 of Lemma 1.12, we can show that
8 Hx() Ōr,Ik
B2 + B3 + Ba
where
B 6A7 diag( Az, . ) E)((a)))))))))))))))) (iag(Ar,,)))
B₂ 2.A Px(x)) diag Ar...) As
B, = -2A (2) diag( A2
B₃ = + (r))
where B2 can be decomposed further
Bg = 2
_ 2 diag((Q)²(z))
E(A.)) (iiag(Ar,)) A
58

---

44g diag(Q (z)
44 (diag( Q (2) (2) diag (I (A x))) diag As, ,
44 diag(Q (2) Az, diag(As, {) 4 diag ((Aa) (Ae, A)
= -44_ diag(Q" (z) Ag,) (i.,,)A.+.,)444 + 4 diag(As,* x)), ( (aagg(A,
where the second step follows from simple algebra, the third step follows from the definition of diag(-), the fourth step follows from simple algebra, and the last step follows from the definition of diag(-). For B3, consider the following:
B3 2 diag((Q (z) E(,,)) diag(Ae,
= 4 diag( (z) - E(A, Az, , diag (A,,,,,)) AA
= 4 (diag (x) - diag ( Ag Ae, ,)) diag( A
_ (iag(("" (z) diag(Ar, A 4 diag (1A)) Az, diag (A,, k))a
-A_ diag(Q (2) A diag(Ae, A² 4 diag(A (aa.g(Aa+.))) ((A()) (Aag((A,, (A,,,,,)) Ag
where the second step follows from simple algebra, the third step follows from the definition of diag(-), , the fourth step follows from simple algebra, and the last step follows from the definition of diag(-).
Finally for B4.
B4 AI( +
I.
diaag (2 (Ag,, o Az,z,))
4 diag (x) An,u.) diag Az,ac))
diag((?" (x) Ae,,, diag(Az,t.) )
+ (A,) diag (A,.. diag( (2,,)))
Eventually, we can show Hessian is C1 .... + Cs. The reason is following:
For the term C1, have
C BI+ B2,2 + B3.2 + B4,5
= (6 + +4+4++0) 6) A diag(Ar,, ()) (E(A)) diag(AA,..))
= 20 AT diag( Az,.. . E( ((A.)) diag( An...)
where B2,2 is the second term of B2, B3,2 is the second term of B3 and B1.5 is the last term of B For the term Cz have
C2 B4.2
For the term C3, we have
C3 B3,1 + Ba33
(4 + 4) A diag(Ar,,.)) diag(Q (er)()
= 8 A diag(Ar,,,)) ( (z),),,)))
For the term Ca, we have
C₁ = B2,1 + Bi,
= _ diag(Q ( (A,, 4A diag(Q Ag k)) diag(Ar, (
SA,__ d (x) diag (Ar, A
For the term Cs. we have
C5 B1.1
59

---

1.6 Gradient and Hessian of F(x)
Lemma L14. Let F(x)) E R be defined as Definition 1.8. Then, we have
Part 1.
@F(=) (r((((()) + Ia) (z dx,
dz,
PPart 2
Ø (x))) + (r((H()) Ia) Ia)-
8
Ork dxy
Dr,azi
0H (x)
0H (x) 0zz
(x)( + 1a)
Proof. Proof of Part 1. We have
(F(x)
0 llgg(do (x) Ia))
0x,
Oz,
tr((H(2) + Id) 0 (H(()) Ta) dz,
tr((H(z) d) -1 (H()) + Old
dz, }
dx,
(x)(H(()) + Ia) -1 (H)))
where the first step follows from the definition of F(x) (see Definition 1.8), the third step follows from the basic calculus rule, and the last step follows from 014 8xy
0.
Proof of Part 2. We have
δ (z)))
dr dx,
(tr((H(z) + Ia)-1 (1((W)
Dak
dr;
0 ((H(x) Ta) Ta)- HH(x)
(H(( Ta) a (H(z)
tr(
Ork
Ox,
Oz
0x,
0 (H (H)) -1 0H(x),
tr (H (x) + -1 á a1l()))
=tr(-(H(z) + Ia)-1
dz
0z
Oxr 0x,
0
0H
x
(x()) + Ta) 0H(x),
= + + 12)
1 tr( + Ta)
Oz, where the first step follows from Part 1, the second step follows from the product rule, and the last step follows from simple algebra.
0x,0zz
Ozk
Using standard algebraic tools in literature [Lee and Sidford, 2019], we can show that, Lemma I.15. We can rewrite as follows
i OF
Òr 20H 6H2 8H3 8H, 8Hs
4G 8G2 8G3 16G4
and
OF Ör dr > 0.
Thus, the function F is convex in
60

---

1.7
Quantity I for Hessian of F
Lemma I.6.. We have
0 aF
20H1 6H2 8H3 8H4 + 8Hs
dr Ör
where
Hi = ((A))
H2 A. diag( Q ( ((.)))
H3 = AT (x)
НW (1)
H₁ AI Q(x))A
Proof. It comes from P.j and the following lemma.
Lemma L.1.. Let H1 Rdxd be Hessian where each entry
(H1) = (x) + 1a 'A diag(Aa,,.s)) diag( (t),
Then, we have
H ((4,)
Proof. We can rewrite (Hi)j, as follows
(Hi), (x) La) A diag( (Ae, diag ( = tr(A.(H( + La) "A_ diag(Au. )) (Ap) diag(An, (A),
(Ap) diag( AA. . (EA)) diag (A,,,,))
=A,,, (A,,
where the first step follows from the Lemma statement, the second step follows from the fact that all of (x) + La) (Aa,,(A,...)). (AA)), and (e, are diagonal matrices, the third step follows from the definition of -,,. (see Definition I.7), and the last step follows from the definition of Σ(Aa)
(see Definition 1.7).
Thus, have
H (AAA,,A)AA
Lemma I.18. Let H2 be defined (H2), tr((H 1a) 11)-1 A. diag (iag(Q) Ar,ss))AN)
Then, we have
H2 A (z)
Proof. We have
(H2), = tr((H(z) + Ia) A diag ( Ar,
tr(E(Am) diag Q 2 (A Asx,))
1TA(( (z) (A,, Az,
" (z))
where the first step follows from the Lemma statement, the second step follows from the definition of ((A,)) (see Definition 1.7), the third step follows from the definition of 1, and the last step follows from the definition of E(A)).
Thus, we have
Hg A (z)g()
61

---

Lemma 1.19. Let H3 be defined as
= tr (t)) - -A diag(Az, (1
Then, we
H3 (=)
Proof. Then, we have
(Ha)j.a = (z((H(x)) + Ia) J diag(A tr((((A,) diag(Ar, A diag Q (A
where the first step follows from the Lemma statement, the second step follows from the definition of E(Aa) (see Definition 1.7), and the last step follows from the property of diag(-).
Thus, we have
H3 = 2
Lemma 1.20. Let H4 be defined
= ()(()))) A diag(A A2) (1
Then, we
H4 = (A
Proof. This proof is very similar to Lemma 1.19, so we omit the details here.
Lemma 1.21. Let Hs be defined as
(Hs),,s 1 (r((H(())++)))))))))))))(()))) Ia) A (Q()))inAAA (4a A.)
Then, we have
H3 = )
Proof. We have
(Hs)j,a = -a) diag( QQ (Q) o A.)
(r({(A,)(Q(z)) diag( A ) diag( A ((T))
Σ (.. diag( A
AT (hac::/adhilllaoo)o))oaoeeoeCCC
AT Eu((A,)(Q..(0))...i. Q(z)))
= A A
Q()))
where the first step follows from the Lemma statement, the second step follows from the definition of Σ(Az) (see Definition L.7), the third step follows from the definition of tr|-], and the last step follows from the fact that E(Az) is diagonal matrix (see Definition I.7).
Thus, we have
H₁ A (Q(=))((.))) Q(=))A,
62

---

L8
Quantity II for Hessian of F
Recall that
0H
2AA diag (1(Aa) (e, Az
Ōzj
+ A (AT))))
Ar
Ox,
2A_ diag (E(Ag)
+ 2AA diag (QQ (I) - A
A diag(2Q (x) - (2x(A))) (A
Lemma 1.22. For the second ilem, we have
② 0F
4G1 + 8G2 8G3 16G4
0z dr
G = (10""
G2 4109
G2 = 444,
G4 (4)
Proof. The proof follows by combining the following lemmas.
Lemma 1.23. Let (G1),,t
(G1); tr((H + Ia)- -1A A (z) + Ia)- diag(Q
Then, we have
G = (==
Proof. We can show
(G1) tr[(H + Ia) diag (z)
A (T)
tr(E(A,) diag(( Ar, Ad Iα) AJ ( An,,k,
tr()(A,) diag 8282 A (AA)) diag(( (2) Ao,,,,
-J. Q (x) (A) Q² (x) Am,,k
Thus, we have
G =
The proofs of G2, G3, GA are similar to G1. So we omit the details here. Lemma I.24. Let (G2),,k
(G2), tr((H -1A- diag (2) -a) A.
Then, we have
G2 = (==44""(334)))3(444441
Lemma 1.25. Let (Ga)j.a
(Ga)j, ur((H + Ia)-A (E(A)) Ia) ' AT (x)
Then, we have
G3 (3=4
Lemma 1.26. Let (G4)j.k (Ga)a. trl(H + Ia) A diag A4. + 1a) -a "A ) A.2
Then, we have
G 111
63

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
The main result of this paper is a class of Dikin walks for log-concave sampling from convex bodies with self-concordant harrier functions. The abstract and introduction reflect the complexity of these walks, and the remainder of the paper is dedicated to prove these assertions.
Justification:
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation as long as is clear that these goals attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification:
The limitations discussed in the last paragraph of Section 4.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification. asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full of assumptions and complete (and correct) proof?
64

---

Answer: [Yes]
Justification: The assumptions are clearly stated in the statement the theorems: convex body contained in ball of radius R. density is log- Lipschitt with parameter L, the walk starts from a w-warm starting point. See Theorem 1.1, 1.2 and 1.3. The proofs of Theo- rem 1.1 and 1.2 results can be found in Section C and D for mixing time, and Section F for per iteration cost. Proofs of Theorem 1.3 can be found in Section G. Guidelines:
The answer NA means that the paper does not include theoretical results. All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that it affects the main claims and/or conclusions the paper (regardless of whether the code and data are provided or not)?
Answer: [NA]
Justification:
The results in this paper are theoretical with no experiments.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurlPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (c.g., to registered users), but it should be possible for other researchers to have some path to reproducing verifying the results.
65

---

5. Open access to data and code
Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, described in supplemental material?
Question:
Answer: [NA]
Justification:
The results this paper are theoretical with no experiments.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (uttps://nips..
wwtpomiissssnoooliooooooomoooeee onPolicy for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https: //nips.
httpa///ggidd////eeuubbibiiiiiii PoodeSummissionPrr
for more details.
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [NA]
Justification:
The results in this paper are theoretical with no experiments.
Guidelines:
The answer means that the paper does not include experiments.
The experimental setting should be presented the core the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment Statistical Significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer:
[NA]
Justification:
The results in this paper are theoretical with no experiments.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi-
dence intervals, statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall with given experimental conditions).
66

---

The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful figures symmetric error bars that would yield results that are out of range (e.g. negative rates).
show in tables or
error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on the puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [NA]
Justification:
The results in this paper are theoretical with no experiments.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary failed experiments that NeurIPS Code of Ethics https htt////eurrpp..peee httppp..c/publll cctppblic//eeeeeecceeceecc htt.hii//ehhisnieeeeiiiileeeeeee
didn't It make it into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the
Answer:
[Yes]
Justification: principles.
We have carefully read Code of Ethics and the paper is adhered to those
Guidelines:
The answer NA means that the authors have reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive socictal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification:
The potential societal impact is discussed in the last paragraph of Section 4.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
67

---

Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper deseribe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: The results in this paper are thcoretical with no experiments.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [NA]
Justification: The results in this paper are theoretical with no experiments.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website). the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms use in package should be provided. For popular datasets, paperswi twwode.hco....c. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
68

---

For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
The results in this paper are theoretical with experiments.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details the dataset/code/moddl as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsoureing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
The results in this paper are theoretical with no experiments.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included the main paper.
According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The results in this paper are theoretical with experiments.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
69


Reviews:
[{"review_id": "a0UZPZuAap", "reviewer": "Reviewer_zXRR", "summary": "In this paper, the authors give improved mixing times for a Markov chain whose goal is to approximately sample from distributions whose densities are proportional to $\\exp(-f(x))$, where $f$ is $L$-Lipschitz and convex, restricted to a convex set defined by the intersection of $n$ hyperplanes that lies in a Euclidean ball of radius $R$. The Markov chain itself is an approximate variant of the Dikin walk that tolerates using spectral approximations to the Hessian instead of the exact Hessian for drawing the next sample.\n\nThe main result is as follows. Given a subroutine to $1/d$-spectrally approximate the Hessian of the barrier that runs in time $C_g$ and given that the barrier itself is $\\nu$-self concordant, the authors give a Markov chain whose iteration complexity to converge to a $\\delta$-TV-approximate distribution to the target is $\\widetilde{O}(\\nu d + d L^2 R^2)\\log(1/\\delta)$. Each step has runtime $C_g + d^{\\omega}$. \n\nThis general framework recovers some of the best known mixing time bounds for special cases of the studied problem, including sampling from the uniform distribution over a convex body. Moreover, for the setting considered, as far as I can tell, this is the first guarantee whose iteration complexity does not depend on the number of constraints defining the polytope (this seems to get replaced by the self-concordance parameter of the barrier, which can be made to be $\\widetilde{O}(d)$ using more sophisticated barriers than the log-barrier).\n\nThe ideas also transfer to giving mixing times for the problem of sampling uniformly from the constraint set of a covering SDP. Again, a key part of the contribution is that the walk can tolerate using spectral approximations to the Hessian of the barrier (in this case, the log det barrier) instead of requiring it exactly. \n\nThe main technical insight is that there are a few properties that one needs the barrier to satisfy in order to give the $O(\\nu d)$ dependence in the mixing time. They consist of a symmetry assumption, a bounded local norm condition, and that the barrier is convex under regularizing it with an identity matrix (see page 6). Under these assumptions, the authors prove that their variant of the Dikin walk satisfies the promised mixing time (the argument is executed in Appendices B, C, D).", "questions": "Can one obtain better mixing times under more structural assumptions on the convex set (e.g. symmetry, some notion of uniform convexity, etc)? Or is this an open question?\n\nDo you think it is possible to get runtime improvements by intelligently recycling Hessian computations? Then, the amortized runtime of finding each approximate Hessian could improve over the stated $C_g + d^{\\omega}$. In particular, self concordance implies some kind of Hessian stability (for any two nearby points, the Hessians at those points are spectrally close) -- so, could you hope that the Hessians of the barrier aren't changing too much between two consecutive steps?", "limitations": "Yes", "rating": 8, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "The problem addressed in this work is an important one to advance our understanding of the geometry of polytopes. The conditions under which mixing times can be obtained are pretty general, and as the authors show, only depend on certain natural properties of the barrier function that is being used. The results also imply significant quantitative improvements over prior work (in particular the result about uniformly sampling from polytopes with an iteration complexity that is independent of the number of constraints when $n$ is much larger than $d$).", "weaknesses": "There are a few questions I have (see below).\n\nBesides that, a very minor comment -- in Section 2.3, the Lewis weight optimization program as written is not actually convex. But, there is a convex formulation for all $p > 0$ given in [LS19].", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We deeply appreciate your thoughtful and positive feedback for our work. We thank you for pointing out the convexity of the Lewis weights program for all $p>0$ as indicated in [LS19], and this is particularly the case for our application when $p$ is $\\mathrm{poly}\\log n$. For your questions,\n\n**Q**: Could one hope to obtain a better mixing time bound on more structural assumptions on the convex sets?\n\n**A**: To the best of our knowledge, the only work that explores this frontier is due to Kook and Vempala [KV24], where they showed that if $f$ satisfies the so-called relatively strongly-convex property, then the mixing time could be independent of $R$.\n\n**Q**: Do you think it’s possible to get further runtime improvements by smartly reusing Hessian computations?\n\n**A**: Yes, we do believe that it’s possible to further reduce per iteration cost by developing a linear algebraic data structure to maintain the inverse Hessian matrix. Similar to many IPM-based algorithms for optimization [CLS19, LSZ19, JSWZ21], since the algorithm takes many steps, the relative change between iterations is small, hence it’s possible to design data structures to amortize the cost. We believe this is an important next step for this line of research, as our work has essentially obtained a good (up to dependence on $L$ and $R$) mixing time.\n\n**References**:\n\n[LS19] Yin Tat Lee and Aaron Sidford. Solving Linear Programs with Sqrt(rank) Linear System Solves. 2019.\n\n[KV24] Yunbum Kook and Santosh Vempala. Gaussian Cooling and Dikin Walks: The Interior-Point Method for Logconcave Sampling. COLT 2024.\n\n[CLS19] Michael Cohen, Yin Tat Lee and Zhao Song. Solving Linear Programs in the Current Matrix Multiplication Time. STOC 2019.\n\n[LSZ19] Yin Tat Lee, Zhao Song and Quiyi Zhang. Solving Empirical Risk Minimization in the Current Matrix Multiplication Time. COLT 2019.\n\n[JSWZ21] Shunhua Jiang, Zhao Song, Omri Weinstein and Hengjie Zhang. Faster Dynamic Matrix Inverse for Faster LPs. STOC 2021.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "8TPvYMzujh", "reviewer": "Reviewer_Cu4d", "summary": "The paper gives a Markov chain Monte Carlo algorithm for estimating the probability of a high dimensional polytope under a log-concave probability distribution. The algorithm improves the mixing time of previous results, while maintaining the best known per-iteration cost. More specifically, in conjunction with a known self-concordant barrier function, they use a cheap spectral approximation of the Hessian matrix required in an iteration, and they modify the Markov chain to converge more quickly despite the approximation.  A second result applies similar ideas to spectrahedrons, and there the improvement over past results is both the asymptotic number of iterations and the asymptotic computational complexity of each iteration.", "questions": "None.", "limitations": "None.", "rating": 7, "confidence": 1, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "This is an important and fundamental problem that has been researched extensively recently, and the result improves upon the best known algorithms.", "weaknesses": "The improvement is somewhat incremental, but nonetheless significant.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank you for your review. We respectfully disagree with your point that the improvement we obtain is incremental, as the prior best algorithms for $d$-dimensional polytopes with $n$ hyperplanes mixes in $nd$ steps, while ours mixes in $d^2$ steps. For the regime where $n\\gg d$, this is a huge improvement, and it is in fact the best known mixing time for uniform sampling over polytopes. Moreover, our framework extends to general convex bodies with well-behaved self-concordant barrier functions. Prior works [MV23] could only handle polytopes.\n\n**Reference**\n\n[MV23] Oren Mangoubi, Nisheeth K. Vishnoi. Sampling from Structured Log-Concave Distributions via a Soft-Threshold Dikin Walk. NeurIPS 2023.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "HzMrYSSl8q", "reviewer": "Reviewer_aoiW", "summary": "The authors propose a sampling framework from a convex body with a barrier.", "questions": "Page 1:\nCan you please motivate the problem a bit better? For a broader CS audience? :)\n\nPage 2:\nCan you please give some intuition about the benefit of using the Hessian in the design of your walk?\n\nPage 3:\nCan you provide further details about the hardness of constructing the Hessian of the universal barrier?\n\nPage 4:\nI do not understand Lines 154 -- 156.\n\nPage 5:\nLine 15:\nWhy is this probability 0.5 \\min\\{\\tau, 1\\}?\n\nPage 6:\nLines 196, 202:\nCan you please give some intuition about the symmetry and convexity properties? :)\n\nPage 7:\nLine 240:\nI do not understand this math display.\n\nCan you please elaborate on Equation (1)?\n\nPage 8:\nLines 285 -- 288:\nI do not understand this sentence :)\n\nCan you please give some details about TensorSRHT?\n\nPage 9:\nCan you please provide a more extensive future work? Thanks!", "limitations": "None.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "The main strength of the paper lies at the novelty and the significance of the results, especially that of Theorem 1.1.", "weaknesses": "See questions below.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate your thoughtful comments and questions, and we would be happy to answer your questions.\n\n**Q**: Can you motivate the problem a bit better, for a broader CS audience?\n\n**A**: Yes. Sampling is a fundamental problem in computer science and machine learning. Consider any constrained optimization problem where the goal is to optimize a function over a constrained set, to run any kind of first-order or second-order algorithm, one has to first construct a good *initial point*. This usually boils down to sample from the constrained set according to particular distributions, which is the central topic studied in this paper. In recent years due to the surge of differentially private machine learning, one could also implement the *exponential mechanism* using log-concave sampling from convex bodies.\n\n**Q**: Can you give some intuition on the benefits of using Hessian in your walk?\n\n**A**: In short, the walk we are using is a variant of ball walk, which could be succinctly described as at the current point $x$, randomly pick a point $y$ in the ball centered at $x$ with radius $\\delta$, and move to $y$ with correct probability. This idea however, does not work so well if the convex body is, say, a long and flat ellipse: we are constrained to pick the ball radius $\\delta$ small, as otherwise we could walk outside of the convex body, hence the walk takes many steps to mix. Computing the Hessian counters the issue by instead computing a maximal volume ellipse around $x$: this enables us to sample points with different scales in different directions, and hence enabling our walk to mix faster.\n\n**Q**: Can you provide further details of the hardness of constructing the Hessian of universal barrier?\n\n**A**: Given a $d$-dimensional convex body $K$, the universal barrier is defined as $\\log\\mathrm{vol}_d (K^\\circ(x))$ where $\\mathrm{vol}_d$ is the $d$-dimensional volume and $K^\\circ(x)$ is the polar set of $K$ with respect to $x$: $K^\\circ(x)=\\\\{ y\\in \\mathbb{R^d}: y^\\top (z-x)\\leq 1, \\forall z\\in K\\\\}$. The main computational obstacle comes from even implementing the zeroth order oracle for the barrier function: estimating the volume of a high dimensional convex body is hard, no deterministic algorithm could estimate the volume to a factor of $(\\frac{d}{\\log d})^d$ [BF86], and randomized algorithms are essentially MCMC types that require generating samples from the convex body. See also [NN94] for a more detailed exposition.\n\n**Q**: I do not understand line 154 - 156. \n\n**A**: Line 154 - 156 says that the cost of computing the Hessian of log-barrier for a spectrahedra is $O(dn^\\omega+n^2 d^{\\omega-1})$, if we consider the regime $n\\geq d$, then the first term $dn^\\omega$ dominates since it has a larger exponent on $n$. If $n\\geq d^{\\frac{3\\omega-4}{\\omega-2}}$ which is roughly $n\\geq d^{8.4}$, then our algorithm is faster than the $dn^\\omega$ time needed by exact computation.\n\n**Q**: Line 15: why is this probability $0.5 \\min\\\\{\\tau, 1 \\\\}$? \n\n**A**: Because we are employing a lazy walk that only moves to the next point with half probability.\n\n**Q**: Can you give more intuition on the symmetry and convex properties?\n\n**A**: The symmetry property is very helpful when one tries to relate the self-concordance parameter with cross ratio distance. Intuitively, we prove that for all barriers of concern, they are also $\\sqrt{\\nu}$-symmetry, and we further prove that the cross ratio distance between any two interior points $u$ and $v$ is at least $\\frac{1}{\\sqrt{\\nu d+dL^2R^2}}$, this establishes our mixing time. The convexity of the regularized barrier is crucial when we try to bound the discrepancy between the log determinants of the current point and the point we are trying to move to. Recall that we define the point to move to is $z=x+\\Phi(x)^{-1/2} \\psi$ where $\\psi$ is an i.i.d. Gaussian vector. Given the convexity, we could bound the log determinant as\n$\\log \\det \\Phi(z) - \\log \\det \\Phi(x) \\geq (z-x)^\\top \\nabla \\log \\det \\Phi(x)$, and the RHS term is precisely the bounded local norm which we could control.\n\n**Q**: Line 240: I do not understand the math display. Can you elaborate more on Equation (1)?\n\n**A**: The point of the math display of line 240 is to show that the probability could be decomposed into a determinant term and a term depending on the local norm of $x$ and $z$, and our subsequent analysis shows how to bound these two terms. For Equation (1), it is the convex program associated with $\\ell_p$ Lewis weights. It involves finding a maximum volume ellipsoid subject to the (soft) constraint that the polytope is contained in the ellipsoid. Here, by ``soft’’, we mean that instead of forcing $a_i^\\top M a_i\\leq 1$ for all $i$ which corresponds to the $\\ell_\\infty$ norm, we instead consider an $\\ell_p$ norm constraint, which could be written as $\\sum_{i=1}^n (a_i^\\top Ma_i)^{p/2}\\leq d$.\n\n**Q**: Line 285 - 288: I do not understand this sentence.\n\n**A**: This sentence compares our sampling algorithm with the IPM-based algorithms for semidefinite optimization [JKLPS20, HJSTZ22]. In those algorithms, the Hessian matrix could be carefully maintained and controlled, further improving the per iteration cost of the algorithm. For us however, while it might still be possible to develop a similar framework to reduce the per iteration cost, the picture is much less clear, and the strategy we employ is a simple method that does not require any sophisticated maintenance of the Hessian matrix.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Meta-review: I would like to thank both the authors and reviewers for the good discussion. The reviewers broadly agree that the paper is interesting and should be accepted. I concur and recommend acceptance.
Decision: Accept (poster)