Para as questões 1 e 2 considerar um valor de semente igual a matricula de um dos componentes da equipe.


Questão 1, Hipercubo Latino - HL

“Você está interessado em gerar uma amostra de valores exponencialmente distribuídos. Gere a amostra com e sem aplicação de técnica de redução de variância. A técnica de redução de variância deverá ser a do Hipercubo Latino. Depois de gerados os valores, peça o histograma e faça o teste Kolmogorov-Smirnov para ambas as amostras e comente. OBS: Você não pode utilizar as funções rexp e randomLHS.”

Método da Transformada Inversa

Para gerar uma variável aleatória que segue uma distribuição Exponencial de paramêtro \(\lambda\) podemos usar o Método da Transformada Inversa, onde, a partir de uma variável uniforme contínua \(U\), onde \(U \sim U(0,1)\) e a função de distribuição acumulada da variável aleatória, \(F_X(x)\), podemos gerar valores aleatórios de \(x\). A técnica usa do fato de que \(F_X(x) = F_X(F^{-1}(U))\), sendo \(F_X(x)\) contínua e monótona não descrecente.

Aplicando a função inversa de \(F_X(x)\) em ambos os lados temos que \(x = F^{-1}(U)\).

Vejamos como ficaria no caso da distribuição \(X \sim Exp(\lambda)\)

Função de Densidade: \[f(x) = \lambda e^{-\lambda x}\] Função de Distribuição: \[F(x) = 1 - e^{-\lambda x}\] Número aleatório: \[x = F^{-1}(u) \Leftrightarrow x = - \frac{\ln(u)}{\lambda}\] Agora usando o software RStudio, podemos aplicar esse resultado na função runif(n) gerando n valores que seguem um variável uniforme contínua (0,1). Para fins de maior esclareciemnto irei descrever passo a passo o que o código está fazendo por meio dos # Comentários

Código em R

Usando R:

set.seed(536837)    # Fixando um valor semente para facilitar a reprodutividade do experimento.
n = 100   # Definindo o número de valores gerados
lambda = 1    # Definindo o parâmetro
u = runif(n)    # Gerando n valores uniformes
X = -log(u)/lambda   # Aplicando o Método da Transformada Inversa

# Gerando o histograma de X

Visualização em R

Redução de Variância, HL

Além disso podemos usar técnicas que permitem reduzir a variância do processo, uma delas é o Hipercubo Latino, onde divide o intervalo de X em n partes iguais e em cada subintervalo faz uma amostragem aleatória.

Vejamos como na pratica isso funciona. Usando a técnica de redução de variância, hipercubo latino.

Código em R

Usando R:

HL <- c()  # Primeiramente defino um vetor vazio para receber os valores de X.
for(i in  1:n){ # Defino uma recursividade para gera os n valores.
  x <- -log((i - runif(1))/n)/lambda  # Usando o Método da Transformada Inversa no intervalo i
  HL[i] <- x  # Armazenando os valores de X, com a técnica de redução de variância, HL.
}

# Gerando o Histograma de X, com a técnica de redução de variância Hipercubo Latino.

Visualização em R

Teste Kolmogorov-Smirnov

Outro aspecto importante é verificar se os métodos estão gerando os valores de maneira satisfatória. O Teste Kolmogorov-Smirnov verifica se há evidência estatística que o vetor aleatório segue uma certa distribuição. Onde a hiposte nula é que a função de distribuição acumulada Fx é igual a alguma função de distribuição.

Usando o RStudio podemos fazer o teste por meio da função ks.test(x,F(x)).

# Testando se o vetor gerado segue uma distribuição Exp(1).
Teste = ks.test(X,pexp,1)
Teste$p.value
## [1] 0.2177719

Para a primeira amostra gerada, o p-valor foi 0,2178 sendo assim, não podemos rejeitar Hipotese nula, ou seja, há evidência estatística que o vetor gerado não segue uma distribuição \(\text{Exp }(\lambda = 1)\).

