os estados 1, 2, ..., m s˜ao transientes e o 0 ´e absorvente. Suponha que o processo assuma apenas valores n˜ao-negativos e que suas taxas de transi¸c˜ao sejam homogˆeneas no tempo. Seja X a v.a. que descreve o tempo gasto pelo processo at´e atingir o estado absorvente. A distribui¸c˜ao de X ´e dita do Tipo Fase, pois cada um dos estados transientes do processo de Markov representa uma das fases.
Denotamos a distribui¸c˜ao do tipo fase por P H(~α, T ), onde ~α ´e o vetor da distribui¸c˜ao inicial das probabilidades, ~α = (αi)i=1,...,m e T ´e a matriz de intensidade de transi¸c˜ao dos
A matriz de intensidade de transi¸c˜ao do processo ser´a na forma T ~t
0 0
,
onde ~t = −T~e com ~e = (1, 1, ..., 1)t. Cada um dos estados transientes representam uma
fase do processo. Vejamos alguns exemplos. Exemplo 1
Suponha que o nosso processo possua apenas uma fase. Podemos represent´a-lo pelo diagrama a seguir:
Para α = 1 e chamando t11 = −β temos, pelas Equa¸c˜oes de Kolmogorov:
P′
10(t) = t10P01(t) + t11P10(t)
P10′ (t) = −βP10(t)
P10(t) = e−βt
Assim, temos que o tempo gasto at´e o processo atingir o estado absorvente 0 tem distribui¸c˜ao exponencial com parˆametro β = −t11. Na nota¸c˜ao, P H(1, −β).
Exemplo 2: Distribui¸c˜ao com m fases
De um modo semelhante ao exemplo anterior, podemos obter uma distribui¸c˜ao com m fases com o diagrama abaixo:
Considere ~α = (1, 0, 0, 0, ..., 0) e ti i+1 = β, i = 1,...,m − 1 e tm = β, e os demais ti j
nulos. Note que este processo possui m fases idˆenticas. Ent˜ao o tempo de absor¸c˜ao, ou seja, o tempo necess´ario pra que o processo chegue no estado absorvente 0 corresponde `a soma dos tempos de transi¸c˜ao entre os estados i, i + 1 partindo do 1 at´e o 0. Assim, a distribui¸c˜ao correspondente ao tempo de absor¸c˜ao ser´a a Gama com parˆametro m (ou Erlang, pois m ´e inteiro) e densidade
βm xm−1
(m − 1)!e
−βx,
que ´e a convolu¸c˜ao de m exponenciais.
Outro exemplo de distribui¸c˜ao do tipo fase ´e a Hiperexponential que pode ser carac- terizada pela existˆencia de duas ou mais fases distintas.
Propriedades de uma distribui¸c˜ao do tipo fase
Seja X uma v.a. cuja distribui¸c˜ao ´e do tipo fase com representa¸c˜ao P H(~α, T ). Ent˜ao 1. A fun¸c˜ao de distribui¸c˜ao de X ´e F (x) = 1 − ~αeT x~1;
2. A fun¸c˜ao geradora de momentos de X ´e MX(s) = E[esX] = ~α(−sI − T )−1~t, onde I
´e a matriz identidade;
3. A densidade de X ´e f (x) = ~αeT x~t, onde ~t = −T~1 e eT = I +P∞
n=1 Tn
n!.
As propriedades acima nos possibilitam trabalhar com as distribui¸c˜oes do tipo fase, uma vez que todas as fun¸c˜oes que caracterizam a distribui¸c˜ao possuem forma fechada. Suas demonstra¸c˜oes ser˜ao omitidas mas podem ser encontradas em Asmussen [1].
Suponha ent˜ao que, no modelo de risco descrito em (1.9), as indeniza¸c˜oes Y tˆem uma distribui¸c˜ao do tipo fase com parˆametros (α, T ) onde
T = −1 00 −2 , e α = (1/2, 1/2) . Seja I = 1 00 1 , ~1 = (1, 1) e t = −T · ~1 = 12 . Temos que a fun¸c˜ao geradora de momentos da vari´avel Y ´e
MY(s) = E[ssY] = α(−sI − T )−1t
Assim,
E[Y ] = d
dsMy(s)|s=0 = 0.75 .
Gostar´ıamos de visualizar a influˆencia da matriz de transi¸c˜ao do processo da taxa de juros nas cotas finais. Para tanto consideramos o espa¸co de estados I = {6%, 8%, 10%} e duas matrizes de transi¸c˜ao diferentes,
P1 = 0 0.9 0.1 0.8 0.2 0 0.9 0.1 0 e P2 = 0.3 0.7 0 0 0.2 0.8 0 0.1 0.9
Calculamos, a seguir, os limitantes superiores para a probabilidade de ru´ına sob as condi¸c˜oes acima descritas. Para o c´alculo das cotas assumimos que s˜ao praticados os valores c = 0.975, θ = 0.1 e x = 5. Note que, nesse caso, C(b) = 0.975 − (1 + 0.1)(1 − b)0.75 = 0.15 + 0, 825b e C(b)E[Z] = 0.15 + 0.825b > 0.75b = bE[Y ], donde os valores satisfazem a condi¸c˜ao de seguran¸ca bE[Y ] < C(b)E[Z].
3.3.2
Desigualdade de Lundberg
Para encontrarmos o coeficiente de Lundberg precisamos garantir a existˆencia de R0 > 0
que satisfaz a equa¸c˜ao
MS(−R0) = 1,
onde MS(s) ´e a fun¸c˜ao geradora de momentos de S := C(b0)Z − b0Y .
Como j´a citado acima, a fun¸c˜ao geradora de momentos de Y ´e dada por MY(s) = E[esY] = α(−sI − T )−1t = −1 2(s − 1) − 1 s − 2 ⇒ MY′ (s) = 2 (2s − 2)2 + 1 (s − 2)2 ⇒ E[Y ] = M ′ Y(0) = 3 4.
Considerando Z com distribui¸c˜ao exponencial com parˆametro 1 e os valores c = 0.975, θ = 0.1 e x = 5 e sabendo que
C(b) = c − (1 + θ)(1 − b)E[Y ]
E[Z] = 0.975 − (1.1)(1 − b)0.75 = 0.15 + 0.825b, podemos calcular E[S], ou seja, E[C(b0)Z − b0Y ]. Temos:
E[C(b0)Z − b0Y ] = C(b0)E[Z] − b0E[Y ] = 0.15 + 0.825b0− 0.75b0 = 0.15 + 0.075b0 > 0.
Note ainda que C(0) = 0.15 > 0, donde bmin = 0 e podemos calcular as cotas superiores
para a probabilidade de ru´ına para todo b ∈ (0, 1] . Ent˜ao, para cada n´ıvel de reten¸c˜ao b ∈ B devemos encontrar uma fun¸c˜ao geradora de momentos de S = C(b)Z − bY .
Com o aux´ılio do software Maple 13 pudemos escrever a fun¸c˜ao geradora de momentos de S em fun¸c˜ao do n´ıvel de reten¸c˜ao b, o que denotamos como segue:
MS(t, b) := −20(3tb + 4)
(6t + 33tb − 40)(tb + 1)(tb + 2) < ∞ para todo t ∈ (−1b,
40
6+33b). No pr´oximo cap´ıtulo detalhamos como a express˜ao acima foi
obtida e porque esta limita-se ao intervalo a valores de s no intervalo dado.
Pelo Lema 2.3.2 podemos, ent˜ao, garantir a existˆencia de um ´unico R0 ≡ R0(b) > 0
para cada n´ıvel de reten¸c˜ao b que satisfaz a equa¸c˜ao MS(−R0, b) = 1.
Com a ajuda do software Maple resolvemos a equa¸c˜ao acima de modo a encontrar R0 = R0(b), chegamos `as seguinte solu¸c˜oes:
0,18 + 59b +p(36 + 396b + 2689b
2)
6(2 + 11b)b ,
18 + 59b −p(36 + 396b + 2689b2)
6(2 + 11b)b .
As duas solu¸c˜oes n˜ao nulas s˜ao positivas para todo b ∈ (0, 1]. Vamos, por ora, adotar a seguinte nomenclatura:
R+(b) :=
18 + 59b +p(36 + 396b + 2689b2)
e R−(b) := 18 + 59b −
p(36 + 396b + 2689b2)
6(2 + 11b)b .
Em Diasparra e Romera [5] os valores apresentados s˜ao provenientes da equa¸c˜ao que resulta no maior valor para R0 a partir de um mesmo n´ıvel de reten¸c˜ao b. ´E f´acil ver
anal´ıtica e graficamente (abaixo) que trata-se da fun¸c˜ao R+(b).
No entanto, quando encontramos uma express˜ao para a fun¸c˜ao geradora de momentos de S (em fun¸c˜ao do n´ıvel de reten¸c˜ao b) ficou claro que tal express˜ao s´o ´e v´alida para t em uma vizinhan¸ca do zero, a saber, para valores de t no intervalo (−1
b, 40
6+33b). Analisando
os valores apresentados em Diasparra e Romera [5] observamos que tais valores n˜ao per- tencem a este intervalo onde a fun¸c˜ao MS(t, b) est´a bem definida. O gr´afico abaixo nos
mostra que os valores encontrados para R0(b) em [5] n˜ao satisfazem a desigualdade
Para este trabalho, ent˜ao,
R0(b) := 18 + 59b −
p(36 + 396b + 2689b2)
6(2 + 11b)b .
O gr´afico a seguir (3.1) mostra a rela¸c˜ao entre a escolha de b e o valor da cota superior de Lundberg para a probabilidade de ru´ına
Figura 3.1: Rela¸c˜ao entre a Cota Superior de Lundberg e o n´ıvel de reten¸c˜ao b Na tabela apresentada no final da se¸c˜ao veremos os valores de cada uma das cotas apre- sentadas calculados para alguns valores particulares de b. Ao longo da se¸c˜ao apresentamos os resultados encontrados para b = 0.5. Temos os seguinte resultados:
R0(0.5) = 0.7732
e
ψ(5) ≤ e−5R0(0.5) = 0.0209 .
Supondo que n˜ao conhecemos a distribui¸c˜ao de Y podemos usar os resultados obtidos no exemplo anterior para obter as seguinte estimativas:
ˆ
e
ψ(5) ≤ e−5 ˆR0(0.5) = 0.0212 .
O gr´afico abaixo, obtido atrav´es do software R, mostra a rela¸c˜ao entre o n´ıvel de reten¸c˜ao b e os valores reais e estimados de R0.
Figura 3.2: R0(b) e ˆR0(b)
Al´em de comparar os valores reais e estimados de R0 podemos tamb´em comparar os
valores das cotas, ou seja, de e−R0x para R
0 real e estimado.
Para a amostra usada na obten¸c˜ao desses gr´aficos podemos afirmar que o resultado da a cota estimada ´e bastante satisfat´orio visto que apresenta-se bastante pr´oxima dos valores reais.
3.3.3
Cota Superior Indutiva
Como vimos no in´ıcio desta se¸c˜ao, a fun¸c˜ao de distribui¸c˜ao de Y ´e dada por F (x) = 1 − ~αeT x~1.
Com o aux´ılio do Maple 13 encontramos a seguinte express˜ao: F (x) = −12e−x− 1 2e −2x + 1 ⇒ ¯F (x) = 1 2e −x+1 2e −2x.
Mostraremos, a seguir, que a distribui¸c˜ao de Y ´e NPUC. ¯ F (y) = 1 2e −y +1 2e −2y e Z ∞ x ¯ F (z)dz = 1 2 Z ∞ x e−z+ e−2zdz, ⇒ ¯F (y) Z ∞ x ¯ F (z)dz = 1 4(e −x−y+ e−x−2y+1 2e −2x−y+ 1 2e −2x−2y) Para x > 0 e y > 0 temos e−x−2y ≤ e−x−y ⇒ 1 4(e −x−y + e−x−2y ) ≤ 12e−x−y e e−2x−y ≤ e−2x−2y ⇒ 1 8e −2x−y+ 1 8e −2x−2y ≤ 1 4e −2x−2y. Ent˜ao ¯ F (y) Z ∞ x ¯ F (z)dz ≤ 12e−x−y+ 1 4e −2x−2y = Z ∞ x+y ¯ F (z)dz, e F ´e NPUC.
Podemos, ent˜ao, usar o resultado do Corol´ario 2.3.6 e calcular a cota superior indutiva segundo a equa¸c˜ao
ψπ(x, i) ≤ (Eπ[eR0bY1])−1Eπ[e−R0x(1+I1)|I
0 = i].
Para x = 5, i = 8%, b = 0.5 e a matriz de transi¸c˜ao P1 obtivemos o seguinte resultado:
ψπ(x, i) ≤ 0.0114.
O valor estimado para a cota n˜ao pode levar em considera¸c˜ao a equa¸c˜ao dada pelo corol´ario, uma vez que a distribui¸c˜ao ´e desconhecida. Portando o c´alculo dessa cota foi feito de modo semelhante ao do exemplo anterior. O valor obtido para b = 0.5 e a matriz P1 foi
ψπ(x, i) ≤ 0.0127.
Repetindo os procedimentos para a matriz P2 e considerando conhecidas a distribui¸c˜oes
obtivemos a seguinte desigualdade
ψπ(x, i) ≤ 0.0101.
No caso em que utilizamos a densidade estimada, obtivemos ψπ(x, i) ≤ 0.0112.
3.3.4
Cota Superior Obtida via Martingales
Assim como no exemplo anterior devemos encontrar, para todo i ∈ I, ρi que satisfa¸ca
Eπ[e−ρiS(1+I1)−1|I
0 = i] = 1
e tomar como R1 o m´ınimo entre eles. Como dito anteriormente, para podermos avaliar
a influˆencia da Matriz de transi¸c˜ao do processo dos juros no resultado da cota superior obtida, faremos esse c´alculo para as duas matrizes j´a apresentadas, P1 e P2.
Com o aux´ılio do software Maple encontramos para a matriz P1
R1(0.5) = 0.8210
e, portanto,
ψ(5) ≤ e−5R1 = 0.0165.
Assumindo desconhecida a distribui¸c˜ao de Y obtivemos, com o aux´ılio do software R, as estimativas
ˆ
R1(0.5) = 0.8180
e
ψ(5) ≤ e−5 ˆR1 = 0.0167.
Os mesmos c´alculos foram feitos para a matriz P2 obtendo os resultados
R1(0.5) = 0.8302 e ψ(5) ≤ e−5R1 = 0.0157 e as estimativas ˆ R1(0.5) = 0.8272 e ψ(5) ≤ e−5 ˆR1 = 0.0160.
Na tabela abaixo, em que L representa a cota superior de Lundberg, I a indutiva e M a cota obtida via Martingales, podemos encontrar o valor real e o estimado das cotas superiores para a probabilidade de ru´ına para alguns valores de b.
b - R0 L I (P1) I - (P2) R1(P1) M(P1) R1(P2) M(P2) 0.5 Real 0.7732 0.0209 0.0114 0.0101 0.8210 0.0165 0.8302 0.0157 0.5 Estimado 0.7703 0.0212 0.0127 0.0112 0.8180 0.0167 0.8272 0.0160 0.75 Real 0.4182 0.1236 0.0818 0.0765 0.4441 0.1085 0.4491 0.1059 0.75 Estimado 0.4133 0.1266 0.0678 0.0635 0.4389 0.1114 0.4438 0.1087 0.85 Real 0.3460 0.1772 0.1226 0.1160 0.3675 0.1592 0.3716 0.1560 0.85 Estimado 0.3412 0.1816 0.1192 0.1129 0.3624 0.1633 0.3664 0.1601 0.95 Real 0.2926 0.2315 0.1656 0.1580 0.3108 0.2114 0.3143 0.2078 0.95 Estimado 0.2880 0.2369 0.1672 0.1597 0.3059 0.2166 0.3093 0.2130 1 Real 0.2709 0.2580 0.1872 0.1792 0.2877 0.2372 0.2910 0.2335 1 Estimado 0.2665 0.2638 0.1916 0.1836 0.2830 0.2429 0.2862 0.2391 Comparando os resultados obtidos neste trabalho com os apresentados no artigo Ine-
qualities for the ruin probability in a controlled dicrete-time risk process para as mesmas
distribui¸c˜oes e sob as mesmas condi¸c˜oes podemos perceber que existe divergˆencia n˜ao s´o nos valores encontrados para R0 mas tamb´em nos valores de R1 e da cota superior indu-
tiva (o que ´e esperado uma vez que a cota indutiva depende de R0). Isso se deve ao fato
de que o equ´ıvoco cometido no artigo de Ros´ario Romera e Maikol Diasparra no sentido de n˜ao observar os intervalos para os quais a fun¸c˜ao geradora de momentos da v.a. S estava bem definida repete-se no c´alculo da cota obtida via Martingales.
Para cada valor de i ∈ I existe, em R, mais de um zero da fun¸c˜ao que chamamos li(r)
que ´e dada pela equa¸c˜ao (2.27), a saber,
li(r) := Eπ[e−r(C(b)Z1−bY1)(1+I1)
−1
|I0 = i] − 1.
´
E importante notar que a fun¸c˜ao li(r) depende do c´alculo da geradora de momentos de
S que chamamos de MS(t), que por sua vez s´o est´a bem definida em uma vizinhan¸ca do
zero para a qual assume valores finitos.
Uma forma simples de verificar a n˜ao validade dos valores apresentados em Diasparra e Romera [5] ´e aplicar os valores encontrados para R1 na equa¸c˜ao (2.26). A desigualdade
li(R1) ≤ 1
n˜ao ´e satisfeita para todo i como demonstramos que deve ser. Logo, os valores apresen- tados em Diasparra e Romera [5] est˜ao incorretos.
Cap´ıtulo 4
Algoritmos de Simula¸c˜ao
4.1
Comandos Maple 13
4.1.1
Exemplo 1
Veja a seguir os comandos usados no software Maple 13 a fim de obter as cotas superiores para a probabilidade de ru´ına do primeiro exemplo apresentado no Cap´ıtulo 3 deste trabalho.
Este software foi utilizado para os c´alculos que consideram conhecidas as distribui¸c˜oes das indeniza¸c˜oes e do tempo entre elas (Y e Z, respectivamente).
O primeiro passo a ser executado para o c´alculo das cotas para a probabilidade de ru´ına ´e a defini¸c˜ao das vari´aveis aleat´orias envolvidas (Y e Z). Em seguida devemos encontrar suas fun¸c˜oes geradoras de momentos e avaliar o intervalo em que est˜ao definidas e s˜ao finitas.
Desigualdade de Lundberg
A seguir definiremos a vari´avel aleat´oria S como sendo S = C(b)Z − bY. Como neste exemplo temos b = 1 e c = 4, S = C(1)Z − Y = cZ − Y = 4Z − Y. Definir esta vari´avel aleat´oria nos permitir´a calcular o valor do coeficiente de Lundberg R0 como sendo a
Note que, pela express˜ao dada acima, MS(s) ´e finita para s ∈ (−12 ,14) e, quando −1
2 < s < 1
4, MS(s) ´e dada pela express˜ao abaixo:
Basta, ent˜ao, resolver a equa¸c˜ao MS(−R0) = 1 para encontrar o valor de R0 e do
Cota superior Indutiva
Para o c´alculo da cota superior Indutiva devemos calcular o valor de β−1 = inft≥0
R∞
t e
R0b0ydF (y)
eR0b0tF (t)¯ .
Este c´alculo foi feito em algumas etapas como segue abaixo:
Com o valor de β j´a conhecido resta calcularmos o valor da cota. O procedure abaixo cumpre esse papel, bastando apenas definir o espa¸co de estados I da Cadeia de Markov dos juros e a segunda linha de sua matriz de transi¸c˜ao pois, uma vez que a taxa de juros inicial desse exemplo ´e de 8%, as probabilidades de transi¸c˜ao de nosso interesse encontram-se nessa linha. Veja os procedimentos e resultados abaixo:
Cota superior obtida via Martingales
Para o c´alculo do limitante superior para a probabilidade de ru´ına obtido via martingales iremos definir a matriz de transi¸c˜ao P do processo da taxa de juros e as fun¸c˜oes dadas em (2.27), a saber,
li(r) := Eπ[e−r(C(b)Z1−bY1)(1+I1)
−1
|I0 = i] − 1.
Vamos restringir o gr´afico ao intervalo onde se encontram as ra´ızes positivas das li’s
e, a seguir, encontrar seus valores.
Sabemos que o valor de R1 ´e o m´ınimo entre os ρ1. Os comandos abaixo nos dar˜ao
o valor de R1 e do limitante superior para a probabilidade de ru´ına completando, assim,
4.1.2
Exemplo 2
Apresentamos a seguir os comandos usados a fim de obter as cotas superiores para a probabilidade de ru´ına do segundo exemplo apresentado no Cap´ıtulo 3 deste trabalho considerando conhecidas as distribui¸c˜oes das indeniza¸c˜oes e do tempo entre elas (Y e Z, respectivamente).
O primeiro passo a ser executado para o c´alculo das cotas para a probabilidade de ru´ına ´e a defini¸c˜ao das vari´aveis aleat´orias envolvidas (Y e Z). Diferentemente do exem- plo anterior, o software Maple 13 n˜ao possui implementada a distribui¸c˜ao do tipo fase que est´a envolvida no processo neste exemplo. Para trabalharmos com tal distribui¸c˜ao devemos definir seus parˆametros e fun¸c˜oes de interesse (que s˜ao conhecidas como vimos anteriormente), e a partir dessas fun¸c˜oes definiremos a v.a. Y .
Ap´os definir os elementos necess´arios para o c´alculo da densidade de Y , podemos achar sua forma sabendo que sua densidade ´e f (x) = ~αeT x~t, onde ~t = −T~1. Uma vez calculada
a densidade de Y podemos definir a v.a. usando o comando Distribution. Em seguida, encontramos a fun¸c˜ao geradora de momentos de Y de acordo com a terceira propriedade das distribui¸c˜oes do tipo fase apresentada no cap´ıtulo anterior.
Como no exemplo anterior, consideramos a distribui¸c˜ao do tempo entre o pagamento de duas indeniza¸c˜oes como uma v.a. Z com distribui¸c˜ao exponencial de parˆametro 1. Abaixo definimos a v.a. Z e encontramos sua fun¸c˜ao geradora de momentos.
Interessados em encontrar as cotas superiores para a probabilidade de ru´ına para todo b ∈ B = (0, 1], devemos definir a fun¸c˜ao C(b) e, em seguida, a v.a. S = C(b)Z − bY . A partir dessa v.a. S = S(b) podemos encontrar uma fun¸c˜ao geradora de momentos MS(s, b)
Se analisarmos a fun¸c˜ao MS(s, b) acima podemos ver que ela converge apenas para
s > −1b e, para esses valores de s, a fun¸c˜ao pode ser reescrita como abaixo:
Novamente analisando a convergˆencia da fun¸c˜ao MS(s, b), podemos mostrar que a
fun¸c˜ao converge para s < 40
Desigualdade de Lundberg
Para encontrar o valor de R0(b) basta resolver a equa¸c˜ao MS(−R0, b) = 1. Fazendo
o c´alculo para b = 0.75 com os comandos abaixo obtemos 3 valores diferentes. Apenas um dos valores positivos encontrados para R0 satisfaz a desigualdade −1b < −R0 < 6+33b40 .
Como b > 0 devemos nos preocupar apenas com a desigualdade −1b < −R0.
Para b = 0.75 temos −1
b = −1.333 < −0.418211. Portanto esse ´e o valor de R0 que
nos interessa. Ou seja, R0(0.75) = 0.418211.
Esse ´e um caso particular. Com os comando abaixo trataremos do caso geral, ou seja, encontraremos uma express˜ao para R0(b) com base no que j´a vimos. A princ´ıpio
encontraremos duas express˜oes para R0. Em seguida mostraremos que uma delas n˜ao
fornece valores correspondentes aos intervalos em que a fun¸c˜ao geradora de momentos de S(b) existe e ´e finita.
O gr´afico abaixo nos permite constatar que a express˜ao R−(b) nos fornece valores que
Em Diasparra e Romera [5] os valores apresentados s˜ao obtidos pela equa¸c˜ao R+(b)
que resulta em valores maiores para R0 que , no entanto, n˜ao se encontram no intervalo
onde a fun¸c˜ao geradora de momentos de S est´a bem definida. Por isso definimos abaixo a fun¸c˜ao R0(b) que nos forneceu os resultados neste trabalho apresentados.
A seguir, calculamos o valor do limitante de Lundberg para a probabilidade de ru´ına para alguns valores de b com X0 = 5 e logo abaixo mostramos graficamente os valores do
Cota superior Indutiva
Como vimos no Cap´ıtulo 3, a distribui¸c˜ao das indeniza¸c˜oes Y neste exemplo ´e NPUC. Sendo assim, podemos usar o resultado da Proposi¸c˜ao 2.3.6 para encontrar a conta superior Indutiva para a probabilidade de ru´ına.
O primeiro passo para o c´alculo da cota ´e definir os elementos envolvidos no mesmo, a saber, o espa¸co de estados I da cadeia de Markov dos juros e as duas matrizes de transi¸c˜ao com rela¸c˜ao as quais faremos o c´alculo com a finalidade de observar sua influencia no valor da cota.
Abaixo denominamos IBP1 a cota indutiva obtida a partir da matriz de transi¸c˜ao P1
e IBP2 a cota indutiva obtida a partir da matriz de transi¸c˜ao P2. Realizamos o c´alculo
O gr´afico abaixo nos permite comparar as cotas j´a calculadas:
Cota superior obtida via Martingales
Assim como no exemplo anterior, a fim de obter o limitante superior para a probabili- dade de ru´ına via martingales definimos as fun¸c˜oes dadas em (2.27), a saber,
li(r) := Eπ[e−r(C(b)Z1−bY1)(1+I1)
−1
|I0 = i] − 1.
Para uma maior simplicidade esse processo foi feito para um valor fixo de b (b = 0.5)1 nos
dois passos a seguir.
1No Apˆendice A apresentamos um algoritmo alternativo que permite-nos obter o valor de R
1(b) para
O valor de b pode ser facilmente substitu´ıdo no segundo argumento da fun¸c˜ao MS
acima, resultando na obten¸c˜ao da cota via martingales para qualquer valor de b ∈ (0, 1] desejado.
Se analisarmos graficamente as fun¸c˜oes li acima definidas podemos perceber que suas
ra´ızes positivas encontram-se bastante pr´oximas.
Restringindo os valores de b a uma vizinhan¸ca das ra´ızes da l′
is podemos visualiza-las
de forma ordenada. Em seguida calculamos seus valores e chamamos de R1 o m´ınimo
Por fim calculamos o valor da cota dado por e−R1x, para x = X
0 = 5.
Abaixo, o valor de R1 encontrado para alguns valores de b.
Analogamente, para encontrar os valores de R1 a partir da matriz P2 basta seguir os
A cota superior obtida via martingales a partir da matriz P2 para b = 0.5 e x = X0 = 5
4.2
Comandos R
Nesta se¸c˜ao inclu´ımos os comandos utilizados no R Software para a obten¸c˜ao dos resultados estimados apresentados neste trabalho. Ao nos citarmos as v.a.’s Y e Z nos referimos a vari´aveis aleat´orias cujas fun¸c˜oes de distribui¸c˜ao s˜ao F e G, respectivamente. Isto ´e, Y tem a mesma distribui¸c˜ao das Yn e Z tem a mesma distribui¸c˜ao das Zn.
Devemos deixar claro que os resultados obtidos pela repeti¸c˜ao dos m´etodos aqui apre- sentados n˜ao ser˜ao iguais aos resultados do Cap´ıtulo 3 deste trabalho pois os mesmos dependem das amostras aleat´orias geradas pelo software.
4.2.1
Exemplo 1
Trataremos do primeiro exemplo apresentado no Cap´ıtulo 3 deste trabalho onde consideramos as indeniza¸c˜oes Yn com distribui¸c˜ao exponencial de parˆametro µ = 12 e o
tempo entre indeniza¸c˜oes Zn com distribui¸c˜ao exponencial de parˆametro λ = 1.
O primeiro passo para a obten¸c˜ao de um estimador ´e gerar amostras das v.a.’s em quest˜ao.
#Gerando amostras de Y e Z de tamanho n n<-1000
dadosy<-rexp(n,1/2) mean(dadosy)
dadosz<-rexp(n,1) mean(dadosz)
Seja S = C(b)Z − bY . Para b = 1 e C(b) = c = 4 temos que S = 4Z − Y . Vamos gerar uma amostra de S de tamanho n.
dadoss<-4*dadosz-dadosy mean(dadoss)
Como vimos anteriormente, a suavidade da fun¸c˜ao de densidade estimada depende tanto do Kernel (ou fun¸c˜ao N´ucleo) utilizado, quanto do parˆametro de suaviza¸c˜ao h, sobre o qual podemos ver mais detalhes em Dias [4] e Silverman [17]. O comando density nos fornecer´a as informa¸c˜oes necess´arias para a obten¸c˜ao do estimador. A seguir, a obten¸c˜ao do parˆametro de suavidade h a partir da amostra de S que geramos anteriormente: hs<-density(dadoss)$bw
Utilizando os comandos abaixo podemos visualizar o histograma da amostra com a qual vamos trabalhar, a densidade real da distribui¸c˜ao da qual obtivemos a amostra e a densidade estimada pelo comando density .
plot.new()#criando um novo plano para os gr´aficos que ser~ao plotados dadosy<-sort(dadosy) #ordenando os dados
hist(dadosy,prob=T,ylim=c(0,0.5))#criando o histograma lines(dadosy,dexp(dadosy,1/2),lwd=2)#densidade real
lines(density(dadosy),xlim=c(0,12),lwd=2,col=2)#densidade estimada Estes comandos geraram a seguinte imagem:
Temos todas as informa¸c˜oes necess´arias para obter a fun¸c˜ao geradora de momentos estimada a partir de uma amostra qualquer.
m<-function(x,dados){ soma1<-0 h<-density(dados)$bw for (i in 1:n) { soma1<-soma1+exp(dados[i]*x)} return((1/n)*exp(1/2*h^2*x^2)*soma1) }
Para encontrar um valor estimado pra R0 basta obter a fun¸c˜ao geradora de momentos
a partir da amostra de S e resolver a equa¸c˜ao
MS(−R0) − 1 = 0,
o que fazemos com os comandos a seguir: g<-function(r){m(-r,dadoss)-1} library(rootSolve)
r0<-uniroot.all(g,c(0.01,1))
Para obter o valor da cota superior basta calcular e−r0x. No nosso exemplo temos
x = 1. O valos estimado da cota superior ´e encontrado, ent˜ao, pelo comando exp(-r0)
Os comandos abaixo nos permitem obter k estimativas para R0 que ser˜ao guardadas
no vetor vetorr0 . O comando final nos fornecer´a a m´edia das estimativas obtidas. k<-1000 vetorr0<-rep(0,k) for (i in 1:k){mdadosy<-rexp(n,1/2) mdadosz<-rexp(n,1) mdadoss<-4*mdadosz-mdadosy mg<-function(r){m(-r,mdadoss)-1} vetorr0[i]<-uniroot.all(mg,c(0.01,1))} mediar0<-mean(vetorr0) mediar0
Passaremos a calcular uma estimativa para a cota superior indutiva para a probabili- dade de ru´ına dada na equa¸c˜ao
ψπ
(x, i) ≤ βX
j∈I
pijEπ[e−R0x(1+j)] = βEπ[e−R0x(1+I1)|I0 = i],
onde β ≡ β(b0) ´e dado por
β−1 = inft≥0
R∞
t e
R0b0ydF (y)
eR0b0tF (t)¯ .
Trabalharemos com a equa¸c˜ao (3.1) segundo a qual R∞ t e R0b0ydF (y) eR0b0tF (t)¯ = MY(R0b0) − Rt −∞e R0b0ydF (y) eR0b0t(1 − F (t)) .
Para tanto precisaremos de express˜oes para a fun¸c˜ao geradora de momentos estimada e fun¸c˜ao de distribui¸c˜ao acumulada estimada de Y , obtidas com os comandos abaixo.
my<-function(x){m(x,dadosy)}
Fun¸c˜ao de distribui¸c˜ao acumulada estimada de Y : pest<-function(x){ hy<-density(dadosy)$bw