Venue: NeurIPS 2022 Conference

Reviews:

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 authors propose a mechanism to produce unbiased samples from expensive posteriors given that it is easy to sample from the prior and possible (although expensive) to evaluate the likelihood. The approach proposes concrete particle dynamics based on NFs that approach the target distribution and corrects them using Metropolis-Hastings. 
strengths_and_weaknesses: Strengths:
- Technically correct
- Clear presentation
- Original dynamics
- Uses the gradient information of the target distribution 
- Produces multiple uncorrelated samples

Weaknesses:
Experimental validation, in terms of:
- Complexity of the models sampled from
- Strength of the baselines
- Ablations
questions: 1. It would seem that a lot of the power of this approach is coming from the strong modeling ability of the NF. I am missing comparisons with simpler approaches that can leverage this power. For instance:

(a) What if we simply run variational inference to optimize the normalizing flow? One can easily compute the gradients with respect to the parameters. (For the SINF case, it might be possible to modify its main algorithm to take an unnormalized log-pdf instead of a set of samples - but the gradients for variational inference can always be approximated stochastically in a NF). Then one can find the closest q(z) within the NF domain and then sample from it and correct using MH. This decouples the number of particles used for learning and the ones used for sampling: More independent samples can be generated on demand after training. And does not require any particle dynamics.

(b) What if we only use SINF and MH? MH produces a new set of samples closer to the posterior, then we use those samples and SINF to update q(z) in closed form. Then we repeat. Again, we have eliminated the separate particle dynamics.

While (a) and (b) might be significantly worse than the proposal of the paper, it'd be reassuring to know that this is the case, and that all the heavy lifting is not being done by the NF being able to learn the energy function. Otherwise, the proposed dynamics don't seem justified enough.

2. The experiments in which the posteriors are calculated are pretty simple. No strong baselines for those models are used, just NUTS which of course is slow.

3. The tested models are _very_ low-dimensional. The high dimensional setting might have very poor performance, this is not known.

4. I expect multimodal likelihoods with broad priors to result in a lack of mixing - I doubt all modes will be discovered even in low dimensional cases, when they are spaced and sharp. Again this is not tested even in a 2-dimensional experiment with many gaussians that are well separated. 


Minor comments:
- Parentheses are used for bibliography citations, should be square brackets, not to be confused with cited equations.
- In the MH update in Algorithm 1, time indices are missing.

[1] Black Box Variational Inference. Rajesh Ranganath et al.



limitations: Addressed in the previous section.
review_title: null

final_score: rating: 6: Weak Accept: Technically solid, moderate-to-high impact paper, with no major concerns with respect to evaluation, resources, reproducibility, ethical considerations.
review_content: summary: The work proposes a Bayesian inference procedure constructed using a deterministic Langevin equation. The resulting algorithm uses a number of particles initialized at random, their dynamics then following the deterministic equation toward its equilibrium distribution, while at regular time intervals a normalizing flow is fit to the collection of particles and used as an independence proposal within Metropolis-Hastings to propose replacing each particle with a new draw. The method is demonstrated on a number of toy examples and standard data sets.

strengths_and_weaknesses: There are some interesting aspects to this work, such as the use of the deterministic Langevin equation, and the empirical results are quite good, especially the ablation studies in the supplementary material that separate out the various components that come together in the final algorithm. There are some gaps, however, that are stopping me from recommending acceptance.

While I think I understand the construction of the algorithm, it is not clear to me that this results in samples distributed according to the target distribution, even up to the approximations introduced (e.g. fitting the normalizing flow to the current stock of particles, time discretization). With that said, the empirical results suggest that the algorithm works as intended numerically; if there is a bias it may be small. I've added some questions for the authors on this topic below and hope they can help me to better understand their work.

It is unclear to me how the results for P-DLMC have been obtained (e.g. in Figure 3). The authors note that these are based "...on the assumption that we can evaluate the likelihood gradient in parallel on N CPU cores" (bottom p.7). Is there actually a parallel implementation here? In my opinion, unless there is a parallel implementation to provide empirical results, this ought to be an item for the discussion, not the experimental section. Otherwise, one is not accounting for the real-world challenges of parallel computing (fewer than N cores, possibly variable execution time between particles leading to load balancing issues, communication overhead in gathering all particles to update the normalizing flow, etc). 

