Faster Sampling from Log-Concave Densities over Polytopes via Efficient Linear Solvers

1. Differences between our randomized estimator for the determinant term, and the estimator of [LLV20]:

Approach of [LLV20]: In [LLV20], the determinant ratio is first estimated by plugging in independent random estimates for the log-determinant with mean $\mathrm{logdet} (H(\theta))- \mathrm{logdet}(H(z))$ into the Taylor expansion for the exponential function $e^t$. This gives an estimate $V$ with mean $E[V] = \frac{\mathrm{det} (H(\theta))}{\mathrm{det}(H(z))}$ (their Lemma 4.3). Roughly speaking, the authors then suggest plugging this estimate $V$ (and another estimate $\hat{V}$ with mean $\frac{\mathrm{det} (H(z))}{\mathrm{det}(H(\theta))}$) into a smooth function $g(t,s)$, to compute an estimate for their Metropolis acceptance rule $\frac{p(z\rightarrow \theta)}{ p(z\rightarrow \theta) + p(\theta\rightarrow z)} = g(\frac{\mathrm{det} (H(\theta))}{\mathrm{det}(H(z))}, \frac{\mathrm{det} (H(z))}{\mathrm{det}(H(\theta))})$. Each time their Markov chain proposes an update $\theta \rightarrow z$, it is accepted with probability $g(V, \hat{V})$.

However, the randomized estimator $g(V, \hat{V})$ used in [LLV20] is not an unbiased estimator for their Metropolis rule $g(\frac{\mathrm{det} (H(\theta))}{\mathrm{det}(H(z))}, \frac{\mathrm{det} (H(z))}{\mathrm{det}(H(\theta))})$, even though $E[V] = \frac{\mathrm{det}(H(\theta))}{\mathrm{det}(H(z))}$ and $E[\hat{V}] = \frac{\mathrm{det} (H(z))}{\mathrm{det}(H(\theta)))}$. This is because, in general, $E[g(V, \hat{V})] \neq g(E[V], E[\hat{V}])$ since $g(V, \hat{V})$ is not a linear function in $V$ or $\hat{V}$. (One can construct a counter-example where the above equality does not hold.) This causes their Markov chain to generate samples from a distribution not equal to the (uniform) target distribution considered in their paper. This problem does not seem to be easily fixable, and requires a different approach.

Our approach to estimating the determinantal term: To obtain an unbiased estimator for the determinantal term, our algorithm bypasses the use of the Taylor expansion for the exponential function $e^t$ used in [LLV20], and instead uses two different infinite series expansions for the sigmoid function (a Taylor series expansion $\mathrm{sigmoid}(t) = \sum_{i=0}^\infty c_i t^i$ centered at $0$ with region of convergence $(-\pi, \pi)$, and a series expansion "at $+ \infty$" with region of convergence $(0, +\infty)$, which is a polynomial in $e^{-t}$) to compute an unbiased estimate with mean equal to a (different) Metropolis acceptance rule. Our Metropolis rule is proportional to $\mathrm{sigmoid}\left(\frac{1}{2}\log \mathrm{det}(\Phi(z)) - \frac{1}{2}\log \mathrm{det}(\Phi(\theta))\right).$ To generate a randomized estimator for this acceptance rule, we compute independent samples $Y_1,\cdots, Y_n$ with mean $\log \mathrm{det} \Phi(\theta) - \log \mathrm{det} \Phi(z)$, and plug these samples into one of the two series expansions, e.g. $\sum_{i=0}^N c_i Y_1…Y_i$, for some small number of terms $N$. To select which series expansion to use, our algorithm chooses a series expansion whose region of convergence contains a ball of radius 1/8 centered at $Y_1$.

Since $Y_1,…,Y_N$ are independent, roughly speaking, the mean of our estimator is proportional to $E[\sum_{i=0}^N c_i Y_1 \cdots Y_i] = \sum_{i=0}^N c_i E[Y_1]\cdots E[Y_i] = \sum_{i=0}^N c_i (\log \mathrm{det} \Phi(\theta) - \log \mathrm{det} \Phi(z))^i$ (or to a simmilar quantity if the other series expansion is selected). We show that, w.h.p., the choice of series expansion selected by our algorithm converges exponentially fast in $N$ when $Y_1,…,Y_N$ are plugged into this expansion (see our proof overview). This implies that, if we choose roughly $N = \log\frac{1}{\gamma}$, our estimator is a (nearly) unbiased estimator with mean within any desired error $\gamma>0$ of the correct value $\mathrm{sigmoid}\left(\frac{1}{2}\log \mathrm{det}(\Phi(z)) - \frac{1}{2}\log \mathrm{det}(\Phi(\theta))\right)$. We then show that this error $\gamma$ is sufficient to guarantee that our Markov chain samples within total variation error $O(\delta)$ from the correct target distribution, if one chooses $\gamma= O(\mathrm{poly}(\delta, 1/d, 1/L))$.

Thank you for your response. I agree with the first assessment, the [LLV20] acceptance rule likely has some issues even for uniform case. A fix for that is appreciated. However, this does not address my major concern of the paper --- from my point of view, this paper pretty straightforwardly combines the last section of [LLV20] and the soft-threshold Dikin walk in [MV23]. I strongly believe to justify the acceptance of this paper, authors should try to add more results in addition to the current $A+B$ result. For example, can you extend soft-threshold Dikin walk to other barriers for a better mixing? For more intricate barriers such as hybrid barrier or Lee-Sidford barrier, can you analogously design this [LLV20]-like fast solvers for them? At its current form, I will keep my score as is.

P.S. (Not for this version), in the next revision, consider comparing with this paper: [KV23]. Authors show that it's possible to do log-concave sampling over polytopes using other barriers such as Lee-Sidford, get better mixing and without dependence on $R$. Moreover, they only require Lipschitz on the density instead of $f$. On the other hand, they do need $f$ to be $\alpha$-relatively strong convex. I think no polynomial dependence on $R$ is essential for Dikin walk, as this (used to be) a key advantage over hit-and-run. Can you extend the [LLV20] machinery to their Dikin walk?

[KV23] Y. Kook and S. Vempala. Gaussian Cooling and Dikin Walks: The Interior-Point Method for Logconcave Sampling. 2023.

Thank you for your response and for contextualizing your work (particularly in your responses to other reviews). Based on your responses, I'm happy to increase my score and confidence. All the best!