BHao_Discussion3

Problem 3.5.17

set.seed(123)

# create simulation function for one season
run_sim = function(days = 90) {
  # create data table to house inputs
  tbl_items = data.table(item = c('oats', 'peas', 'beans', 'barley'),
                  cost = c(1.05, 3.17, 1.99, 0.95),
                  price = c(1.29, 3.76, 2.23, 1.65),
                  minq = c(0, 0, 0, 0),
                  maxq = c(10, 8, 14, 11))
  tbl_items = tbl_items[, profit := price - cost]

  # create data frame to house results
  tbl_sim = data.frame(day = rep(0, days), total_revenue = rep(0, days), 
                       total_cost = rep(0, days), total_profit = rep(0, days))

  for (i in 1:days) {
    qty = NULL
    revenue = 0
    cost = 0
    profit = 0
    for (r in 1:nrow(tbl_items)) {
      qty = sample(tbl_items[r, minq]:tbl_items[r, maxq], 1)
      revenue = revenue + tbl_items[r, price] * qty
      cost = cost + tbl_items[r, cost] * qty
      profit = profit + tbl_items[r, profit] * qty
    }
    tbl_sim[i,] = c(i, revenue, cost, profit)
  }
  return(data.table(tbl_sim))
}

# create simulation function for x seasons
run_sim_season = function(seasons = 1000, days = 90) {
  tbl_daily = data.frame(day = rep(0, seasons*days), total_revenue = rep(0, seasons*days), 
                         total_cost = rep(0, seasons*days), total_profit = rep(0, seasons*days))
  tbl_season = data.frame(season = rep(0, seasons), total_revenue = rep(0, seasons), 
                          total_cost = rep(0, seasons), total_profit = rep(0, seasons))
  for (i in 1:seasons) {
    ## REPLACE LOGIC BELOW WITH A SIMPLE ARRAY AND THEN RETURN A DATA.TABLE 
    sim = run_sim(days = days)
    tbl_daily[(i*days-days+1):(i*days),] = sim
    tbl_season[i,] = c(season = i, sim[, lapply(.SD, sum), .SDcols = 2:4])
  }
  return(list(tbl_daily = data.table(tbl_daily), tbl_season = data.table(tbl_season)))
}

# run simulation x times and aggregate results
sim_results = run_sim_season(seasons = 1000, days = 90)

summary(sim_results$tbl_daily)

##       day       total_revenue     total_cost     total_profit   
##  Min.   : 1.0   Min.   : 0.00   Min.   : 0.00   Min.   : 0.000  
##  1st Qu.:23.0   1st Qu.:35.34   1st Qu.:28.02   1st Qu.: 6.800  
##  Median :45.5   Median :46.13   Median :37.05   Median : 9.090  
##  Mean   :45.5   Mean   :46.17   Mean   :37.08   Mean   : 9.089  
##  3rd Qu.:68.0   3rd Qu.:57.09   3rd Qu.:46.13   3rd Qu.:11.360  
##  Max.   :90.0   Max.   :92.35   Max.   :74.17   Max.   :18.180

summary(sim_results$tbl_season)

##      season       total_revenue    total_cost    total_profit  
##  Min.   :   1.0   Min.   :3612   Min.   :2910   Min.   :702.1  
##  1st Qu.: 250.8   1st Qu.:4056   1st Qu.:3255   1st Qu.:797.8  
##  Median : 500.5   Median :4157   Median :3339   Median :818.8  
##  Mean   : 500.5   Mean   :4156   Mean   :3338   Mean   :818.0  
##  3rd Qu.: 750.2   3rd Qu.:4251   3rd Qu.:3417   3rd Qu.:838.6  
##  Max.   :1000.0   Max.   :4610   Max.   :3741   Max.   :905.8

sim_results$tbl_daily %>% select(-day) %>%
  gather('metric', 'value') %>%
  mutate('metric' = factor(metric, levels = c('total_revenue', 'total_cost', 'total_profit'))) %>%
  ggplot(aes(x = value)) +
  geom_histogram() +
  facet_grid(. ~ metric, scales = 'free')

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

sim_results$tbl_season %>% select(-season) %>%
  gather('metric', 'value') %>%
  mutate('metric' = factor(metric, levels = c('total_revenue', 'total_cost', 'total_profit'))) %>%
  ggplot(aes(x = value)) +
  geom_histogram() +
  facet_grid(. ~ metric, scales = 'free')

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