How to estimate spatial models with the use of the likelihood function André Braz Golgher1 Abstract This paper is part of a series that discusses introductory concepts of spatial econometrics 2. The texts were written in Portuguese and intend to present this field of study to students at upper undergraduate to graduate levels in Economy and in Regional Sciences. The main objective of this third text is to discuss how to apply the likelihood function to estimate different spatial models. First, I present this function in a general perspective and obtain the maximum likelihood estimator for the normal distribution. Then, I show some of the proprieties of this estimator that makes it very popular as an estimation method, in particular for spatial models. Afterwards, I develop the mathematical expressions of the maximum likelihood estimators for the OLS, spatial error and spatial lag models. Next, I present some illustrative simulations using Matlab that address the concepts discussed in the text. Finally, I obtain the covariance matrix for the spatial lag model as an example for the other spatial models. Key words: spatial models; likelihood function; maximum likelihood estimators. 1 Associate Professor at the Economics Department at the Cedeplar/FACE/UFMG, visiting scholar at the Regional Research Institute (RRI) at the West Virginia University (WVU) and visiting scholar at the Carolina Population Center (CPC) at the University of North Carolina (UNC) in Chapel Hill. 2 A large part of the references selected for these texts comes from the readings assigned by professor Donald Lacombe in his course ARE 693L Spatial Econometrics (spring 2012) at the WVU (http://community.wvu.edu/~djl041/teaching.html). 1 Estimando modelos espaciais com a função de verossimilhança André Braz Golgher3 Resumo Esse texto faz parte de uma série que apresenta pontos introdutórios da econometria espacial. O objetivo principal deste terceiro texto é discutir como utilizar a função de verossimilhança para estimar diferentes modelos espaciais. Inicialmente, apresenta-se essa função e discute-se como obter o estimador de máxima verossimilhança utilizando a distribuição normal como exemplo. Em seguida, apresentam-se algumas das propriedades do estimador, que o tornam muito popular em diversos tipos de estimação, em particular na estimação dos modelos espaciais. Depois, desenvolvem-se as expressões matemáticas que são utilizadas na obtenção dos estimadores de máxima verossimilhança para os modelos de mínimos quadrados ordinários, de lag espacial e de erro espacial. Posteriormente, apresentam-se simulações feitas no Matlab que ilustram a aplicação dos conceitos discutidos no texto. Por fim, apresenta-se a matriz de covariância para o modelo de lag espacial como exemplo para os demais modelos espaciais. Palavras chave: modelos espaciais; função de verossimilhança; estimador de máxima verossimilhança. 3 Professor do Cedeplar/FACE/UFMG, pesquisador visitante do Regional Research Institute (RRI) da West Virginia University (WVU) e pesquisador visitante do Carolina Population Center (CPC) na University of North Carolina (UNC) em Chapel Hill. 2 1 - Introdução Esse texto faz parte de uma série que apresenta a econometria espacial em pontos introdutórios4. No primeiro texto, “Introdução à Econometria Espacial”, foram discutidos alguns conceitos introdutórios sobre a econometria espacial, onde foram apresentados alguns dos modelos espaciais, incluindo motivações teóricas para o uso desses modelos. No segundo texto da série, “Interpretando os coeficientes dos modelos espaciais”, discutiu-se como interpretar os coeficientes obtidos nos diversos modelos, e também como calcular os efeitos diretos, indiretos e totais de diversas ordens. Neste terceiro texto, discute-se como utilizar a função de verossimilhança para estimar diferentes modelos espaciais. Inicialmente, apresenta-se a função de verossimilhança e discute-se o estimador de máxima verossimilhança, utilizando a distribuição normal como exemplo. Em seguida, apresentam-se algumas das propriedades do estimador, que o tornam muito popular em diversos tipos de estimação. Depois, desenvolve-se passo-a-passo as expressões que são utilizadas na obtenção dos estimadores de máxima verossimilhança de diferentes modelos espaciais. Uma vez discutidos os conceitos teóricos envolvidos na estimação dos modelos espaciais com estimadores de máxima verossimilhança, apresenta-se simulações feitas no Matlab que exemplificam a aplicação desses conceitos. Por fim, apresenta-se a matriz de covariância para o modelo de lag espacial como exemplo para os demais modelos espaciais. Segundo Elhorst (2010), são três os principais métodos desenvolvidos para a estimação de modelos espaciais: método da máxima verossimilhança (ML), variáveis instrumentais ou o método de momentos generalizados (IV/GMM), e a abordagem Bayesiana com Monte Carlo via cadeias de Markov (MCMC). Segundo Anselin (1988), o uso da função de máxima verossimilhança era a abordagem mais familiar para a estimação e teste de hipóteses de modelos espaciais na época de publicação de seu influente livro. Greene (2003) afirmou que os modelos Bayesianos estavam se tornando muito populares em aplicações da econometria, mas que o método da máxima verossimilhança permanecia sendo o preferido na maioria dos campos de estudo. Neste texto discutimos o método 4 Grande parte do material citado aqui foi selecionado da ementa do curso ARE 693L Spatial Econometrics (spring 2012) (http://community.wvu.edu/~djl041/teaching.html) ministrado por Donald Lacombe do RRI da WVU. 3 mais usado que é o de máxima verossimilhança. Para uma discussão sobre o método Bayesiano ver LeSage e Pace (2009). Uma das razões principais do desenvolvimento do método com variáveis instrumentais ou o método de momentos generalizados (IV/GMM) foi para facilitar os cálculos computacionais na estimação dos modelos, que podem ser demandantes devido a matrizes de ordem n x n presentes no processo de obtenção dos estimadores. Segundo Elhorst (2010), LeSage e Pace (2009) mostram evidências conclusivas que essas limitações computacionais são coisas do passado. Para uma discussão sobre esse método ver a extensiva série de publicações de Kelejian e Prucha (Kelejian, 2011) . Este texto foi dividido em dez seções, incluindo essa introdução. Na seção 2 é apresentado o estimador de máxima verossimilhança, utilizando como exemplo a distribuição normal. Em seguida, discutem-se algumas das propriedades desse estimador. A seção 4 apresenta os modelos espaciais e discute algumas das similaridades entre eles que facilitam a estimação dos modelos. Nas seções 5, 6 e 7 são desenvolvidos respectivamente os estimadores de máxima verossimilhança para o MQO, para o modelo de lag espacial e para o modelo de erro espacial. A seção 8 mostra simulações ilustrativas no Matlab sobre os conceitos discutidos nas seções anteriores. Na seção seguinte, obtém-se a matriz de covariância do modelo de lag espacial, onde, como veremos, os cálculos são bastante extensos. A última seção conclui o texto. Além dessas seções são incluídos três apêndices econométricos no fim do texto, como forma de tornar a discussão dos conceitos do corpo do texto mais completa. 2-O estimador de máxima verossimilhança Nesta seção apresentamos alguns pontos relacionados com a função de verossimilhança e com a utilização desta para a obtenção de estimadores. Inicialmente, apresenta-se a função de densidade de probabilidade (pdf), utilizando a distribuição normal como exemplo. Em seguida, obtém-se a função de verossimilhança e derivam-se os estimadores de máxima verossimilhança para essa distribuição. O objetivo é dar ao leitor uma ideia inicial desse método de estimação, que será posteriormente utilizada na obtenção das funções de verossimilhança dos modelos espaciais. 4 As funções de densidades de probabilidade, f (x) , de distribuições contínuas, como a normal, b têm as seguintes propriedades: f ( x) 0 , Pr ob(a x b) f ( x)dx, e 0 Pr ob(a x b) 1 a (Greene, 2003). No caso particular da distribuição normal com média e variância 2 , N ( , 2 ) , essa função tem o seguinte formato: f ( y; , ) 1 2 (2 2 ) 1 e 1 2 2 y 2 . 2 Se tivermos duas observações independentes, y1 e y 2 , obtemos a densidade conjunta destas multiplicando as respectivas expressões: 1 1 2 y1 2 2 y2 2 1 1 . 2 2 F ( y1 , y2 ; , ) e e (2 2 ) 12 (2 2 ) 12 2 Note que para distribuição normal temos dois parâmetros, e 2 . Para simplificar e generalizar a notação, definimos ( , 2 ) , como o vetor composto por esses parâmetros e y ( y1 , y2 ) como o vetor com as observações. Reescrevemos a expressão acima como: 1 2 yi 2 1 . 2 F ( y; ) e 1 2 i 1 ( 2 ) 2 2 Seguindo esse raciocínio, a densidade conjunta de n observações independentes e identicamente distribuídas, que no caso acima é a distribuição normal, é dada pelo produto das densidades individuais, e é denominada função de verossimilhança: n L( y ) f ( yi ). i 1 Ao contrário da expressão da pdf, f ( yi ), onde temos os dados condicionados aos parâmetros ( yi ) , note que na função de verossimilhança, L( y), temos ( y ) , ou seja, os parâmetros são condicionados aos dados. Assim, o objetivo do uso da função de verossimilhança (FV) eh, a 5 partir dos dados empíricos, estimar os parâmetros da função. Estima-se os parâmetros que maximizam essa função, isto eh, maximizam a probabilidade de ocorrência dos dados empíricos. Seguindo a discussão com a distribuição normal, a FV para n observações é escrita como: n L( , y ) i 1 1 2 (2 ) 2 1 e 1 2 2 yi 2 n 1 (2 ) 2 2 n e 1 2 2 yi 2 2 i 1 1 (2 ) 2 n e 1 2 2 n yi 2 i 1 2 Para acharmos o máximo dessa função, podemos derivar parcialmente essa expressão com relação aos dois parâmetros, e 2 , igualando cada uma das derivadas a zero. Entretanto, em geral, é mais simples trabalhar com o logaritmo dessa função. Note que o máximo da FV é o mesmo ponto que o máximo do logaritmo dessa função, pois o logaritmo é uma função monotonicamente crescente. Vejamos. 1 2 yi 2 n n 1 2 i 1 2 2 ln (2 2 ) 2 1 yi 2 ln L( , y ) l ( , y ) ln e n 2 2 i 1 (2 2 ) 2 1 n n n yi 2 . ln 2 ln 2 2 2 2 2 i 1 n Uma vez obtido o logaritmo da FV, derivamos a expressão, inicialmente com relação a , obtemos o estimador de máxima verossimilhança (MV) para esse parâmetro: l ( , 2 y ) n n n 1 n 2 yi 2 12 yi 0 2 ln 2 2 ln 2 2 i 1 i 1 yi 0 n i 1 ny n 0 MV y O estimador de MV para MV é a média entre as observações. Em seguida, derivamos com relação à 2 , e obtemos o estimador para a variância: 6 l ( , 2 y ) 2 n n 1 n 2 2 ln 2 ln 2 yi 2 2 2 2 i 1 n 1 2 2 2 2 n n yi 2 1 n 1 yi y 2 0 2 2 2 2 i 1 i 1 n 2 MV 1 n 2 yi y 0 2 i 1 1 n 2 yi y n i 1 Como discutido, a FV é muito utilizada na obtenção de estimadores em diversos campos de aplicação e, em particular, também são aplicados aos modelos espaciais. Esse fato é devido as interessantes propriedades desses estimadores (ver Greene 2003 para uma discussão mais detalhada). A próxima seção discute algumas delas. 3 – Propriedades do estimador de MV Os estimadores de MV apresentam algumas características que os tornam muito populares em diversos tipos de aplicações, inclusive na estimação dos modelos espaciais. Os estimadores de MV são consistentes e assintoticamente eficientes, com matriz de covariância assintótica 1 2 ln L . Além disso, os definida no limite inferior de Cramér-Rao: Var (ˆ) [ I ( )]1 E ' estimadores são assintoticamente normalmente distribuídos: ˆ ~ a ~ N[ , I ( ) 1 ] (ver apêndice econométrico 1 para maiores detalhes sobre esses conceitos). Para que esses estimadores tenham essas propriedades, existem certas condições formais que devem ser satisfeitas. Para o caso particular dos modelos espaciais, essas foram descritas em estudos realizados na década de 80 (Anselin, 1988). Dentre outras, essas condições formais são: o logaritmo da FV deve existir para os valores dos parâmetros sendo considerados; essa função deve ser continuamente derivável até a segunda ou a terceira ordens; e as derivadas parciais devem assumir valores finitos. 7 Como vimos na seção anterior quando foram obtidos os estimadores de MV da distribuição normal, todas essas condições são satisfeitas para essa distribuição. Uma vez apresentadas as propriedades desses estimadores, continua-se esse exemplo, obtendo a matriz de covariância dos estimadores, que deve ser positivamente definida. Obtivemos os estimadores na seção anterior com o uso das derivadas de primeira ordem: l 1 n yi 2 i 1 1 n l 1 2 n y y i 2 i 1 2 2 2 As covariâncias dos estimadores da distribuição normal são obtidas a partir do limite de CramérRao, com a seguinte matriz: 2 ln L 2 2 ln L E 2 I ( ) E ln L ' 2 2 ln L 2 2 ln L ( 2 ) 2 Para determinar matriz de covariância, calculamos as derivadas de segunda ordem: 2l 1 n n 2 yi 2 . 2 i 1 1 n 2l 2l l 1 n 2 2 yi . y i 2 2 2 2 2 i 1 ( ) i 1 2l 2 2 ( ) 2 1 1 n 2 2 n 2 yi y 2 i 1 1 n 1 1 1 n 2 2 n y y y y 2 2 2 i i 2 2 2 i 1 2( ) 2 i 1 n 1 n 2 y y . i 3 2 2 2( ) 2 i 1 8 Em seguida, obtemos os valores esperados dessas expressões, lembrando que E yi e que E ( yi ) 2 n 2 : 2l n n E 2 E 2 2 1 n 2l 1 n E E 2 2 yi 4 E yi 0 2 i 1 ( ) i 1 n 1 n n 2l yi y 2 n 1 E yi y 2 E E 3 3 2 2 2 2 2 2 2( ) 2 i 1 2( ) 2 i 1 ( ) n 1 n (n 2 ) 3 2 2 2 2 2( ) 2 2( ) Obtemos assim a matriz de covariância com o limite de inferior de Cramér-Rao: 2 ln L 2 [ I ( )] E 2 ln L 2 2 n 1 [ I ( )] 0 2 ln L n 2 0 2 2 ln L n 4 0 ( 2 ) 2 2 4 2 n 0 Além dessas propriedades discutidas, outra propriedade do estimador de MV é que ele é invariante a uma mudança por uma função injetiva e continuamente derivável. Ou seja, se temos um estimador ˆ para um parâmetro e temos uma função injetiva e continuamente derivável f (x), como uma reta f ( x) ax b, então o estimador de f ( ) será f (ˆ) . Essa propriedade será utilizada na obtenção dos estimadores para os modelos econométricos. Uma vez obtidos os estimadores de MV para a distribuição normal, utilizaremos esses mesmos conceitos para discutir os estimados de MV dos modelos espaciais. A próxima seção apresenta 9 esses modelos e discute algumas regularidades apresentadas por eles, que simplificarão a discussão subsequente. 4-Modelos espaciais Antes de discutirmos como estimar cada um dos modelos espaciais, eles são analisados em pares, pois isso vai simplificar em muito a discussão futura. A figura 1, similar a figure 1 em Elhorst (2010), mostra a relação entre diversos modelos espaciais. Como mostra essa figura, o 7 - modelo de lag de X (SLX) é descrito pela seguinte equação: Y X WX , ~ N (0, 2 I n ). Note que podemos reescrever esse modelo com Z [ X ,WX ] e [ , ] da seguinte forma: y Z , ~ N (0, 2 I n ). Essa expressão é exatamente igual ao 8 - MQO. Assim, os modelos são estimados de forma similar depois de feitas essas mudanças no modelo SLX e, assim, não abordaremos este ultimo em separado o modelo SLX. De forma similar, o 4 - modelo de erro espacial de Durbin (SDEM) é descrito pela seguinte equação: Y X WX u, u Wu , ~ N (0, 2 I n ) Fazendo a mesma transformação acima, obtemos a equação do 6 - modelo de erro espacial. Assim, também não abordaremos o modelo SDEM em separado. Por sua vez, o 2 - modelo espacial de Durbin (SDM) pode ser escrito como o 5 - modelo de lag espacial com essas mesmas transformações (LeSage e Pace, 2009) e esse primeiro também não é discutido em separado. Fazendo uso dessa mesma transformação, 8 - modelo de Manski pode ser escrito como o 3 - modelo de Kelejian-Prucha. Assim, esse primeiro também não é discutido. 10 Em todos os modelos mostrados na figura 1 assumimos que ~ N (0, 2 I n ) . Assim, utilizamos a distribuição normal como base para toda a discussão subsequente. Seguimos a discussão de como estimar os modelos fazendo uso da FV inicialmente com o MQO, e em seguida com o modelo de lag espacial e o modelo de erro espacial. O modelo de Kelejian-Prucha engloba esses dois últimos modelos. Assim, pode-se estender a discussão apresentada sobre a metodologia de estimação dos modelos de lag espacial e de erro espacial na obtenção das expressões dos estimadores de MV desse primeiro. Assim, para não estender em demasiado este texto, esse modelo não é discutido em separado. 11 Figura 1 – Modelos espaciais 1 - Modelo de Manski Se Se Se 0 0 0 4 - Modelo erro espacial de Durbin (SDEM) 3 - Modelo de Kelejian-Prucha (SAC) 2 - Modelo espacial de Durbin (SDM) 0 Se 0 Se Se 0 Se Se 0 0 7 – Modelo de lag de X (SLX) 0 6 – Modelo de erro espacial (SEM) 5 – Modelo de lag espacial (SAR) Se Se Se 0 0 0 8 - MQO 12 5-Estimadores de MV para o MQO Essa seção apresenta os estimadores de MV para o MQO, que, apesar de não ser um modelo espacial, serve de base de comparação para os modelos espaciais, em particular para os modelos de lag espacial e de erro espacial. Além disso, como vimos na seção anterior, os estimadores do modelo de lag de X (SLX) são similares aos do modelo de MQO. O apêndice econométrico 3 apresenta algumas das hipóteses do MQO, além de mostrar o processo de obtenção do estimador de mínimos quadrados ordinários. Essa apresentação é incluída neste texto para que seja feita uma comparação entre este estimador e o estimador de MV. O MQO tem como equação a seguinte expressão: yi X i i , i ~ N (0, 2 ) Note que o erro tem distribuição normal, com média zero e variância 2 . Assim, a função de densidade de probabilidade para o erro tem o seguinte formato: f 1 (2 2 ) 1 e 1 x 2 2 2 2 1 (2 2 ) 1 e 1 i 0 2 2 2 2 1 (2 2 ) 1 e 1 i 2 2 2 2 A função de verossimilhança é escrita como: L 1 (2 2 ) n e 2 n 1 2 2 i 2 i 1 1 (2 2 ) n e 1 ' 2 2 2 Entretanto, no modelo MQO desejamos obter o estimador para e 2 a partir das variáveis dependentes e independentes, obtidas empiricamente. Assim, reescrevemos a expressão acima com: i yi X i O Jacobiano da transformação de n variáveis por m funções (ver apêndice econométrico 4) é dado por: 13 x1 y1 J abs ... xn y1 x1 y n ... ... . xn ... y n ... Assim, para o caso particular do MQO, o Jacobiano da transformação de i para yi tem o seguinte valor: J i 1. yi Incorporando essas expressões na FV do erro, obtemos a FV do MQO: L 1e 1 (2 2 ) n 1 2 2 ( y X )'( y X ) 2 1 (2 2 ) n e 1 2 2 ( y X )'( y X ) 2 Tomando o logaritmo da função, utilizando a relação para matrizes ( A B)' A' B' , e a relação para vetores X 'Y Y ' X , temos: n n 1 l ln 2 ln 2 2 ( y X )' ( y X ) 2 2 2 n n 1 ln 2 ln 2 2 y' y ( X )' y y' X ( X )' ( X ) 2 2 2 n n 1 ln 2 ln 2 2 y' y 2( X )' y ( X )' ( X ). 2 2 2 Derivando com relação ao vetor , obtemos o estimador para esses parâmetros: l 1 2 2 X ' y 2 X ' X 0 2 X ' y X ' X ˆ ( X ' X ) 1 X ' y MQO Note que esse estimador é exatamente igual ao obtido para o estimador de mínimos quadrados ordinários, discutido no apêndice econométrico 3. 14 Derivando com relação a 2 , obtemos o estimador para a variância: l n 1 ( y X )' ( y X ) 0 2 2 2 2 2 2( ) 1 n 2 ˆ MQO ( y X )' ( y X ) e' e . n Por fim, obtemos a matriz de covariância utilizando o limite de Cramér-Rao. l 1 1 X'X 2 2 X ' y 2 X ' X 2 2 X ' X 2 ' 2 2 2l X'X E 2 ' l 1 2 X ' y 2 X ' X 12 2 2 X ' y 2 X ' X 2 2 2 ' 2 2( ) 1 1 2l E E 2 2 X ' y X ' X 2 2 E X ' y X ' X 0 2 ( ) ( ) l n 1 ( y X )' ( y X ) 2 2 2 2 2 2 2 2( ) n n 1 1 2 3 ( y X )' ( y X ) 2 3 ' 2 2 2 2 2( ) ( ) 2( ) ( ) n 2l 1 n 1 2 3 ' 2 3 E ' E E 2 2 2 2 2 2 ( ) 2( ) ( ) 2( ) ( ) n n 1 2 3 (n 2 ) 2 2 2 2 2( ) ( ) 2( ) 15 2 ln L 2 1 [ I ( )] E 2 ln L 2 ' 2 ln L 2 2 ln L ( 2 ) 2 1 2 X ' X 0 4 2 n 0 Como discutido na seção anterior, esse procedimento descrito aqui também é utilizado para o modelo SLX. 6 – Estimadores de MV para o modelo de lag espacial Segue a discussão de como obter os estimadores para o modelo de lag espacial e, como vimos, também para o modelo espacial de Durbin. A matriz de covariância será obtida posteriormente na seção nove do texto, devido a extensão dos cálculos. A discussão abaixo é baseada em Doreian (1981). Como vimos, o modelo de lag espacial é expresso pela seguinte equação: Y WY X Seguindo os mesmos procedimentos da seção anterior, isolamos o termo do erro, obtemos o Jacobiano da transformação de i para yi , e substituindo os resultados na FV do erro: y Wy X ( I W ) y X J i I W . yi Assim, a FV do modelo de lag espacial é dada por: ( y Wy X )' ( y Wy X ) L( , , ; y, X ) (2 2 ) n / 2 I W exp 2 2 Tomando o logaritmo da função e simplificando a expressão, temos: 16 ( y Wy X )' ( y Wy X ) l ( , 2 , ; y, X ) ln{(2 2 ) n / 2 I W exp } 2 2 ( y Wy X )' ( y Wy X ) ln( 2 2 ) n / 2 ln I W ln[exp ] 2 2 n 1 ( y Wy X )' ( y Wy X ) ln( 2 ) ln( 2 ) ln I W 2 2 2 Derivando com relação a 2 , obtemos o estimador da variância: 1 l n ( y Wy X )' ( y Wy X ) 0 2 2 2 2 2 2( ) 1 1 2 n 2 ( y Wy X )' ( y Wy X ) 0 2 1 n 2 ˆ SAR ( y Wy X )' ( y Wy X ) Entretanto, não podemos estimar a variância, pois não sabemos os valores de nem de , que estão presentes na expressão. Assim, devemos utilizar alguma expressão para estimarmos esses parâmetros para depois estimarmos ˆ 2 . Assim, derivando o log da verossimilhança com relação a , temos: l 1 1 2 X ' ( y Wy X ) 2 X ' (( I W ) y X ) 0 Para prosseguir o raciocínio, fazemos uma mudança de variável com z ( I W ) y. Daí ficamos com: X ' ( z X ) 0 Manipulando essa expressão, obtemos o estimador do modelo: ˆSAR ( X ' X ) 1 X ' z 17 Note que para estimarmos ̂ SAR devemos estimar primeiro o parâmetro , para depois substituirmos em z. Substituindo a mudança de variável z ( I W ) y e o resultado para ̂ SAR na expressão do estimador da variância, lembrando que ( AB)' B' A' , ficamos com: 1 n 2 ˆ SAR ( z X )' ( z X ) 1 n 2 ˆ SAR ( z X ( X ' X ) 1 X ' z )' ( z X ( X ' X ) 1 X ' z ) 1 ( z ' z ' ( X ( X ' X ) 1 X ' )' )( z X ( X ' X ) 1 X ' z ) n 1 ( z ' z 2 z ' ( X ( X ' X ) 1 X ' ) z z ' X ( X ' X ) 1 X ' X ( X ' X ) 1 X ' z ) n 1 ( z ' z 2 z ' ( X ( X ' X ) 1 X ' ) z ) z ' X ( X ' X ) 1 IX ' z ) n 1 1 z ' ( I ( X ( X ' X ) 1 X ' ) z z ' Mz, n n onde M I ( X ( X ' X ) 1 X '. Note que agora o estimador da variância depende somente de . Ou seja, neste ponto da 2 ( ). Mas discussão os dois estimadores já discutidos são função desse parâmetro: ˆSAR ( ) e ˆ SAR ainda não sabemos o valor de . Assim, devemos estimar o parâmetro para depois estimarmos os demais via essas expressões. Retomamos a expressão do log da FV. n 1 l ( ˆ ( ),ˆ 2 ( ), ; y, X ) ln( 2 ) ln(ˆ 2 ) ln I W ( y Wy Xˆ )' ( y Wy Xˆ ) 2 ˆ 2 2 1 2 Sabendo que ˆ SAR ( y Wy Xˆ )' ( y Wy Xˆ ) , temos: n 18 n ln( 2 ) ln ˆ 2 ) ln I W 2 n 1 ( y Wy Xˆ )' ( y Wy Xˆ ) ˆ ( y Wy X )' ( y Wy Xˆ ) 2 l n 1 ln(2 ) n ln 1 ( y Wy Xˆ )' ( y Wy Xˆ ) ln I W 2 2 n Um ponto chave na obtenção do estimador de MV é o logaritmo do Jacobiano. Para simplificar essa relação devemos reescrevê-lo. Note que podemos calcular o determinante de uma matriz A como o produto de seus autovalores: A i . i Daí, temos: W i , onde i são os autovalores da matriz de peso. i I W ( i ) , i I W (1 i ) . i Tomando o logaritmo, temos: ln I W ln (1 i ) ln(1 i ). i i Note que o logaritmo só é definido para valores positivos. Assumindo que os autovalores são reais, dessa relação, temos para todo i: (1 i ) 0 1 i Se i é positivo: 1 i . Como essa relação é valida para todo i, daí temos que 1 max , onde max é o maior autovalor real positivo. 19 Se i é negativo: 1 min 1 i . Como essa relação também é valida para todo i, dai temos que , onde min é o menor autovalor real negativo. 1 1 , Assim obtemos o intervalo para os valores de , onde a FV é definida: min max . Se a matriz de peso for normalizada temos os seguintes valores: min 1 e max 1 . Esse mesmo raciocínio é valido para o parâmetro de correlação espacial do modelo de erro espacial. Substituindo a relação acima na função de log de FV e substituindo K n 1 ln(2 ) 2 constante, temos: l ˆ ˆ n 1 ln(2 ) n ln ( y Wy X )' ( y Wy X ) ln I W 2 2 n n l K ln ˆ 2 ln(1 i ) 2 i Sabemos que: 1 1 1 n n n 1 1 ( y ' y 'W ' ) M ( y Wy ) ( y ' My 2 y ' MWy 2 y 'W ' MWy) n n ˆ 2 z ' Mz (( I W ) y )' M ( I W ) y ( y Wy)' M ( y Wy) Substituindo essa expressão na função do log da FV, obtemos a seguinte equação: n l K ln ˆ 2 ln(1 i ) 2 i n 1 K ln ( y ' My 2 y ' MWy 2 y 'W ' MWy) ln(1 i ) 2 n i 20 Simplificando essa expressão, suprimindo os valores constantes e reescrevendo, temos: 2 g ( ; y, X ) ln y' My 2 y' MWy 2 y'W ' MWy ln(1 i ) n i Devemos obter o mínimo dessa função para obtermos o valor do parâmetro que maximiza a FV. Em vez de resolvermos essa expressão analiticamente, ela é minimizada por procedimentos numéricos com a utilização de métodos computacionais. Isto é, variam-se os valores de passoa-passo em pequenos degraus em um intervalo pré-definido e obtém-se o mínimo da função. Normalmente se utiliza o intervalo [0,1). Os exemplos de correlação positiva são muito mais numerosos que os de correlação negativa, o que justifica essa escolha. Valores negativos para os parâmetros espaciais são observados em alguns trabalhos empíricos, como na competição entre empresas em um mesmo espaço urbano ou na disputa por espaço geográfico entre espécies vegetais ou animais (Griffth e Arbia, 2010) . Ou seja, na prática os parâmetros e não são incluídos na otimização da função de verossimilhança. Isso facilita em muito os procedimentos de maximização. Para uma discussão mais aprofundada ver LeSage e Pace (2009). Uma vez obtido o parâmetro que minimiza a expressão acima, são estimados os demais parâmetros, ˆ ( ) e ˆ 2 ( ), e o processo de obtenção dos parâmetros se encerra. A seção seguinte discute o modelo de erro espacial de forma similar ao apresentado aqui. 7- Estimadores de MV do modelo de erro espacial Seguindo um raciocínio similar ao descrito acima obtemos os estimadores do modelo de erro espacial (SEM). Como vimos, o modelo tem a seguinte equação: Y X u u Wu , ~ N (0, 2 I n ). 21 Em seguida, isolamos : Y X ( I W ) 1 ( I W )Y ( I W ) X ( I W )(Y X ). A partir dessa expressão, calculamos o Jacobiano da transformação de para Y : J ( I W ). Y De posse dessa informação, obtemos a FV do modelo: [( I W )(Y X )]'[( I W )(Y X )] L( , , ; y, X ) (2 2 ) n / 2 I W exp 2 2 [(Y X )' ( I W )' ][( I W )(Y X )] L( , , ; y, X ) (2 2 ) n / 2 I W exp 2 2 Tomando o logaritmo dessa função, temos: 1 n l ( , , ; y, X ) ln( 2 2 ) ln I W (Y X )' ( I W )' ( I W )(Y X ) 2 2 2 Derivando com relação a , obtemos o estimador ̂ SEM : l 1 (2)( X )' ( I W )' ( I W )(Y X ) 0 2 2 X ' ( I W )' ( I W )(Y X ) 0 [( I W ) X ]' ( I W )(Y X ) 0 [( I W ) X ]' ( I W )Y [( I W ) X ]' ( I W ) X [( I W ) X ]'[( I W )Y ] [( I W ) X ]'[( I W ) X ] Fazemos as seguintes mudanças de variáveis: 22 Y * ( I W )Y X * ( I W )Y Ficamos com: ( X *)' Y * ( X *)' ( X *) ˆ (( X *)' ( X *))1 ( X *)' Y * SEM Note que o estimador acima, assim como observado no modelo de lag espacial, depende do valor da variável de correlação espacial, que no caso do modelo de erro espacial é o valor de , que ainda não conhecemos. Prosseguindo, derivamos com relação à 2 e obtemos o estimador da variância: l n 1 (Y X )' ( I W )' ( I W )(Y X ) 0 2 2 2 2 2 2( ) 1 1 2 n 2 (Y X )' ( I W )' ( I W )(Y X ) 0 2 1 n ˆ SEM (( I W )(Y Xˆ ))' ( I W )(Y Xˆ ) 2 Note que esse estimador também depende do valor de . Ou seja, ambos estimadores, ̂ SEM e ˆ SEM , dependem de . Assim devemos estima-lo e em seguida obtemos os demais. 2 Retornamos a expressão do log da FV. n 1 l (ˆ 2 ( ), ; y, X ) ln( 2ˆ 2 ) ln I W (Y X )' ( I W )' ( I W )(Y X ) 2 2 2ˆ Para simplificar essa relação devemos reescrever o logaritmo do determinante como: ln I W ln(1 i ), onde i são os autovalores da matriz de peso (utiliza-se essa i nomenclatura para autovalores para não confusão com o parâmetro espacial ). 23 Substituindo a expressão do estimador da variância, ficamos com: n 1 l ( ˆ ,ˆ 2 , ; y, X ) ln( 2ˆ 2 ) ln(1 i ) (Y Xˆ )' ( I W )' ( I W )(Y Xˆ ) 2 ˆ 2 2 i n n ln( 2 ) ln ˆ 2 ) ln(1 i ) 2 2 i Simplificando, suprimindo as constantes e reescrevendo, obtemos a função que será minimizada por procedimentos numéricos: 2 g ( ) ln ˆ 2 ln(1 i ) n i Na próxima seção todo o procedimento discutido para os modelos MQO, de lag espacial e de erro espacial é implementado em simulações ilustrativas no Matlab. 8 – Simulação ilustrativa Este seção apresenta simulações ilustrativas feitas no Matlab que trata dos conceitos discutidos nas seções anteriores. Como ponto de partida, temos o DGP do modelo de Kelejian-Prucha: Y WY X u u Wu , ~ N (0, 2 I n ) Nas simulações estipulamos valores arbitrários para os termos desse DGP. Esses valores são baseados em Florax et al (2003) e em Anselin et al (1996). Maiores detalhes das simulações discutidas por esses autores são apresentados no próximo texto dessa serie. As simulações foram feitas com 100 observações, como sugerido em Mur e Angulo (2009) como número mínimo para uma estimativa confiável. As observações foram definidas como regiões, que foram dispostas em uma estrutura quadrada regular, tipo tabuleiro de xadrez, com cada lado de dimensão 10. Assim temos 10 x 10 = 100 observações. Um ponto importante de toda essa análise é referente a relação entre o tamanho da amostra e dificuldades computacionais de estimação. Note que para uma amostra com n observações, a 24 matriz de peso é n x n. Assim, uma amostra maior aumenta a acurácia da estimação, mas implica em maior dificuldade computacional, pois as matrizes a serem manipuladas no processo de estimação tem n2 elementos. Note, porém, que, de forma geral, as matrizes de peso têm muitos elementos nulos, sendo que a proporção destes tende a aumentar quando o tamanho da amostra aumenta. O uso de técnicas de matrizes esparsas, não utilizada aqui, pode ser usado para facilitar a estimação dos parâmetros (Pace, 1997). A matriz de peso foi definida como de contiguidade por torre, pelos lados, como uma torre se movimento no xadrez, como mostra o diagrama abaixo. Note que a contiguidade em uma estrutura regular pode ser definida de outras maneiras, como por exemplo, pelo movimento da rainha no xadrez, que incluem também as regiões que se tocam apenas no vértice. A matriz foi normalizada na linha. Uma vez definida a matriz de peso, o próximo passo é obter a matriz X das variáveis independentes. Assume-se que temos apenas uma variável exógena, sendo que o valor de cada uma das 100 observações foi obtido a partir de uma distribuição uniforme com valores entre 0 e 10: X ~ U[0,10] . Os valores são mantidos os mesmos em todas as simulações. Para os parâmetros temos 0 1 , referente ao intercepto, e 1 1 , referente a variável exógena X. Esses valores foram incorporados ao DGP do modelo de Kelejian-Prucha. Foram analisados três casos particulares desse modelo: o modelo MQO com 0 e 0 ; o modelo de lag espacial com 1 1 e 0 ; e o modelo de erro espacial com 0 e . 2 2 25 Em seguida, determinam-se os erros estocásticos. Esses foram obtidos a partir de uma distribuição normal com média zero, e assumindo 2 1 . Foram gerados 100 valores a partir dessa distribuição para cada uma das simulações. No total são 100 simulações para o MQO, 20 para o modelo de lag espacial e 20 para o modelo de erro espacial. Foram assim gerados 1040 vetores com 100 valores aleatórios para o erro estocástico. De posse de toda essa informação, obtêm-se com o DGP do modelo de Kelejian-Prucha, os valores da variável dependente para esses modelos, também em 1040 vetores com 100 valores cada um. Ou seja, sabemos quais são os parâmetros que geraram os valores da variável dependente para cada um dos modelos. Os resultados são discutidos em separado para o MQO, o modelo de lag espacial e o modelo de erro espacial. 8.1 – MQO Com a informação descrita acima, geramos os 100 valores da variável dependente e os 100 valores da variável independente para cada uma das 1000 simulações. De posse dessa informação, estima-se o modelo MQO 1000 vezes utilizando as expressões: ˆMQO ( X ' X ) 1 X 'Y 2 ˆ MQO e' e . n Compara-se, assim, o modelo real que gerou os dados com o modelo estimado. Como vimos, o valor real para os parâmetros são 0 1 , 1 1 e 2 1 . O valor médio das 1000 estimações para esses parâmetros foi respectivamente de 1,0103, 0,9977 e 0,9749. Ou seja, os valores estimados são muito próximos dos reais. Os histogramas a seguir mostram a distribuição de cada um dos parâmetros, todos aproximadamente normais, com a média descrita acima. Note que a dispersão dos valores de 1 , que é o coeficiente de maior interesse empírico, era menor que para os demais parâmetros. Na 26 grande maioria das simulações, os valores ficaram entre 0,95 e 1,05, variações relativamente pequenas em torno do valor real. Porém, note que os valores absolutos dessas dispersões dependem dos valores arbitrários incluídos nas estimações. Histograma 1 – Dispersão dos valores de 0 250 200 150 100 50 0 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 Histograma 2 – Dispersão dos valores de 2 1 300 250 200 150 100 50 0 0.85 0.9 0.95 1 1.05 1.1 1.15 1.2 27 Histograma 3 – Dispersão dos valores de 2 300 250 200 150 100 50 0 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 8.2 – Modelo de lag espacial Também com a informação descrita acima, foram obtidos os 100 valores da variável dependente e da variável independente para o modelo de lag espacial. Estimou-se o modelo empírico 20 vezes. Note que para cada uma dessas estimativas devemos fazer um procedimento que demanda mais do computador do que o MQO. Essa é a razão desse número menor de simulações. O primeiro passo é obter o valor de que minimiza a expressão abaixo: 2 g ( ;Y , X ) ln Y ' MY 2 Y ' MWY 2Y 'W ' MWY ln(1 i ), n i onde M I ( X ( X ' X ) 1 X '. Assim, inicialmente, obtém-se a matriz M, e estimam-se os autovalores da matriz de peso. Substituem-se esses resultados em conjunto com o vetor Y e a matriz W na expressão acima. Esses valores serão constantes para todos os valores de . Qual valor de minimiza essa expressão? Isso é obtido numericamente. No caso da simulação aqui, variou-se o parâmetro de 0,05 até 0,95 em acréscimos de 0,05, em um total de 19 valores. Note que poderíamos ter feito de 0 a 0,95 em acréscimos de 0,01. A única diferença seria o tempo maior do computador. Com cada um dos valores de , calculou-se o valor da função g. 28 Esse mesmo procedimento foi feito 20 vezes. Note que dependendo do erro gerado na simulação, temos um valor de variável dependente distinto, o que leva a estimações de distintas. Lembrando que o valor real foi estipulado como 1 , o diagrama abaixo mostra os resultados 2 obtidos nas 20 simulações. Na maioria das simulações, o valor de que minimizava a expressão acima ficou em torno desse valor real. O valor médio das 20 simulações foi 0,485. Gráfico 1 – Valores de g ( ; Y , X ) para diferentes valores de 6.5 6 5.5 5 4.5 4 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Uma vez obtido o valor de que minimiza a expressão acima, estima-se os demais parâmetros do modelo a partir das expressões abaixo: ˆSAR ( X ' X ) 1 X ' ( I W )Y 1 n 2 ˆ SAR ( I W )Y ' M ( I W )Y A média dos valores das 20 replicações para 0 , 1 e 2 foram respectivamente 1,0318, 0,9981 e 0,9729, também muito próximos dos valores reais. 29 8.3 – Modelo de erro espacial De forma similar, de posse da informação descrita acima, foram obtidos os valores para a variável dependente e independente do modelo de erro espacial e estimaram-se os modelos empíricos. O primeiro passo foi estimar os 100 autovalores da matriz de peso, W, que são expressos aqui como i . Depois, deve-se obter o valor de que minimiza a expressão abaixo, sendo n o número de observações: 2 g ( ) ln ˆ 2 ln(1 i ), n i 2 1 onde ˆ (( I W )(Y Xˆ ))' ( I W )(Y Xˆ ), n ' ˆ (( X *) ( X *))1 ( X *)' Y *, X * ( I W ) X , e Y * ( I W )Y Inicialmente, substituem-se os valores encontrados para os autovalores, e os vetores X, Y e a matriz W nesta expressão. Assim, o único parâmetro a determinar é . Varia-se o valor do parâmetro numericamente e calcula-se o valor da função g para cada valor de . No caso da simulação aqui, variou-se o parâmetro de 0,05 até 0,95 em acréscimos de 0,05, em um total de 19 valores. O diagrama abaixo mostra os resultados, onde na maioria das simulações, o valor ficou em torno de 1 , o valor real. O valor médio das 20 simulações foi de 0,482. 2 30 Gráfico 2 – Valores de g (; Y , X ) para diferentes valores de 6 5.8 5.6 5.4 5.2 5 4.8 4.6 4.4 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 De posse do valor de , utilizamos as expressões acima para obter os valores dos demais parâmetros. As médias dos valores das 20 replicações para 0 , 1 e 2 foram respectivamente 0,9970, 0,9954 e 1,0009, também muito próximos dos valores reais. 9 - Matriz de covariância Nas seções 3 e 5 foram obtidas as matrizes de covariância respectivamente para a distribuição normal e para o MQO. Devido à extensão dos cálculos para a obtenção dessa matriz para os modelos espaciais, a obtenção da matriz de covariância desses modelos será feita apenas para o modelo de lag espacial. As passagens são feitas passo-a-passo e, assim, procedimentos análogos podem ser adaptados para os demais modelos. Como vimos, o estimador de máxima verossimilhança tem com propriedade ser assintoticamente eficiente com covariância assintótica positivamente definida no limite inferior de Cramér-Rao: 1 2 ln L . Var (ˆ) [ I ( )]1 E ' 31 Para o modelo de lag espacial, essa matriz de covariância tem o seguinte formato: 2 ln L 2 2 2 ' ln L [ I ( )]1 E 2 ' 2 ln L 2 ' 2 ln L 2 ' 2 ln L ' 2 ln L ' 2 ln L 2 ' 2 ln L ' 2 ln L ' 1 Então partimos do log da FV do modelo: l n 1 ( y Wy X )' ( y Wy X ), ln( 2 ) ln( 2 ) ln(1 i ) 2 2 2 i e obtemos as derivadas de primeira ordem. As duas primeiras já foram feitas anteriormente: l 1 2 2 X ' y 2 X ' X 2 1 l n ( y Wy X )' ( y Wy X ) 2 2 2 2 2 2( ) A terceira é realizada aqui pela primeira vez: n 1 ( y Wy X )' ( y Wy X ) l ln( 2 ) ln( 2 ) ln(1 i ) 2 2 2 i i i 1 2 (Wy)' ( y Wy X ) (1 i ) Em seguida, obtemos as derivadas de segunda ordem de todas as expressões acima, lembrando que a matriz é simétrica. Além disso, estimamos os valores esperados de cada uma das expressões. 1 2 ln L l n ( y Wy X )' ( y Wy X ) 2 2 2 2 2 2 2 2 ' 2 2( ) 2 n ( y Wy X )' ( y Wy X ) 2 2 2( ) 2( 2 ) 3 32 Obtemos o valor esperado dessa expressão: n n 2 2 E ( y Wy X )' ( y Wy X ) E ' 2 2 2 2 2( 2 ) 3 2( 2 ) 3 2( ) 2( ) n n 2 2 n E ' n 2 2 2 2 3 4 2 3 4 2( ) 2( ) 2( ) 2( ) 2( ) 1 2 ln L l n ( y Wy X )' ( y Wy X ) 2 2 2 2 2 2 2( ) 1 1 2 2 (Wy )' ( y Wy X ) 2 2 ' Wy ( ) ( ) Para simplificar a notação na obtenção do valor esperado definimos a matriz B W I W . 1 Utilizamos as seguintes expressões: E[ X ' ] E[ X ' ] 0; para uma matriz 1 x 1, tr (a11) a11; e temos a seguinte relação para o traço de uma matriz, tr ( AB) tr ( BA). O valor esperado da derivada acima tem o seguinte valor: 1 1 1 1 E 2 2 ' Wy 2 2 E ' Wy 2 2 E ' W I W ( X ) ( ) ( ) ( ) 1 1 1 2 2 E ' B ( X )) 2 2 E( ' B ) 2 2 E tr ( ' B ) ( ) ( ) ( ) 1 1 1 2 2 tr BE ' 2 2 tr BI 2 2 tr B ( ) ( ) 1 2 ln L l n ( y Wy X )' ( y Wy X ) 2 ' 2 2 2 2( 2 ) 2 1 1 2 2 ( X )' ( y Wy X ) 2 2 X ' ( ) ( ) 33 Obtemos o valor esperado dessa expressão: 1 E 2 2 X ' 0 ( ) Utilizando a notação i (i ) 2 para simplificar, obtemos a próxima expressão: (1 i ) 2 i 2 ln L l 1 2 (Wy )' ( y Wy X ) ' i (1 i ) (i ) 2 1 1 2 (Wy )' (Wy ) 2 ( y 'W 'Wy 2 i (1 i ) 1 2 (( I W ) 1 ( X ))'W ' ( B)( X )) 1 2 ( ' X ' ' ) B' B( X ) 1 2 [( ' X ' B' BX 2 ' B' BX ' B' B ] 1 E{ 2 [ ' X ' B ' BX 2 ' B' BX ' B ' B ]} 1 2 E[ ' X ' B' BX 2 ' B ' BX ' B' B ] 1 2 [ ' X ' B' BX E[2 ' B' BX ] E[ ' B' B ]] 1 2 [ ' X ' B' BX 0 E[tr ' B' B ]] 1 2 [ ' X ' B' BX trB ' BE[ ' ]] 1 2 [ ' X ' B' BX trB ' B ( I 2 )] 1 2 [ ' X ' B' BX ] trB ' B 34 1 2 ln L l 1 2 ( X )' ( y Wy X ) 2 X ' (Wy) ' 1 1 1 E 2 X ' Wy 2 EX ' Wy 2 EX ' B( X ) 1 2 X ' BX Finalmente: 1 2 ln L l 1 2 ( X )' ( y Wy X ) 2 X ' X ' 1 X'X E 2 X ' X 2 Substituindo todas essas expressões, obtemos a matriz de covariância: n 1 0 2 tr B 4 2( ) 1 1 1 [ I ( )]1 2 tr B 2 [( ' X ' B' BX ] trB' B 2 X ' BX 1 X'X 0 2 X ' BX 2 1 10 – Conclusão Esse texto faz parte de uma série que apresenta a econometria espacial em pontos introdutórios. No primeiro texto da serie foram discutidos alguns conceitos introdutórios sobre a econometria espacial. No segundo texto, discutiu-se como interpretar os coeficientes obtidos nos diversos modelos. Este terceiro texto da serie, discute como utilizar a função de verossimilhança para estimar alguns dos modelos espaciais comumente utilizados em análises empíricas. Conceitos básicos da função de verossimilhança foram discutidos, e os estimadores de máxima 35 verossimilhança para os modelos de lag espacial e de erro espacial foram derivados. A metodologia de estimação de vários dos demais modelos, a exceção do modelo de KelejianPrucha, é similar ao desses dois primeiros. Simulações no Matlab exemplificam a aplicação dos conceitos discutidos no texto. Como vimos, o modelo de Kelejian-Prucha engloba os modelos de lag espacial e de erro espacial. Assim, pode-se obter a metodologia de estimação desse primeiro diretamente a partir do apresentado para esses dois últimos. Assim, para que este texto não ficasse demasiadamente extenso, esse modelo não foi discutido em separado. Porém, existe uma característica na estimação desse modelo que deve ser mencionada, uma vez que ela não aparece nos demais modelos espaciais citados. Como vimos na estimação do modelo de lag espacial, o valor de é obtido por procedimentos numéricos. Ou seja, varia-se o valor desse parâmetro em pequenos acréscimos com o uso do computador, e verifica-se qual valor minimiza uma função em particular. De forma similar, o valor de também é obtido por procedimentos numéricos na estimação dos modelos de erro espacial. Nos dois casos temos análises unidimensionais: apenas um parâmetro, ou , varia de cada vez e temos um vetor com os valores de uma função especifica que está sendo minimizada. No modelo de KelejianPrucha devemos variar esses dois parâmetros de forma simultânea. Ou seja, em vez de um vetor com o valor da função sendo minimizada, obtemos uma matriz, e verificamos para qual par de valores de e que a função é minimizada. Assim, conceitualmente as estimações dos três modelos são similares, mas computacionalmente estimar o modelo de Kelejian-Prucha é mais demandante. O objetivo principal deste texto foi fornecer uma base para que o leitor possa entender de forma mais sólida como aplicar o estimador de máxima verossimilhança em estudos com modelos espaciais. Assim, caso seja necessário a utilização de modelos espaciais distintos dos normalmente utilizados, o próprio leitor pode desenvolver a sintaxe de estimação. 36 Referências Anselin, L. (1988) Spatial Econometrics: Methods and Models. Kluwer Academic, Dordrecht. Anselin, L., Bera, A., Florax, R. e Yoon, M. (1996) Simple diagnostic tests for spatial dependence. Regional Science and Urban Economics 26 (1), 77-104. Doreian, P. (1981) Estimating Linear Models with Spatially Distributed Data. Sociological Methodology 12, 359-388. Elhorst, J. (2010) Applied spatial econometrics: raising the bar. Spatial Economic Analysis 5 (1). Florax, R., Folmer, H. e Rey, S. (2003) Specification searches in spatial econometrics: the relevance of Hendry's methodology. Regional Science and Urban Economics 33 (5), 557-579. Greene, W. (2003) Econometric Analysis. Prentice Hall, 5º edition. Griffith, D. e Arbia, G. (2010) Detecting negative spatial autocorrelation in georeferenced random variables. International Journal of Geographical Information Science 24 (3), 417-437 Kelejian, H. (2011) Spatial models in Econometrics. Draft for the Spatial Econometric Advanced Institute. LeSage, J e Pace, R. (2009) Introduction to spatial econometrics. Taylor & Francis Group, Boca Raton. Mur, J. e Ângulo, A. (2009) Model selection strategies in a spatial setting: some additional results. Regional Science and Urban Economics 39, 200-213 Pace, K. (1997) Performing Large Spatial Regressions and Autoregressions. Economic Letters 54 (3), 283–291. Simon, C. e Blume, L.(2004) Matemática para Economistas. Bookman, 1ª edição. 37 Apêndice econométrico 1 Esse apêndice discute os seguintes pontos estimadores consistentes, matriz de covariância, matriz positivamente definida, estimadores eficientes e limite inferior de Cramér-Rao. Para uma discussão mais abrangente sobre esses temas ver Greene (2003) e Simon e Blume (2004). Estimador consistente Uma variável aleatória xn converge em probabilidade para uma constante se, lim Pr ob xn 0 para todo 0. Ou seja, a probabilidade de se observar uma n estimação muito diferente de é muito pequena, tendendo para zero, quando o tamanho da amostra tende a infinito. No caso de um estimador não enviesado, com E[ˆn ] , temos, lim Pr ob ˆn 0, ou n p lim ˆn , onde n é o tamanho da amostra. Um estimador assim é consistente. Matriz de covariância Dada uma distribuição multivariada de n variáveis, xi , temos como valor esperado de cada uma delas o seguinte vetor: 1 E[ x1 ] 2 E[ x2 ] E[ x]. ... ... E[ x ] n n A matriz de covariância é dada por: 38 ( x1 1 )( x1 1 ) ( x1 1 )( x2 2 ) ... ( x 2 )( x1 1 ) E[( x )( x )' ] 2 ... ... ( x )( x ) ... n 1 1 n ... ( x1 1 )( xn n ) ... ... . ... ... ... ( xn n )( xn n ) Ou de forma mais sintética: 11 12 ... Var ( x) 21 ... ... n1 ... ... 1n ... ... . ... ... ... nn Matriz positivamente definida Uma matriz A é positivamente definida se z' Az 0 para todo z 0 . Estimador eficiente Tome dois estimadores não enviesados ˆ1 e ˆ2 de : E[ˆ1 ] e E[ˆ2 ] . O estimador ˆ1 é mais eficiente que ˆ2 se a variância do primeiro for menor que a do segundo, Var (ˆ1 ) Var (ˆ2 ). Quando temos um vetor de parâmetros, o estimador ˆ1 é mais eficiente que ˆ2 se a matriz obtida pela subtração das respectivas matrizes de covariâncias, Var (ˆ2 ) Var (ˆ1 ) , for positivamente definida. 39 Limite inferior de Cramér-Rao Assumindo certas condições formais, a variância de um estimador não enviesado será pelo 1 1 ln L( ) 2 2 ln L( ) E menos tão grande quanto: [ I ( )] E . ' 1 Assim, um estimador não enviesado terá uma variância igual ou maior que esse limite. Um estimador com variação definida por esse limite é, portanto, eficiente, como é o caso dos estimadores de MV. Apêndice econométrico 2 – Função injetiva Uma função diz-se injetiva (ou injetora) se, para quaisquer valores de x1 e x2 , tais que x1 x2 , pertencentes ao domínio da função, tivermos f ( x1 ) f ( x2 ). Apêndice econométrico 3 - Estimador mínimos quadrados ordinários Aqui são incluídos alguns pontos introdutórios sobre a estimação de modelos MQO via mínimos quadrados ordinários, que justamente dão nome ao modelo, que foram julgados úteis para o entendimento das estimações utilizando a função de verossimilhança. A apresentação aqui se baseia em Greene (2003). O modelo clássico de regressão linear múltipla tem alguns pressupostos. Seguem alguns deles. 1) Existe uma linearidade na relação entre a variável dependente e as independentes, y X ; 2) A matriz das variáveis explicativas tem posto máximo, isto é, essas variáveis são independentes uma das outras, e, assim, não podem ser escritas como uma combinação linear das demais; 3) As variáveis explicativas são exógenas, ou seja, E[ X ] 0; 40 4) Os erros têm como características a homocedasticidade e a inexistência de autocorrelação, E[ ' X ] 2 I ; e 5) Os erros são distribuídos de forma normal, X ~ N (0, 2 I ) . Nos modelos MQO estimasse os parâmetros minimizando a soma dos quadrados dos resíduos. Para obtermos o estimador para os parâmetros do modelo, partimos do pressuposto 1): y X . Isolamos o erro e manipulamos a expressão: y X ' y X ' y X . Lembrando que ( A B)' A' B' e que ( AB)' B' A' , obtemos a seguinte expressão: S ( ) ' y X ' y X y'( X )' y X y' y ( X )' y y' ( X ) ( X )' ( X ) Lembrando-se da relação para vetores, X 'Y Y ' X , temos: S ( ) y' y 2( X )' y ' X ' X y' y 2 ' X ' y ' X ' X Derivando com relação a , obtemos a condição de primeira ordem para minimizar a expressão: S ( ) ( y' y 2 ' X ' y ' X ' X ) 2 X ' y 2 X ' X 0 Daí, temos: X ' X X ' y Como X tem posto máximo, como descrito no segundo pressuposto, a matriz X ' X tem inversa, e obtemos o estimador de mínimos quadrados ordinários: ˆMQO ( X ' X ) 1 X ' y . 41 Apêndice econométrico 4 – O Jacobiano Como vimos, o Jacobiano é incluído nas funções de verossimilhança. Esse apêndice econométrico apresenta alguns conceitos referentes a este tema, inicialmente pela definição e, em seguida, com a utilização na FV. A apresentação aqui se baseia em Greene (2003). A4.1 - Definição de Jacobiano de uma transformação A matriz Jacobiana é a matriz formada pelas derivadas parciais de primeira ordem de uma função vetorial. O valor absoluto do determinante dessa matriz, como veremos, é utilizado quando são feitas mudanças de função que define uma densidade de probabilidade conjunta, como a função de verossimilhança. Seguem exemplos. Partimos de uma função injetiva e de sua inversa: f ( x) y x f 1 ( y) Derivando essa função inversa, temos: dx d f 1 ( y ) . dy dy Se f (x) é linear, por exemplo, y ax b , temos: x y b a a J dx 1 dy a A função acima foi definida como f : . No caso de n funções e n variáveis, podemos escrever o Jacobiano como: 42 x1 y x 1 J ... y ' xn y 1 x1 yn ... ... xn ... yn ... Segue um exemplo para um sistema de equações lineares n x n, onde podemos escrever o sistema como: Y AX Se a inversa da matriz existir, temos: X A1Y , onde J A1 Exemplificando com um sistema 2 x 2, ficamos com: y1 a11 a12 x1 y 2 a21 a22 x2 1 x1 a11 a12 y1 x2 a21 a22 y2 x a11 a12 J y' a21 a22 1 O valor absoluto do determinante da matriz acima é o Jacobiano da transformação: abs J . Uma vez definido o Jacobiano da transformação, segue uma explicação do por que do uso deste na FV. 43 A4.2 – O uso do Jacobiano na função de verossimilhança Tome uma densidade de probabilidade contínua, definida por f x (x) . Assim, temos: b Pr ob( x b) f ( x)dx. x Seja y g (x) uma transformação monotônica de x e, portanto, temos x g 1 ( y) . Qual é a densidade de probabilidade de y? Em outras palavras, qual é a função f y ( y ) , tal que: a Pr ob( y a) f y ( y)dy ? Note que: f x ( x) f x ( g 1 ( y)) e que dx [ g 1 ( y )]' dy Substituindo esses termos na integral acima com g(a) = b, e, ainda, tomando apenas valores positivos para não termos pdf com valores negativos, temos: a Pr ob( y a) f x ( g 1 ( y )) [ g 1 ( y )]' dy Assim, f y ( y) f x ( g 1 ( y)) [ g 1 ( y)]' f x ( g 1 ( y)) J Ou seja, no caso da densidade conjunta, temos: Ly ( ; y) Lx ( g 1 ( ; y)) J 44