|
1 | 1 | ---
|
2 |
| -title: Bayesian Inference on a Pendulum using Turing.jl |
| 2 | +title: Bayesian Inference on a Pendulum using DiffEqBayes.jl |
3 | 3 | author: Vaibhav Dixit
|
4 | 4 | ---
|
5 | 5 |
|
6 | 6 | ### Set up simple pendulum problem
|
7 | 7 |
|
8 | 8 | ```julia
|
9 |
| -using DiffEqBayes, OrdinaryDiffEq, RecursiveArrayTools, Distributions, Plots, StatsPlots |
| 9 | +using DiffEqBayes, OrdinaryDiffEq, RecursiveArrayTools, Distributions, Plots, StatsPlots, BenchmarkTools, TransformVariables |
10 | 10 | ```
|
11 | 11 |
|
12 | 12 | Let's define our simple pendulum problem. Here our pendulum has a drag term `ω`
|
@@ -76,7 +76,7 @@ length of the pendulum L is probably around 3.0:
|
76 | 76 | priors = [Uniform(0.1,3.0), Normal(3.0,1.0)]
|
77 | 77 | ```
|
78 | 78 |
|
79 |
| -Finally let's run the estimation routine from DiffEqBayes.jl using the Turing.jl backend |
| 79 | +Finally let's run the estimation routine from DiffEqBayes.jl with the Turing.jl backend to check if we indeed recover the parameters! |
80 | 80 |
|
81 | 81 | ```julia
|
82 | 82 | bayesian_result = turing_inference(prob1,Tsit5(),t,data,priors;num_samples=10_000,
|
@@ -104,3 +104,20 @@ plot(bayesian_result, colordim = :parameter)
|
104 | 104 | Notice that after awhile these chains converge to a "fuzzy line", meaning it
|
105 | 105 | found the area with the most likelihood and then starts to sample around there,
|
106 | 106 | which builds a posterior distribution around the true mean.
|
| 107 | + |
| 108 | +DiffEqBayes.jl allows the choice of using Stan.jl, Turing.jl and DynamicHMC.jl for MCMC, you can also use ApproxBayes.jl for Approximate Bayesian computation algorithms. |
| 109 | +Let's compare the timings across the different MCMC backends. We'll stick with the default arguments and 10,000 samples in each since there is a lot of room for micro-optimization |
| 110 | +specific to each package and algorithm combinations, you might want to do your own experiments for specific problems to get better understanding of the performance. |
| 111 | + |
| 112 | +```julia |
| 113 | +@btime bayesian_result = turing_inference(prob1,Tsit5(),t,data,priors;syms = [:omega,:L],num_samples=10_000) |
| 114 | +``` |
| 115 | + |
| 116 | +```julia |
| 117 | +@btime bayesian_result = stan_inference(prob1,t,data,priors;num_samples=10_000) |
| 118 | +``` |
| 119 | + |
| 120 | +```julia |
| 121 | +@btime bayesian_result = dynamichmc_inference(prob1,Tsit5(),t,data,priors,as(Vector, asℝ₊, 1)) |
| 122 | +``` |
| 123 | + |
0 commit comments