post sobre estimação de médias e método dos momentos generalizado

.Rproj.user/shared/notebooks/paths

@@ -1,3 +1,7 @@
-C:/Users/pedro/azul/content/post/2022-06-19-note-taking/index.pt-br.Rmd="BAED8761"
-C:/Users/pedro/azul/content/post/first-class-deliverables/post.Rmd="52A11494"
-C:/Users/pedro/azul/public/post/2021-06-23-Categorias-espaco-institucional/post.Rmd="C42B8134"
+/Users/danielcoutinho/Documents/GitHub/azul/.Rprofile="AEDAA334"
+/Users/danielcoutinho/Documents/GitHub/azul/config.toml="7AFA2F44"
+/Users/danielcoutinho/Documents/GitHub/azul/content/post/2021-10-10-estimador-james-stein-e-admissibilidade/index.pt-br.Rmd="6BFF49B4"
+/Users/danielcoutinho/Documents/GitHub/azul/content/post/2022-01-06-tratamento-sequencial/codigo.R="DF5348D9"
+/Users/danielcoutinho/Documents/GitHub/azul/content/post/2022-01-09-mqo-mqo-mqo/index.pt-br.Rmd="6D8749EA"
+/Users/danielcoutinho/Documents/GitHub/azul/content/post/2023-05-18-estima-o-da-m-dia-m-todo-generalizado-dos-momentos-e-cramer-rao/index.pt-br.Rmd="64CE4BCA"
+/Users/danielcoutinho/Documents/GitHub/azul/content/post/first-class-deliverables/post.Rmd="EE7DC204"

config.toml

@@ -1,6 +1,6 @@
 baseURL = "https://azul.netlify.app/"
 title = "AZUL"
-copyright = "Copyright © 2008–2020, Pedro Cavalcante & Daniel Coutinho; all rights reserved."
+copyright = "Copyright © 2008–202?, Pedro Cavalcante & Daniel Coutinho; all rights reserved."
 paginate = 2
 languageCode = "pt-br"
 DefaultContentLanguage = "pt-br"

content/post/2023-05-18-estima-o-da-m-dia-m-todo-generalizado-dos-momentos-e-cramer-rao/index.pt-br.Rmd

