Venue: NeurIPS 2022 Conference

Reviews:

final_score: rating: 8: Strong Accept: Technically strong paper, with novel ideas, excellent impact on at least one area, or high-to-excellent impact on multiple areas, with excellent evaluation, resources, and reproducibility, and no unaddressed ethical considerations.
review_content: summary: This paper gives a method for sampling log-concave distributions supported on compact convex sets, based on Riemannian Hamiltonian Monte Carlo (RHMC). Equality constraints are enforced via a Lagrange multiplier and inequality constraints are enforced via a carefully chosen Riemannian metric constructed from a self-concordant barrier function. Algorithms for exploiting sparsity for matrix inversion enable sampling from high-dimensional distributions. Tests show substantial improvement in sampling times in comparison to existing constrained samplers. Formal guarantees are provided for convergence of the associated integration scheme the convergence to the stationary distribution is bounded.
strengths_and_weaknesses: Strengths:
This is a nice paper. It gives strong practical results for an important, but largely overlooked problem and gives a thorough theoretical analysis. The paper is very well-organized, with helpful guides to the structure of the appendix. 

Weaknesses:
Probably the biggest weakness is that some of the assumptions are not clear from the problem setup. In particular, it is just stated that the constraint set, $\mathcal{K}$, must be convex. However, it really must be compact. Many of the arguments implicitly depend on this. These restrictions should definitely be clarified in the problem setup, and ideally also in the abstract.

If $\mathcal{K}$ is not compact, the cross-ratio distance from Definition 5 is not well-defined. It might be possible to extend this definition to the case in which $\mathcal{K}$ does not contain a line by using appropriate limits. However, in this case the constant $G$ from Lemma 17 could be zero. (Consider the case of the non-negative orthant in $\mathbb{R}^2$.)

If $\mathcal{K}$ is not compact and the $f$ is linear, then the corresponding density may not exist, rendering arguments about stationarity from B.3 irrelevant.

There are also some small typos and unclear arguments. Unfortunately, I did not have the time to step through every  argument in detail. But, I did catch a few issues:

- In line 575, it should read $MNM=M$
- In many places, "dynamic" is used where "dynamics" would be more common. For example, in 603 "the HMC dynamic" would normally be "the HMC dynamics"
- Similarly, "the dynamic on x" would normally be "the dynamics of x", in Lemma 4.
- $P$ gets redefined in line 601, and clashes with the definition from Lemma 2 in a confusing way.
- In line 798, we should have $\|x\|_A = \sqrt{x^\top A x}$
- Definition 5 implies that $\mathcal{K}$ must be compact. This should be stated more explicitly in the problem setup.   
- I think there is a typo in the first inequality after line 967. I feel like this is based on a lower bound on the right side of the equation after line 964. So, the second side should have an integrand of $\mathcal{T}_x(S^c)$.
- I think the factor on the second inequality below line 967 should be $\frac{1}{40}$ and not $2/5$. This is because the bound seems to follow from lower bounding the integrands by $0.05=\frac{1}{20}$ for both terms. I don't think this really matters for the rest of the proof. 
- The last inequality below line 968 should be strict.

questions: Would it be possible to extend this to noncompact sets? How about non-convex sets and non-log-concave distributions? What are the challenges?
limitations: The authors address many technical limitations, but as described above, some fairly serious technical limitations are not adequately discussed without digging into the technicalities. 

The authors indicate not potential negative societal impact, and I agree.
review_title: null

final_score: rating: 8: Strong Accept: Technically strong paper, with novel ideas, excellent impact on at least one area, or high-to-excellent impact on multiple areas, with excellent evaluation, resources, and reproducibility, and no unaddressed ethical considerations.
review_content: summary: The paper develops a constrained version of Riemannian Hamiltonian Monte Carlo, maintaining both sparsity and constraints. The paper shows that ill-conditioned, non-smooth, constrained distributions in very high dimensions, upwards of 100,000, can be sampled efficiently. The work is supported by several theoretical results on the correctness and computational efficiency, as well as bounds on the mixing time of the CRHMC and its discretizations.
strengths_and_weaknesses: ### Strengths

