3- Modelos Geocêntricos¶

Distribuição Normal¶

fonte: wikipedia

fonte: Distribuição Gausiana

Aula - Statistical Rethinking Winter 2019 Lecture 03

Regressão Linear¶

São simples golems estatísticos.

Modelo de média e variância normalmente (Gaussiano) distribuído.
Média como uma combinação aditiva dos pesos das varáveis que a compõem.
Variância é constante.

Gerando um processo de flutuação normal¶

Através da flutuação normal (sobe 1 ou desce 1) iremos construir um processo que tenha como resultante natural o surgimento de um comportamento cuja distribuição podemos descrever como a normal.

# Necessário para desbloquear o asyncio no Jupyter

# Fonte: https://pystan.readthedocs.io/en/latest/faq.html
# Fonte: https://github.com/microsoft/playwright-python/issues/178

# Docs Pystan: https://pystan.readthedocs.io/en/latest/index.html

# -------------
# Instalar nest_asyncio: pip install nest_asyncio - Estára no requirements.txt
# Versões 
# python==3.8.0
# numpy==1.21.1
# pystan==3.2.0
# nest_asyncio==1.5.1
# -------------

# Rodar esse comando antes de import a stan (pystan versão 3.x)
import nest_asyncio
nest_asyncio.apply()

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats

# Add fonts to matplotlib to run xkcd

from matplotlib import font_manager

font_dirs = ["fonts/"]  # The path to the custom font file.
font_files = font_manager.findSystemFonts(fontpaths=font_dirs)

for font_file in font_files:
    font_manager.fontManager.addfont(font_file)
    
# plt.xkcd()   

# Construindo um passeio conforme a distribuição Bernoulli(0.5) em +1 ou -1
n_linhas = 500  # Quantas linhas iremos simular?
n_passos = 30   # Quantidade de passos que iremos simular?

corte_1 = 3   # Indica o primeiro corte no passeio
corte_2 = 14  # Indica o segundo corte no passeio
corte_3 = 26  # Indica o terceiro corte no passeio

aleatorizar_grafico = True  # True para ver a simulação um pouco mais aleatorizada visualmente. False para ver sem efeito.

passeios = []  # Armazenar todos os passeios gerados num array

# Parâmetros da simulação
plt.figure(figsize=(17, 9))
plt.rcParams['axes.facecolor'] = 'lightgray'  # Alterando a cor de fundo, para ficar mais elegante. 

for linha in range(n_linhas):  # Quantas linhas iremos plotar
    
    passeio = [0]  # Inicializando todos os pontos com o valor 0 (conforme dito acima!)
    novo_ponto = 0  # variável de controle 
    
    for passo in range(n_passos):  # Quantos passos iremos simular?
        
        if aleatorizar_grafico:
            # +1 (subida) ou -1 (descida) - Com efeito de aleatorizacão (apenas para melhorar a visualização).
            novo_ponto = novo_ponto + np.random.uniform(1, -1) 
        else:    
            # +1 (subida) ou -1 (descida) - Sem efeito de aleatorizacão - "Pixelizado"
            novo_ponto = novo_ponto + np.random.choice([1, -1]) 
        
        passeio.append(novo_ponto)
    
    passeios.append(passeio)  # Armazendo todas as linhas geradas
    plt.plot(passeios[linha], color='darkblue', alpha=0.5, linewidth=0.1)  # Plotando a linha gerada
    
# Retas horizontais dos cortes
plt.vlines(corte_1, -15, 15, color='darkred', ls='--',  linewidth=3)
plt.vlines(corte_2, -15, 15, color='orange', ls='--', linewidth=3)
plt.vlines(corte_3, -15, 15, color='darkgreen', ls='--', linewidth=3)

# Configurando infos do gráfico    
plt.title("Passeio Aleatório {+1, -1} com p(0.5)")
plt.xlabel('Número de passos')
plt.ylabel('Valor')

plt.grid(axis='y', ls='--', color='white')

plt.show()

passeios = np.array(passeios).T  # Transpondo os vetores de passeios para ajuste na plotagem.

# ------------
# Gerando os histogramas de todas as linhas para alguns passos específicos ao longo do passeio.
fig, [ax1, ax2, ax3] = plt.subplots(1, 3, figsize=(17, 6))

ax1.hist(passeios[corte_1], color='darkred', rwidth=0.8, density=True)
ax1.grid(axis='y', ls='--', color='white')
ax1.set_title("Histograma \n Normalmente distribuído \n  passo " + str(corte_1))
ax1.set_xlabel('Valor')
ax1.set_xlim(-20, 20)
ax1.set_ylim(0, 0.25)

