# \[ANN-RFC\] KissABC.jl- Approximate Bayesian Computation

**URL:** https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668
**Category:** Package Announcements
**Tags:** announcement, bayesian-inference
**Created:** [June 3, 2020, 1:25pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668 "2020-06-03T13:25:12Z")
**Posts on this page:** 20
**Page:** 2

<div class="post-metadata">

### Author: ![mcreel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcreel/32/30088_2.png) [@mcreel](https://discourse.julialang.org/u/mcreel)
#### Post date: [June 10, 2020, 6:04pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/22 "2020-06-10T18:04:25Z")

</div>

Yes, I know, but, asymptotically, the quantiles of the posterior will define correct confidence intervals under the frequentist interpretation (see [An MCMC approach to classical estimation - ScienceDirect](https://www.sciencedirect.com/science/article/pii/S0304407603001003?casa_token=M0j7AvCXgY8AAAAA:xATDbr571oeyFlQfWr5jd-h2iiDUo7OtUs3gpE3GifT2e0FWlR9Vixl0o6wf8NV-VXxPCJllMvo), theorem 3, for example).

The issue of the variance being larger than what one gets with full sample inference, and thus, leading to broader confidence intervals, is different than the confidence intervals having correct coverage.

So, it seems to me that it would be desirable to be able to find a tuning method so that what I’m loosely calling confidence intervals would contain the true parameter values the appropriate percentage of the time, at least when the sample is large enough so that asymptotic results should be a good approximation.

---

<div class="post-metadata">

### Author: ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)
#### Post date: [June 10, 2020, 6:10pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/23 "2020-06-10T18:10:02Z")

</div>

I don’t have access to that paper unfortunately, but in general it’s easy to get Bayesian models that have no frequentist interpretation (though in your case that’s not the issue since you’re literally using an RNG). And how far into the asymptotic range you have to go can be variable based on the model.

It looks like you’re using 5000 samples in your model. Does it get more like a calibrated confidence interval if you move to 10000 or 100000 ?

---

<div class="post-metadata">

### Author: ![francesco.alemanno](https://avatars.discourse-cdn.com/v4/letter/f/e8c25b/32.png) [@francesco.alemanno](https://discourse.julialang.org/u/francesco.alemanno)
#### Post date: [June 10, 2020, 6:55pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/24 "2020-06-10T18:55:48Z")

</div>

At 5000 samples It should be already damn close to MLE prediction, Indeed now, with an improved distance function the intervals are much tighter, and even now It Is not optimal, ABC methods require a lot of tuning to yield exact bayesian posteriors  
@mcreel you can improve the results of ABC via regression adjustement, for now i do not Plan on introducing such methods in KissABC since they are not general purpose, but they can be pretty effective as a post processing procedure on the results of and ABC inference algorithm

---

<div class="post-metadata">

### Author: ![mcreel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcreel/32/30088_2.png) [@mcreel](https://discourse.julialang.org/u/mcreel)
#### Post date: [June 10, 2020, 7:00pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/25 "2020-06-10T19:00:02Z")

</div>

I think the main issue is that the distance measure is not very good for the purpose of arriving to good confidence intervals. The KissABC.jl documentation notes “now we need a distance function to compare datasets, this is possibly the worst distance we could use, but it will work out anyway”. The problem, I believe, is that the distance measure is giving equal weight to statistics that probably have very different variances. The tail quantiles very likely have much larger variances than do the central quantiles. With improper weighting of statistics, that Theorem 3 I cited above will not hold. So, I believe that with other choices of distance measure, the issue could be resolved. But, possibly, there may be ways of tuning the algorithm to get around the problem even with the chosen distance measure. That’s what my question is focused on.

I want to note that the distance measure is only an example, it’s not part of the algorithms that the package provides.

---

<div class="post-metadata">

### Author: ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)
#### Post date: [June 10, 2020, 7:08pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/26 "2020-06-10T19:08:10Z")

</div>

That makes sense, yes using all the quantiles from 0:0.01:1 as equally important means you will need to go deep into high sample size before the distance metric itself is asymptotically anything… The 1% and 99% quantiles are high variance for example.

In essence fitting an entire distribution is an infinite dimensional problem, and there doesn’t exist an asymptotic sample number when the dimensionality is infinite. Of course in this case you’ve restricted it to a 100 dimensional finite dimension approximation, and eventually you’d get to an asymptotic range, but I imagine for that model the asymptotic range is well above 5000 samples.

---

<div class="post-metadata">

### Author: ![francesco.alemanno](https://avatars.discourse-cdn.com/v4/letter/f/e8c25b/32.png) [@francesco.alemanno](https://discourse.julialang.org/u/francesco.alemanno)
#### Post date: [June 15, 2020, 11:13am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/27 "2020-06-15T11:13:33Z")

</div>

Just released KissABC 1.3, now supporting Kernel Differential Evolution, which is way more sample efficient than hard threshold ABC!

For those interested look for the method `KABCDE`, and as always let me know if you find any issues, i’ll do my best to solve them almost in real time 🙂

---

<div class="post-metadata">

### Author: ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)
#### Post date: [June 15, 2020, 12:41pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/28 "2020-06-15T12:41:22Z")

</div>

awesome! I am looking forward to trying this out at some point.

---

<div class="post-metadata">

### Author: ![robsmith11](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/robsmith11/32/29641_2.png) [@robsmith11](https://discourse.julialang.org/u/robsmith11)
#### Post date: [June 15, 2020, 2:36pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/29 "2020-06-15T14:36:58Z")

</div>

I gave it a try with my same toy mixture model.

- Compared with version 1.0, I get slightly tighter posteriors for ABCDE with `earlystop=true` and similar run times.
- But with the new default setting (`earlystop=false`), it seems to never converge, despite it reaching “completion = 1.0” quickly. I’m not sure if there’s a termination bug or it really takes that much longer to run (I killed it after it was still running 10x longer than with earlystop).
- KABCDE seems to work okay, but with the same settings it runs 2.5x longer than ABCDE and results in a posteriors that are 5-10x wider.

---

<div class="post-metadata">

### Author: ![francesco.alemanno](https://avatars.discourse-cdn.com/v4/letter/f/e8c25b/32.png) [@francesco.alemanno](https://discourse.julialang.org/u/francesco.alemanno)
#### Post date: [June 15, 2020, 2:49pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/30 "2020-06-15T14:49:40Z")

</div>

Hey, thank you for trying It out, the new changes are due to the fact that some hard problems require refinement steps even after reaching the target tolerance, so `earlystop` defaults to `false`, the reason the algorithm Is running longer Is that the parameter `generations` defaults to `500`,

I think you should tune It according to your problem, by the way the methods `KABCDE` uses the epsilon value in a completely different way, so ti obtain comparable results between kernelised and non kernelised methods the epsilon must be chosen accordingly (and you must consider the weights when using the kernelised methods in order to have proper CI’s)

The past performance can be easily attained and surpassed with some care, I simply didnt want to use defaults that could lead to wrong results

---

<div class="post-metadata">

### Author: ![robsmith11](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/robsmith11/32/29641_2.png) [@robsmith11](https://discourse.julialang.org/u/robsmith11)
#### Post date: [June 15, 2020, 4:42pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/31 "2020-06-15T16:42:03Z")

</div>

Thanks! After letting ABCDE run a bit longer (needed about 14x longer without earlystop), I got similar results, so I guess my problem just isn’t difficult enough as you said.

As for KABCDE, yes, I forgot to consider the weights! After running 2.5x as long as ABCD with earlystop, I now get CI that are only 1.5x wider. So still not quite as good. Is there anything you would recommend tuning? Or do you only expect it to be more sample efficient than ABCD without earlystop?

---

<div class="post-metadata">

### Author: ![francesco.alemanno](https://avatars.discourse-cdn.com/v4/letter/f/e8c25b/32.png) [@francesco.alemanno](https://discourse.julialang.org/u/francesco.alemanno)
#### Post date: [June 15, 2020, 7:10pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/32 "2020-06-15T19:10:23Z")

</div>

@robsmith11

i deeply care about improving performance as much as possible, so i experimented a bit with your problem

```julia
function sim((u1, p1), params; n=10^6, raw=false)
 u2 = (1.0 - u1*p1)/(1.0 - p1)
 x = randexp(n) .* ifelse.(rand(n) .< p1, u1, u2)
 raw && return x
 [std(x), median(x)]
end

function dist(s, s0)
 sqrt(sum(((s .- s0)./s).^2))
end
plan=ABCplan(Factored(Uniform(0,1), Uniform(0.5,1)), sim, [2.2, 0.4], dist)

res,del,conv=ABCDE(plan, 0.01, nparticles=100,generations=30,parallel=true)

using Statistics
function getCI(x::Vector{<:Number})
    quantile(x,[0.25,0.5,0.75])
end
function getCI(x::Vector{<:Tuple})
    [getCI(getindex.(x,i)) for i in 1:length(x[1])]
end

```

results

```julia
via getCI(res)
240 generations:
 [0.48958933397111065, 0.4924062224370781, 0.49559446402487584]
 [0.879783065265908, 0.8816472031816496, 0.8835803050367947]
120 generations:
 [0.4893221894893949, 0.49278449533673585, 0.494863093578758]
 [0.8795982153875357, 0.8816146345915951, 0.8829915018185673]
60 generations:
 [0.4887524655164148, 0.49234470862896673, 0.49502567359353133]
 [0.8796953221457162, 0.8814094610516047, 0.8833899764047788]
30 generations:
 [0.48893164848681747, 0.49163740259340305, 0.4944261757524125]
 [0.8795333045815598, 0.8809134893725196, 0.8829819098348083]
early stop:
 [0.49006664933267297, 0.49313860531909304, 0.49497013116625105]
 [0.8804136291097875, 0.8819843728641816, 0.8834306754737902]

```

i see very good convergence on this problem, even 30 generations lead reasonable CI’s, so on `ABCDE` i do not see a clear cut performance regression, surely without enforced early stopping there is a bit of tuning involved, what i did see instead is that `KABCDE` is not performing very well on this problem, i will think of something to improve it’s performance.

is there anything i’ve missed?

---

<div class="post-metadata">

### Author: ![StevenSiew](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevensiew/32/218393_2.png) [@StevenSiew](https://discourse.julialang.org/u/StevenSiew)
#### Post date: [June 16, 2020, 2:53am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/33 "2020-06-16T02:53:40Z")

</div>

Can someone please explain the purpose of the argument **param**  
it does not seems to do anything!

```julia
sim((μ,σ), param) = randn(100) .* σ .+ μ

```

---

<div class="post-metadata">

### Author: ![pjentsch0](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pjentsch0/32/15148_2.png) [@pjentsch0](https://discourse.julialang.org/u/pjentsch0)
#### Post date: [June 16, 2020, 3:37am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/34 "2020-06-16T03:37:30Z")

</div>

`param` is used to pass any parameters that your `sim` function might need, but which you don’t want to do SMC on.

---

<div class="post-metadata">

### Author: ![robsmith11](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/robsmith11/32/29641_2.png) [@robsmith11](https://discourse.julialang.org/u/robsmith11)
#### Post date: [June 16, 2020, 4:33am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/35 "2020-06-16T04:33:12Z")

</div>

Thanks. Agreed, by cutting the number of generations ABCDE can be made very fast with good results even without early stopping.

---

<div class="post-metadata">

### Author: ![StevenSiew](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevensiew/32/218393_2.png) [@StevenSiew](https://discourse.julialang.org/u/StevenSiew)
#### Post date: [June 16, 2020, 7:31am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/36 "2020-06-16T07:31:19Z")

</div>

```julia
tdata=randn(1000).*0.04.+2

sim((μ,σ), param) = randn(100) .* σ .+ μ

```

The next thing which I do not understand is what the TESTDATA tdata has 1000 elements in an array but the sim function only generate an array with only 100 elements in it.

Why 100 elements? Why not 1000 elements to make it a parity with the TESTDATA? Or why not just 10 elements or 1\_000\_000 elements? It’s all very very mysterious and confusing.

---

<div class="post-metadata">

### Author: ![francesco.alemanno](https://avatars.discourse-cdn.com/v4/letter/f/e8c25b/32.png) [@francesco.alemanno](https://discourse.julialang.org/u/francesco.alemanno)
#### Post date: [June 16, 2020, 7:36am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/37 "2020-06-16T07:36:46Z")

</div>

That is a mistake It should be the same, i Will fix It when i get a chance

---

<div class="post-metadata">

### Author: ![robsmith11](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/robsmith11/32/29641_2.png) [@robsmith11](https://discourse.julialang.org/u/robsmith11)
#### Post date: [June 16, 2020, 9:31am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/38 "2020-06-16T09:31:22Z")

</div>

Is there a reason why the number of samples in the simulated and test data should match?

I usually just adjust the number simulated samples to keep the noise of my chosen statistic within a desired range while minimizing the simulation cost. Depending on the amount of test data available and the specifics of the statistic, this may me mean taking more or fewer simulation samples.

---

<div class="post-metadata">

### Author: ![StevenSiew](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevensiew/32/218393_2.png) [@StevenSiew](https://discourse.julialang.org/u/StevenSiew)
#### Post date: [June 16, 2020, 10:40am UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/39 "2020-06-16T10:40:58Z")

</div>

I found a youtube video on

ABC  
Approximate Bayesian Computation

[The ABC’s of ABC (Approximate Bayesian Computation)](https://youtu.be/MsgdXDXXP_0)

---

<div class="post-metadata">

### Author: ![francesco.alemanno](https://avatars.discourse-cdn.com/v4/letter/f/e8c25b/32.png) [@francesco.alemanno](https://discourse.julialang.org/u/francesco.alemanno)
#### Post date: [June 16, 2020, 1:03pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/40 "2020-06-16T13:03:45Z")

</div>

In theory they must be the same, in practice it is very convenient to have them different

---

<div class="post-metadata">

### Author: ![pjentsch0](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pjentsch0/32/15148_2.png) [@pjentsch0](https://discourse.julialang.org/u/pjentsch0)
#### Post date: [July 5, 2020, 9:07pm UTC](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668/41 "2020-07-05T21:07:16Z")

</div>

I’m having trouble getting similar results with KissABC 2.0 compared to v1.3, will there be more documentation coming for the various new parameters?

Does the method `AIS()` correspond to `ABCDE` in some way?

Thanks

[Previous page](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668.md?page=1)

[Next page](https://discourse.julialang.org/t/ann-rfc-kissabc-jl-approximate-bayesian-computation/40668.md?page=3)
