Vou começar lendo os dois arquivos.

Estou fechando a análise. Li o enunciado e o seu script e, como o R não está instalado aqui, reproduzi a Tarefa 1 em Python para ter os números reais do seu seguro. Eles devem bater com os do seu R até a 4ª ou 5ª casa. Só o p-valor do bootstrap muda com a semente.

---

O cenário

Sua carteira tem duas fontes de incerteza: quantos sinistros acontecem no ano ($N$) e quanto custa cada um ($X$). O que a seguradora paga no ano é $S = X_1 + \dots + X_N$. O Exercício 1 pergunta só isto: que “personalidade” têm $N$ e $X$?

Um detalhe que o enunciado não diz: $87+79+\dots+103 = 1.038$. Ou seja, os 1.038 sinistros são todos os sinistros dos 12 anos, não uma amostra separada.

---

1a. Estatísticas-resumo, histogramas e dispersão

🟢 Versão fácil

	Severidade (R$)	Frequência (nº/ano)
Média	4.243,13	86,5
Desvio-padrão	8.543,89	13,30
CV	2,01	0,154
Mediana	1.416,53	81,5
q90 / q95 / q99	9.957 / 16.363 / 48.017	102,4 / 107,5 / 111,9
Mín. / Máx.	83,76 / 87.732	65 / 113
Assimetria	5,07	0,52
Curtose	35,75	2,60

A história da cafeteria. O sinistro “típico” custa cerca de R$ 1,4 mil (a mediana), mas a média é R$ 4,2 mil. É como o salário de uma firma pequena em que o dono ganha 100 mil: a média sobe e a mediana não. O CV (desvio ÷ média) diz o mesmo: na severidade o desvio é o dobro da média; na frequência é só 15%. Quase toda a incerteza do seu seguro vem do tamanho dos sinistros, não da quantidade.

Quem paga a conta? Os 10% maiores sinistros (104 deles) respondem por 55,1% do valor total. Os 10 maiores, por 14,7%. O q99 (48 mil) é cerca de 34 vezes a mediana.

Frequência: razão variância/média = 2,05. Se os sinistros “caíssem do céu” a uma taxa fixa (Poisson), a variância seria igual à média, ou seja, a razão seria ≈ 1. Aqui é o dobro, e o teste dá p = 0,021, abaixo de 5%. Isso sugere sobredispersão: a taxa de sinistros provavelmente muda de ano para ano. Um indício concreto: numa Poisson(86,5), um ano com 113 ou mais sinistros tem probabilidade de 0,36%, e um ano com 65 ou menos, de 0,96%. Os dois aparecem nos 12 anos. Mas são só 12 anos, então a evidência é moderada, não esmagadora.

Mistério: o excesso médio $e(u)$. A pergunta é: “dado que o sinistro já passou de $u$, quanto passa em média?”

$u$	mediana	q75	q90	q95	q98
$e(u)$	6.409	9.109	13.394	17.759	20.739

Numa Exponencial (sem memória), $e(u)$ seria constante. Aqui ele triplica: quanto maior o sinistro, maior o “resto” esperado. Essa é a assinatura de cauda pesada.

🔵 Versão rigorosa

Variância amostral, CV e momentos (fórmulas em LaTeX):

$$s^2=\frac{1}{n-1}\sum_{i=1}^n (x_i-\bar x)^2,\qquad CV=\frac{s}{\bar x}$$

$$m_k=\frac1n\sum_{i=1}^n (x_i-\bar x)^k,\qquad \gamma_1=\frac{m_3}{m_2^{3/2}},\qquad \gamma_2=\frac{m_4}{m_2^{2}}$$

• Por que $n-1$ em $s^2$? Porque $\bar x$ foi estimada com os mesmos dados, o que “gasta” um grau de liberdade. Sem a correção, $s^2$ subestimaria a variância em média.
• Por que dividir por $m_2^{k/2}$? Para a medida ficar sem unidade e invariante a escala ($X\to cX$ não muda $\gamma_1,\gamma_2$). Na Normal, $\gamma_1=0$ e $\gamma_2=3$.
• O script usa $1/n$ nos momentos e $n-1$ no sd. Com $n=1038$ a diferença é desprezível, mas com $n=12$ (frequência) vale saber.
• Quantil empírico: o R usa o tipo 7, interpolando linearmente entre estatísticas de ordem em $h=(n-1)p+1$.