If it is a model rather than actual implementation, have the authors simply divided through by N? Or have they used a model where only some p portion of the sequential program can be parallelized, thus using a scaling factor such as (1 - p + p/N)? There is a comment about the assumption that evaluating the likelihood gradient takes about a minute, while the cost of the normalizing flow is in the seconds, but even that limits the maximum speedup to much less than N (say 60 seconds likelihood gradient parallelized down to essentially zero, 1 second for normalizing flow, that still limits speedup to factor of 61, i.e. 61 seconds down to 1 second, even when N >> 61). These are some of the challenges of scaling these sort of population/ensemble algorithms in practice (SMC has similar challenges).

questions: * Before the Metropolis-Hastings step, aren't the particles only approximately distributed according to the target distribution, given time discretization to get them there, and applying the Metropolis-Hastings kernel will not leave them properly distributed, either? This is not the case with HMC or MALA, say, where the time-discretized move is itself the proposal (over the current state of the Markov chain), and the Metropolis-Hastings rule accepts or rejects it. Here it seems there is a particle drawn from the time discretization, and a proposal drawn from the normalizing flow fit to the collection of those particles, and we accept or reject the swap. Perhaps the authors could better explain how this yields samples from the target distribution, or clarify the claims (e.g. first paragraph of conclusion claims "To make it asymptotically unbiased we add MH step...").

* Independence proposals rely on the proposal being a good approximation of the target distribution in order to have a reasonable acceptance rate, and so for the Markov chain to mix well. I can accept that a normalizing flow can achieve a good approximation with sufficient size, but I do not see any empirical results on this part of the algorithm, e.g. what is the average acceptance rate of this step in the examples?

* Bottom p.7 there is a comment that, "...we find there is a tradeoff between number of particles and number of iterations, such that starting with a larger number of particle does not always imply a higher overall computational cost." I can see how this might be the case, but can also seem the limits of it, as in e.g. the extreme case of reducing to one particle. A plot showing this tradeoff would be a good addition.

* As above, specific numerical results on the swap moves (e.g. acceptance rate) would be a good addition.

* As above, some clarity on the presentation of the parallel version.

limitations: No concerns.
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: The authors of this work propose using the deterministic “dual” of the langevin equation for simulating from a target posterior. This dual presents difficulties since it requires estimating the score of the FPK equations solutions, to do this the authors propose using normalizing flows, furthermore , the authors “precondition” the flows in latent space via simulating what seems to be latent Gaussian dynamics at equilibrium, this “preconditioning” improves/stabilizes the overall procedure.
strengths_and_weaknesses: Pros:
1. A sound methodology for simulating the deterministic “dual” of the Langevin equation using normalizing flows.
2. Sound and successful experimental procedure.

Cons:

1. The overall guarantees of the method (theoretically) are not clear and there's no acknowledgement if the method is mostly heuristic. 
2. Furthermore, intuition behind the method's preconditioning step is also not presented clearly. There are some potential conceptual/theoretical issues illustrated by the reviewer in point 3 (with a potential contradiction).
3. Lack of motivation when moving from langevin to its stochastic counterpart.

questions: Comments/Suggestions:

