Archive for randomisation

random variate generation with [finite] guarantees

Posted in Books, Statistics, University life with tags , , , , , , , , on August 17, 2025 by xi'an

I came across this paper by Feras A. Saad and Wonyeol Lee, to appear in Proc. ACM Program. Lang. It is calling for a finite precision assessment of random numbers generators. Rather than the “fictitious infinite-precision (Real-RAM)” model. With the following illustration

“Mironov (2012) demonstrates that floating-point effects in the Laplace random variate generator from existing software libraries can entirely destroy the real-world privacy guarantees of algorithms”

Their solution is to resort to finite precision computation of the CDF of a target distribution, and then to apply the inverse CDF transform to a chain of random bits. I did not go through any of the technical (gory) details of the implementation, presented as an optimised version of the original Knuth and Yao method, but the author compute the cdf of standard distributions from

“the GNU Scientific Library (GSL) by reusing high-quality CDF implementations. The built-in GSL Gaussian generators often have complex implementations spanning hundreds of lines of code, and each specify different output distributions which are all intractable to estimate. Indeed, any GSL random variate generator that makes just two (or more) calls to uniform is already intractable to analyze” 

and claim faster execution times, larger ranges of output, and a minimal overhead for extended-accuracy generators. I wonder if an MCMC study is under production towards handling intractable CDFs.

differentially private distributed Bayesian linear regression with MCMC

Posted in Books, pictures, Statistics, University life with tags , , , , , , , , , on August 30, 2023 by xi'an

An ICML 2023 paper by Barıs¸ Alparslan, Sinan Yıldırım¸ and Ilker Birbil that (re)addresses the issue of privacy when running a Bayesian regression analysis. Resorting to the common notion of differential privacy, imposing a limited variability if a single observation is modified, and a Gaussian randomisation of the observations.

“A differentially private algorithm constrains the difference between the probability distributions of the output values obtained from neighbouring data sets”

In the super classical setup of simple Normal linear regression, y=Xθ+σε. Summary statistics are chosen as

S=X’X and z=X’y,

(why the separation?) then randomised. (Keeping Ŝ definite positive? Not necessarily, it appear.) Inspired directly from Dwork & al. (2014). The authors still manage to spend an entire column in (re)deriving the conditional Normal distribution of z conditional on S and (θ,σ)… Which is later exploited for integrating z out in the MCMC algorithm.

“some important differences between our work and that of Bernstein & Sheldon (2019) [stem] from the choice of summary statistics and the consequent hierarchical structure used for modelling linear regression [and]lead to significant differences in the inference methods as well as significant computational advantages [O(d³) vs. O(d⁶)]”

In a distributed setting several agents are handling their own data and keep their privacy by the same mechanishttps://www.slideshare.net/xianblog/discussion-of-icml23pdfm [as in the top graph from the paper]. On principle, a Bayesian analysis of the resulting hierarchical model should directly consider the posterior on the global parameter by considering the distributions of the randomised pairs (ẑ,Ŝ). The elephant in the room is the distribution of the regressors, which is customarily unknown and not accounted for in a traditional Bayesian analysis. It is needed here due to the division in S and z, plus the randomisation step that calls for the posterior distribution of S given Ŝ. Elephant that is exfiltrated by either assuming Normality or substituting Ŝ for S without accounting for the noise! Definitely not exactly Bayesian. Another column is spent on the Metropolis-within-Gibbs simulation of the posterior…

Overall, I remain reserved about this approach, since it does not follow a clear Bayesian pathway and in particular does not incorporate privacy as part of the Bayesian decision analysis.

a football post?!

Posted in Statistics with tags , , , , , , , , , , , , , , on June 22, 2022 by xi'an

I am not interested in football, neither as a player (a primary school trauma when I was the last being picked!) or as a fan, contrary to my dad (who was a football referee in his youth) and my kids, but Gareth Roberts (University of Warwick) and Jeff Rosenthal wrote a paper on football draws for the (FIFA) World Cup, infamously playing in Qatar by the end of the year, which Gareth presented in a Warwick seminar.

For this tournament, there are 32 teams, first playing against opponent teams supposedly drawn from a uniform distribution over all draw assignments, within 8 groups of 4 teams, with constraints like 1-2 EU teams per group, 0-1 from the other regions. As done at the moment and on TV, the tournament is filled one team at time by drawing from Pot 1, then Pot 2, then Pot 3, & Pot 4. &tc.. Applying the constraints one draw at a time, conditional on the past draws and the constraints, rather obviously creates non-uniformity! Uniformity would be achievable by rejection sampling (with a success probability of 1/540!) But this is not televisesque enough…

A debiasing solution is found by using several balls for each team in the right proportion, correcting for the sequential draws. Still impractical when requiring 10¹⁴ balls…!