@@ -0,0 +1,166 @@
+---
+title: Estimação da Média, Método Generalizado dos Momentos e Cramer Rao
+author: Daniel Coutinho
+date: '2023-05-18'
+slug: estimacao-da-media-metodo-dos-momentos
+categories:
+  - Estatística
+tags:
+  - Cramer Rao
+  - Método Generalizado dos Momentos
+images: []
+authors: ["danielc"]
+output:
+  blogdown::html_page:
+    pandoc_args: 
+      [
+      "--lua-filter=../script_number_and_braces.lua"
+      ]
+---
+
+Faz muito tempo que eu não escrevo aqui. O que o primeiro ano do doutorado não faz! Eu trago um post que apresenta um problema que é extremamente simples e com uma solução não óbvia - pelo menos para mim. Eu vou apresentar o problema, discutir o método generalizado dos momentos (_Generalized Method of Moments_) e resolver o problema. Na sequência, eu vou discutir uma extensão do problema, e usar o limite inferior de Crámer Rao. 
+
+O modelo é bem simples: nós temos duas variáveis aleatórias, $X$ e $Y$. $Y$ tem média 0 e $X$ tem uma média desconhecida, $\mu$. Nós também sabemos que $X$ e $Y$ tem variâncias finitas, $\sigma_x^2$ e $\sigma_y^2$, respectivamente; e covariância $\sigma_{xy}$. Nós queremos estimar a média de $X$, $\mu$. 
+
+O estimador mais natural é a média amostral, $\bar{X} = \frac{1}{n} \sum_{i=1}^n X_i$, que tem variância $\sigma_x^2/n$. A pergunta aqui é: existe um estimador melhor que esse? A intuição é simples, e assuma que a covariância entre $X$ e $Y$ é positiva: sabemos que a média de $Y$ é zero. Suponha que nós observamos uma amostra e calculamos a média amostral de $Y$. Se a média é maior que zero, nós podemos usar a correlação positiva para afirmar que a nossa média amostral de $X$ é maior do que o valor esperado de $X$ e obter um estimador que é mais preciso que a média amostral de $X$ usando a informação obtida com a média amostral de $Y$. 
+
+Vamos formalizar isso usando o método generalizados dos momentos.
+
+# O Método Generalizado dos Momentos
+
+Vamos começar com o que a gente já (?) sabe: o método dos momentos é uma técnica de estimação que nos diz para pegar qualquer momento teórico e substituir ele por um momento amostral. Você quer estimar a média, $E[X]$? O método dos momentos nos diz para usar a média amostral, $\frac{1}{n} \sum_{i=1}^n X_i$. E a variância, $E[X^2] - E[X]^2$? Bom, use $\frac{1}{n} \sum_{i=1}^n X_i^2 - \bar{X}^2$ - lembre que $\bar{X} = \frac{1}{n}\sum_{i=1}^n X_i$. 
+
+Para motivar o método generalizado dos momentos, vamos considerar que $X$ tem uma distribuição Poisson, que tem um único parâmetro $\lambda$. A coisa legal da Poisson é que a média _e_ a variância são iguais a $\lambda$. Você pode se sentir tentado a usar a média ou a variância para estimar o parâmetro. Mas usar um ou outro parece jogar informação fora e talvez uma combinação dos dois seja o ideal? 
+
+O método generalizado dos momentos faz exatamente isso. Para combinar os momentos, nós minimizamos uma forma quadrática. Me mantendo no exemplo da Poisson, nós poderíamos fazer:
+
+$$
+\min_{\lambda} \left(\frac{1}{n} \sum_{i=1}^n X_i - \lambda\right)^2 + \left(\frac{1}{n} \sum_{i=1}^n (X_i - \lambda)^2 - \lambda\right)^2
+$$
+Eu posso ser mais geral que isso e colocar alguns pesos! Eu posso colocar até pesos que fazem o termo entre os dois (ou mais) momentos! Eu vou escrever isso em forma matricial: seja $g(X,\theta)$ a nossa função dos momentos e $W$ uma matriz de pesos (ela tem que ser positiva definida - isso vai garantir que a nossa função quadrática tem um mínimo e não um máximo). O que o método generalizado dos momentos faz é:
+
+$$
+\min_{\theta} g(X,\theta)^T W g(X,\theta)
+$$
+O exemplo com os momentos da Poisson usa $W = I$ e:
+
+$$
+g(X,\theta) = \begin{bmatrix}
+\frac{1}{n} \sum_{i=1}^n X_i - \lambda\\
+\frac{1}{n} \sum_{i=1}^n (X_i - \lambda)^2 - \lambda
+\end{bmatrix}
+$$
+
+O método generalizado dos momentos é muito próximo de variáveis instrumentais: variáveis instrumentais impõe uma condição de momento, $E[Zu] = 0$. Se nós temos mais variáveis instrumentais do que regressores, nós temos sobreidentificação. A solução usual é mínimos quadrados em dois estágios, que é o método generalizado dos momentos $G = Z^T(Y - X\beta)$ e $W = (Z^T Z)^{-1}$. 
+
+# De volta ao problema
+
+Nosso problema é achar um estimador para a média de $X$ quando nós temos uma amostra iid $\{X_i,Y_i\}_{i=1}^n$ tal que a média de $Y_i$ é 0 e nós variâncias finitas e covariância entre $X$ e $Y$. 
+
+Nós vamos usar o método generalizado dos momentos. Os momentos que vamos usar são as médias amostrais, então nós temos:
+
+$$
+\min_{\mu} \begin{bmatrix} \bar{X} - \mu & \bar{Y} \end{bmatrix} W \begin{bmatrix} \bar{X} - \mu \\ \bar{Y} \end{bmatrix} \tag{1}
+$$
+
+E os pesos? Se eu usar $W = I$, o problema é resolvido com $\hat{\mu} = \bar{X}$. A matriz "ótima" de pesos - significando que minimiza a variância do estimador - é a inversa da variância dos momentos. O bom desse exemplo é que a variância dos momentos é fácil de calcular: a variância da média amostral é só a variância da variável aleatória dividida pelo tamanho da amostra. Você obtém a covarância usando a mesma ideia (e o fato de que a amostra é iid):
+
+$$
+W = \begin{bmatrix} \frac{\sigma_x^2}{n} & \frac{\sigma_{xy}}{n}\\
+\frac{\sigma_{xy}}{n} & \frac{\sigma_y^2}{n} \end{bmatrix}^{-1}
+$$
+
+Matriz $2 \times 2$ [é inversível na mão](https://pt.wikipedia.org/wiki/Matriz_inversa#Invers%C3%A3o_de_matrizes_2%C3%972) e é fácil o suficiente que nem eu sou capaz de errar:
+
+$$
+W = \frac{n}{\sigma_x^2 \sigma_y^2 - \sigma_{xy}^2}
+\begin{bmatrix}
+\frac{\sigma_{y}^2}{n} & -\frac{\sigma_{xy}}{n} \\
+- \frac{\sigma_{xy}}{n} & \frac{\sigma_x^2}{n}
+\end{bmatrix} \tag{2}
+$$
+
+O determinante na divisão vai multiplicar a expressão toda, e como ele é positivo, não faz nenhuma diferença pro problema em mãos. Nós iremos minimizar:
+
+$$
+\min_{\mu} \sigma_y^2(\bar{X} - \mu)^2 - 2\sigma_{xy} (\bar{X} - \mu)\bar{Y} + \sigma_x^2 \bar{Y}^2
+$$
+
+Isso é obtido usando a equação 2 na equação 1 e fazendo todas as multiplicações de matriz. Agora, resta a única coisa que um economista precisa saber: derivar e igualar a zero!:
+
+$$
+-2\sigma_y^2 (\bar{X} - \mu) + 2\sigma_{xy} \bar{Y} = 0 \therefore \sigma_y^2 (\bar{X} - \mu) - \sigma_{xy} \bar{Y}\\
+\hat{\mu} = \bar{X} - \frac{\sigma_{xy}}{\sigma_y^2} \bar{Y}
+$$
+
+Verificando a nossa intuição lá em cima: suponha que a covariância entre $X$ e $Y$ é positiva. Se a média amostral de $Y$ for maior que zero, então a média amostral de $X$ vai ser acima do valor esperado $\mu$. Sem muito esforço, a gente pode verificar que este estimador tem uma variância mais baixa que a média amostral:
+
+$$
+\text{Var}(\hat{\mu}) = \text{Var}\left(\bar{X} - \frac{\sigma_{xy}}{\sigma_y^2} \bar{Y}\right) = \text{Var}(\bar{X}) + \frac{\sigma_{xy}^2}{(\sigma_y^2)^2} \text{Var}(\bar{Y}) - 2\frac{\sigma_{xy}}{\sigma_y^2} \text{Cov}(\bar{X},\bar{Y}) = \\
+= \frac{\sigma_{x}^2}{n} + \frac{\sigma_{xy}^2}{(\sigma_y^2)^{\cancel{2}}} \frac{\cancel{\sigma_y^2}}{n} - 2\frac{\sigma_{xy}}{\sigma_y^2}\sigma_{xy} = \frac{\sigma_x^2}{n} - \frac{\sigma_{xy}^2}{n\sigma_y^2}
+$$
+
+É claro que $\frac{\sigma_{xy}^2}{\sigma_y^2}$ é positivo, então a gente está subtraindo alguma coisa da variância da média amostral de $X$ e o nosso estimador tem uma variância menor que a variância da média amostral.
+
+A estrutura que eu impus no modelo é mínima. O fascinante aqui é que ter uma variável com média conhecida que é correlacionada com a variável cuja média nós desconhecemos permite obter uma estimativa melhor do que a média amostral. 
+
+Esse estimador é não viesado:
+
+$$
+E[\hat{\mu}] = E[\bar{X} - \sigma_{xy}/\sigma_y^2 \bar{Y}] = \mu - \frac{\sigma_{xy}}{\sigma_y^2} 0 = \mu
+$$
+
+# Limite de Crámer Rao
+
+Se você teve um excelente curso de estatística, você já ouviu falar do limite de Crámer Rao. A ideia é bem simples: entre todos os estimadores não viesados, qual alcança a menor variância? O limite de Cramer Rao nos diz que, se $\hat{\theta}$ é um estimador não viesado de $\theta$ (para uma amostra i.i.d.) e $f(X|\theta)$ é a densidade de $X$:
+
+$$
+\text{Var}(\hat{\theta}) \geq \left(nE\left[\left(\frac{\partial \log(f(X|\theta))}{\partial \theta}\right)^2\right]\right)^{-1}
+$$
+
+Por exemplo, pro caso da normal univariada com média $\mu$ e $\sigma^2$, nenhum estimador (não viesado) consegue obter uma variância menor que $\frac{\sigma^2}{n}$. Isso abarca o exemplo acima, mas o univariado está fazendo todo o trabalho aqui. 
+
+Eu não impus nenhuma distribuição no exemplo acima: eu só impus uma média finita e segundo momento finito. Será que existe algum estimador melhor do que o do método generalizado dos momentos no caso bivariado? Crámer Rao pode nos dar uma resposta, mas como depende da distribuição, a gente vai ter que impor uma distribuição. A mais conveniente é a normal. A normal bivariada tem densidade:
+
+$$
+f(X,Y|\theta) = \frac{1}{\sqrt{2\pi \sigma_x^2 \sigma_y^2 (1 - \rho^2)}} \exp \left(-\frac{(X - \mu)^2}{2\sigma_x^2 (1 - \rho^2)} + \frac{\rho(X-\mu)Y}{\sigma_x \sigma_y (1 - \rho^2)} - \frac{Y^2}{2 \sigma_y^2(1-\rho^2)}\right)
+$$
+
+Eu to usando $\theta = (\mu,\sigma_{x}^2, \sigma_y^2, \rho)$ e $\rho$ é a correlação entre $X$ e $Y$, logo $\rho = \frac{\sigma_{xy}}{\sigma_x \sigma_y}$. Vamos tirar o log:
+
+$$
+\log(f(X,Y|\theta)) = -\frac{1}{2}\log(2\pi \sigma_x^2 \sigma_y^2 (1 - \rho^2)) -\frac{(X - \mu)^2}{2\sigma_x^2 (1 - \rho^2)} + \frac{\rho(X-\mu)Y}{\sigma_x \sigma_y (1 - \rho^2)} - \frac{Y^2}{2 \sigma_y^2(1-\rho^2)}
+$$
+
+Vamos tirar a primeira derivada:
+
+$$
+\frac{\partial \log(f(X,Y|\theta))}{\partial \mu} = \frac{X - \mu}{\sigma_x^2 (1 - \rho^2)} - \frac{\rho Y}{\sigma_x \sigma_y (1 - \rho^2)}
+$$
+
+Eu poderia calcular o quadrado disso e obter uma expressão que depende de $(X-\mu)^2$, $Y^2$ e $(X - \mu)Y$. Honestamente, isso dá muito trabalho! Eu digo, sem fornecer nenhuma prova, que:
+
+$$
+E\left[\left(\frac{\partial \log(f(X|\theta))}{\partial \theta}\right)^2\right] = - E\left[\frac{\partial^2 \log(f(X|\theta))}{\partial \theta \partial \theta^T}\right]
+$$
+
+Obtendo a segunda derivada:
+
+$$
+\frac{\partial^2 \log(f(X,Y|\theta))}{\partial \mu^2} = -\frac{1}{\sigma_x^2(1-\rho^2)}
+$$
+
+Então, o limite inferior de Cramer Rao é:
+
+$$
+\left(-nE\left[\frac{\partial^2 f(X|\theta)}{\partial \mu^2}\right]\right)^{-1} = \frac{\sigma_x^2 (1 - \rho^2)}{n}
+$$
+
+Agora, use a definição de $\rho = \frac{\sigma_{xy}}{\sigma_x \sigma_y}$, para obter:
+
+$$
+\left(-nE\left[\frac{\partial^2 f(X|\theta)}{\partial \mu^2}\right]\right)^{-1} = \frac{\sigma_x^2 (1 - \rho^2)}{n} = \frac{\sigma_x^2}{n} - \frac{1}{n} \frac{\sigma_x^2 \sigma_{xy}^2}{\sigma_x^2 \sigma_y^2} = \\
+= \frac{\sigma_x^2}{n} - \frac{\sigma_{xy}^2}{n\sigma_y^2}
+$$
+
+Essa é exatamente a variância obtida usando o método dos momentos. O nosso estimador obtido pelo método generalizado dos momentos é o mais eficiente se a distribuição conjunta das variáveis é normal.