MonteCarloMeasurements.jl: Nonlinear Propagation of Arbitrary Multivariate Distributions by means of Method Overloading
This manuscript outlines a software package that facilitates working with probability distributions by means of Monte-Carlo methods, in a way that allows for propagation of multivariate probability distributions through arbitrary functions. We provid…
Authors: Fredrik Bagge Carlson
Abstract type P articles R eal 1 Introduction 𝑦 = 𝑓 (𝑥) = 𝐴𝑥 + 𝑏 𝑥 ∼ 𝑁 (𝜇, Σ) 𝑦 𝑦 ∼ 𝑁 (𝐴𝜇 + 𝑏, 𝐴Σ𝐴 T ) 𝑦 2 MonteCarloMeasurements.jl 𝑥 +,-,*,/,sin,exp +,-,*,/,sin,exp 𝑥 P articles{T,N} 𝑁 𝑇 P articles 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 𝑓 (𝑥) 𝑓 (𝑥) 𝑓 R eal P articles 𝑝 P articles 𝑝 𝑝 2.1 Basic Usage julia> using MonteCarloMeasurements , Plots , Distributions julia> a = π ± 0 . 1 # Construct Gaussian uncertain parameters using ± (\pm) P art500(3.142 ± 0.1) julia> b = 2 ∓ 0 . 1 # ∓ (\mp) creates StaticParticles (with StaticArrays) SP art100(2.0 ± 0.1) julia> std ( a ) # Ask about statistical properties 0.09997062445203879 julia> sin ( a ) # Use them like any real number P art500(1.255e-16 ± 0.0995) julia> sin ( a ) / cos ( a ) - tan ( a ) # Self-correlation is naturally handled P art500(0.0) julia> plot ( a ) # Plot them julia> b = sin .( 1:0 . 1:5 ) . ± 0 . 1 ; # Create multivariate uncertain numbers julia> plot ( b ) # Vectors of particles can be plotted julia> c = Particles ( 500 , Poisson ( 3 .)) # Create uncertain numbers distributed according to a given distribution ↪ P art500(2.896 ± 1.71) Complex{P articles} 2.2 Interaction wit h the Julia Ecosystem Normal(p) Normal 2.3 P er formance Adv antages 2.3.1 Shared computations function least_squares ( A , y ) Q = qr ( A ) return Q \ y end 𝑦 𝐴 𝑁 𝑁 𝐴 𝑦 𝑁 2.3.2 Optimal usage of SIMD instr uctions 2.4 Systematic Sampling 1/𝑁 𝑁 2.5 Multivariate Distr ibutions P articles P articles(N, d::Distribution) 𝑑 𝑝 𝑦 = 𝐴𝑝 𝑦 𝐴Σ𝐴 T julia> p = [ 1 ± 1 , 5 ± 2 ] # Create a vector of uncorrelated particles 2-element Array{Particles{Float64,500},1}: 1.0 ± 1.0 5.0 ± 2.0 julia> A = randn ( 2 , 2 ); # Create a random matrix julia> y = A * p # Transform particle vector 2-element Array{Particles{Float64,500},1}: -8.04 ± 3.1 2.4 ± 1.5 julia> cov ( y ) # Covariance of posterior multivariate particles 2×2 Array{Float64,2}: 9.61166 -3.59812 -3.59812 2.16701 julia> A * Diagonal ([ 1^2 , 2^2 ]) * A ' # Theoretical posterior covariance 2×2 Array{Float64,2}: 9.4791 -3.53535 -3.53535 2.15126 2.5.1 Sigma Points and t he Unscented T ransform sigma points 2.6 Limitations function negsquare ( x ) x > 0 ? x ^2 : - x ^2 end p = 0 ± 1 negsquare(p) 𝑥 > 0 negsquare P articles 3 Usage Examples and Benchmarks 3.1 Uncertainty Propagation—ODE simulation 𝑢 1 = 𝜃 𝑢 2 = − 𝑔 𝐿 sin (𝜃) 𝑔 = 9.7 9 ± 0.02 𝐿 = 1.00 ± 0.01 𝑢 0 = [0 ± 0, 𝜋 /3 ± 0.02] 𝜇 ± 𝜎 ∼ 𝑁 (𝜇, 𝜎 2 ) 𝑔 𝐿 0.0 0.1 0.2 0.3 0.4 0.5 0.00 0.25 0.50 0.75 1.00 Time [s] Linear MCM 195 196 197 198 199 200 - 3 - 2 - 1 0 1 2 3 Time [s] Linear MCM MCM Σ Float32 Linear MCM MCM Σ N aive MC 𝑢 0 2𝑛 + 1 𝑛 = 3 0 2 4 6 8 Σ Algorit hm 1 ± P articles{Float64,500} using MonteCarloMeasurements , OrdinaryDi ff Eq , Plots g = 9 . 79 ± 0 . 02 # Gravitational constant L = 1 . 00 ± 0 . 01 # Length of the pendulum u ₀ = [ 0 . 0 ± 0 . 0 , π / 3 . 0 ± 0 . 02 ] # Initial speed and initial angle function pendulum ( u , u , p , t ) θ , θ = u [ 1 ], u [ 2 ] u [ 1 ] = θ u [ 2 ] = - ( g / L ) sin ( θ ) end prob = ODEProblem ( pendulum , u ₀ , ( 0 . 0 , 2 . 0 )) sol = solve ( prob , T sit5 (), reltol = 1 e -6 ) plot ( sol ) ± P articles 3.2 MCMC poster ior P articles P articles P articles 3.3 Robust Optimization cost' cost Algorit hm 2 𝑐𝑥 + 𝑑𝑦 > 10 ∀ 𝑐, 𝑑 ∈ 𝑃 using MonteCarloMeasurements , Optim , Zygote const c = 1 ∓ 0 . 1 # These are the uncertain parameters const d = 1 ∓ 0 . 1 # These are the uncertain parameters function cost ( pars ) x , y = pars - ( 3 x +2 y ) + 10000* ( maximum ( c * x + d * y ) > 10 ) end pars = [ 1 ., 1 ] # Initial guess res = Optim . optimize ( cost , cost ' , pars , BFGS (), inplace = false ) 4 Concluding Remarks Ref erences Advances in Engineering Software SIAM review CoRR http : / / arxiv . org / abs / 1810 . 07951 JuliaSt ats/Distributions.jl: a Julia pac kage f or pr oba- bility distributions and associated functions 10.5281/zenodo.2647458 IEEE T ransactions on automatic contr ol Journal of Open Sour ce Softwar e 10.21105/joss.00615 The Journal of Open Researc h Softw are 10.5334/jors.151
Original Paper
Loading high-quality paper...
Comments & Academic Discussion
Loading comments...
Leave a Comment