Building your first model
Installing MCHammer
Install the package as usual using Pkg.
using Pkg
Pkg.("MCHammer")
If you need to install direct, we recommend using ']' to go in the native Pkg manager.
(v1.1) pkg> add https://github.com/etorkia/MCHammer.jl
Loading MCHammer
To load the MCHammer package
using MCHammer
Getting your environment setup for modeling
In order to build your first model, you will need to get a few more packages installed:
Distributions.jl : To build a simulation, you need distributions as inputs. Julia offers univariate and multivariate distributions covering most needs.
StatsBase.jl and Statistics.jl : These packages provide all the functions to analyze results and build models.
To load the support packages:
julia> using Distributions, Statistics, StatsBase, DataFrames
Building a simple example
EVERY MONTE CARLO MODEL HAS 3 COMPONENTS
- Inputs: Ranges or Single Values
- A Model: Set of mathematical relationships f(x)
- Outputs: The variable(s) of interest you want to analyze
Main Distributions for most modeling situations
Though the most used distributions are cite below, Julia's Distributions package has an impressive array of options. Please check out the complete library of distributions at Distributions.jl
- Normal()
- LogNormal()
- Triangular()
- Uniform()
- Beta()
- Exponential()
- Gamma()
- Weibull()
- Poisson()
- Binomial()
- Bernoulli()
In order to define a simulated input you need to use the rand function. By assigning a variable name, you can generate any simulated vector you want.
using Distributions
input_variable = rand(Normal(0,1),100)
100-element Array{Float64,1}: 0.49396120494898144 -2.4870283129165918 0.10925582834072056 -1.6518022744960785 0.403174916030235 0.5986110827594519 -0.6228945272767029 -0.22163629640554214 1.0655641482279776 -0.7986035728253849 ⋮ 1.0473636603577106 -0.5371622435326651 -1.2775729760328973 0.3324112256295076 1.2043968012107904 1.127021306414491 -0.3452499296043534 -1.0180713406985824 0.15525964099054929
Creating a simple model
A model is either a visual or mathematical representation of a situation or system. The easiest example of a model is
PROFIT = REVENUE - EXPENSES
Let's create a simple simulation model with 1000 trials with the following inputs:
Setup environment and inputs
using Distributions, StatsBase, DataFrames, MCHammer
n_trials = 1000
Revenue = rand(TriangularDist(2500000,4000000,3000000), n_trials)
Expenses = rand(TriangularDist(1400000,3000000,2000000), n_trials)
1000-element Array{Float64,1}: 1.8648501671239818e6 2.2251974334942773e6 2.6722253150278404e6 2.357605171301006e6 2.4028904622442266e6 2.5241743834304223e6 2.20253588791027e6 2.660244382520066e6 2.7900800843616487e6 2.6811719637771677e6 ⋮ 1.5918528020011084e6 1.9596502667604194e6 1.970887427773498e6 1.6600644350288461e6 1.7100699218384176e6 2.1496188478392866e6 2.022554128384966e6 1.859481542314854e6 2.0060707304267227e6
Define a Model and Outputs
# The Model
Profit = Revenue - Expenses
#Trial Results : the Profit vector (OUTPUT)
Profit
1000-element Array{Float64,1}: 1.5016031612712722e6 604659.7581555275 337059.57111664256 1.5798958104149648e6 997416.2776403851 412889.7564387056 487290.56057081604 952514.2037599627 616392.6264003362 561434.7237059474 ⋮ 1.3502740072905058e6 1.5832496495287719e6 1.3222483619034346e6 1.7869750731157272e6 1.6329631947467462e6 945383.6533520473 1.7502791890084809e6 1.106063015354592e6 1.5319127140207249e6
Analyzing the results in Julia
# `fractiles()` allows you to get the percentiles at various increments.
fractiles(Profit)
11×2 Array{Any,2}: "P0.0" -245007.0 "P10.0" 412380.0 "P20.0" 606819.0 "P30.0" 759795.0 "P40.0" 892232.0 "P50.0" 1.01079e6 "P60.0" 1.13439e6 "P70.0" 1.24538e6 "P80.0" 1.43235e6 "P90.0" 1.63644e6 "P100.0" 2.24382e6
density_chrt(Profit)
Sensitivity Analysis
First we need to create a sensitivity table with hcat() using both the input and output vectors.
#Construct the sensitivity input table by consolidating all the relevant
#inputs and outputs.
s_table = hcat(Profit, Revenue, Expenses)
#We then need to convert to a DataFrame and add names
s_table = DataFrame(s_table)
names!(s_table, [:Profit, :Revenue, :Expenses])
#To produce a sensitivity tornado chart, we need to select the output against
#which the inputs are measured for effect.
sensitivity_chrt(s_table, 1, 3)