vpl_sim: Simulação com o Método de Monte Carlo

View source: R/involutivo.R

vpl_simR Documentation

Simulação com o Método de Monte Carlo

Description

Simulação com o Método de Monte Carlo

Usage

vpl_sim(
  Nsim,
  ranges,
  variables,
  distribution = "uniform",
  params,
  dependencia = diag(length(ranges))
)

Arguments

Nsim

Número de simulações

ranges

Intervalos de variação de cada variável

variables

Variáveis utilizadas para o computo do VPL

distribution

Distribuição a priori a ser utilizada para gerar as variáveis

params

parâmtros a serem utilizados pelas distribuições a priori.

dependencia

matriz de covariância entre as variáveis

Examples

## Simulacao de Monte Carlo com distribuição uniforme e dependencia total

set.seed(1)
ranges <- list(vgv = c(min = 0.9, max = 1.1),
               cc = c(min = 0.9, max = 1.1),
               bdi_i = c(min = 0.9, max = 1.1),
               bdi_c = c(min = 0.9, max = 1.1)
)
variables <- list(vgv = vgv, wv = wv, cc = cc, wc = wc, bdi_i = bdi_i,
                  bdi_c = bdi_c, cor = cor, tma = tma)
dependencia100 <- matrix(data = c(1, -1, -1, -1,
                                  -1, 1, 1, 1,
                                  -1, 1, 1, 1,
                                  -1, 1, 1, 1),
                       nrow = 4, byrow = TRUE,
                       dimnames = list(names(ranges), names(ranges)))

Nsim <- 500

vpl_unif100 <- vpl_sim(Nsim, ranges = ranges, variables = variables,
                   distribution = "uniform", dependencia = dependencia100)
m_unif100 <- mean(vpl_unif100$vpl)
std_unif100 <- sd(vpl_unif100$vpl)
hist(vpl_unif100$vpl, freq = FALSE)
curve(dnorm(x, mean = m_unif100, sd = std_unif100),
        col = "darkblue", lwd = 2, add = TRUE, yaxt = "n")
summary(vpl_unif100$vpl)

# Its possible to compute the probabilities based in the simulations:
mean(vpl_unif100$vpl < 0.85*m_unif100) # probability that the VPL is less
than 85% of the mean

# Or compute the probabilities based on the normal curve with mean and sd
equals to that of the simulation.
pnorm(0.85*m_unif100, mean = m_unif100, sd = std_unif100)

## Simulacao de Monte Carlo com distribuição uniforme e dependencia 50%

dependencia50 <- matrix(data = c(1, -.5, -.5, -.5,
                                -.5, 1, .5, .5,
                                -.5, .5, 1, .5,
                                -.5, .5, .5, 1),
                         nrow = 4, byrow = TRUE,
                         dimnames = list(names(ranges), names(ranges)))

vpl_unif50 <- vpl_sim(Nsim, ranges = ranges, variables = variables,
                   distribution = "uniform", dependencia = dependencia50)
m_unif50 <- mean(vpl_unif50$vpl)
std_unif50 <- sd(vpl_unif50$vpl)
hist(vpl_unif50$vpl, freq = FALSE)
curve(dnorm(x, mean = m_unif50, sd = std_unif50),
        col = "darkblue", lwd = 2, add = TRUE, yaxt = "n")
summary(vpl_unif50$vpl)

## Simulacao de Monte Carlo com distribuição uniforme e variáveis totalmente
independentes

dependencia0 <- diag(4)
dimnames(dependencia0) <- list(names(ranges), names(ranges))

vpl_unif <- vpl_sim(Nsim, ranges = ranges, variables = variables,
                   distribution = "uniform", dependencia = dependencia0)
m_unif <- mean(vpl_unif$vpl)
std_unif <- sd(vpl_unif$vpl)
hist(vpl_unif$vpl, freq = FALSE)
curve(dnorm(x, mean = m_unif, sd = std_unif),
          col = "darkblue", lwd = 2, add = TRUE, yaxt = "n")
summary(vpl_unif$vpl)

## Simulacao de Monte Carlo com distribuição beta e dependencia 100%
params <- list(vgv = c(shape1 = 2, shape2 = 2),
               cc = c(shape1 = 2, shape2 = 2),
               bdi_i = c(shape1 = 2, shape2 = 2),
               bdi_c = c(shape1 = 2, shape2 = 2)
)
vpl_beta <- vpl_sim(Nsim, ranges = ranges, variables = variables,
                     distribution = "beta", params = params,
                     dependencia = dependencia100)

m_beta <- mean(vpl_beta$vpl)
std_beta <- sd(vpl_beta$vpl)
ggplot(as.data.frame(vpl_beta), aes(vpl_beta)) +
 geom_histogram(aes(y =..density..), col="red", fill="green", alpha=.2) +
  stat_function(
    fun = dnorm,
    args = list(mean = m_beta, sd = std_beta)
  )

lfpdroubi/appraiseR documentation built on April 14, 2024, 10:27 p.m.