/*
** Lesson 8.2: Logit Model of Economic Education
** Greene [1999], Example 19.1
** See also Spector and Mazzeo [1980]
*/
use gpe2;
output file=gpe\output8.2 reset;

n=33;
load data[n,4]=gpe\grade.txt;
gpa=data[2:n,1];
tuce=data[2:n,2];
psi=data[2:n,3];
grade=data[2:n,4];
z=gpa~tuce~psi~ones(rows(grade),1);

call reset;

@ logit model: estimation @
_nlopt=2;  @ using component log-likelihood @
_method=4;
_iter=50;
_b={0.5,0.0,0.5,0};

call estimate(&logitf,grade~z);
/*
_deriv=&logitf1;
call estimate(&logitf,grade~z);
*/

@ logit model: interpretation @
p=1./(1+exp(-z*__b));
print " Probability                    Slopes";;
print p~(p.*(1-p).*__b[1:rows(__b)-1]');

end;

/* log-likelihood function of logit model */
proc logitf(x,b);
    local k,z,f;
    k=rows(b);
    z=x[.,2:k+1]*b;
    f=1./(1+exp(-z));    @ same as: f=exp(z)./(1+exp(z)); @
                         @ logistic distribution function @
    retp(x[.,1].*ln(f)+(1-x[.,1]).*ln(1-f));
endp;

/* 1st derivatives of log-likelihood function of logit model */
proc logitf1(x,b);
    local z,k,f;
    k=rows(b);
    z=x[.,2:k+1]*b;
    f=1./(1+exp(-z));    @ same as: f=exp(z)./(1+exp(z)); @
                         @ logistic distribution function @
    retp((x[.,1]-f).*x[.,2:k+1]);
endp;