Na aula anterior vimos que as ferramentas que aprendemos em LEGO II nem sempre são as melhores dependendo do tipo de problema que queremos enfrentar. Se trabalhamos com variáveis dependentes qualitativas, que assumem categorias de valores e não valores contínuos, as premissas por trás da aplicação do método dos Mínimos Quadrados Ordinários são violadas. Embora o modelo linear sirva como aproximação exitem métodos que fornecem melhores estimadores.

Na primeira aula vimos um desses métodos, o método da Máxima Verossimilhança onde observamos o quão ‘verossímel’ é uma série de valores hipóteticos dos parâmetros que queremos estimar na geração dos dados que temos. O valor mais ‘verossímil’, ou seja com a máxima verossimilhança, é o valor que mais se aproxima de um valor verdadeiro do parâmetro, condicional aos nossos dados.

Na aula passada introduzi o modelo logit que utiliza a distribuição logística para estimar os parâmetros que relacionam variáveis explicativas a uma variavél dependente dicotômica, isto é, que assume apenas dois valores, em geral 0 e 1. Hoje vamos ver como podemos estimar esses valores usando a máxima verossimilhança.

A primeira coisa a fazer é atribuir um modelo de distribuição de probabilidade às nossas observações, isto é, atribuir probabilidades a todos os valores que uma variável aleatória \(Y_{i}\) podem assumir. Assim teremos o componente estocástico de nosso modelo. No caso do referendo no Chile os valores possíveis seriam o ‘No’ e o ‘Sí’, que assumiriam os valores 0 e 1 com a seguinte probabilidade:

\[\begin{align} &Pr(Y_{i}=1)=\pi \\ &Pr(Y_{i}=0)=1-\pi \\ &Pr(Y\neq0,1)=0 \end{align}\]

\(\pi\) é o parâmetro da distribuição. O intervalo de valores que o parâmetro pode assumir é o espaço paramétrico simbolizado por \(\Theta\). No caso do referendo e de qualquer variável binária, esse intervalo vai de 0 a 1, isto é, \(\pi \in \Theta = [0,1]\).

Mais conveniente do que listar as probabilidades de todos os valores que \(Y_{i}\) pode assumir é usar alguma equação que ‘resuma’ essa lista. Essas equações são derivadas de certas premissas como, no caso de variáveis dependentes dicotômicas, a premissa de que ela só pode assumir dois valores com probabilidade não zero, isto é, a variável \(Y_{i}\) deve ter dois resultados mutuamente exclusivos, \(y_{i} = 0\) ou \(y_{i} = 1\) que são exaustivos. Essa premissa deve se apoiar em princípios substantivos, isto é, em nossa teoria. Se estamos examinando a votação para um cargo majoritário e o eleitor só pode escolher um de dois candidatos, então as premissas do modelo estatístico correspondem a de nossa teoria. Se o eleitor tivesse mais de um voto ou se pudesse escolher mais de um candidato, essas premissas não se aplicariam e terímaos que procurar um modelo que atendesse melhor esses requisitos teóricos. Isto é: os modelos estatístico que utilizamos devem ser sempre derivado de nosso modelo teórico substantivo.

No caso de variáveis binárias o modelo estatístico que melhor se aproxima de um processo como o da votação em um referendo é a distribuição de Bernoulli cuja equação resume o caráter exclusivo e exaustivo dos valores possíveis de \(Y_{i}\):

\[\begin{equation} Y_{i} \sim f_{bern}(y_{i}|\pi) \begin{cases} \pi^{y_{i}}(1-\pi)^{1-y_{i}}, & \text{para }y_{i} = 0,1 \\ 0, & \text{do contraŕio} \end{cases} \end{equation}\]

Se fixamos o valor do parâmetro \(\pi\) conseguimos, com a fórmula acima, calcular \(Pr(Y_{i} = 1)\). Mas, como vimos na primeira aula, o que nos interessa é conhecer o valor do parâmetro \(\pi\) dado observações de \(Y_{i}\). Para isso, partindo também de outra premissa do modelo estatístico, a de independência de \(Y_{i}\) e \(Y_{j}\) para todo \(i \neq j\) obtemos o modelo completo:

\(Pr(Y|\pi) = \prod_{i = 1}^{n}\pi_{i}^{y_{i}}(1-\pi_{i})^{1-y_{i}}\)

O modelo logit que vimos nos dá o componente sistemático de nosso modelo estatístico e é formalizado por:

\(\frac{\pi}{1-\pi} = exp(\alpha + \beta_{1}X_{i1} + \ldots + \beta_{k}X_{ik})\)

Que modificado nos dá:

\(\pi = \frac{1}{1 + \text{exp}(-x{i}\beta)}\)

Onde \(\beta\) é um vetor com os parâmetros do modelo. Substituindo no componente estocástico (“reparameterizando”) temos:

\[ Pr(Y|\beta) = \prod_{i = 1}^{n}\big(\frac{1}{1 + exp(-x_{i}\beta)}\big)^{y_{i}}\big(1-\frac{1}{1 + exp(-x_{i}\beta)}\big)^{1-y_{i}} \\ = \prod_{i = 1}^{n}[1 + exp(-x_{i}\beta)]^{-y_{i}} [1 + exp(x_{i}\beta)]^{-(1-y_{i})} \]

Podermos tornar essa equação mais ‘tratável’ utilizando o logaritmo e o fato de que a verossimilhança \(\mathcal{L}(\beta|y)\) é proporcional à equação acima. Obtemos:

\[\begin{align} \mathcal{L}(\tilde\beta|y) &= \sum_{i = 1}^{n}-y_{i}\text{ln}[1 + exp(-x_{i}\tilde\beta)]-(1-y_{i})\text{ln}[1 + exp(x_{i}\beta)]\\ &= \sum_{i = 1}^{n}y_{i}x_{i}\beta - log(1-exp(x_{i}\beta)) \end{align}\]

Para achar o valor que maximiza o logaritmo da verossimilhança acima diferenciamos a equação e igualamos a derivada a zero e resolvemos para \(\tilde\beta\).

A derivada parcial de \(\mathcal{L}(\tilde\beta|y)\) comrelação a \(\beta\) é

\[ \frac{\partial \mathcal{L}\beta}{\partial \beta} = \sum_{i = 1}^{n} x_{i}(y_{i} - \frac{1}{1 + exp(-x_{i}\beta)})\]

(In)felizmente não há como resolver isso com métodos analíticos e temos que recorrer a métodos numéricos (computacionais) para estimar \(\tilde\beta\).

Se observarmos o comando glm() que utilzamos no R para estimar os parâmetros de um modelo logit veremos a seguinte passagem:

##------------- THE Iteratively Reweighting L.S. iteration -----------
#        for (iter in 1L:control$maxit) {
#            good <- weights > 0
#            varmu <- variance(mu)[good]
#            if (anyNA(varmu))
#                stop("NAs in V(mu)")
#            if (any(varmu == 0))
#                stop("0s in V(mu)")
#            mu.eta.val <- mu.eta(eta)
#            if (any(is.na(mu.eta.val[good])))
#                stop("NAs in d(mu)/d(eta)")
# (...)

Essa passagem nos diz que o R usa um algoritmo chamado Iterative Reweight Least Square (IRLS) para estimar o parâmetro \(\beta\). O que esse algoritmo faz é determinar uma objetivo (na verdade uma loss function) e, por tentativa e erro ir se aproximando deste objetivo. O objetivo, em geral, é minimizar uma função quadrática do erro, a discrepância entre os valores observados e os valores preditos de \(Y\). Se fizermos os gráficos da distribuição dessa função com relação aos valores preditos veremos um “U” ou um a“vale”. O valor de \(\beta\) que minimiza a função se encontra no fundo deste “vale”. Para descobrir o fundo do “vale” além da derivada acima usamos a segunda derivada ou, no caso de um vetor \(\beta\), uma matriz de segundas derivadas conhecida como a matriz Hessiana:

\[ \frac{\partial^2 \mathcal{L}\beta}{\partial \beta \partial \beta^{T}} = -\sum_{i = 1}^{n} x_{i}x_{i}^{T} \Big(\frac{1}{1 + exp(-x_{i}\beta)}\Big) - \Big(1 -\frac{1}{1 + exp(-x_{i}\beta)}\Big)\]

No caso do IRLS o objetivo é minimizar a distância entre o valor \(y_{i}\) observado e o estimado \(X\beta\), quase como no MQO:

\[\underset{\boldsymbol \beta}{ \operatorname{arg\,min} } \big\| \mathbf y - X \boldsymbol \beta \|_p = \underset{\boldsymbol \beta}{ \operatorname{arg\,min} } \sum_{i=1}^n \left| y_i - X_i \boldsymbol\beta \right|^p \]

O algoritmo testa vários valores de \(\beta\) atualizando esse valor do seguinte modo:

\[ \beta_{novo} = \beta_{antigo} - \left( \frac{\partial^2 \mathcal{L}\beta}{\partial \beta \partial \beta^{T}} \right)^{-1}\frac{\partial \mathcal{L}\beta}{\partial \beta} \]

e vai ‘aprendendo’ com os testes. Escrevendo em formato de matrizes temos:

\[\begin{align} \frac{\partial \mathcal{L}\beta}{\partial \beta} &= \boldsymbol{X}^{T}(\boldsymbol{y - p}) \\ \\ \frac{\partial^2 \mathcal{L}\beta}{\partial \beta \partial \beta^{T}} &= \boldsymbol{X^{T}WX} \end{align}\]

Desse modo a atualização de \(\beta\) é dada por

\[ \boldsymbol\beta^{(t+1)} = (X^{\rm T} W^{(t)} X)^{-1} X^{\rm T} W^{(t)} \mathbf{y} \]

Para isso inicializamos uma matriz de pesos \(W\) que vai sendo atualizada iterativamente. A atualização é dada por:

\[ w_i^{(t)} = \big|y_i - X_i \boldsymbol \beta ^{(t)} \big|^{p-2} \]

Isto tudo quer dizer que aos poucos, por tentativa e erro, procura-se um parâmetro que minimiza o erro absoluto e assim vai se aproximando do parâmetro com a máxima verossimilhança. É como se aos poucos fomos descendo o “vale” aprendendo com os erros e corrigindo nossas estimativas. Veremos outros algoritmos como esse nas aulas sobre simulação e aprendizado de máquina.

A mesma estratégia está por trás da estimativa de parâmetros de outros modelos estatísticos. Derivamos a distribuição, ou seja o componente estocástico, com base em nossa teoria e utilizamos metodos numéricos para estimar o componente sistemático. Quando utlizamos funções como a glm() no R ela nos pede uma fórumla e uma família. A fórmula é definida como uma descrição do modelo a ser ajustado. A família é uma descrição da distribuição, equivalente ao componente estocástico. A função também pede um ‘link’, que equivale a forma funcional de nosso componente sistemático.

No exemplo do pleibiscito no Chile, nossso modelo teórico procurava prever o voto ‘Sim’ ou ‘Não’ dadas certas características do eleitor. Como vimos acima esse modelo nos levou à escolha de uma distribuição de bernoulli. Essa distribuição faz parte da família da distribuição binomial, com o logit como link (a forma funcional que assume o parâmetro que estamos interessados, \(\pi\)). No R escrevemos:

# glm(formula, family = binomial(link = "logit"), ...)

No caso da nomeação de mulheres para compor o ministério nos governos Dilma e Bolsonaro, que vimos na primeira aula, a distribuição é da família binomial, com o link logit, mas onde, diferente do caso do pleibiscito, temos apenas uma soma: um número de ‘sucessos’ (ministras) em certo número de ‘tentativas’ (ministérios) ou seja, a soma de uma série de processos de Bernoulli (uma sequência de ‘Sins’ e ‘Nãos’), cada um independente do outro e com a mesma probabilidade de ocorrência \(\pi\).