- Significance:
    - The paper considers an important problem in machine learning that is of great interest to the broader machine learning community.
    - The proposed framework can address sampling problems that are deemed infeasible previously. The paper shows that ill-conditioned, non-smooth, constrained distributions in very high dimensions, upwards of 100,000, can be sampled efficiently.
    - The experiments also demonstrate the approach's superiority in comparison to existing methods in terms of rate of convergence and total sampling time.
- Quality and Clarity:
    - The paper is well-written and well-organized.
    - The theoretical analyses of the work are rigorous and of very high quality. The experiments are extensive and support the findings of the theoretical part.

### Weaknesses

- Quality:
    - I believe the proof of Theorem 8 (Appendix, on the correctness of the discretized version of CRHMC) is incorrect. In principle, the detailed-balance conditions need to be verified on two arbitrary measurable sets. This could be done by verifying the condition pointwise (as done in Theorem 8) if the state space is discrete or the proposal distribution is absolutely continuous. HMC proposals, however, contain point masses, and the detailed-balance condition can be verified pointwise only if the proposal is volume-preserving. This can be seen in the proof of Theorem 6 for the ideal version of RCHMC but is absent in Theorem 8. 

      It is worth noting that this a not a minor point since the rest of the theoretical analyses (e.g., bounds on mixing time) relies on the convergence of discretized CRHMC obtained through this Theorem.
    - This brings a broader point about the choice of integrator for RCHMC, namely, why the implicit midpoint is preferred over the leapfrog integrator (which is known to be volume-preserving). The implicit midpoint is described in the paper with no context or references and very few details about its properties are spelled out. I think such an important component of the algorithm deserves some details treatment on its background and properties.
    - A minor point: the paper seems to be reluctant to state that RCHMC is ergodic or converges into the formal statement of the theorems: Theorem 6 states that RCHMC satisfies detailed balance; Theorem 8 states that RCHMC has a stationary distribution. Mathematically, those two statements are much weaker than ergodicity and are not strong enough to support later discussions such as Line 816.
- Presentation:
    - To continue on the point above, while the development of CRHMC follows well-founded rationales, I find it very hard to separate the contributions of the current construction from pre-existing methods in the field, simply because the paper doesn’t provide enough contexts. The construction of the constraint term in the Hamiltonian $\lambda(x, v)^T c(x)$ is rigorous, but is it an original concept, or does it draw inspiration from pre-existing works? The barrier function is a well-known tool in constrained HMC, but was introduced with no proper references. Similarly, several aspects of efficient computation and simplifications in Section 3 seems to base on are pre-existing ideas. I believe a more thorough review of the key ideas would strengthen the manuscript greatly.
questions: (see above)
limitations: None noted.
review_title: null

final_score: rating: 5: Borderline accept: Technically solid paper where reasons to accept outweigh reasons to reject, e.g., limited evaluation. Please use sparingly.
review_content: summary: The authors propose a method to sample from a distribution using Riemannian Hamiltonian Monte Carlo while respecting a set of constraints. The efficiency of the method relies on a combination of tricks from linear algebra, while the proposed Riemannian metric is based on the Jacobian of the constraints and a barrier function. In the experiments it is demonstrated the behavior of the sampler in different scenarios.
strengths_and_weaknesses: - Originality: As the authors mention their approach is based on a previous work, which they extend such that to make it much more efficient computationaly. I think that the related work is cited and discussed very briefly in the intro. I am not an expert in the topic of sampling, but the proposed approach seems like a sufficient contribution.

- Quality: The technical part of the paper seems ok, but I have not checked it carefully. The claims are supported by theoretical analysis and experimentally validated.

- Clarity: I think that the paper is ok written, but the Section 2 seems a bit too technical without the necessary explanations to help the reader. For example, a lot of linear algebra results are used, but I have the feeling that the gist and intuition are missing. I believe that this makes the paper hard to access from non-expert readers like me. (Please see questions)

