Simulación por el Método Montecarlo

Tal como se ha señalado en los textos previos el método Montecarlo es un método fundamentalmente estocástico. Es exactamente lo opuesto a un método deterministico ( o si lo prefieres que se expresa con un modelo matemático).

Generlidades

La familia de los métodos montecarlo son especialmente interesantes cuando tenemos un modelo matemático determinístico. Sin embargo pueden ser usados aún en casos en los que no tengamos el modelo matemático y esto es la principal ventaja de los métodos estocásticos frente a los determinísticos. Utilizar Montecarlo es declarar que no tenemos en realidad la más mínima idea de como se relacional las variables del problema y apenas podemos o queremos tener un marco de referencia para ver que riesgo asumimos al decidir.

En este texto analizaremos tres problemas típicos en los que se puede partir del modelo determínistico matemático y enriqueceremos el análisis con Montecarlo. La técina es utilizada por gran cantidad de paquetes de software comercial, entre ellos los más notables son Crystal Ball de Oracle y @Risk de Palizade Software.

Los casos que veremos tienen que ver con el modelo de Wilson para lote económico, análisis del VAN / TIR y Programación por Camino Crítico. Por una cuestión de practiciad repasaremos el modelo determinístico sólo en el caso del modelo de Wilson y dejamos al lector el repaso de los otros modelos para tener un punto de referencia con el la simulación por el método Montecarlo.

Flujo de Trabajo

Para trabajar en este método es necesario contar con un dataset (conjunto de datos históricos) o en su defecto con un modelo formal estocástico que nos permita generar ese dataset.

El mecanismo de trabajo práctico implica un enriquecimiento de la información que nos ofrece Montecarlo respecto a métodos como los que nos aporta un modelo formal (determinístico).

A efectos de poder comprender mejor esto que se expresa iniciaremos con el estudio de un modelo determinístico, para luego pasar a tomar ese modelo y someterlo al análisis ampliado del Método Montecarlo.

Seguiremos la secuencia o flujo de trabajo que se detalla:

A- Modelar un diagrama de influencia

B- Obtener Datos Históticos de las Variables

C- Si no es posible armar el dataset buscar en la bibliografía modelos que contemplen las variables

D- Si no es posible B y/o C Modelar el problema

E- Validar el modelo (alimentar con datos históricos para obtener valores del presente)

F- Alimentar el Modelo con datos del presente y predecir el futuro

G- Simular con Montecarlo (alimentar con datos del presente alterados por ruido estadístico) y obtener un resultado estadístico

H- Determinar el pronóstico, el nivel de confianza y establecer el riesgo del negocio o problema simulado con Montecarlo

Diagrama de Influencia

En sí mismo esta técnica ya puede considerarse un modelado del problema. Es un modelado muy simplificado, pero la mayor parte de las veces que se realiza, sobre todo si se hace en una reunión, los resultados pueden ser más que sorprendentes.

En la técnica de simulación por Dinámica de Sistemas hay un software desarrollado por el MIT que realiza este tipo de diagrama. El software se llama Ventana Simulator y puede descargarse desde:

https://vensim.com

Diagrama de Influencia

Este diagrama se construye a partir del hexágono ROJO en el que se tiene a la variable que se pretende controlar u optimizar. Dentro de los cuadriláteros hay algunos que representan ratios o índice de tasas, algunos controlables por la empresa y otros independientes (tal como los impuestos). Finalmente las elipses representan variables que influyen directa o indirectamente sobre el hexágono y los colores refieren a la categorización de externas (independientes, inciertas, exógenas, o no controlables); e internas (dependientes o en las que es posible influir en parte o sobre las que tenemos certidumbre).

Modelos determinísticos

Para ofrecer un ejemplo de los pasos del plujo de trabajo abordaremos un caso de gestión de inventarios.

El Caso de la Freight Forwarder

Un caso de emprendedora exitoso en la región vitı́cola de Mednoza - Argentina está ligado a una egresada de la carrera de ingenierı́a industrial que a comienzos de la década de 1990. Ante la incertidumbre de las crisis económicas decidió, sin capital alguno, lanzarse como free lance en el negocio de mover containers con vino a los puertos de Valparaı́so y San Antonio - Chile. En los años 2000 esto se conoció como el negocio de Freight Forwarder en el que todos los activos de trabajo son subcontratado.

Su empresa actuaba como intgradora de servicios entre cada eslabon de la cadena de suministros, entre ellos transportistas, despachantes de aduana, bodegas, navieras y clientes de brokers en el extranjero.