# Testando se o vetor gerado, como a téc. de redução HL, segue uma distribuição Exp(1).
Teste = ks.test(HL,pexp,1)
Teste$p.value
## [1] 1

Para a segunda amostra gerada, o p-valor foi 1 sendo assim, não podemos rejeitar Hipotese nula, ou seja, há evidência estatística que o vetor gerado não segue uma distribuição \(\text{Exp }(\lambda = 1)\).


Questão 2, Amostragem Descritiva

“Você está interessado em gerar uma amostra de valores exponencialmente distribuídos. Gere a amostra com e sem aplicação de técnica de redução de variância. A técnica de redução de variância deverá ser a da Amostragem Descritiva. Depois de gerados os valores, peça o histograma e faça o teste Kolmogorov-Smirnov para ambas as amostras e comente. OBS: Você não pode utilizar a função rexp.”

Método da Transformada Inversa

Aplicando a função inversa de \(F_X(x)\) em ambos os lados temos que \(x = F^{-1}(U)\).

Vejamos como ficaria no caso da distribuição \(X \sim Exp(\lambda)\)

Função de Densidade: \[f(x) = \lambda e^{-\lambda x}\] Função de Distribuição: \[F(x) = 1 - e^{-\lambda x}\] Número aleatório: \[x = F^{-1}(u) \Leftrightarrow x = - \frac{\ln(u)}{\lambda}\]

Agora usando o software RStudio, podemos aplicar esse resultado na função, runif(n) gerando n valores que seguem um variável uniforme contínua (0,1).

Código em R

Usando R:

set.seed(536837)  # Fixando um valor semente para facilitar a reprodutividade do experimento.
n = 100  # Definindo o número de valores gerados
lambda = 1  # Definindo o parâmetro
u = runif(n)  # Gerando n valores uniformes
X = -log(u)/lambda  # Aplicando o Método da Transformada Inversa

# Gerando o histograma de X

Visualização R

Redução de Variância, AD

Outra técnica que permitem reduzir a variância do processo, é Amostragem Descritiva, onde divide o intervalo de X em n partes iguais e em cada subintervalo faz uma amostragem do valor médio do intervalo.

Vejamos como na prática isso funciona. Usando a técnica de redução de variância, amostragem descritiva.

Código em R

Usando R:

AD <- c()  # Primeiramente defino um vetor vazio para receber os valores de X.
for(i in  1:n){ # Defino uma recursividade para gera os n valores.
  y <- -log((i - 0.5)/n)/lambda  # Usando o Método da Transformada Inversa no intervalo i
  AD[i] <- y  # Armazenando os valores de X, com a técnica de redução de variância, AD.
}

# Gerando o Histograma de X, com a técnica de redução de variância Amostragem Descritiva.

Visualização R

Teste Kolmogorov-Smirnov

Outro aspecto importante é verificar se os métodos estão gerando de maneira satisfatória. O Teste Kolmogorov-Smirnov verifica se há evidência estatística que o vetor aleatório segue uma certa distribuição. Onde a hiposte nula é que a função de distribuição acumulada \(F_X\) é igual a alguma função de distribuição.

Usando o RStudio podemos fazer o teste por meio da função ks.test(x,F(x)).

# Testando se o vetor gerado segue uma distribuição Exp(1).
Teste = ks.test(X,pexp,1)
Teste$p.value
## [1] 0.2177719

Para a primeira amostra gerada, o p-valor foi 0,2178 sendo assim, não podemos rejeitar Hipotese nula, ou seja, há evidência estatística que o vetor gerado não segue uma distribuição \(\text{Exp }(\lambda = 1)\).

# Testando se o vetor gerado, como a téc. de redução HL, segue uma distribuição Exp(1).
Teste = ks.test(HL,pexp,1)
Teste$p.value
## [1] 1

