The formula implemented in (https://github.com/mathnet/mathnet-numerics/blob/master/src/Numerics/Distributions/NegativeBinomial.cs) on line 260 is incorrect: ### var lambda = Gamma.SampleUnchecked(rnd, r, p); should be: ### var lambda = Gamma.SampleUnchecked(rnd, r, p/(1-p)); See e.g., https://mathworld.wolfram.com/NegativeBinomialDistribution.html or https://timothy-barry.github.io/posts/2020-06-16-gamma-poisson-nb/ or https://gregorygundersen.com/blog/2019/09/16/poisson-gamma-nb/ or https://probabilityandstats.wordpress.com/tag/poisson-gamma-mixture/