\[P(y) = \binom{N}{y} \pi^y(1-\pi)^{N-y} = \frac{N!}{y!(N-y)!}\pi^y(1-\pi)^{N-y}\] \[ \pi = \text{logit}^{-1}(\boldsymbol{X \beta})\]

Substituindo a componente sistemático no componente estocástico temos a verossimilhança de cada observação. Para todas as observações temos:

\[ \begin{align} \mathcal{L}(\pi\mid y, N) &= \prod_{i = 1}^{n}\binom{N}{y} (\frac{1}{1 + e^{-X \beta}})^y(1-\frac{1}{1 + e^{-X \beta}})^{N-y} \\ \text{ln}\mathcal{L}(\pi\mid y, N) &= \sum_{i = 1}^{n}y_{i}\text{ ln}\left(\frac{1}{1 + e^{-X \beta}}\right) + (N-y)\text{ ln}\left(1- \frac{1}{1 + e^{-X \beta}}\right) \\ &= \sum_{i = 1}^{n} -y_{i}\text{ ln}\left(1 + e^{-X \beta}\right) - (N-y)\text{ ln}\left(1 + e^{-X \beta}\right) \end{align} \]

Métodos númericos permitem encontrar o máximo a partir da derivada do log da verossimilhança e assim estimar os coeficientes de \(\boldsymbol \beta\). A interpretação destes coeficientes é a mesma que no caso da variável binária. Com já examinamos um exemplo em aula anterior não vamos nos deter mais tempo aqui.

Modelos de contagem - Poisson

Vamos supor agora que nossos dados nos dão uma contagem sem um limite, não mais o número de ‘sucessos’ em um número definido de ‘tentativas’, mas um número de eventos em certo intervalo de tempo, como a ocorrência de golpes de estado, guerras, retweets etc. onde o número de ocorrências não tem um limite natural, isto é, \(Y_{i} \in (0,1,2,\ldots)\). Isso seria equivalente a se estimar o número de ocorrências de eventos que acontecem completamente ao acaso a uma taxa de \(\lambda\) por unidade de tempo.

Uma maneira de estimar essa taxa, o parâmetro de interesse, é partir da premissa de que no momento \(t = 0\), ainda não temos ocorrências, então \(N (0) = 0\). Agora dividimos a meia-linha \([0, \infty)\) em pequenos subintervalos de comprimento \(\delta\), conforme mostrado na Figura 11.

Poisson-coinflip
Dividindo a semi-linha \([0, \infty)\) em pequenos intervalos de comprimento \(\delta\).

Cada subintervalo corresponde a um intervalo de tempo de comprimento \(\delta\). Assim, os intervalos são \((0, \delta]\), \((\delta, 2 \delta]\), \((2 \delta, 3 \delta]\), \(\cdots\). De maneira geral, o \(k\) th intervalo é \(\big((k-1) \delta, k \delta \big]\). Vamos supor que em cada intervalo de tempo, lançamos uma moeda para a qual \(P(H) = p = \lambda \delta\). Se a moeda cair cara, dizemos que temos uma ocorrência nesse subintervalo. Caso contrário, dizemos que não temos nenhuma ocorrência nesse intervalo. A Figura 2 mostra esse processo. Aqui, temos uma ocorrência no tempo \(t = k \delta\), se o lançamento da \(k\)ésima moeda resultar em cara.

Poisson-coinflip-2
Processo Poisson como o limite de um processo de Bernoulli.

Agora, seja \(N(t)\) definido como o número de ocorrências (número de caras) do tempo \(0\) ao tempo \(t\). Existem \(n \approx \frac{t}{\delta}\) slots de tempo no intervalo \((0, t]\). Assim, \(N(t)\) é o número de caras em \(n\) lançamentos de moeda. Conclui-se que \(N(t) \sim Binomial(n, \pi)\). Observe que aqui \(p = \lambda \delta\) e lembre que a média de uma distribuição binomial é dada por \(n\pi\) então

\[ \begin{align*} n\pi &= n \lambda \delta\\ &=\frac{t}{\delta} \cdot \lambda \delta\\ &=\lambda t. \end{align*} \]