Para a segunda amostra gerada, o p-valor foi 1 sendo assim, não podemos rejeitar Hipotese nula, ou seja, há evidência estatística que o vetor gerado não segue uma distribuição \(\text{Exp }(\lambda = 1)\).


Questão 3, Monte Carlo

“Para a modelagem de um contrato de seguro discreto e temporário por 5 anos através de simulação monte carlo:

  a) Simule as mortes considerando distribuição poisson;

Em R, iremos usar a função rpois, com \(\lambda = \text{n . q}_X\), onde \(\text{q}_X\) é probabilidade de morte e \(\text{n}\) o número de apólice.

  b) Idade de contrato 50 anos;

Em R, idade é definida na variável idade.atual.

  c) Tábua AT83;

Importando a tábua AT83 em R, read.table.

  d) Benefício de R$10.000;

Em R, beneficio.

  e) Juros 6%;

Em R, taxa.juros

  f) Quantidade de apólices: 1000;

Em R, n.apolice.

  g) Simule o VP financeiro da obrigação, considerando 1000 simulações;

Em R, VP.

  h) Calcule o prêmio individual."

Em R, premio.ind.

Monte Carlo

A Simulação de Monte Carlo, também conhecida como Método de Monte Carlo ou uma simulação de probabilidade múltipla, é uma técnica matemática, que é usada para estimar os possíveis resultados de um evento incerto. Sabendo disso, vamos usar para definir o valores para o seguro discreto e temporário, onde durante 5 anos o segurado paga a seguradora e em caso de morte a seguradora paga aos beneficiario o valor que lhe é devido e se após os 5 anos não houver morte a seguradora não paga nada. Para simular as mortes vamos usar como base a Distribuição de Probabilidade Poisson, onde a taxa de morte por ano é 1.

Seguro discreto e temporário por 5 anos - Poisson

Código em R

Importando a tábua AT83, em R.
qx <- read.table('tabua.txt',header=F,dec=',')[,1] 
Definindo parâmetros
n.sim <- 1000 # Número de simulações
n.apolice <- 1000 # Número de apólices
idade.atual <- 50 # Idade no início do contrato
temp <- 5 # Temporaridade do contrato
beneficio <- 10000 # Valor a ser pago aos beneficiários
taxa.juros <- 0.06 # Taxa de juros anual
Aplicação Monte Carlo
S <- rep(0,n.sim) # Valor da obrigação
loss <- c() # Valor a ser pago em caso de morte
v <- 1/(1+taxa.juros) # Desconto
    for (i in 1:n.sim){
        n.apol <- n.apolice
        for (j in idade.atual:(idade.atual+temp-1)){
            mortes <- rpois(1,n.apol*qx[j])  #  Número de mortes
            loss[i] <- mortes*beneficio*(v^(j-50+1)) # VPF da perda na data do contrato
            n.apol <- n.apol-mortes # Número de apólice restante após 1 ano de contrato 
            S[i] <- S[i] + loss[i] # Valor Presente financeiro da obrigação.
        }
    }
VP <- mean(S) # Valor presente médio da Obrigação. 
premio.ind <- VP/n.apolice # Premio a ser pago por cada individuo.
## Valor Presente financeiro médio da obrigação após  1000 simulações. VP
## [1] 423529.9
## Valor presente médio do prêmio individual após  1000 simulações. VP/1000
## [1] 423.5299

Visualização R


Questão 4, Bootstrap

“Através da Simulação Bootstrap não paramétrica, construa um intervalo de confiança (α = 5%) para a média do VP financeiro da obrigação calculado na questão anterior. Você deve considerar uma amostra aleatória com 20 elementos. A quantidade de reamostragens deve ser igual a 100. Você pode escolher simular no R ou Excel.”

Bootstrap

O método Bootstrap não paramétrico, é uma técnica estatística de reamostragem tem como finalidade obter informações de características da distribuição de alguma variável aleatória. A partir de uma amostra aleatória é feito um reamostragem, com repetição, do elementos presente na amostra original, sendo assim não existe limite para o número de amostras que podem ser geradas. Para isso vou usar a função runif, para determinar qual elemento da amostra original foi sorteada para compor a nova amostra.

