/* jagmnl1.prg does maximum likelihood for binary logit regression
on Jaggia data using Jaggia data. It is more advanced than jagmle1.prg
TO RUN THIS PROGRAM ....
(1) You need an implementation of Gauss that has Maxlik (a Gauss add-on)
(2) You need my program SELDATA.SRC in the same directory (to select data)
(3) You need my program GHMNL.SRC in the same directory (for gradient etc.)
You can get seldata.src and ghmnl.src from my web site.
Program jagmnl1.prg
1. Reads in a flat ascii data set.
2. Creates new data from transformation of original data
This includes making discrete chocie variables for dependent.
3. Converts the data to a Gauss data set, with names.
4. Selects a subset of the data to be used in the OLS regression
(This uses my separate procedure seldata.src)
5. Calls my procedure ghmnl.src to get logl, gradient and hessian
6. Calls the Gauss add-on maxlik to perform maximum likelihood.
The original data are from Sanjiv Jaggia and Satish Thosar, 1993,
"Multiple Bids as a Consequence of Target Management Resistance"
Review of Quantitative Finance and Accounting, 447-457.
The data in the ascii file jaggia.asc are on 126 firms, 12 variables per firm
For description of the variables see jagols1.prg
In addition we need to create a constant and Size Squared
*/
/******* Initial setup: libraries used etc. and output file ********
*/
New; @ clear workspace @
library maxlik; @ So Gauss knows we use procedures in maxlik @
maxset; @ A maxlik command @
disable; @ Allows missing values in calculation @
outfile = "jagmnl1.out"; @ all output goes to this file @
output file = ^outfile reset; @ delete earlier contents @
print "File jagmnl1.out is output from Gauss program jagmnl1.src";
print "Date is:";
print date;
print;
/******* Include additional Gauss code in separate files **********
seldata.src: Gauss program to select the subset of data
and variable names used in analysis.
This includes procedures opendat and readdat used below.
ghmnl.src: Gauss program with gradient and hessian for logit and mnl3
This includes procedures logli, loggrd, loghess
and mnl3li, mnl3grd, mnl3hess used below.
I assume you have these files in the same directory as this program.
*/
#include seldata.src;
#include ghmnl.src;
@ or if e.g. in directory gauss then include c:\gauss\seldata.src @
/******* Read in the data and create variables **********
Here we first read the data in as a flat ascii file.
The data are ordered by observation, with four variables
per line so that the first observation (on 12 variables)
takes up the first three lines etc.
Variables in order are DOCNO WEEKS NUMBIDS TAKEOVER BIDPREM INSTHOLD
SIZE LEGLREST REALREST FINREST REGULATN WHTKNGHT
Then we create a Gauss data set with names provided.
*/
load data1[126,12] = jaggia.asc; @ Data read in to matrix @
one = ones(126,1); @ A column vector of ones @
sizesq = data1[.,7].*data1[.,7]; @ A column vector = seventh column squared @
numbids = data1[.,3];
d2numbid = 1*(numbids.>=1); @ 0,1 indicator = 1 if numbids > 0 @
d3numbid = 1*(numbids.>=1) + 1*(numbids.>=2); @ 0, 1, 2 indicator @
data2 = data1~one~sizesq~d2numbid~d3numbid;
@ New variables added to existing data @
/******* Create Gauss data set with names *********
The created Gauss data set will be called jaggiatr.dat
and the associated names will be in file jaggiatr.dht
Note that an alternative is to use the Gauss utility program ATOG386.EXE
(where ATOG means ascii to Gauss). See Gauss Systems and Graphics manual.
*/
dataset = "jaggiatr"; @ or e.g. "c:\\gauss\\myex\\jaggiatr" @
let vnames = DOCNO WEEKS NUMBIDS TAKEOVER BIDPREM INSTHOLD SIZE LEGLREST
REALREST FINREST REGULATN WHTKNGHT ONE SIZESQ D2NUMBID D3NUMBID;
call saved(data2,dataset,vnames);
/******** Extract the subset of the Gauss data set to be used *********
*/
data ="jaggiatr"; @ or e.g. "c:\\gauss\\myex\\jaggiatr" @
print;
call dstat(data,0); @ check that data okay @
print;
let YVAR = D2NUMBID;
let XVAR = ONE REALREST BIDPREM SIZE;
/* Set data size dimensions and some other variables
*/
m = 1; @ Number of dependent variables @
k = rows(xvar); @ Number of regressors @
nobs = 126; @ Maximum number of rows to read in dataset @
/* Call opendat, rad, readdat in seldata.src
This gives
x - nobs*k matrix of explanatory variables
y - nobs*m matrix of dependent variables
indy - Index for the dependent variables
indx - Index for the exogenous variables
f1 - File handle for the open data set
n - Number of observations in the data set
nror - number of times data set needs to be read to be complete
rest - number of obervations in the time data is read
The data are read into memory one time.
Nobs must be equal to number of observations in the data set
*/
{f1,n,indy,indx} = opendat(data,yvar,xvar,&indices,m);
{nror, rest } = rad(nobs,n);
{ y,x } = Readdat(f1,1,nobs,indy,indx);
/******** Estimation by maxlik version 4 *********
The log-likelihood, gradient and hessian procedures
used appear in ghmnl.src which was included earlier
*/
@ Some setting up not needed for maxlik @
z = y~x;
clear y,x;
parnames = XVAR;
output on;
" BINARY LOGIT MAXIMUM LIKELIHOOD using Maxlik and ghmnl1.src";
" Dataset ";;Data;
" Nr of observations ";;n;
" Nr of times data is read ";;nror;
" Nr of obs. read last chew ";;rest;
" Dependent variable ";;$yvar;
" Parameters (Variables) names ";
$parnames';
"";
output off;
/* Analytical derivatives and hessian.
Can delete if instead use numerical derivatives
*/
_max_HessProc = &loghess; @ analytical hessian provided in hgpoiss.src @
_max_GradProc = &loggrd; @ analytical gradient provided in hgpoiss.src @
/* Setting various parameters for maxlik.
These are generally ones that over-ride the default values in maxlik.
For default values and descriptions see the maxlik manual
or the Gauss code in maxlik.src
*/
_max_Parnames = xvar; @ Parameter names which here = variable names @
_mlditer = 40; @ Not sure @
_max_Algorithm = 4; @ Newton Raphson = 4 BFGS = 2 (def) @
_max_MaxIters = 250; @ maximum number of iterations @
_mlmtime = 1e+5; @ Not sure 1e+5 = default @
_max_CovPar = 3; @ Method to compute covariance matrix @
@ 3 = heteroskedastic-consistent estimate @
_max_GradTol = 0.00001; @ Convergence tolerance for gradient of coeffs @
__output = 0; @ Minimal output. Use for Unix without X-windows @
start = (0.01*rndn(k,1)); @ start values for coeffs here randomly chosen @
@ probably gives the dim of parameter vector @
/* Call maxlik.
Arguments are
data = the name of the Gauss data set (or matrix if data in matrix)
yvar|xvar = the variables to be selected from the Gauss data set
loglli = procedure to compute logit log-likelihood (in ghmnl.src)
start = starting values of the parameters
*/
output on;
{blogit,phi0,grdlogit,covlogit,retcode}
= maxlik(data,yvar|xvar,&loglli,start);
/* Print the results from maxlik */
call maxprt(blogit,phi0,grdlogit,covlogit,retcode);
output off;
f1 = close(f1);