Simulación de Sistemas - Tidy Simulation - Newspaper problem
Daniel Chavarría
9/10/2019
Comment is in spanish at the moment, let me know if you have any questions in twitter @TetricAttack! Have fun!
Enunciado - Problema del vendedor de periódicos:
Don Ricardo, compra periódicos a $1.300 cada uno y los vende a $1.600 la unidad. Al finalizar cada día la empresa de periódicos le paga a Don Ricardo $750 por cada periódico que no haya vendido. La demanda diaria de periódicos tiene una función de probabilidad como la siguiente:
Demanda diaria de periódicos
| (x) | Probabilidad p(x) |
|---|---|
| 25 | 0,12 |
| 30 | 0,08 |
| 35 | 0,3 |
| 40 | 0,25 |
| 45 | 0,1 |
| 50 | 0,1 |
| 55 | 0,05 |
Se desea determinar por simulación el número óptimo de periódicos que el vendedor debe ordenar cada día. Para ello realice la actividad entregable propuesta a continuación. Actividad entregable:
Debe simular 1000 días la venta de periódicos y responder las siguientes preguntas:
- ¿Cuántos periódicos debe ordenar diariamente Don Ricardo de tal forma que se maximice la utilidad esperada con la venta de los periódicos?
- Calcule estadística descriptiva de la demanda diaria de periódicos como: mínimo, máximo, media, mediana y desviación estándar e interprételas.
- Calcule un intervalo de confianza con una confianza del 95% para el promedio de las utilidades diarias de la venta de periódicos.
library(tidyverse)
prob_ventas <- tibble(ventas = c(25, 30, 35, 40, 45, 50, 55),
probabilidad = c(0.12, 0.08, 0.3, 0.25, 0.1, 0.1, 0.05)) %>%
mutate(acumulada = cumsum(probabilidad))En este caso:
Ganancias totales = Ingresos de venta - Costo de compra - Utilidad perdida por exceso de demanda + Ingresos precio reducido
Inicialmente exploremos el funcionamiento de la función sample()
sample(prob_ventas$ventas, prob = prob_ventas$probabilidad, size = 10000, replace = TRUE) %>%
hist(freq = F)
Las probabilidades parecen corresponder con la probabilidad de cada una, por lo que procedemos con la simulación.
Pregunta 1. ¿Cuántos periódicos ordenar?
k <- 10000 #Número de días a simular
simulacion_dias <- function(k, periodico) {
tibble(dia = 1:k,
demanda = sample(prob_ventas$ventas, prob = prob_ventas$probabilidad, replace = TRUE, size = k),
ingresos = case_when(periodico > demanda ~ demanda * 1600,
TRUE ~ periodico * 1600),
costos = periodico * 1300,
ganancia_perdida = case_when(periodico > demanda ~ (periodico - demanda) * (1600 - 1300),
TRUE ~ 0),
ingreso_remate = (periodico - demanda) * 750,
ganancia = ingresos + ingreso_remate - costos - ganancia_perdida) -> sim_dias
return(sim_dias)
}
simulacion_dias(k, 40) #Prueba de 1000 días simulados y 70 periódicos ordenados## # A tibble: 10,000 x 7
## dia demanda ingresos costos ganancia_perdida ingreso_remate ganancia
## <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 25 40000 52000 4500 11250 -5250
## 2 2 35 56000 52000 1500 3750 6250
## 3 3 35 56000 52000 1500 3750 6250
## 4 4 50 64000 52000 0 -7500 4500
## 5 5 35 56000 52000 1500 3750 6250
## 6 6 40 64000 52000 0 0 12000
## 7 7 40 64000 52000 0 0 12000
## 8 8 30 48000 52000 3000 7500 500
## 9 9 35 56000 52000 1500 3750 6250
## 10 10 40 64000 52000 0 0 12000
## # ... with 9,990 more rowsLa función simulacion_dias entrega una tabla completa de k días simulados con el número de periódicos ordenados. Ahora, sólo nos falta variar el número de periódicos y calcular las ganancias para cada una de estas tablas.
vector_ganancia <- vector()
for (i in 0:80) {
simulacion_dias(k, i) %>%
summarize(total = sum(ganancia)) %>%
as.numeric() -> ganancia_dia
vector_ganancia <- append(vector_ganancia, ganancia_dia)
}
resultados_sim <- tibble(periodicos_ordenados = 0:80, ganancia_total = vector_ganancia)
resultados_sim## # A tibble: 81 x 2
## periodicos_ordenados ganancia_total
## <int> <dbl>
## 1 0 -286953750
## 2 1 -276390000
## 3 2 -265053750
## 4 3 -255138750
## 5 4 -244275000
## 6 5 -232515000
## 7 6 -222836250
## 8 7 -212538750
## 9 8 -202608750
## 10 9 -191606250
## # ... with 71 more rowsEn este momento tenemos un resultado de simulación de 1000 días para el número de periódicos ordenados entre 0 y 80, veamos los resultados gráficamente.
resultados_sim %>%
ggplot(aes(x = periodicos_ordenados, y = ganancia_total)) +
geom_col() +
scale_y_continuous(labels = scales::dollar_format()) +
labs(title = "Número de periódicos ordenados y ganancia total asociada",
subtitle = "k = 10000 simulaciones",
x = "Número de periódicos ordenados",
y = "Ganancia total")
De esta gráfica concluimos que los números de periódicos a ordenar que nos interesa analizar se encuentran entre 30 y 50, con un máximo al parecer en 40, veámoslo nuevamente:
Con esta gráfica podemos responder la pregunta 1 - Don Ricardo debería ordenar 40 periódicos.
Pregunta 2. Estadística descriptiva de la demanda diaria
Para responder esta pregunta, debemos ver una de las tablas individuales para un número fijo de periódicos, analicemos el caso de 40 períodicos que corresponde a las ganancias máximas. mínimo, máximo, media, mediana y desviación estándar
set.seed(1995)
sim <- simulacion_dias(k, 40)
sim %>%
summarize(mínimo = min(demanda), máximo = max(demanda), mediana = median(demanda),
media = mean(demanda), desv_est = sd(demanda)) -> estadistica_descriptiva
estadistica_descriptiva## # A tibble: 1 x 5
## mínimo máximo mediana media desv_est
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 25 55 40 38.2 7.99Para el caso de 40 periódicos ordenados, se obtuvo una demanda diaria entre 25 y 55, que corresponden a los valores máximos y mínimos de la f.d.p. y una mediana de 40 periódicos ordenados, que significa que la mitad de los datos se encuentran por encima y por debajo respectivamente. Se concluye además que el promedio de periódicos ordenados en las 1000 simulaciones es de 38.1625 (a pesar de no ser un valor viable ya que las demandas son discretas) con una desviación estandar de 7.9885841 periódicos asociada que habla de una variabilidad no muy alta.
Pregunta 3. Calcule un intervalo de confianza con una confianza del 95% para el promedio de las utilidades diarias de la venta de periódicos.
library(broom)
sim$ganancia %>%
t.test(conf.level = 0.95) %>% tidy() -> conf_intervalo
conf_intervalo## # A tibble: 1 x 8
## estimate statistic p.value parameter conf.low conf.high method
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <chr>
## 1 5544. 103. 0 9999 5439. 5650. One S~
## # ... with 1 more variable: alternative <chr>Con una confianza del 95% para el caso de 40 periódicos, se estima que la media de la ganancia estará entre 5439.012068 y 5649.837932 pesos por día con una media de 5544.425.