ax2.hist(passeios[corte_2], color='darkorange', rwidth=0.8, density=True)
ax2.grid(axis='y', ls='--', color='white')
ax2.set_title("Histograma \n Normalmente distribuído \n  passo " + str(corte_2))
ax2.set_xlabel('Valor')
ax2.set_xlim(-20, 20)
ax2.set_ylim(0, 0.25)

ax3.hist(passeios[corte_3], color='darkgreen', rwidth=0.8, density=True)
ax3.grid(axis='y', ls='--', color='white')
ax3.set_title("Histograma \n Normalmente distribuído \n passo " + str(corte_3))
ax3.set_xlabel('Valor')
ax3.set_xlim(-20, 20)      
ax3.set_ylim(0, 0.25)
    
plt.show()

Todas as flutuações geradas se equilibriam entre si, gerando assim uma curva simétrica. Porém cada uma delas individualmente, e não necessariamente, se equilibram entre si. Isso é a \(Normal\)!

Porque a Normal?¶

O interessante, e frustrante ao mesmo tempo, é que mesmos nós sabendo como cada uma das curvas são geradas, não temos acesso a como cada uma delas individualmente faz o seu percurso, pois temos muito pouca informação sobre como é o processo gerador dessas curvas.

Assim, quando vimos um comportamento de um evento que é normalmente distribuído, por exemplo, não sabemos o que tem la dentro! Não temos um acesso intuitivo da informação contida no comportamento individual de cada passeio.

Mas sabendo que o processo como um todo tem o comportamento, aproximadamente, normalmente distribuído, podemos usar essa informação para falsificar ou não hipóteses de um estudo particular.

\[ x \sim normal(\mu, \sigma) \]

E isso é muito legal, pois podemos descrever o processo como um todo, com apenas por dois números, pois toda a informação do funcionamento da maquinaria subjacente do processo foi eliminada. E tudo que foi preservado foi a média (\(\mu\)) e o desvio padrão (\(\sigma\)), por isso que precisamos apenas desses dois números para descrever uma distribuição normal.

Note

Essa é toda a informação que conseguimos extrair da maquinaria subjacente, a Natureza!

O mais terrível é que não podemos saber o funcionamento do processo a partir de um simples histograma! Temos que realmente fazer ciência para isso, temos de cavar mais fundo, medir coisas mais difíceis para só assim descobrir mais algumas coisas nas profundezas do mecanismo de geração do processo.

Isso é verdade não apenas para a \(Normal\), mas para muitas outras distribuições! Todas as distribuições de máxima entropia tem essa mesma propriedade, que muitos processos diferentes possuem a mesma distribuição de frequências.

Perspectiva Ontológica¶

Ontologia é o ramo da filosofia que estuda a natureza do ser, da existência e da própria realidade. Vamos observar a distribuição \(Normal\) sobre uma perspectiva ontológia:

Processo que adiciona flutuações amortecidas como resultado.
As flutuações amortecidas se aproximam de uma gaussiana.
Mas nenhuma informação sobrou do processo gerador, apenas a média e a variância.
Não podemos inferir o processo a partir dessa distribuição.

Se quisermos construir um modelo visando responder de modo mais conservador possível, no qual tudo que estiver disposto a dizer sobre algum conjunto de medidas, como medida de alturas, é que eles têm variância finita e podemos usar a distribuição gaussiana.

Mesmo se eles forem distorcidos ou alguma outras coisas, a distribuição gaussiana cobrirá uma faixa mais ampla de valores do que qualquer outra com a mesma média e a mesma variância!

Essa é a distribuição mais é conservadora que podemos assumir. Qualquer outra distribuição será mais restrita por isso terá mais informação incorporada, por isso a opção mais conservadora que podemos propor é a gaussiana, no qual tudo que precisamos para ela é a média a variância. (Isso será demostrado mais para frente no curso.)

Perspectiva Epistemologica¶

Epistemologia, também conhecida como a Teoria do Conhecimento, é o ramo da filosofia que estuda como o ser humano ou a própria ciência adquire e justifica seus conhecimentos. Vamos observar a \(Normal\) sobre essa perspectiva:

Conhecemos apenas a média e a variância.
A menos surpreendente e mais conservadora (máxima entropia) distribuição gaussiana.
É a distribuição natural da máxima entropia.

Modelos Lineares¶