- Significance: From the claims and the empirical results it seems that the method indeed solves the problems better (more efficient) than previous approaches. It is quite likely that this work would be useful when sampling under constraints.
questions: In particular, it is not clear to me: 
1. why the metric should not be invertible
2. how the associated null spaces are related to the problem
3. why the metric is then constructed as in Eq. 2.6 
4. how the function $\phi$ is defined in practice.

In general, I find hard to access the Section 2.1. I think some discussion about the geometrical intuition/properties of these steps together with some figures could make the paper much more accessible. I believe that this part can be improved.
limitations: The authors discuss some limitations of the proposed method.
review_title: null

final_score: rating: 7: Accept: Technically solid paper, with high impact on at least one sub-area, or moderate-to-high impact on more than one areas, with good-to-excellent evaluation, resources, reproducibility, and no unaddressed ethical considerations.
review_content: summary: The manuscript proposes a MCMC sampler based on the Constrained HMC method for sampling from a log log-concave density $\pi(x) \propto \exp[-f(x)]$ constrained on a convex subset $x \in K \subset \mathbb{R}^d$ and a constraint $c(x)=0$. The main innovations are computational tricks to scale the method to high-dimensional settings and exploit the sparsity pattern of constraints of the type $Ax=b$ for a sparse matrix $A$.
strengths_and_weaknesses: The proposed method seems to deliver state-of-the-art results when evaluated on a range of challenging settings (the numerical studies are well executed). Furthermore, the authors have incorporated their method within a popular Bioinformatics library.

I have a number of questions, as well as remarks on the presentation of the paper, that are listed in the next section.

questions: The manuscript proposes a MCMC sampler based on the Constrained-HMC method for sampling from a log-concave density $\pi(x) \propto \exp[-f(x)]$ and constrained to stay within a convex subset $x \in K \subset \mathbb{R}^d$ and and an implicit manifold defined through the constraint $c(x)=0$. The main innovations are computational tricks to scale the method to high-dimensional settings and exploit the sparsity pattern of constraints of the type $Ax=b$ for a sparse matrix $A$.

The proposed method seems to deliver state-of-the-art results when evaluated in a range of challenging settings (the numerical studies are well executed). Furthermore, the authors have incorporated their method within a popular Bioinformatics library.


1. I am not finding the discussion of sparsity well-motivated or explained. The authors mention that sparsity is preserved, but it takes quite a bit of time for the reader to realize that (I think) it is related to the sparsity pattern of A within the constraint Ax=b.

2. The authors are not clear within the main paper (although it is somehow clarified in the Supp Material) whether they are considering MC methods with or without an accept reject. For example, several of the cited papers on (variation of) the Langevin dynamics do not have an accept-reject step (and are consequently only approximately sampling from the posterior).

3. It would be best to make it clearer whether the paper concentrates on constraints of the type Ax=b, or more general implicit manifolds. It seems like all the experiments are with linear constraints.

4. It was not clear to me whether the idea of using barrier functions within an HMC scheme is new or not. If it is not, it would be best to clarify. Something that is not very well explained is the following. Suppose that one runs a standard HMC algorithm (with explicit leapfrog integrator) in $\mathbb{R}^d$ constrained to a convex set $K \subset \mathbb{R}^d$ (i.e. without an equality constraint $c(x)=0$). For a barrier function $\phi$ and positive definite mass matrix $M(x) = \nabla^2 \phi(x)$, there does not seem anything that prevents the Markov chain to step out of the convex set $K$ because of discretization error. Is it different when using the implicit midpoint method for integrating the Hamiltonian dynamics? Can the authors clarify this? 

5. Do the authors use the assumption that the target density is log-concave? where?

6. For general constraint $c(x)=0$, the CHMC requires a "projection" step to bring back the Markov chain back to the manifold $c(x)=0$. In practice (i.e. without further assumptions on the function $c$), this projection step may fail. It is the reason why the "reversibility checks"  have been developed (eg. see [1]).

7. There is a discussion that describes how having a position-dependent mass matrix can help take the geometry of the distribution into account for improved mixing properties (i.e. Riemaniann-HMC/etc..). But, if I am not mistaken, the authors are simply using as mass matrix the projection of $\nabla^2 \phi(x)$ on the tangent space to the manifold $c(x)=0$, which is *not* taking the geometry of the target distribution (eg. anisotropy of scales) into account (i.e. the target distribution $\pi(dx)$ is not used to design the mass matrix)

