Primeira vez aqui? Seja bem vindo e cheque o FAQ!
x

Você pode me apresentar uma simulação Monte Carlo para explorar a diferença da estimação do estimador de mínimos quadrados e variável instrumental na presença de endogeneidade?

0 votos
119 visitas

1 Resposta

0 votos
respondida Ago 27, 2015 por danielcajueiro (5,261 pontos)  

Considere que você deseja estimar os parâmetros da equação (1) abaixo que pertence ao sistema de equações simultâneas:

\(y_1=\alpha_1 + \beta_1 y_2 + \gamma_1 z_1 + u_1\; (1)\)

\(y_2= \alpha_2 + \beta_2 y_1 + \gamma_2 z_2 + u_2\; (2),\)

onde \(y_1\) e \(y_2\) são endógenas e \(z_1\) e \(z_2\) são exógenas.

Podemos reescrever esse modelo da seguinte forma (isolando os valores de \(y_1\) e \(y_2\)):

\(y_1 (1-\beta_1 \beta_2)=\alpha_1 + \beta_1 \alpha_2 + \beta_1 \gamma_2 z_2 + \gamma_1 z_1 + \beta_1 u_2 + u_1\; (3) \)

\(y_2 (1-\beta_1 \beta_2)=\alpha_2 + \beta_2 \alpha_1 + \beta_2 \gamma_1 z_1 + \gamma_2 z_2 + \beta_2 u_1 + u_2 \; (4)\)

Vamos proceder fazendo simulações Monte Carlo. Com o conhecimento dos valores dos \(\alpha\)s, dos \(\beta\)s e dos \(\gamma\)s e de amostras das variáveis \(z\)s e \(u\)s geramos várias amostras de variáveis \(y\)s usando as equações (3) e (4) acima.

No código abaixo, implementamos dois estimadores:

a) Métodos de mínimos quadrados (que é viesado) supondo que estamos estimando a Equação (1) independente do conhecimento da segunda.

b) Estimador de variáveis instrumentais a partir do método dos momentos e supondo que

\(E[u_1]=0\; (5)\)

\(E[z_1 u_1]=0\; (6)\)

\(E[z_2 u_1]=0\; (7).\)

Note que para a Equação (7) ser válida (ou seja, dar informação adicional as equações (5) e (6)), nós precisamos que o coeficiente \(\delta_2\) na Equação (8) abaixo seja significante

\[y_2=\delta_0 + \delta_1 z_1 + \delta_2 z_2 + \epsilon.\]

Veja abaixo também que esse coeficiente ser significante implica que \(z_2\) consegue recuperar a casualidade de \(y_2\) para \(y_1\).

A figura abaixo apresenta as estimativas do parâmetro \(\beta_1\) na Equação (1), onde a linha vermelha apresenta o valor correto desse parâmetro. Note que o estimador (a) [figura da esquerda] é viesado e o estimador (b) [figura da direita]não é:

A imagem será apresentada aqui.

O código usado para implementar essa figura está apresentado a seguir:

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

plt.figure(num=None, figsize=(8, 4), dpi=80, facecolor='w', edgecolor='k')

def geracaoDados(n):
    alpha1=0.5
    alpha2=0.7    
    beta1=-1
    beta2=1
    gamma1=-0.5
    gamma2=2

    z1=np.random.normal(2,1,n)
    z2=np.random.normal(2,1,n)
    u1=np.random.normal(0,1,n)
    u2=np.random.normal(0,1,n)

    y1=((alpha1+beta1*alpha2)*np.ones(n)+beta1*gamma2*z2+gamma1*z1+beta1*u2+u1)/(1-beta1*beta2)   
    y2=((alpha2+beta2*alpha1)*np.ones(n)+beta2*gamma1*z1+gamma2*z2+beta2*u1+u2)/(1-beta1*beta2)   

    return y1,y2,z1,z2


def estimacaoVarInstrumental(y1,y2,z1,z2):
    sumZ1=np.sum(z1)
    sumZ2=np.sum(z2)
    sumY1=np.sum(y1)
    sumY2=np.sum(y2)
    sumY1Z1=np.dot(y1,z1)
    sumY2Z1=np.dot(y2,z1)
    sumY1Z2=np.dot(y1,z2)
    sumY2Z2=np.dot(y2,z2)
    sumZ1Z1=np.dot(z1,z1)
    sumZ1Z2=np.dot(z1,z2)
    A=np.array([[n,sumY2,sumZ1],[sumZ1,sumY2Z1,sumZ1Z1],[sumZ2,sumY2Z2,sumZ1Z2]])
    b=np.array([sumY1,sumY1Z1,sumY1Z2])
    beta=np.linalg.solve(A, b)
    return beta

def estimacaoOls(x,y):
    hatBeta = np.linalg.lstsq(x,y)[0]
    return hatBeta

if __name__ == '__main__':
    numeroAmostras=1000
    n=50 # numero de observacoes
    # Estimacao OLS
    vetorOLS=np.empty([numeroAmostras])
    for i in range(0,numeroAmostras):
        [y1,y2,z1,z2]=geracaoDados(n)
        x=np.vstack([np.ones(n),y2,z1]).T
        hatBetaOls=estimacaoOls(x,y1)
        vetorOLS[i]=hatBetaOls[1]
    axes=plt.subplot(121)    
    plt.setp(axes, xticks=[-1.3,-1.1,-0.9,-0.7,-0.5,-0.3,-0.1,0.1],yticks=[20,40,60,80,100,120,140])    
    plt.hold(True)
    plt.plot([-1, -1],[0,100],'r-')
    plt.hist(vetorOLS,bins=30)
    plt.title('OLS: Amostra de tamanho 1000')

    # Estimacao Variavel Instrumental
    vetorInst=np.empty([numeroAmostras])
    for i in range(0,numeroAmostras):
        [y1,y2,z1,z2]=geracaoDados(n)
        hatBetaInst=estimacaoVarInstrumental(y1,y2,z1,z2)
        vetorInst[i]=hatBetaInst[1]
    axes=plt.subplot(122)    
    plt.hold(True)
    plt.plot([-1, -1],[0,100],'r-')
    plt.hist(vetorInst,bins=30)
    plt.title('I: Amostra de tamanho 1000')
...