The fun in their paper is that the problem can be formulated as a particle filter, estimating the right probabilities by randomising the number of balls [hidden randomness] and estimating the probability for team j to be included by a few thousands draws. With some stratified sampling on the side to minimise randomness. Removing the need for the (intractable?) distribution is thus achieved by retrospective sampling, as in pseudo-marginal MCMC. Alternatively, one could swap pairs of teams by a simplistic MCMC algorithm, with no worry about stationarity and the possibility of on-screen draws. (Jeff devised a Java applet to simulate an actual draw.) Obviously, it is still a far stretch that this proposal will be implemented for the next World Cup. If so, I will watch it!

multilevel linear models, Gibbs samplers, and multigrid decompositions

Posted in Books, Statistics, University life with tags , , , , , , , , , , , , , on October 22, 2021 by xi'an

A paper by Giacommo Zanella (formerly Warwick) and Gareth Roberts (Warwick) is about to appear in Bayesian Analysis and (still) open for discussion. It examines in great details the convergence properties of several Gibbs versions of the same hierarchical posterior for an ANOVA type linear model. Although this may sound like an old-timer opinion, I find it good to have Gibbs sampling back on track! And to have further attention to diagnose convergence! Also, even after all these years (!), it is always a surprise  for me to (re-)realise that different versions of Gibbs samplings may hugely differ in convergence properties.

At first, intuitively, I thought the options (1,0) (c) and (0,1) (d) should be similarly performing. But one is “more” hierarchical than the other. While the results exhibiting a theoretical ordering of these choices are impressive, I would suggest pursuing an random exploration of the various parameterisations in order to handle cases where an analytical ordering proves impossible. It would most likely produce a superior performance, as hinted at by Figure 4. (This alternative happens to be briefly mentioned in the Conclusion section.) The notion of choosing the optimal parameterisation at each step is indeed somewhat unrealistic in that the optimality zones exhibited in Figure 4 are unknown in a more general model than the Gaussian ANOVA model. Especially with a high number of parameters, parameterisations, and recombinations in the model (Section 7).

An idle question is about the extension to a more general hierarchical model where recentring is not feasible because of the non-linear nature of the parameters. Even though Gaussianity may not be such a restriction in that other exponential (if artificial) families keeping the ANOVA structure should work as well.

Theorem 1 is quite impressive and wide ranging. It also reminded (old) me of the interleaving properties and data augmentation versions of the early-day Gibbs. More to the point and to the current era, it offers more possibilities for coupling, parallelism, and increasing convergence. And for fighting dimension curses.

“in this context, imposing identifiability always improves the convergence properties of the Gibbs Sampler”

Another idle thought of mine is to wonder whether or not there is a limited number of reparameterisations. I think that by creating unidentifiable decompositions of (some) parameters, eg, μ=μ¹+μ²+.., one can unrestrictedly multiply the number of parameterisations. Instead of imposing hard identifiability constraints as in Section 4.2, my intuition was that this de-identification would increase the mixing behaviour but this somewhat clashes with the above (rigorous) statement from the authors. So I am proven wrong there!

Unless I missed something, I also wonder at different possible implementations of HMC depending on different parameterisations and whether or not the impact of parameterisation has been studied for HMC. (Which may be linked with Remark 2?)

back to the Bernoulli factory

Posted in Books, Statistics, University life with tags , , , on April 7, 2020 by xi'an

“The results show that the proposed algorithm is asymptotically optimal for the mentioned subclass of functions, in the sense that for any other fast algorithm E[N] grows at least as fast with p.”

Murray Pollock (Warwick U. for a wee more days!) pointed out to me this paper of Luis Mendo on a Bernoulli factory algorithm that estimates functions [of p] that can be expressed as power series [of p]. Essentially functions f(p) such that f(0)=0 and f(1)=1. The big difference with earlier algorithms, as far as I can tell, is that the approach involves a randomised stopping rule that involves, on top of the unlimited sequence of Bernoulli B(p) variates a second sequence of Uniform variates, which sounds to me like a change of paradigm, given the much higher degree of freedom brought by Uniform variates (as opposed to Bernoulli variates with an unknown value of p). Although there exists a non-randomised version in the paper. The proposed algorithm is as follows, using a sequence of d’s issued from the power series coefficients:

1. Set i=1.
2. Take one input X[i].
3. Produce U[i] uniform on (0,1). Let V[i]=1 if U[i]<d[i] and V[i]=0 otherwise.
If V[i] or X[i] are equal to 1, output X[i] and finish.
Else increase i and go back to step2.

As the author mentions, this happens to be a particular case of the reverse-time martingale approach of Łatuszynski, Kosmidis, Papaspiliopoulos and Roberts (Warwick connection as well!). With an average number of steps equal to f(p)/p, surprisingly simple, and somewhat of an optimal rate. While the functions f(p) are somewhat restricted, this is nice work