Econ 455W/555: Monte Carlo Simulation Example

Setup

# Load required packages
# Load required packages
library(magrittr) # for data manipulation
library(ggplot2) # for data visualization

# Define input parameters
num_users <- rnorm(1000, mean = 5000, sd = 1000)
revenue_per_user <- runif(1000, min = 0, max = 5000)
operating_costs <- rnorm(1000, mean = 10000, sd = 200)
development_cost <- 2000000
discount_rate <- 0.06

Calculate net present value (NPV) for each sample

npv <- ((revenue_per_user * num_users) - operating_costs) / ((1 + discount_rate) ^ 1) +
  ((revenue_per_user * num_users) - operating_costs) / ((1 + discount_rate) ^ 2) +
  ((revenue_per_user * num_users) - operating_costs) / ((1 + discount_rate) ^ 3) +
  ((revenue_per_user * num_users) - operating_costs) / ((1 + discount_rate) ^ 4) +
  ((revenue_per_user * num_users) - operating_costs) / ((1 + discount_rate) ^ 5) - development_cost

Analyze the results

expected_value <- mean(npv)
uncertainty <- sd(npv)
confidence_interval <- c(expected_value - 1.96 * uncertainty / sqrt(length(npv)), expected_value + 1.96 * uncertainty / sqrt(length(npv)))

expected_value

## [1] 50978827

confidence_interval

## [1] 48965659 52991994

Plot the results

ggplot(data.frame(npv), aes(x = npv)) +
  geom_histogram(aes(y = ..density..), fill = "lightblue", color = "black") +
  geom_vline(xintercept = expected_value, color = "red", linetype = "dashed") +
  geom_vline(xintercept = 0) +
  geom_vline(xintercept = confidence_interval[1], color = "green", linetype = "dashed") +
  geom_vline(xintercept = confidence_interval[2], color = "green", linetype = "dashed") +
  labs(x = "Net Present Value", y = "Density", title = "Monte Carlo Simulation Results") +
    theme(legend.text = element_text(size = 15, color = "gray30"), legend.title = element_text(size = 16)) +
  theme(axis.text.x = element_text(size = 14, angle = 0, color = "gray30"),
        axis.text.y = element_text(size = 14, angle = 0, color = "gray30"),  
        axis.title.x = element_text(size = 16, angle = 0, color = "gray30"),
        axis.title.y = element_text(size = 16, angle = 90, color = "gray30")) +
  theme_minimal() + 
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank()) 

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.