Índice de dispersão e por que $\chi^2$.

$$I=\frac{s^2}{\bar n},\qquad (n-1),I=\frac{\sum_i (N_i-\bar N)^2}{\bar N};\overset{H_0:\text{Poisson}}{\approx};\chi^2_{n-1}$$

Para $\lambda$ não muito pequeno, $N_i\approx \mathcal N(\lambda,\lambda)$. Então $\sum\big((N_i-\bar N)/\sqrt\lambda\big)^2\sim\chi^2_{n-1}$, e trocamos $\lambda$ por $\bar N$. O teste é unilateral (só nos interessa $I$ grande demais).

Excesso médio e sua versão empírica:

$$e(u)=E[X-u\mid X>u],\qquad \hat e(u)=\frac{\sum_i (x_i-u)\mathbb 1{x_i>u}}{\sum_i \mathbb 1{x_i>u}}$$

Exponencial: $e(u)=1/\lambda$ (constante). Pareto: linear crescente. Lognormal: crescente e suave. Perto do q98 poucos pontos sobram, então $\hat e$ fica ruidoso.

---

1b. Estimação por máxima verossimilhança (MV)

🟢 Versão fácil

Ideia. Moeda com 7 caras em 10 lançamentos: o $p$ que torna esse resultado “menos surpreendente” é $p=0{,}7$. A MV faz o mesmo: cada distribuição é um molde de bolo com botões (parâmetros), e a MV gira os botões até os 1.038 sinistros observados ficarem o mais plausíveis possível.

Seus candidatos (resultados do ajuste):

Distribuição	Parâmetros estimados	Personalidade
Exponencial	taxa = 0,000236	“sem memória”, CV sempre = 1 (o seu é 2,0)
Gama	forma 0,643; taxa 0,000151	flexível, mas cauda exponencial
Lognormal	meanlog 7,4015; sdlog 1,3221	cauda pesada: $\ln X$ é Normal
Weibull	forma 0,721; escala 3.242	forma < 1: sinistros grandes mais comuns que na Exp.
Pareto	forma 1,682; escala 3.187	cauda pesadíssima (potência)

Leitura do seu caso: $e^{7{,}4015}\approx 1.638$ é a mediana do modelo lognormal (a observada é 1.417).

Frequência:

• Poisson: $\hat\lambda = 86{,}5$ (a própria média, como no exemplo da moeda).
• Binomial Negativa: $\hat\mu=86{,}5$ e $\hat r=102{,}0$. Variância do modelo = $86{,}5+86{,}5^2/102 = 159{,}8$, contra 86,5 da Poisson e 177 da amostra.

História da Binomial Negativa. Há anos de tempestade e anos calmos, e cada ano tem sua taxa $\Lambda$. Dentro do ano, sinistros são Poisson($\Lambda$). Se $\Lambda$ varia como uma Gama, o resultado é a Binomial Negativa.

🔵 Versão rigorosa

Princípio.

$$\hat\theta=\arg\max_\theta,\ell(\theta),\qquad \ell(\theta)=\sum_{i=1}^n\ln f(x_i;\theta)$$

O log transforma o produto de números minúsculos em soma estável, com o mesmo ponto de máximo.

Exponencial. $\ell=n\ln\lambda-\lambda\sum x_i;\Rightarrow;\frac n\lambda-\sum x_i=0;\Rightarrow;\hat\lambda=1/\bar x$.

Lognormal. Com $y=\ln x$:

$$\ell=-\sum\ln x_i-\frac n2\ln(2\pi\sigma^2)-\frac{1}{2\sigma^2}\sum(\ln x_i-\mu)^2$$

$$\hat\mu=\frac1n\sum\ln x_i,\qquad \hat\sigma^2=\frac1n\sum(\ln x_i-\hat\mu)^2$$

