What is a Monte Carlo Simulation? Well, think about it as a computation process that utilized random numbers to derive an outcome(s). So instead of having fixed inputs, probability distributions are assigned to some or all of the inputs. This will generate a probability
distribution for the output after the simulation is ran.
Here is an example. A firm that sells product X under a pure/perfect competition market* wants to know the probability distribution for the profit of this product and the probability that the firm will loss money when marketing it.
The equation for the profit is: TP = TR - TC = (Q*P) - (Q*VC+FC)
Assumption:
* The Quantity Demanded (Q) flucturates between 8,000 and 12,000 units and is uniformly distributed.
* Variable Cost (VC) is normally distributed (with mean = 7, Sd = 2) truncated on both sides (with a minimum of 7 / 2 and a maximum of 10).
* Market Price (P) is normally distributed (with mean = 10, Sd = 3) truncated on the left-hand side (with a minimum of 1).
* Fixed Cost (FC) is $5,000.
Also See Log-Normal Probability Distribution
* Under perfect competition market, the firm does not have the influence to affect the price of this product - the firm takes the market price as a given, dP/dQ = 0.
Outcome:
The average profit for this investment is $29,546 as shown on cells G25 after 50,000 iteration is ran. The probability that the profit of the investment turns out to be negative (loss money) is 22.28% as shown on cell C24. The probability distribution of the profit > X is display on column F and G. For example, there is 65% of chance that the profit will be greater than $12,481. The probability distribution is quite normal as shown on the figure. The mean is also very close to the median. This is due to the probability distribution that we assigned to the variables.
************************************************************************************************
************************************************************************************************
Option Explicit
Option Base 1
'**********************************************************************************
'* Run Monte Carlo Simulation *
'**********************************************************************************
Sub MonteCarlo()
Dim Iteration As Long, i As Long
Dim Q As Double, P As Double, TR As Double
Dim VC As Double, FC As Double, TC As Double
Dim SdVC As Double, MeanVC As Double, SdP As Double, MeanP As Double
Dim MinQ As Double, MaxQ As Double, AverageTP As Double, SumTP As Double
Dim ProfitX As Double, CountNo As Double
Iteration = Range("C3").Value
FC = Range("C7").Value
MinQ = Range("C13").Value
MaxQ = Range("C14").Value
MeanVC = Range("C15").Value
SdVC = Range("C16").Value
MeanP = Range("C17").Value
SdP = Range("C18").Value
ProfitX = Range("B24").Value
ReDim TP(Iteration) As Double
SumTP = 0
CountNo = 0
For i = 1 To Iteration: Cells(12, 3) = i
VC = Truncate_Normal_VC(MeanVC, SdVC, MeanVC / 2, MeanP)
P = Truncate_Normal_P(MeanP, SdP, 1)
Q = Int((MaxQ - MinQ + 1) * Rnd + MinQ)
TC = FC + VC * Q
TR = P * Q
TP(i) = TR - TC
'Comment out the following will make the simulation run faster
Cells(5, 3) = Q
Cells(6, 3) = P
Cells(8, 3) = VC
Cells(9, 3) = TC
Cells(10, 3) = TR
Cells(11, 3) = TP(i)
If TP(i) > ProfitX Then CountNo = CountNo + 1
SumTP = SumTP + TP(i)
Next i
AverageTP = SumTP / Iteration
Cells(25, 7) = AverageTP
Cells(24, 3) = 1 - CountNo / Iteration
Call Sort(Iteration, TP)
Call Hist(Iteration, 40, TP(1), TP(Iteration), TP)
For i = 1 To 20
Cells(i + 3, 6) = 1 - (0.05 * i)
Cells(i + 3, 7) = TP(Int(Iteration / 20 * i))
Next i
Cells(3, 6) = "Close to 100%"
Cells(13, 6) = "Median = 50%"
Cells(23, 6) = "Close to 0%"
Cells(3, 7) = TP(1)
End Sub
'**********************************************************************************
'* Return random numbers from a Truncated Normal Distribution for VC *
'**********************************************************************************
Function Truncate_Normal_VC(MeanX, SdX, leftLimit, RightLimit)
Dim x As Double
Dim fac As Double, r As Double, V1 As Double, V2 As Double
5 x = gauss * SdX + MeanX
If RightLimit < x Or x < leftLimit Then GoTo 5
Truncate_Normal_VC = x
End Function
'**********************************************************************************
'* Return random numbers from a Truncated Normal Distribution for Price *
'**********************************************************************************
Function Truncate_Normal_P(MeanX, SdX, leftLimit)
Dim x As Double
Dim fac As Double, r As Double, V1 As Double, V2 As Double
5 x = gauss * SdX + MeanX
If x < leftLimit Then GoTo 5
Truncate_Normal_P = x
End Function
'**********************************************************************************
'* Return random numbers from Standard Normal Distribution *
'**********************************************************************************
Function gauss()
Dim fac As Double, r As Double, V1 As Double, V2 As Double
10 V1 = 2 * Rnd - 1
V2 = 2 * Rnd - 1
r = V1 ^ 2 + V2 ^ 2
If (r >= 1) Then GoTo 10
fac = Sqr(-2 * Log(r) / r)
gauss = V2 * fac
End Function
'***********************************************************************************
'* Sort the numbers *
'***********************************************************************************
Sub Sort(n As Variant, arr() As Double)
Dim Temp As Double
Dim i As Long
Dim j As Long
For j = 2 To n
Temp = arr(j)
For i = j - 1 To 1 Step -1
If (arr(i) <= Temp) Then GoTo 10
arr(i + 1) = arr(i)
Next i
i = 0
10 arr(i + 1) = Temp
Next j
End Sub
'**********************************************************************************
'* Construct Historgram Distribution *
'**********************************************************************************
Sub Hist(n As Variant, M As Long, Start As Double, Right As Double, arr() As Double)
Dim i As Long, j As Long, Find As Long
Dim Length As Double
ReDim breaks(M) As Single
ReDim freq(M) As Single
For i = 1 To M
freq(i) = 0
Next i
Length = (Right - Start) / M
For i = 1 To M
breaks(i) = Start + Length * i
Next i
For i = 1 To n
If (arr(i) <= breaks(1)) Then freq(1) = freq(1) + 1
If (arr(i) >= breaks(M - 1)) Then freq(M) = freq(M) + 1
For j = 2 To M - 1
If (arr(i) > breaks(j - 1) And arr(i) <= breaks(j)) Then freq(j) = freq(j) + 1
Next j
Next i
For i = 1 To M
Cells(i + 1, 9) = breaks(i)
Cells(i + 1, 10) = freq(i)
Next i
End Sub