Portanto, se ignorarmos a constante \(t\) temos a igualdade \(\pi = \frac{\lambda}{n}\).

Como acabmos de ver acima a distribuição binomial é dada por

\[P(y) = \binom{N}{y} \pi^y(1-\pi)^{N-y} = \frac{N!}{y!(N-y)!}\pi^y(1-\pi)^{N-y}\] Se substituirmos \(\pi\) por \(\lambda/n\) nesta equação temos

\[P(Y = y) = \frac{N!}{y!(N-y)!}\left(\frac{\lambda}{n} \right)^y \left(1-\frac{\lambda}{n} \right)^{N-y}\] Como em um processo de Poisson podemos dividir o intervalo de tempo \(t\) em tantos subintervalos \(\delta\) quanto quisermos isso equivale a aproximar \(n\) (o númerode tentativas da binomial) ao infinito.

\[\lim_{n \to \infty}P(Y = y) = \lim_{n \to \infty}\frac{N!}{y!(N-y)!}\left(\frac{\lambda}{n} \right)^y \left(1-\frac{\lambda}{n} \right)^{N-y}\] Depois de algum manbo-jambo algébrico2 obtemos:

\[\left( \frac{\lambda^{y}}{y!} \right)\lim_{n \to \infty}\frac{N!}{(N-y)!}\left(\frac{1}{n^{y}} \right)^y \left(1-\frac{\lambda}{n} \right)^{-y} = \left( \frac{\lambda^{y}}{y!} \right)(1)(e^{-\lambda})(1) = \left( \frac{\lambda^{y}e^{-\lambda}}{y!} \right)\]

Aqui usamos a definição de \(e\):

\[ e = \lim_{x \to \infty} \left(1 + \frac{1}{x} \right)^x \\ e^{\delta} = \lim_{x \to \infty} \left(1 + \frac{1}{x} \right)^{x\delta} \] Definimos \(x = -\frac{n}{\lambda}\) e fazemos

\[\lim_{n \to \infty}\left(1-\frac{\lambda}{n} \right)^{n} = \lim_{n \to \infty}\left(1-\frac{1}{x} \right)^{x(-\lambda)} = e^{-y}\]

Portanto, como \(\delta \rightarrow 0\) e, portanto \(N \rightarrow \infty\), pelo teorema central do limite a distribuição de \(N(t)\) converge para uma distribuição de Poisson com taxa \(\lambda\). De forma mais geral, podemos argumentar que o número de ocorrências em qualquer intervalo de comprimento \(\tau\) segue uma distribuição \(Poisson(\lambda \tau)\) quando \(\delta \rightarrow 0\).

O Processo Poisson pode ser descrito como:

Seja \(\lambda> 0\) fixo. O processo de contagem \(\{N(t), t \in [0, \infty)\}\) é chamado de processo de Poisson com taxa s \(\lambda\) se todas as seguintes condições são válidas:
  1. \(N (0) = 0\);
  2. \(N (t)\) tem incrementos independentes;
  3. o número de chegadas em qualquer intervalo de comprimento \(\tau> 0\) tem distribuição \(Poisson(\lambda \tau)\).

A fórmula geral para o Poisson é escrita como:

\[ \begin{align} \nonumber f(y_{i}\mid\lambda) = \left\{ \begin{array}{l l} \frac{e^{-\lambda}\lambda^{y_{i}}}{y_{i}!} &\quad \text{para } \lambda > 0 \text{ e } y_{i} = 0,1,\ldots,\\ & \quad \\ 0 & \quad \text{do contrario} \end{array} \right. \end{align} \]

Essa fórmula, que vocês encontram nos textos de estatística, é equivalente à derivação da binomial e, mais uma vez, nada mais é do que uma maneira de contarmos eventos.

A estimativa do parâmetro \(\lambda\) pode ser obtida por Máxima Verossimilhança. Para isso partimos da distribuição conjunta das probabilidades de todos os valores de \(\lambda\).

\(Pr(Y\mid\lambda) = \prod_{i = 1}^{n}\frac{e^{-\lambda}\lambda_{i}^{y_{i}}}{y_{i}!}\)

Para derivar o componente sistemático observamos que \(E(Y_i)=\lambda\) e, como \(Y_i\) é sempre positivo sem um lmite superior e tem, por definição, o limite inferior = 0 uma função linear não seria uma boa aproximação, por assumir valores negativos e uma função logit também não, por ter limite superior = 1. Precisamos então de uma função da classe das funções de valores positivos (entre 0 e \(\infty\)). É natural também pensar que o efeito \(\beta\) é decrescente, isto é, seria mais provável passarmos de 1 para 2 ocorrências do que de 20 para 21. Levando isso em conta chegamos à família de funções exponenciais e o componente sistemático pode ser escrito como

\(E(Y_i)=\lambda = \text{exp}(x_{i}\beta)\)

Para derivar o log da verossimilhança, substituimos essa função (componente sistemático) na distribuição de probabilidade (componente estocástico) e tiramos o log \[ \begin{align} \text{ln}\mathcal{L}(\tilde\beta \mid y) &= \sum_{i = 1}^{n}{(x_{i}\tilde\beta)y_{i} - \text{exp} (x_{i}\tilde\beta) - \text{ln }y_i} \\ &= \sum_{i = 1}^{n}{(x_{i}\tilde\beta)y_{i} - \text{exp} (x_{i}\tilde\beta)} \end{align} \]

Podemos elmiminar \(- \text{ln }y_i\) por não ser função de \(\tilde\beta\). Para estimar \(\tilde\beta\) usamos algoritmos como o IRLS na função glm().

[Exemplo de Gelman e Hill - Poisson, exposure & overdispersion]

VD categórica com mais de uma categoria

Uma extensão do modelo binomial se torna necessária quando temos mais do que duas categorias. No caso do plebiscito, por exemplo, a variável dependente original trazia a categoria dos indecisos. Repare que as três categorias, “sí”, “no” e “indeciso”, não possui uma ordem entre elas, o que caracateriza um modelo logit multinomial. Neste modelo calculamos qual a probabilidade do voto do eleitor \(i\) cair em alguma das categorias. Seja \(\pi_{ij}\) a probabilidade da \(iª\) observação cair na \(jª\) categoria, isto é, \(\pi_{ij} \equiv \text{Pr}(Y_{i} = {j})\) para \(j = 1,\dots,m\) categorias:

\[\begin{align} \pi_{ij} &= \frac{\text{exp}(\gamma_{0j}+\gamma_{1j}X_{i1}+\dots+\gamma_{kj}X_{ik})}{ 1 + \sum_{l=1}^{m-1}\text{exp}(\gamma_{0j}+\gamma_{1j}X_{i1}+\dots+\gamma_{kj}X_{ik})} \text{ para } j = 1,\dots,m-1\\ \pi_{im} &= 1 - \sum_{l=1}^{m-1}\pi_{ij} \end{align} \]

Com alguma manipulação algébrica temos

\[ \text{log}\frac{\pi_{ij}}{\pi_{im}} = \gamma_{0j}+\gamma_{1j}X_{i1}+\dots+\gamma_{kj}X_{ik} \text{ para } j = 1,\dots,m-1\] Vemos que a última categoria se torna um tipo de base, e \(\pi_{im}\) entra no denominador do cálculo da probabilidade nas outras categorias. O coeficiente, portanto, nos dá o log da chance de se pertencer a categoria \(j\) contra a categoria base \(m\) [como fica a equação acima no caso de termos apenas duas categorias?]

[Exemplo Fox e Weisberg 1]

Muitas vezes temos mais do que duas categorias exaustivas e mutuamente exclusivas. Avaliação de governo em geral envolvem as categoras “Ótimo”, “Bom”, “Regular”, “Ruim” e “Péssimo”. Neste caso temos 5 categorias exaustivas e mutuamente exclusivas e que, ainda, mantém uma ordem entre si. Podemos pensar que essas cinco categorias são partições de um a’variável latente’ \(z\). Essa variável, não observada (ex: ‘conformismo’), é função linear de nossas variáveis independentes. Conforme \(z\) vai assumindo diferentes valores muda-se a categoria de \(Y\). Podemos estimar pontos latentes a partir do qual ocorre a transição entre categorias:

\[ \begin{align} \nonumber y_{i} = \left\{ \begin{array}{l l} 1 \text{ se } z_{i} < c_{1.5}\\ 2 \text{ se } z_{i} \in(c_{1.5},c_{2.5})\\ 3 \text{ se } z_{i} \in(c_{2.5},c_{3.5})\\ 4 \text{ se } z_{i} \in(c_{3.5},c_{4.5})\\ 5 \text{ se } z_{i} > c_{4.5}\\ \end{array} \right. \\ \\ z_{i} \sim \text{logistic}(x_{i}, \sigma^2) \end{align} \]

Onde \(c_{(.)}\) é o valor dos pontos de transição na escala de \(x\) e \(\sigma\) nos dá o quão forte é a transição. A distribuição de probabilidade acumulada de \(y\) é dada por:

\[ \begin{align} \text{Pr}(y \leq j \mid \boldsymbol{x}) &= \text{Pr}(z \leq c_{j} \mid \boldsymbol{x}) \\ &= \text{Pr}(\alpha + \beta_{1}x_{1} + \dots + \beta_{k}x_{k} + \epsilon \leq c_{j} \mid \boldsymbol{x}) \\ &= \text{Pr}(\epsilon \leq c_{j} - \alpha - \beta_{1}x_{1} - \dots - \beta_{k}x_{k} \mid \boldsymbol{x}) \end{align} \]

Para \(j = 1,2 \dots, q-1\), isto é, a quebra de \(z\) em \(q\) categorias acrescenta \(q-1\) parâmetros ao modelo. Assumindo que \(\epsilon\) segue uma distribuição logística então temos o modelo logit ordinal:

\[ \begin{align} \text{logit}[\text{Pr}(y >j\mid \boldsymbol{x})] &= \text{log} \frac{\text{Pr}(y >j\mid \boldsymbol{x})}{\text{Pr}(y \leq j\mid \boldsymbol{x})} \\ &= (\alpha - c_{j}) - \beta_{1}x_{1} - \dots - \beta_{k}x_{k} \end{align} \]

Repare que a única diferença na equação acima a única difernça entre os diferentes \(j\)’s está no intercepto \((\alpha - c_{j})\), portanto as curvas em um gráfico desta função serão paralelas para cada \(j\). Isso quer dizer que os logits, ou as chances, são proporcionais. Por isso alguns chama esse modelo de modelo de chance proporcinal ou, em inglês, proportional-odds model. Ao utilizar um programa para estima-lo ele retornará o coeficiente de \(x\) e os pontos de transição \(c_{(.)}\).

[Exemplo Fox e Weisberg 2]

Alternativamente podemos estimar esse modelo com logits separados por categoria, se forem poucas, ou por MQO, se forem muitas. O problema do MQO é que ele considera os intervalos entre as categorias como sendo do mesmo tamanho, o que, por se tratar de variáveis qualitativas nem sempre é verdade.

Outros tipos de VD

Algumas vezes parte de nossos dados é qualitativa e parte quantitativa. Algumas vezes registramos apenas valores críticos, como a quantidade de votos apenas de candidatos eleitos. Para os que não foram eleitos temos apenas essa informção: “Não eleito”. Um modelo para este tipo de distribuição é o Tobit que primeiro estima a probabildade de ser eleito e então combina com a distribuição de votos.

Como podem ver as possibilidades são inúmeras. Além desses temos modelos de duração, de tempo de espera etc. Além dos ínumeros modelos não paramétricos que examinaremos quando estudarmos aprendizado de máquina.

  1. Sigo aqui o livro online Introduction to Probability, Statistics, and Random Processes de Hossein Pishro-Nik.↩

  2. Os detalhes da derivação da Poisson a partir da binomial pode ser encontrada na página de Andrew Chamberlain (https://medium.com/@andrew.chamberlain/deriving-the-poisson-distribution-from-the-binomial-distribution-840cc1668239)↩