Para as questões 1 e 2 considerar um valor de semente igual a matricula de um dos componentes da equipe.
“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.”
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
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
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.
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.
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)\).
“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.”
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).
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
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.
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.
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)\).
“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.
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
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
“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.”
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
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)
}
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