“Modelo lineares generalizados”: teste t, regressão simples, regressão multipla, ANOVA, ANCOVA, MANOVA, MANCOVA, etc, etc, etc. Todos esses modelos são a mesma coisa, são todos modelos os lineares.

A seguir, iremos construir modelos lineares manualmente e do zero para entendermos como é o processo de pensamento de construção de hipóteses e, também, iremos construir um gráfico no qual iremos visualizar a nossa incerteza do sistema.

Warning

Iremos aprender estratégias e não procedimentos!

Acordando a linguagem que iremos trabalhar¶

Relembrando o primeiro modelo que nós fizemos anteriormente, o lançamento do globo, nós tinhamos:

\[ w \sim Binomial(N, p) \]

\[ p \sim Uniform(0, 1) \]

\(w\): Resultado
\(\sim\): “É distribuído como…”
\(Binomial( N, p)\): É a função de distribuição de probabilidade dos dados (verossimilhança, ou likelihood no inglês)
\(p\): É o parâmetro que iremos estimar
\(Uniform(0, 1)\): Nossa distribuição à priori

No exemplo do lançamento do globo, o que estávamos interessados em saber era a proporção de água na superfície (\(p\)) da Terra. O processo que usamos para coletar os dados, a amostragem, automaticamente nos propõem, de modo quase impositivo, a utilização da estrutura Binomial.

A estrutura que distribui os dados binominalmente, nada mais é do que a contagem do número de vezes que é possível acontecer tal fato, dado que a proporção \(p\) é fixa para todas às vezes \(N\) que iremos retirar uma amostra.

Essa é a linguagem que iremos usar nesse curso. A notação matemática padrão, na qual consiste em uma maneira de comunicar a todos os nossos colegas de trabalho o que foi pensado, quais foram as suposições, quais foram as hipóteses (isto é, as sugestões) e, também, qual a estrutura foi proposta para um possível caminho para se entender o problema.

Note

A linguagem matemática será a forma de comunicar ao mundo sua forma de pensar!

Linguagem da modelagem¶

Do mesmo modo que escrevemos o modelo acima iremos escrever todos os outros modelos, inclusive, a programação usando a linguagem probabilística Stan, o qual tem a sua sintaxe próxima a esse formato.

Para qualquer modelo que iremos elaborar, seja uma regressão linear simples ou modelos mais elaborados, temos que criar uma lista com todas as variáveis que irão participar desse modelo. Algumas coisas que iremos observar nos dados, como a contagem do número de vezes que o globo foi lançado. E também haverá outras coisas que não vamos poder observar, como a inclinação de uma reta de uma modelo linear simples, ou proporção de água no globo, pois elas não são entidades observáveis.

Por isso temos que listar todas as variáveis e então defini-las.

Listando as variáveis:

\[ y_i \]

\[ \mu_i \]

\[ \beta \]

\[ \sigma \]

\[ x_i \]

Um modelo de regressão, assim como qualquer outro modelo, será escrito da mesma forma que vimos anteriormente, porém esse modelo terá muito mais símbolos pois temos muito mais variáveis participando da explicação, mas é a mesma coisa, apenas precisamos definir a cada um desses símbolos.

O motor dos modelos de regressão linear é a segunda linha da definição abaixo (a seguir iremos construir esse modelo desde o início):

\[ y_i \sim Normal(\mu_i, \sigma) \]

\[ \mu_i = \beta x_i \]

\[ \beta \sim Normal(0, 10) \]

\[ \sigma \sim Exponencial(1) \]

\[ x_i \sim Normal(0, 1) \]

Essa segunda linha, \( \mu_i = \beta x_i \), geralmente a parte que é mais confusa para se entender (iremos ver um exemplo a seguir), mas ela significa que média da distribuição normal de cima (média do \(y_i\)) é geralmente definida por uma equação, o que define em termos de alguma outra variável que observamos, \(x_i\).

O \(x_i\) é uma variável que ajuda a explicar o comportamento do \(y_i\), ou seja, \(x_i\) é uma variável explicativa.

Mas repare que \(x_i\) também tem uma distribuição, nós normalmente não nos preocupamos em definir uma distribuição para as variáveis explicativas, pois não iremos prevê-las, porém existe uma grande vantagem em fazer essas suposições de distribuição, pois podemos fazer coisas muito legais com esse fato, iremos ver essas coisas mais à frente, tais como medidas de erros e também dados faltantes.

Assim, como todas as variáveis têm sua própria definição de distribuição, se você não sabe algo sobre algumas delas podemos colocar dentro do modelo e assim iremos ganhar automaticamente um poder inferencial do que estávamos perdendo antes.