O termo $-\sum\ln x_i$ (o jacobiano) não afeta $\hat\mu,\hat\sigma$, mas afeta o AIC. Para comparar famílias diferentes, todas as verossimilhanças têm de estar na mesma escala (a original). O script faz isso.

Gama. $\ell=n\alpha\ln\beta-n\ln\Gamma(\alpha)+(\alpha-1)\sum\ln x_i-\beta\sum x_i$. Derivando em $\beta$: $\hat\beta=\hat\alpha/\bar x$. Substituindo na derivada em $\alpha$:

$$\ln\hat\alpha-\psi(\hat\alpha)=\ln\bar x-\overline{\ln x}$$

Não tem forma fechada ($\psi$ é a digama), então precisa de otimização numérica.

Pareto (do actuar). É a Pareto II (Lomax), não a Pareto clássica:

$$f(x)=\frac{\alpha\theta^\alpha}{(x+\theta)^{\alpha+1}},\quad x>0,\qquad E[X]=\frac{\theta}{\alpha-1};(\alpha>1)$$

Com $\hat\alpha=1{,}68$: a média é finita, mas a variância é infinita (exige $\alpha>2$). Isso seria um problema sério na Tarefa 2.

O truque da escala (ESC=1000). O otimizador usa derivadas por diferenças finitas, e uma taxa de ~$2\times10^{-4}$ pode ser “pisada” para valores negativos. Ajusta-se $Y=X/1000$ e converte-se de volta:

$$\text{taxa}_X=\text{taxa}_Y/c,\quad \text{escala}_X=c\cdot\text{escala}_Y,\quad \hat\mu_X=\hat\mu_Y+\ln c,$$

com formas e $\sigma$ inalterados. A log-verossimilhança é recalculada na escala original ($\ell_X=\ell_Y-n\ln c$).

Poisson. $\ell=\sum(-\lambda+k_i\ln\lambda-\ln k_i!)\Rightarrow\hat\lambda=\bar n$.

Binomial Negativa (parametrização $r,\mu$; $p=\frac{r}{r+\mu}$):

$$P(N=k)=\binom{k+r-1}{k}p^r(1-p)^k,\qquad E[N]=\mu,\quad Var(N)=\mu+\frac{\mu^2}{r}$$

Por que é Poisson-Gama: com $\Lambda\sim\text{Gama}(r,\text{taxa }\beta=r/\mu)$,

$$P(N=k)=\int_0^\infty \frac{e^{-\lambda}\lambda^k}{k!},\frac{\beta^r\lambda^{r-1}e^{-\beta\lambda}}{\Gamma(r)},d\lambda=\frac{\Gamma(k+r)}{k!,\Gamma(r)}\Big(\frac{\beta}{1+\beta}\Big)^r\Big(\frac1{1+\beta}\Big)^k$$

A variância vem da lei da variância total:

$$Var(N)=E[\Lambda]+Var(\Lambda)=\mu+\frac{\mu^2}{r}$$

Ou seja, a variância da Poisson mais a variância extra causada pela taxa oscilar. Com $r\to\infty$ volta a Poisson. Nessa parametrização, $\hat\mu=\bar n$ e só resta otimizar $r$ (verossimilhança perfilada, no script via optimize em $\ln r$).

---

1c. Escolha do melhor ajuste

🟢 Versão fácil

É um campeonato com três juízes.

Dist.	logLik	AIC	KS: D	KS: p	Peso AIC
Exponencial	−9.708,5	19.418,9	0,225	≈ 0	≈ 0
Gama	−9.628,0	19.260,0	0,131	5·10⁻¹⁶	≈ 0
Lognormal	−9.445,4	18.894,8	0,048	0,016	≈ 1
Weibull	−9.571,2	19.146,4	0,099	2·10⁻⁹	≈ 0
Pareto	−9.488,9	18.981,8	0,071	6·10⁻⁵	≈ 0

