/*
** Lesson 7.3: Hypothesis Testing 
** CES Production Function: b[4]=1/b[3]
** Judge, et. al. [1988], Chapter 12
*/
use gpe2;
output file=gpe\output7.3 reset;
load x[30,3]=gpe\judge.txt;

call reset;
_nlopt=1; @ MAXLIK @
_method=5;
_iter=100;
_tol=1.0e-5;
_conv=1;
_jacob=0;

/* Lagrangian Multiplier Test */
@ based on constrained estimation @
_b={1.0,0.5,-1.0};
call estimate(&rfc,x);
b1=__b|(1/__b[3]);  @ original parametrization, b[4]=1/b[3] @
ll1=__ll;           @ log-likelihood @
e=rf(x,b1);         @ estimated errors @
s2=meanc(e^2);      @ estimated error variance @
g=gradp2(&rf,x,b1); @ gradient of error function @
lm=(e'g)*invpd(g'g)*(g'e)/s2;

/* Wald Test */
@ based on unconstrained estimation @
_b={1.0,0.25,-1.0,-1.0};
call estimate(&rf,x);
b2=__b;             @ estimated parameters @
vb2=__vb;           @ estimated var-cov. of parameters @
ll2=__ll;           @ log-likelihood @
w=eqc(b2)'*invpd(gradp(&eqc,b2)*vb2*gradp(&eqc,b2)')*eqc(b2);

/* Likelihood Ratio Test */
lr=-2*(ll1-ll2);

print "Wald Test = " w;
print "Lagrangian Multiplier Test = " lm;
print "Likelihood Ratio Test = " lr;
end;

proc rf(data,b); @ unconstrained residual function @
    local l,k,q,e;
    l=data[.,1];
    k=data[.,2];
    q=data[.,3];
    e=ln(q)-b[1]-b[4]*ln(b[2]*l^b[3]+(1-b[2])*k^b[3]);
    retp(e);
endp;

proc rfc(data,b); @ constrained residual function @
    local l,k,q,e;
    l=data[.,1];
    k=data[.,2];
    q=data[.,3];
    e=ln(q)-b[1]-(1/b[3])*ln(b[2]*l^b[3]+(1-b[2])*k^b[3]);
    retp(e);
endp;

proc eqc(b);  @ constraint function @
    retp(b[3]*b[4]-1);
endp;