Un factor clave de éxito ( KPI ) era el tamaño de lote óptimo a transportar en cada envı́o. Es bien conocido el modelo de Wilson que data de la década de 1920, pero que hoy resulta casi una garantı́a al fracaso si se utiliza en contextos tan turbulentos como los que nos toca vivir. A pesar de ello el modelo aún se estudia y no es descartado, pues sirve de base para ser mejorado. El modelo tiene tres variables independientes y una dependiente, a saber:

Variables Independientes

* Demanda Media Anual D (en Containers por año)

* Costo de almacenamiento unitario Ca ($ por mes)

* Costo de emisión o gestión Ce ($ por orden)

Variable Dependiente

* Tamaño de Lote Mensual EOQ

La recomendación del Método Montecarlo es salir a buscar en la bibliografía y en los papers modelos determinísticos que vinculen a estas variables. Cuanto más actuales sean los modelos mayor certeza de los resultados obtenidos tendremos.

Modelo Matemático de Wilson

\[ {EOQ} = \sqrt {(\frac {2*D * Ce} {Ca} )}\]

El caso de la Bodega Torrenti operada por el Freight Forwarder

Esta bodega tiene tres destinos principales de exportación.

  • Estados Unidos demanda 500 TEUs al año
  • Brasil demanda 53 TEUs al año
  • Cadandá demand 28 TEUs al año

  • Costo de Almacenamiento U$D 2.500

  • Costo de Gestión U$D 750

Modelado numérico en R-CRAN

Diagrama Modelo Base

D1 <-500
D2 <-52
D3 <-28

D<- D1+D2+D3

Ca <- 2500
Ce <- 750

EOQ <- sqrt(2*D*Ce/Ca)

EOQ
## [1] 18.65476

La conclusión final es que la bodega debería enviar 18 containers por mes para operar en un punto que minimice sus costos logísticos totales.

Lote Económico resuelto por Montecarlo

Este método no concibe a las variables independientes como fijas, sino que las considera variables aleatorias e introduce Ruido en el modelo matemático. En esta parte deberíaos analizar el comportamento de cada variable y asociarle una distribución (función paramétrica) de probabilidad. Por sencillez supondremos en este ejemplo que las varaibles siguen una distribución normal.
En rigor deberíamos estudiar que distribución siguen los datos con los que alimentaremos el modelo. <bR R-Cran tiene por defecto definidas las siguientes distribuciones de probabilidad.

  • Normal
  • Logística
  • Uniforme
  • Gamma
  • Lognormal
  • Weibull
  • Cauchy
  • Exponencial
  • Chi-cuadrado
  • F
  • T-Student.


Puede explorar el paquete rriskDistributions que le dará una idea más clara. Generaremos una muestra al azar de 47 valores con media 40 y desvio standard 5

#    library(rriskDistributions)
#    muestra=rnorm(47,40,5)
#    muestra
#    fit.cont(muestra)

Estocatización de Variables

  • Estados Unidos demanda promedio 500 TEUs al año con una varianza de 90

  • Cadandá demanda promedio 28 TEUs con una varianza de 6 al año

  • Brasil demanda promedio 93 TEUs al año con una varianza de 8


  • Costo de Almacenamiento U$D 2.500 Varianza 1.000

  • Costo de Gestión U$D 750 varianza 250

Cálculo del EOQ Montecarlo

#Genero tres variables aleatorias normales con 
#una cantidad de muestras m igual a tres veces #el tamaño de la media 

m1 <- 500*3
D1 <- rnorm(m1,500,90)
m2 <- 93*3
D2 <- rnorm(m2,93,8)
m3 <- 28*3
D3 <- rnorm(m3,28,6)
D <- c(D1,D2,D3)
# m_d muestras de demanda
m_d <- length(D)

#Generamos las variables aleatorias de costos igual a la cantidad de muestras de demanda que tenemos.

Ca <- rnorm(m_d,2500,1000)

Ce <- rnorm (m_d,750,250)

# Preparo área de graficos 3 filas 1 columna
#par(mfcol=c(1,3))

#Grafico 1
plot(density(D), 
xlab="Demanda TEU", #Change the x-axis label
ylab="Densidad", #y-axis label
main="Muestras de Demanda 5 años")#Main title

#Grafico 2
plot(density(Ca),
xlab="Costo Almacenamiento", #x-axis label
ylab="% eventos registrados", # y label
main="Población Costos Almacenamiento")     

