Simulation Modeling Functions

Overview

Though most of your modeling can be realized in raw Julia, some of the most important features in a Monte-Carlo simulation package have to do with analyzing and applying correlation in models. MCHammer's correlation approach is based on "Ronald L. Iman & W. J. Conover (1982) A distribution-free approach to inducing rank correlation among input variables, Communications in Statistics - Simulation and Computation"

The simulation and correlation functions are designed to quickly obtain risk and decision analysis metrics such as moments, percentiles and risk over time.

Risk Events and Conditional Distributions

Risk Events allow you to model a joint distribution accounting for the probability of it occurring and the conditional impact when it does. The process simulates the Probability x Impact formula correctly.

MCHammer.RiskEvent — Function

RiskEvent(Prob, Distribution, Trials; seed=0)

Risk Events are defined as conditional distributions that will inflate the 0.

The RiskEvent() allows you to conditionally sample any distribution for as many trials as defined. Prob is the conditional probability of sampling the impact Distribution Distribution is any univariate distribution from Distributions.jl Trials is the number of total iterations. seed allows you to set a seed for the RiskEvent. Left blank or set to 0, the seed is set to random.

# Simulate a conditional risk event with a 30% chance of occurring and an impact
# that is distributed along a standard Normal. 10 trials are generated and about
# 3 should have non-zero outcomes.


RiskEvent(0.3, Normal(0,1), 10)

10-element Array{Float64,1}:
  0.0
  0.9045264902106298
  2.0679171427126892
  1.547808345914827
 -0.4315039376200465
  0.0
 -1.5737165904503339
 -0.0
  0.0
  0.0

Correlation and Covariance

MCHammer.cormat — Function

cormat(ArrayName, RankOrder=1)

Cormat calculates a symetric correlation matrix using both PPMC and Rank Order. Rank Order is default because this is what it used in the Iman-Conover method for correlating of simulated variables.

RankOrder = 1 calculates the Spearman rank order correlation used in MCHammer (this argument is optional and defaults to Spearman)

RankOrder = 0 calculates the Pearson Product Moment Correlation

Random.seed!(1) #hide
test = rand(Normal(),1000,5)
cormat(test)

# output

5×5 Array{Float64,2}:
  1.0         0.045012   0.00247197  -0.0455839   0.0126131
  0.045012    1.0        0.0534       0.0449149   0.0219751
  0.00247197  0.0534     1.0          0.0194396   0.0504692
 -0.0455839   0.0449149  0.0194396    1.0        -0.0301272
  0.0126131   0.0219751  0.0504692   -0.0301272   1.0

MCHammer.covmat — Function

covmat(ArrayName)

Calculates the covariance matrix.

covmat(test)

# output

5×5 Array{Any,2}:
  1.00069    0.0286309  0.0102903  -0.0356815   0.0132867
  0.0286309  1.09233    0.0598539   0.0536936   0.0216141
  0.0102903  0.0598539  1.07241     0.0140108   0.0583517
 -0.0356815  0.0536936  0.0140108   0.942015   -0.0240827
  0.0132867  0.0216141  0.0583517  -0.0240827   1.02767

Correlating simulation variables

MCHammer.corvar — Function

corvar(ArrayName, n_trials, correl_matrix)

The corvar function correlates simulation inputs unsing the Iman Conover Method. Your array must contain >2 simulated inputs. Remember to hcat() your inputs into tables reflecting your input correlation matrices.

n_trials: is the number of trials in the simulation. This must be consistent.

correl_matrix: must be defined as a Square Positive Definite correlation matrix. This can be calculated from histroical data using cormat() function.

Random.seed!(1)
n_trials = 1000
sample_data = [rand(LogNormal(0, 0.5),n_trials) rand(Normal(3,2),n_trials) rand(Gamma(1, 0.5),n_trials) rand(LogNormal(0, 0.5),n_trials) rand(Normal(3,2),n_trials) rand(Gamma(1, 0.5),n_trials)]

test_cmatrix = [1 0 0 0 0 0; 0 1 0 0 0 0; 0 0 1 0 0 0;0 0 0 1 0.75 -0.7; 0 0 0 0.75 1 -0.95; 0 0 0 -0.7 -0.95 1 ]

Random.seed!(1)
cormat(corvar(sample_data, n_trials, test_cmatrix))

# output

6×6 Array{Float64,2}:
  1.0          0.045012    0.00247197  -0.0455839  -0.0138308   0.0112554
  0.045012     1.0         0.0534       0.0449149   0.0592791  -0.0355262
  0.00247197   0.0534      1.0          0.0194396   0.0532426  -0.0468971
 -0.0455839    0.0449149   0.0194396    1.0         0.719585   -0.662708
 -0.0138308    0.0592791   0.0532426    0.719585    1.0        -0.939008
  0.0112554   -0.0355262  -0.0468971   -0.662708   -0.939008    1.0

Analyzing Simulation Results

MCHammer.GetCertainty — Function

GetCertainty(ArrayName, x, AboveBelow=0)

This function returns the percentage of trials Above (1) or Below(0) a target value of x.

Random.seed!(1)
test = rand(Normal(),1000)

GetCertainty(test, 0, 1)

# output

0.502

MCHammer.fractiles — Function

fractiles(ArrayName, Increment=0.1)

The fractiles function calculates percentiles at equal increments. The default optional argument for Increments is 0.1 for deciles but can be set to anything such as 0.05 for quintiles or 0.01 for percentiles.

fractiles(test)

# output

11×2 Array{Any,2}:
 "P0.0"    -3.882
 "P10.0"   -1.34325
 "P20.0"   -0.860748
 "P30.0"   -0.526089
 "P40.0"   -0.274806
 "P50.0"    0.00474446
 "P60.0"    0.218685
 "P70.0"    0.472055
 "P80.0"    0.800966
 "P90.0"    1.2504
 "P100.0"   3.12432

Misc functions

MCHammer.cmd — Function

cmd(x)

Shell /Dos Command wrapper to run batch and shell commands in script. This is used to process SQL from the command line or perform system level operation in a script using a command prompt.

MCHammer.truncate_digit — Function

function  truncate_digit(num, digits=2)

Truncation algorithim to remove decimals (ported by anonymous author from Maple) e.g.

  0.066 = 0.06
  0.063 = 0.06

Result_1 = truncate_digit(0.667)
Result_2 = truncate_digit(0.661)
Result_1 == Result_2

# output
true