8. No discussion of how to tune the leapfrog step size and momentum refresh rate (i.e. coefficient $\beta$)?

9. The authors reports the ESS, but do not specify the ESS of what functional of the trajectory is considered



**Minor**:
1. if space is an issue, I do not think that it is necessary to spend so much time explaining to NeurIPS readers that sampling is important and difficult. The introduction can be shortened, if needed.


**[Answers by the authors acknowledged 4/8/2022]**
Thank you for the careful answers that have clarified my misunderstandings. It is indeed a very solid piece of work!

**References:**
1. "Monte Carlo on manifolds: sampling densities and integrating functions" by Emilio Zappa, Miranda Holmes-Cerfon, Jonathan Goodman
limitations: OK.
review_title: null

final_score: rating: 4: Borderline reject: Technically solid paper where reasons to reject, e.g., limited evaluation, outweigh reasons to accept, e.g., good evaluation. Please use sparingly.
review_content: summary: This paper provides a variation of the Constraint Riemannian Hamiltonian Monte Carlo (RHMC). Constraint variations of HMC deal with  sampling from a target distribution that by means of some constraints, is (implicitly) defined on a manifold. In the present work it is assumed that the constraint has the form c(x) = Ax − b. Such assumptions has led to fast computations and as it is mentioned, the proposed algorithm can be run on target spaces with dimensionality as high as 100K. 
strengths_and_weaknesses: * It seems to me that the contribution of the presented work is incremental. The idea of constraint HMC is already proposed in a more general setting by: 
    Brubaker et al, "A Family of MCMC Methods on Implicitly Defined Manifolds" 
where Constrained Riemann Manifold HMC (CRHMC) is mentioned as a specific case. The present work correctly gives this credit to Brubaker et al. (in Line 135) but then claims Brubaker's CRHMC "requires eigenvalue decomposition and is not efficient for large problems". As far as I understand, Brubaker's work does not necessarily require an eigenvalue decomposition and the present work is indeed an instance of Brubaker's work where extra assumptions/restrictions mentioned in Lines 192 to 220, has led to efficient computations. 

* If I have understood correctly and the contribution is efficient computations under extra assumptions/restrictions, then the robustness of this work could possibly be demonstrated by extensive empirical experiments. Nonetheless, the provided experimental results seem minimal to me (the proposed algorithm is only compared against one benchmark) which is suitable for sampling from (constraint) uniform distributions only. I will return to this in the "Questions" section. 

* In general the paper's writing could be improved and notations/terminology could be defined more rigorously.   

questions: * The proposed work is only compared against (2 implementations of) CHRR algorithm. I suggest that the comparison be with other variations of constraint HMC (e.g. Brubaker's constraint baseline HMC, constraint Metropolis-Hastings) and/or more recent implementations of CHMC algorithms such as Graham's "Manifold Markov chain Monte Carlo methods for Bayesian inference in diffusion models".

* I think CHRR is for sampling from (constraint) uniform distributions. Are the presented experimental target densities uniform? If no, how CHRR is used and if yes, what is the usage of Riemannian HMC (since gradients are always 0) ?


limitations: 
Not applicable. This work which is theoretical.
review_title: null

Meta-review: metareview: The focus of the submission is practically efficient sampling of non-smooth constrained log-concave distributions in high dimension as formulated in (1). In order to tackle this task, the authors design a constrained variant of the RHMC (Riemannian Hamiltonian Monte Carlo) technique, relying on self-concordant barrier functions. They demonstrate the practical efficiency of the proposed method (achieving significantly improved sampling time) compared to existing constrained samplers, with additional theoretical results in the supplement.

Sampling in high dimension is a fundamental problem of data science with various successful applications; designing new methods in the area is of clear interest for the NeurIPS community. As it was assessed by the reviewers, the authors deliver important new tools in this context. As it was also noted, the manuscript could be made slightly more accessible to wider audience.
confidence: Certain
Decision: Accept