1. The deterministic Langevin Equation is well known and used in the diffusion score matching literature [1] (see Eq X) and should be acknowledged appropriately and maybe briefly discussed. Notice the velocity term coincides with the score. 
2. Equation 2 is broken U(x) should be U(x(t)) in order to be consistent, however this is also wrong, the FPK equation is defined pointwise for all x \in R^d  that is its defined as q(x,t) or q_t(x) , it's not only defined pathwise on x(t), this notation is inaccurate. The f(x(t)) notation should be reserved for the deterministic langevin equation, but the FPK equation is defined pointwise. Across the whole paper the usage of q(x), q(x(t)) and x vs x(t) in general is quite inconsistent.
3. From flows setup you would expect p(z|y) to be Gaussian but this is never made clear, instead you refer to \pi(z)  being Gaussian which was never explicitly defined in the context of p(x|y) ? let's assume that in your notation \pi(z) = p(z|y) = N(0,1) then indeed the score of q_t(z) is infact -z (i.e. V(z)=z) however this implies that z(t) evolves entirely according to linear Gaussian dynamics , not only that but if the score is indeed -z this implies that q_t(z) = N(z; 0,1) which is at equilibrium at does not depend on time whatsoever. What is the point of simulating gaussian noise at equilibrium in latent space ?  Formally I cannot see how this can obtain finite time estimates of q_t(x)  (t< \infty) 
4. To clarify the previous point a bit further if indeed V(z)=z then this implies directly that -ln q_t = || z||^2/2 which implies directly that q_t is an independent Gaussian where q_t(z) solves the FPK equation for the density implied by the ODE dynamics in Eq 8. Remember U(z) and V(z) cannot be chosen independently, given U(z)  then grad(-U(z)) implies an FPK equation  where V(z) depends on the solution to that equation. One can't just go  oh let V(x) = ln Gaussian without first proving that q_t(z) = Gaussian.
5. Continuing the line of reasoning dz_t =grad(-U(z))dt +dW_t should have as its steady state N(0,1)  if z_0 where samples from \pi(z)=p(z|y) then in fact V(z) = z can easily be derived, but again you'd be just simulating an OU process at equilibrium … What is the advantage of this? Can you formally prove an advantage of doing this? Or at least motivate it conceptually ?


Questions on methodology:

1. In your algorithm you don't have any pseudocode for the z(t) dynamics ? Why is this ? I found this quite confusing.
2. For estimating q_t(x) you just say run normalizing flows  on x(t) ? so here you just do density estimation with normalizing flows correct ? (you learn flows in the standard way) . It really feels this is the essence to the method, all that is being done is running normalizing flows to estimate q_t(x) and then V(x) and then some adhoc stabilizing noise procedure is done with a lack of motivation. 
3. How is the preconditioning used then ? you do some extra gaussian dynamics in latent space with the learned flows to stabilize things a bit ?  (Please clarify and give more motivation)
4. I can't seem to find an ablation for the effect of preconditioning/simulating the Gaussian, what effect does the adhoc preconditioning routine have? You can still use NF to fit the score without this step.

limitations: Overall there are clear counter points in this work that lack detailed discussion:

1. When going from Langevin to the deterministic version, more approximation is required, in  sense it is a harder problem, whilst this allows for a hybrid version between ULA / normalizing flows, the authors do not motivate well why one would want to pursue such a problem (which requires more difficult approximation schemes).
2. The whole derivation of the latent dynamics in the flow/gaussian space seems lacking and also is motivated. At the end of the day the equations you have seem to simulate Gaussian dynamics at equilibrium, this is either not particularly useful or wrong (if somehow they are out of equilibrium in practice). 



[1] Song, Y., Sohl-Dickstein, J., Kingma, D.P., Kumar, A., Ermon, S. and Poole, B., 2020. Score-based generative modeling through stochastic differential equations. arXiv preprint arXiv:2011.13456.
review_title: null

Meta-review: metareview: This paper proses an inference method that combines gradient ascent and normalizing flows. The idea is that one could, in principle, simulate the deterministic Fokker-Planck equation, but this would require access to the density of the evolving approximating density, which is intractable. Thus, the paper proposes to maintain a set of particles and update a normalizing flow to approximate this density. The resulting procedure is deterministic, with an accuracy that depends on the number of particles and the power of the normalizing flow.

Reviewers agreed this was an interesting approach and the experimental results are promising (albeit fairly low-dimensional). However, there were a few apparent weaknesses: Firstly, there was a lack of clarity about the theoretical guarantees. The authors state that this is all clear in the manuscript, but readers would undoubtedly benefit from a much more centralized/explicit description of what approximations are involved, and under what guarantees the method is claimed to work, which could be put in a single place. In addition, in trying to understand exactly what was done in the experiments, it is difficult to understand several of the details. Algorithm 1 is very helpful in this regard, but it would be beneficial to have a self-contained elaboration of all of the points (perhaps even in an appendix). Finally, the experimental results are all relatively low-dimensional. Often particle methods do not scale gracefully to higher dimensions. (I do not see this as a huge flaw because the method could be useful even if it does not scale, but reader would benefit from evidence on this point either way.)

In the end, however, most of the above issues are issues of clarity and I am willing to trust the authors to fix these before final submission. The paper appears to present a novel idea and the community would benefit from seeing it and discussing it.
confidence: Less certain
Decision: Accept