Rincón de Práctica Stata

Contenido

Rincón de Práctica Stata#

MLE usando Stata#

Ejemplo 1.#

Sea $\epsilon_i\sim i.i.d.\, N(\mu,\sigma^2)$.

Planear la función log-verosimilitud.
Realizar en el computador simulación de 100 datos con $\mathbb{N}(\mu_0,\sigma_0^2)=\mathbb{N}(1.7,3^2)$ y estimación ML.

Respuesta propuesta. Como vimos, la log-likelihood function estará dada por

\[\ell(\boldsymbol{\theta})=-\frac{n\cdot log(2\pi)}{2}-\frac{n\cdot log(\sigma^2)}{2}-\sum_{i=1}^N{\frac{(\epsilon_i-\mu)^2}{2\cdot \sigma^2}}\]

Así, una propuesta a la parte computacional en Stata sería la siguiente:

import pandas as pd
import ipystata

%%stata
  * 1. Simular datos 
  set    seed 123
  set    obs  100
  scalar mu    = 1.7
  scalar sigma = 3
  
  gen eps = rnormal(mu,sigma)
  sum eps

  * 2. Escribir la función log-likelihood
  cap prog drop my_mle
  pro def my_mle
	  args lnf mu sigma
	  quie replace `lnf' = -0.5*ln(2*_pi*`sigma'^2)-($ML_y1 - `mu')^2/(2*`sigma'^2)
  end
  
  * 3. Optimizar la log-likelihood
  ml model lf my_mle (mu: eps=)(sigma: eps=)
  ml max

Number of observations (_N) was 0, now 100.

    Variable |        Obs        Mean    Std. dev.       Min        Max
-------------+---------------------------------------------------------
         eps |        100    1.516204    3.095611  -12.11381    8.02918

  1.   args lnf mu sigma
  2.   quie replace `lnf' = -0.5*ln(2*_pi*`sigma'^2)-($ML_y1 - `mu')^2/(2*`sigma'^2)
  3.   end

initial:       log likelihood =     -<inf>  (could not be evaluated)
feasible:      log likelihood = -2126.5089
rescale:       log likelihood = -279.44898
rescale eq:    log likelihood = -260.90153
Iteration 0:   log likelihood = -260.90153  
Iteration 1:   log likelihood = -254.80682  
Iteration 2:   log likelihood = -254.39028  
Iteration 3:   log likelihood = -254.38986  
Iteration 4:   log likelihood = -254.38986  

                                                           Number of obs = 100
                                                           Wald chi2(0)  =   .
Log likelihood = -254.38986                                Prob > chi2   =   .

------------------------------------------------------------------------------
             | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
mu           |
       _cons |   1.516204   .3080094     4.92   0.000     .9125165    2.119891
-------------+----------------------------------------------------------------
sigma        |
       _cons |   3.080094   .2177955    14.14   0.000     2.653223    3.506965
------------------------------------------------------------------------------

Ejemplo 2.#

Bajo normalidad, MLE y OLS dan el mismo vector estimado de parámetros para el modelo lineal.

Sea,

\[Y=\boldsymbol{X}\boldsymbol{\beta}+u \hspace{0.5cm};\hspace{0.5cm}u_i\sim N(0,\sigma^2)\]

Planear la función log-verosimilitud
Encontrar $\hat{\boldsymbol{\theta}}_{MLE}=(\hat{\boldsymbol{\beta}},\hat{\sigma}^2)$
Encontrar la matriz de información
Realizar en el computador simulación de 100 datos para $y_i=\beta_0+\beta_1x_{1i}+\beta_2x_{2i}+u_i$, donde $x_1\sim U[0,10]$, $x_2\sim\chi^2_1$, $u_i\sim\mathbb{N}(0,\sigma^2=3^2)$ y $\boldsymbol{\beta}'=(\beta_0,\beta_1,\beta_2)'=(1,0.75,3.2)$. Realizar estimación ML de los parámetros y comparar con OLS.

Respuesta sugerida.

Como vimos, para $\boldsymbol{x}=(1,X_1,x_2)'$ y $u_i=y_i-\boldsymbol{x}_i'\boldsymbol{\beta}$ con $\boldsymbol{\beta}'=(\beta_0,\beta_1,\beta_2)'$ se tiene que

\[\ell(\boldsymbol{\theta})=-\frac{n\cdot log(2\pi\sigma^2)}{2}-\frac{(\boldsymbol{y}-\boldsymbol{X}\boldsymbol{\beta})'(\boldsymbol{y}-\boldsymbol{X}\boldsymbol{\beta})}{2\sigma^2}\]

y a partir de las condiciones de primer orden (F.O.C.) evaluadas en el estimador ML, se obtiene: $$\hat{\boldsymbol{\beta}}_{MLE}=(X'X)^{-1}(X'Y)\hspace{0.5cm},\hspace{0.5cm}\hat{\sigma}^2_{MLE}=(Y-X\hat{\beta})'(Y-X\hat{\beta})\cdot n^{-1}$$

Notar que la expresión para $\beta$ corresponde a la misma que se obtiene por OLS.

Además, después de resolver, se obtiene la siguiente matriz de información

\[\begin{split}I(\boldsymbol{\beta}_0)= \left[\begin{array}{cc} \frac{1}{\sigma_0^2}(X'X) & 0\\ 0 & \frac{N}{2\sigma_0^4} \end{array}\right]\end{split}\]

Así, una propuesta para el código de programación en Stata es la siguiente:

%%stata
  * 1. Simular datos 
    set seed 123
    set obs  100

    scalar sigma = 3

    gen x1 = runiform(0,10)
    gen x2 = rchi2(1)
    gen u  = rnormal(0,sigma)

    gen y = 1 + 0.75*x1 + 3.2*x2 + u

  * 2. Escribir la función log-likelihood
    cap prog drop my_ols
    pro def my_ols
	    args lnf Xb sigma
	    quie replace `lnf' = ln(normalden($ML_y1,`Xb',`sigma'))
    end

  * 3. Optimizar la log-likelihood
    ml model lf my_ols (xb: y = x1 x2)(sigma:)
    ml max
  
  * 4. Comparación con OLS
    reg y x1 x2

Number of observations (_N) was 100, now 100.

  1.     args lnf Xb sigma
  2.     quie replace `lnf' = ln(normalden($ML_y1,`Xb',`sigma'))
  3.     end

initial:       log likelihood =     -<inf>  (could not be evaluated)
feasible:      log likelihood = -822.51142
rescale:       log likelihood =  -348.3095
rescale eq:    log likelihood =  -348.3095
Iteration 0:   log likelihood =  -348.3095  (not concave)
Iteration 1:   log likelihood = -305.30761  
Iteration 2:   log likelihood = -272.72806  
Iteration 3:   log likelihood = -251.69736  
Iteration 4:   log likelihood = -251.33146  
Iteration 5:   log likelihood = -251.33126  
Iteration 6:   log likelihood = -251.33126  

                                                        Number of obs =    100
                                                        Wald chi2(2)  = 445.22
Log likelihood = -251.33126                             Prob > chi2   = 0.0000

------------------------------------------------------------------------------
           y | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
xb           |
          x1 |   .8781539   .1068187     8.22   0.000     .6687931    1.087515
          x2 |   3.229557   .1698408    19.02   0.000     2.896675    3.562439
       _cons |   .2802955   .6265632     0.45   0.655    -.9477458    1.508337
-------------+----------------------------------------------------------------
sigma        |
       _cons |   2.987312   .2112349    14.14   0.000       2.5733    3.401325
------------------------------------------------------------------------------

      Source |       SS           df       MS      Number of obs   =       100
-------------+----------------------------------   F(2, 97)        =    215.93
       Model |  3973.17594         2  1986.58797   Prob > F        =    0.0000
    Residual |  892.403429        97  9.20003535   R-squared       =    0.8166
-------------+----------------------------------   Adj R-squared   =    0.8128
       Total |  4865.57937        99  49.1472664   Root MSE        =    3.0332

------------------------------------------------------------------------------
           y | Coefficient  Std. err.      t    P>|t|     [95% conf. interval]
-------------+----------------------------------------------------------------
          x1 |   .8781539   .1084579     8.10   0.000     .6628949    1.093413
          x2 |   3.229557   .1724472    18.73   0.000     2.887297    3.571817
       _cons |   .2802955   .6361785     0.44   0.660    -.9823427    1.542934
------------------------------------------------------------------------------