#Grafico 3
plot(density(Ce),
xlab="Costo Gestión Inventario", #x-axis label
ylab="% eventos registrados", # y label
main="Población Costos Emisión")     

Cálculo de muestras de lote Económico EOQ

En este caso realizaremos las mismas operaciones que antes, pero las variables independientes ahora serán matrices y el resultado será un conjunto de valores con su propia distribución de probabilidades.

EOQ <- sqrt(2*D*Ce/Ca)
## Warning in sqrt(2 * D * Ce/Ca): NaNs produced
#Media del Lote Económico
mean(EOQ, na.rm=TRUE)
## [1] 16.43874
#Varianza del Lote Económico
var(EOQ, na.rm=TRUE)
## [1] 110.7694
#Desvio Estandard del Lote Económico
sd(EOQ, na.rm=TRUE)
## [1] 10.52471
plot(density(EOQ,na.rm=TRUE))

Modelo Montecarlo

Hallaremos los valores máximos y mínimos de las variables dependientes e independientes

#Índice Mínimos
which.min(Ca)
## [1] 36
which.min(Ce)
## [1] 507
which.min(EOQ)
## [1] 1830

Valores Mínimos

Ca[which.min(Ca)]
## [1] -1317.145
Ce[which.min(Ce)]
## [1] -326.3086
EOQ[which.min(EOQ)]
## [1] 1.322674
#Máximos
which.max(Ca)
## [1] 959
which.max(Ce)
## [1] 1754
#Atención excuiremos los NA
which.max(EOQ)
## [1] 322

Riesgo del Modelo

Generación de Dataset de Inventario Simulado

Se ha visto como a partir del modelo numérico pudimos hacer el trabajo de simulación. Si tuviésemos datos históricos con solamente calcular los parámetros \(\mu\) y \(\sigma\) podríamos reporducir el experimento y la simuación En algunas oportunidades las estadísticas como las de kaggle, CEPAL, Banco Mundial, etc. nos entregan los parámetros para construir el historial. Veremos ahora cómo generar el dataset a partir de los parámetros.

Para poder realizar este tipo de simulación es necesario que tengamos un historial de los costos de gestión de un inventario. En este ejercicio estamos pensando en un historial de datos que varios años. Por caso plantemos un ejemplo que abarca desde 1950 a 2021.

Si el caso es utilizar datos reales, todos los métodos de simulación requieren que se haga un trabajo de tratamiento y depuración de la información. Tendremos que realizar cosas y tomar decisiones que implican:

    1. Remover datos de circunstancias excepcionales o capturados erróneamente.
    1. Revisar los NA
    1. Decidir como se completaran los NA (¿Serán reemplazados por el promedio anual?, ¿Se interpolará?, ¿Se tomará una aproximación basada en el mismo dato de la misma semana del años pasado incrementado por la tendencia?, etc)
    1. Reralizar una actualización a precios constantes si corresponde
    1. Realizar entrevistas y reuniones con expertos para validar
    1. Comparar con datos de cámaras o estaddísticas nacionales Inicia en 1950, finaliza en Dic 2021
    1. Repasar la consistencia metódica de los datos

Demanda

Constituiremos el primer dataset que representa la demanda de containers adquiridos por la empresa

D1 <- round(rnorm(426,300,30))
D2 <- round(rnorm(426,200,50))
D <- c(D1,D2)
hist(D)

Costo de Emisión

Ce1 <- rnorm (426,30,4)
Ce2 <- rnorm (426,25,4)
Ce <- c(Ce1,Ce2)
plot(density(Ce))

Costo de Almacenamiento

Ca1 <- rnorm(426,72,12)
Ca2 <- rnorm(426,95,8)
Ca <- c(Ca1,Ca2)
plot(density(Ca))

Recuperación de datos desde la Web

Inventario Simulado

Si no deseas generar el datasets, puedes cargarlo con load(“Wilson.Rda”)

https://themys.sid.uncu.edu.ar/ /rpalma/R-cran/Inventario_Simulado/Wilson.Rda

load ("Wilson.Rda")
library(psych)
pairs.panels(Wilson)

Lote Económico

EOQ <- sqrt(2*D*Ca/Ce)
hist(EOQ)

Este método que hemos usado utilizado nos permite crear un dataset semejante al que podríamos obtener si capturamos los datos de la operación de un negocio formal.