• Juiz 1, AIC (nota com desconto por complexidade): a Lognormal vence com 87 pontos de vantagem sobre a Pareto. Na escala do AIC, diferenças acima de 10 já são “sem chance”.
• Juiz 2, KS (a curva teórica anda colada à empírica?): todas reprovam a 5%, inclusive a Lognormal (p = 0,016). Com $n=1038$ o teste é tão sensível que acusa até desvios pequenos. A Lognormal é a “menos errada”, não a “perfeita”. Isso deve constar no relatório.
• Juiz 3, QQ-plot em escala log. Ele compara cada sinistro ordenado com o valor que o modelo esperaria naquela posição. Nos extremos, o maior observado é 87.732 e o modelo esperaria ≈ 128.748. O menor observado é 83,76 e o modelo esperaria ≈ 20,85.

Frequência:

	logLik	AIC	Peso
Poisson	−48,81	99,61	39%
Bin. Negativa	−47,37	98,73	61%

A Binomial Negativa tem o menor AIC, mas a diferença é de 0,88. Abaixo de 2 é praticamente empate técnico. Como $\hat r=102$ é grande, ela quase é uma Poisson. O teste de dispersão (p = 0,021) a favorece, e ela é a escolha conservadora: mais variância em $S$.

Escolha final: severidade Lognormal, frequência Binomial Negativa.

Mistério: o buraco no canto esquerdo. A Lognormal prevê sinistros minúsculos (≈ R$ 20) que não aparecem: o menor real é R$ 83,76. Pode haver um piso, como uma franquia ou um valor mínimo de abertura de sinistro. É uma hipótese minha, que os dados não confirmam.

Mistério: a MV não preserva a média. O modelo lognormal tem $E[X]=3.926$, e a amostra tem 4.243, ou seja, 7,5% a menos. Isso vem do fato de a MV casar a média dos logs, não a média dos valores. Na Tarefa 2 isso se transmite: $E[S]=86{,}5\times3.926\approx339{,}6$ mil, contra ≈ 367 mil se usássemos a média empírica.

🔵 Versão rigorosa

AIC, BIC e peso de Akaike.

$$AIC=2k-2\hat\ell,\qquad BIC=k\ln n-2\hat\ell,\qquad w_i=\frac{e^{-\Delta_i/2}}{\sum_j e^{-\Delta_j/2}},\quad\Delta_i=AIC_i-AIC_{\min}$$

Por que a penalidade $2k$? O $\hat\ell$ é otimista porque os mesmos dados escolhem $\hat\theta$ e avaliam o ajuste. O AIC vem de estimar a divergência de Kullback-Leibler esperada entre o modelo e a verdade. O viés do $\hat\ell$ é ≈ $k$, então $\hat\ell-k$ é quase não viesado, e multiplicar por $-2$ dá o AIC. O BIC penaliza mais ($\ln n>2$ quando $n>7$) porque vem de uma aproximação bayesiana da verossimilhança marginal. Aqui, AIC e BIC concordam.

Teste KS.

$$D=\sup_x|F_n(x)-F(x)|,\qquad F_n(x)=\frac1n\sum\mathbb 1{x_i\le x}$$

O valor crítico aproximado a 5% é $1{,}36/\sqrt n=0{,}042$, e o nosso $D=0{,}048$ já passa dele. Mas o p-valor do ks.test pressupõe parâmetros fixos. Como estimamos com os mesmos dados, a $\hat F$ fica artificialmente perto de $F_n$ e o p sai alto demais. A correção (Lilliefors via bootstrap paramétrico) é: simular 1.000 amostras do modelo ajustado, reestimar os parâmetros em cada uma, recalcular $D$ e comparar. Resultado: p ≈ 0,001 (o mínimo possível com 1.000 réplicas). A rejeição é ainda mais forte do que o p = 0,016 ingênuo, então ela é robusta.

QQ-plot. Plota-se $\big(\ln F^{-1}(p_i),\ \ln x_{(i)}\big)$ com $p_i=(i-\tfrac12)/n$ (o ppoints). Se o modelo serve, os pontos ficam na diagonal. A escala log amplia as caudas, que é o que importa para o risco.

Mistério adicional: assimetria. A assimetria teórica da Lognormal ajustada é

$$\gamma_1=(e^{\sigma^2}+2)\sqrt{e^{\sigma^2}-1}\approx 16{,}9,$$