Simulação Bootstrap – Método não paramétrico, Amostra escolhida.
amostra <- S[1:20]  # Selecionei os 20 resultados obtidos aleatóriamente na questão anterior
amostra # Apresentando os dados
##  [1] 367749.3 413058.6 477715.8 510376.8 340849.5 491678.6 391448.4 469356.7
##  [9] 409696.7 375618.3 350564.0 342169.4 325815.4 455019.2 446114.1 373085.8
## [17] 479759.5 395569.2 452421.5 488924.6

Codigo em R

Definindo uma função que recebe a amostra original e o número de reamostragens.
bootstrap <- function (Amostra, B, alpha = 0.05){
    amostra <- Amostra  # Definindo amostra Original
    size <- length(amostra)  # Verificando o tamanho da amostra
    mediasB <- rep(0,B) # Criando o vetor da médias de cada amostra
    amostraB <- rep(0,size) # Reamostragem
      for (i in 1:B){
        for (j in 1:size){
            # Sorteando qual posição do elemento escolhido na amostra original 
            y <- trunc(runif(1)*size) + 1
            # Atribuindo o elemento sorteado na reamostragem. 
            amostraB[j] <- amostra[y] # Reamostragem gerada
        }
      mediasB[i] <- mean(amostraB) # Média da reamostragem i.
      }
    ordem <- sort(mediasB) # Ordenando o vetor das médias
    q1 <- quantile(ordem, alpha/2) # Definindo quantil 0.025, com alfa = 0.05
    q2 <- quantile(ordem, 1-alpha/2) # Definindo o quantil 0.975.
    
    # Mostrando os resultados
    cat("\n Média da Amostra inicial.\n")
    print(mean(amostra))
    cat("\n Médias das ",B,"reamostras Bootstrap.\n")
    print(sort(mediasB))
    icb <- round(c(q1,q2),2)
    cat("\n Intervalo Bootstrap com ",(1-alpha)*100,"% de confiança, usando",B,"reamostras.\n")
    print(icb)
}

Visualização em R

x = bootstrap(amostra, 100)
## 
##  Média da Amostra inicial.
## [1] 417849.6
## 
##  Médias das  100 reamostras Bootstrap.
##   [1] 391538.7 391956.8 392530.8 395122.7 397757.4 398401.5 398980.1 399115.7
##   [9] 400258.6 400343.6 400435.3 400961.0 401415.7 401521.9 402009.8 402690.9
##  [17] 402755.9 402801.6 402913.9 404380.4 406421.1 407609.9 407863.6 407919.6
##  [25] 408620.9 409175.7 409431.3 409747.1 410478.6 410732.9 410734.3 410742.8
##  [33] 411068.8 412536.1 413604.4 413847.6 413960.4 414545.4 415489.5 415855.9
##  [41] 416748.0 417231.1 417295.9 417318.9 418246.1 418348.2 418407.8 418409.0
##  [49] 418588.4 418697.7 418924.6 419105.4 419407.6 420155.9 421168.4 421667.2
##  [57] 422652.8 422751.6 422923.0 423234.3 423359.4 423633.6 423697.4 423736.4
##  [65] 424598.1 424815.9 425453.0 426045.9 426371.3 426445.3 426562.9 426633.7
##  [73] 426774.3 427198.1 427460.8 427593.9 427598.7 427672.6 427981.2 428070.1
##  [81] 428071.5 428630.0 430190.9 431029.0 431107.7 431231.5 431315.8 431813.2
##  [89] 431925.5 432230.8 434079.5 434283.2 434656.3 435326.6 436525.0 436648.3
##  [97] 437267.5 440046.2 440133.3 452558.9
## 
##  Intervalo Bootstrap com  95 % de confiança, usando 100 reamostras.
##     2.5%    97.5% 
## 393762.0 438726.3