Tengamos en cuenta que todos los pasos anteriores no son importantes para el trabajo que vamos a realizar. De hecho la formula del lote económico de Wilson (EOQ) fue creado en 1912. Tiene más de 100 años y en la actualidad , como modelo, garantiza el fracaso en los cotos de gestión de inventario.

Discutiremos en clase más sobre este tema.

Guardaremos este archivo y a partir de el realizaremos una simulación por el método Montecarlo.

Grabar Archivo

Guardaremos en archivo en la carpeta que nos indica el comando getwd(). Lo haremos utilizando los formatos csv (que ya conocemos) y el formato propio de R-Cran denominado .Rda

getwd()
## [1] "/media/rpalma/Datos/2023/Posgrado/MBA/3er_Clase/Inventario_Simulado"
Wilson <- cbind(D,Ca,Ce,EOQ)
save (Wilson, file="Wilson.Rda")
write.csv(Wilson,file="Wilson.csv")

Análisis exploratorio

Realizaremos un análisis exploratorio de cada componente del modelo. Nuestro modelo tiene una variable Dependiente llamada EOQ, que representa la cantidad de contenedores pedidos para satisfacer la demanda. y tres variavles independientes:

  • Costo de Almacenamiento Ca
  • Costo de sostenimiento o gestión del inventario Ce
  • Demanda anual D

Costo de almacenamiento

summary(Ca)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   37.84   71.58   85.87   83.15   95.45  115.29

Valores máximos y mínimos

which.max(Ca)
## [1] 521
which.min(Ca)
## [1] 403

¿Que valor tiene el costo mínimo?

Ca[3]
## [1] 69.25299

¿Qué valores están por arriba del tercer cuanti?

which(Ca> 90) 
##   [1]   6  14  32  34  55  98 118 121 139 140 169 196 197 198 199 230 233 344
##  [19] 365 391 414 420 427 428 430 432 433 434 435 436 437 439 440 441 443 445
##  [37] 446 447 449 450 451 453 456 459 461 462 464 466 468 469 470 471 473 474
##  [55] 475 476 477 478 479 481 483 484 485 487 488 490 491 492 493 495 496 497
##  [73] 498 499 500 501 502 503 505 508 510 511 512 514 516 517 518 519 521 522
##  [91] 524 525 526 527 528 530 532 533 534 535 536 537 538 540 541 543 545 546
## [109] 547 548 549 551 552 555 557 558 559 561 562 563 564 565 566 567 568 569
## [127] 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 586 587 588
## [145] 589 590 591 592 593 594 595 596 598 599 600 601 602 603 606 607 609 610
## [163] 611 612 614 616 618 620 621 622 623 625 626 627 629 631 632 635 636 638
## [181] 639 640 641 642 644 646 647 648 649 655 656 657 659 661 663 664 665 666
## [199] 668 669 670 672 673 674 675 678 681 683 686 687 688 690 691 692 694 695
## [217] 696 697 698 699 700 701 702 703 704 705 706 708 709 710 711 712 715 716
## [235] 717 720 723 724 725 726 729 730 731 732 734 736 738 739 740 741 742 744
## [253] 745 747 750 751 753 754 756 757 759 761 763 764 766 767 768 769 770 771
## [271] 773 774 776 778 779 780 782 783 784 785 786 787 788 789 790 791 792 793
## [289] 794 795 796 797 798 800 801 802 803 804 805 806 807 809 811 812 813 814
## [307] 815 816 817 818 819 821 822 823 824 827 828 829 830 831 832 833 835 837
## [325] 838 839 840 841 842 843 845 846 847 848 850 851

Histogramas

Evolución de costos de almacenamiento

plot(sort(Ca), main ="costos ordenados", type="l")

DataExplorer

DataExplorer simplifica y automatiza el proceso EDA y la generación de informes. El paquete escanea automáticamente cada variable que realiza el perfilado de datos y ofrece varias funciones útiles para generar diferentes gráficos en características discretas y continuas.

Es importante destacar que el paquete permite generar un informe HTML completo sin esfuerzo al invocar la útil función create_report en un conjunto de datos (por ejemplo, create_report(df)). Es posible pasar argumentos adicionales, como una variable de respuesta y, para agregar varios análisis bivariados al informe:

Busca en la carpeta en la que has guardado este archivo. La biblioteca habrá creado un archivo llamado Reporte.html Las líneas están comentadas por que este comando no se puede correr dentro de un archivo Mackdown

# library(DataExplorer)
# Wilson %>% create_report(output_file = "Report")