enquanto a amostral é 5,07. Momentos de ordem alta de caudas pesadas são dominados por valores raríssimos que 1.038 observações quase nunca capturam. Em parte, isso explica a diferença, mas também é um sinal de que a Lognormal pode ser “pesada demais” em alguma região. Confirmação visual: o modelo prevê 51,6% do valor nos 10% maiores (observado: 55,1%) e 15,8% nos 1% maiores (observado: 14,7%), então a concentração está razoável.

Frequência, justificativa. Com $\Delta AIC=0{,}88<2$, o AIC não distingue os modelos. Os argumentos que desempatam a favor da Binomial Negativa são (i) o teste de dispersão ($p=0{,}021$) e (ii) o fato de que, se a taxa realmente variar entre anos e usarmos Poisson, subestimamos a variância de $S$ (o erro é para o lado perigoso). O contra-argumento é a amostra de apenas 12 anos, e vale ressalvar isso.

---

📖 Glossário

Símbolo / variável no R	Definição e para que serve	Valor no seu caso
sin ($X_i$)	Valor individual de cada sinistro	1.038 obs.; média 4.243,13
fr ($N$)	Nº de sinistros em cada ano	12 anos; média 86,5
$S$	Sinistro agregado anual $=X_1+\dots+X_N$	Tarefas 2 em diante
premio ($P$), capital ($C$)	Prêmio anual da carteira e capital inicial	400.558,54 e 255.459,09 (Tarefas 3 e 4)
media ($\bar x$), dp ($s$)	Média e desvio-padrão amostral ($n-1$)	ver tabela 1a
cv	$s/\bar x$: dispersão relativa, sem unidade	2,01 (sev.) / 0,154 (freq.)
assimetria ($\gamma_1$)	Grau de assimetria; > 0 = cauda direita longa	5,07 / 0,52
curtose ($\gamma_2$)	Peso das caudas; Normal = 3	35,75 / 2,60
q90, q95, q99	Valores abaixo dos quais ficam 90%, 95%, 99% das obs.	sev.: 9.957 / 16.363 / 48.017
razao ($I$)	Variância ÷ média da frequência; Poisson ≈ 1	2,046
p_disp	p-valor de $(n-1)I\sim\chi^2_{n-1}$ sob Poisson	0,0207
me ($e(u)$)	Excesso médio acima de $u$; crescente = cauda pesada	6.409 → 20.739
$\theta$, par	Vetor de parâmetros de um modelo	por distribuição
shape ($\alpha$), rate ($\beta$), scale	Forma; taxa ($=1/\text{escala}$); escala	ver tabela 1b
meanlog ($\mu$), sdlog ($\sigma$)	Média e desvio de $\ln X$ na Lognormal	7,4015 e 1,3221
ESC	Fator (1000) para estabilizar o otimizador	1000
loglik ($\hat\ell$)	Log-verossimilhança máxima; maior = ajusta melhor	−9.445,4 (lognormal)
$k$	Nº de parâmetros do modelo	1 ou 2
aic, bic	Critérios com penalidade por complexidade; menor = melhor	18.894,8 / 18.904,7 (lognormal)
dAIC ($\Delta_i$), peso_AIC ($w_i$)	Diferença para o melhor AIC e probabilidade relativa de ser o melhor	lognormal: 0 e ≈ 1
KS_p, $D$	Distância máxima entre FDA empírica e ajustada, e seu p-valor	0,048; 0,016
p_boot	p-valor KS por bootstrap paramétrico (Lilliefors)	≈ 0,001
lambda_hat ($\hat\lambda$)	Taxa estimada da Poisson	86,5
size_hat ($\hat r$), mu_hat ($\hat\mu$)	Parâmetros da Bin. Negativa (dispersão e média)	102,02 e 86,5
melhor_sev, melhor_freq	Modelos escolhidos, usados nas tarefas seguintes	lnorm e nbinom

---

Resumo do Exercício 1: a severidade é assimétrica, com cauda pesada, e a Lognormal é a melhor candidata, mas é rejeitada pelo KS e não preserva a média. A frequência tem leve sobredispersão, e a Binomial Negativa é uma escolha defensável, embora o AIC quase empate. Quando quiser, seguimos para o Exercício 2 (distribuição de $S$).