-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathAula 3 - Estágio Supervisionado II.Rmd
152 lines (108 loc) · 3.67 KB
/
Aula 3 - Estágio Supervisionado II.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
---
title: "Aula 3 - Estágio Supervisionado II"
author: "Filipe Duarte"
date: "8/14/2019"
output:
html_document: default
pdf_document: default
---
# Aula 3
Nesta aula veremos como estimar os parâmetros das distribuições pelo método da máxima verossimilhança. Em seguida, aprenderemos como selecionar o modelo com o melhor ajuste para a amostra.
## Método da Máxima Verossimilhança
O método da máxima verossimilhança é um método para estimar os parâmetros de um modelo probabilístico. Dessa maneira, dada uma amostra e dado um modelo probabilístico, a estimativa por máxima verossimilhança estima os valores para os diferentes parâmetros do modelo.
Considerando $X$ uma variável aleatória e sendo $f(x,\theta)$ a função de densidade, onde o parâmetro $\theta$ é desconhecido. Dessa forma, dada uma amostra aleatória simples de $X$,de tamanho $n$, $X_1, X_2, ..., X_n$, e sejam $x_1, x_2, ..., x_n$ os valores observados, a função de verossimilhança $L$ define-se por:
$$
L(\theta; x_1,..., x_n) = f(x_1;\theta)*...*f(x_n;\theta).
$$
Dessa forma, dizemos que $\hat{\theta}$ é um estimador de máxima verossimilhança (EMV) para $\theta = \theta_0$.
## fitditsrplus
O pacote `fitdistrplus` estima os parâmetros dos modelos pela máxima verossimilhança. Portanto, faremos uso dessa técnica para estimar os parâmetros dos modelos para os dados dos sinistros.
Carregando biblioteca:
```{r}
library(fitdistrplus)
```
Carregando dados:
```{r}
sinistro <- read.csv("sinistros.csv", header = TRUE, sep=",")
```
Visualizando os dados:
```{r}
head(sinistro)
```
Criando um vetor apenas com os valores dos sinistros:
```{r}
sinistro <- sinistro$x
```
Estatísticas descritivas:
```{r}
# summary
summary(sinistro)
# desvio padrão
sd(sinistro)
```
Histograma
```{r}
# histograma
hist(sinistro)
```
# Ajuste das distribuições
```{r}
descdist(sinistro)
```
Vamos testar para gamma, weibull e lognormal:
```{r}
fit <- fitdist(sinistro, "gamma")
fitdist(sinistro, "weibull")
fitdist(sinistro, "lnorm")
```
Vamos plotar o histograma com a curva gerada a partir da estimação:
```{r}
hist(sinistro, pch=20, breaks = 20, prob = TRUE, main = "")
curve(dgamma(x, fit$estimate[1], fit$estimate[2]), add = TRUE, col = "red")
```
Vamos fazer o histograma e a distribuição acumulada:
```{r}
plotdist(sinistro, histo = TRUE, demp = TRUE)
```
Agora vamos estimar novamente para criar as 3 funções
```{r}
fit_w <- fitdist(sinistro, "weibull")
fit_g <- fitdist(sinistro, "gamma")
fit_ln <- fitdist(sinistro, "lnorm")
```
```{r}
summary(fit_w)
summary(fit_g)
summary(fit_ln)
par(mfrow=c(2,2))
plot.legend <- c("Weibull", "lognormal", "gamma")
denscomp(list(fit_w, fit_g, fit_ln), legendtext = plot.legend)
cdfcomp (list(fit_w, fit_g, fit_ln), legendtext = plot.legend)
qqcomp (list(fit_w, fit_g, fit_ln), legendtext = plot.legend)
ppcomp (list(fit_w, fit_g, fit_ln), legendtext = plot.legend)
```
Testando outras distribuições:
Carregar outro pacote:
```{r}
library(actuar)
```
Testando a loglogistica e pareto:
```{r}
fit_ll <- fitdist(sinistro, "llogis", start = list(shape = 1, scale = 500))
fit_P <- fitdist(sinistro, "pareto", start = list(shape = 1, scale = 500))
```
```{r}
summary(fit_ll)
summary(fit_P)
```
```{r}
par(mfrow=c(2,2))
plot.legend <- c("lognormal", "loglogistic", "Pareto")
denscomp(list(fit_ln, fit_ll, fit_P), legendtext = plot.legend)
cdfcomp (list(fit_ln, fit_ll, fit_P), legendtext = plot.legend)
qqcomp (list(fit_ln, fit_ll, fit_P), legendtext = plot.legend)
ppcomp (list(fit_ln, fit_ll, fit_P), legendtext = plot.legend)
```
```{r}
gofstat(list(fit_w, fit_g, fit_ln, fit_ll, fit_P), fitnames = c("Weibull","Gama","LogNormal", "Loglogística", "Pareto"))
```