The RATS Software Forum

I want to modify a probit model to make it non-standard so that I can't use DDV, so I began programming using MAXIMIZE to maximize the likelihood function. The problem is that I can't even get MAXIMIZE to replicate or even come close to the output of the probit example given in the WINRATS example program probit.rpf. I'm even using the example from page RM-298 of the reference manual which shows how to estimate a probit model using MAXIMIZE. I must be doing something wrong but I can't see it. I'm attaching the program and pasting it below. Any suggestions would be much appreciated.

* THE FIRST PART IS JUST AN EXTRACT FROM THE WINRATS EXAMPLE PROG probit.rpf
* RATS User's Guide, Example from Section 12.2

* PROBIT.RPF
* Example of estimation of probit and logit models
*
open data probit.dat
data(org=obs) 1 95 public1_2 public3_4 public5 private $
years teacher loginc logproptax yesvm
*
* Linear probability model
*
linreg yesvm
# constant public1_2 public3_4 public5 private years teacher loginc logproptax
*
ddv(dist=logit) yesvm
# constant public1_2 public3_4 public5 private years teacher loginc logproptax
ddv(dist=probit) yesvm
# constant public1_2 public3_4 public5 private years teacher loginc logproptax

* NEXT PART JUST REDEFINES THE VARIABLES

set x2 = public1_2
set x3 = public3_4
set x4 = public5
set x5 = private
set x6 = years
set x7 = teacher
set x8 = loginc
set x9 = logproptax
ddv(dist=probit) yesvm
# constant x2 x3 x4 x5 x6 x7 x8 x9

** NEXT PART SETS UP THE PROBIT MODEL FOR MAXIMIZE AND IS ADAPTED FROM THE EXAMPLE GIVEN ON PAGE RM-298 OF THE REFERENCE MANUAL

compute b1=%beta(1)
compute b2=%beta(2)
compute b3=%beta(3)
compute b4=%beta(4)
compute b5=%beta(5)
compute b6=%beta(6)
compute b7=%beta(7)
compute b8=%beta(8)
compute b9=%beta(9)
comp success=0.0
do i = 1, 95
comp success = yesvm(i) + success
end do i
comp reps = 95.0

frml zfrml = b1+b2*x2+b3*x3+b4*x4+b5*x5+b6*x6+b7*x7+b8*x8+b9*x9
frml probit = (z=zfrml(t)) , $
success*log(%cdf(z))+(reps-success)*log(1-%cdf(z))
nonlin b1 b2 b3 b4 b5 b6 b7 b8 b9
maximize(iters=500,SUBITERS=200,method=bfgs) probit

**** MAXIMIZE DOES NOT APPEAR TO WORK - WHAT AM I DOING WRONG?

That example has REPS and SUCCESSES as separate data series (i.e. they're different for each individual). You're generating a single set of values across all data. If you're actually trying to work out how to do that with that data set, you would want

set success = yesvm
set reps = 1

replacing your current calculation of SUCCESS. Then the rest is the same. The point estimates will be the same--the standard errors won't be, and will be quite wrong because you're starting with the converged estimates. If you use zeros as guess values, the standard errors will end up being similar, but not the same, since probit can do full 2nd derivatives, while BFGS can only approximate them.

Note, BTW, that your original calculation of the total number of successes can be done with a single SSTATS instruction

SSTATS / YESVM>>SUCCESS

And that works perfectly and solves my problem - Thank You, Tom!