In replicating some work done with package GAUSS,
I wish to maximize a log-likelihood function. Naturally,
I want to use analytic derivatives. However, "sneq"
(the name given to function) does not appear in the
derivative table.
Any suggestions for how to use deriv() with a multiline
function?
Thanks,
Bruce
ps - the particular function appears below. It is not mine; I translated it
as best I could from the GAUSS program. I am well aware that the function
could be much better written and that bounding the likelihood to prevent zeroes
is illegitimate - this is replication work, and I have to try to stick to the
original
code.
sneq <- function(sed,d0,cm1,cm2,cm3,cm4,cm5,cm6,cm7,cm8,cm9,cm10,cm11,cm12,
cm13,cm14,cm15,cm16,sep,betar0,betad0,bgov,bpre,bjim,bvinc,bvblk,bin,bedu,
tm1,tm1,tm3,tm4,tm5,tm6,tm7,tm8,tm9,tm10,outst,rho,zoe0,winst,TRIEHEAT,
JULYEC2,PRESINC,VP4,HOME4X,ADAACA,DEV3,DEV2,LEGIS,ZECON41,DSOUTH4,SOUTH4,
SOUTH64B,WEST76,NCENT72,NEWENG64,TIEBREAK,VOTPOP,GOVRACE,RAIN,PCAPREAL,
EDUCAT,JIMCROW,PCBLCK,PCIN,OUTST,DUMMY88,DUMMY84,DUMMY80,DUMMY76,DUMMY72,DUMMY68,
DUMMY64,DUMMY60,DUMMY56,DUMMY52){
DBAR <- d0 + cm1*TRIHEAT + cm2*JULYEC2 + cm3*PRESINC + cm4*VP4+cm5*HOME4X+
cm6*ADAACA + cm7*DEV3 + cm8*DEV2 + cm9*LEGIS + cm10*ZECON41 + cm11*DSOUTH64 +
cm12*SOUTH4 + cm13*SOUTH64B + cm14*WEST76 + cm15*NCENT72 + cm16*NEWENG64
gama <- exp(betar0-betad0)
alpha <- gama/(gama+1)
DELTA <- gama/( (1/DOOV) - (1-gama) )
DBAR <- DBAR*(1-a) + DELTA*a
Ef <- (TIEBREAK + winst)*rho*alpha*(1-alpha)*pnorm( (alpha-DBAR)/sed)/sed -
zoe0*VOTPOP^zoe1
TM <- tm1*DUMMY88 + tm2*DUMMY84 + tm3*DUMMY80 + tm4*DUMMY76 + tm5*DUMMY72 +
tm6*DUMMY68 + tm7*DUMMY64 +
tm8*DUMMY60 + tm9*DUMMY56 + tm10*DUMMY52
Xb <- betad0 + bgov*GOVRACE + bpre*RAIN + bvinc*PCAPREAL + bedu*EDUCAT +
bjim*JIMCROW + bvblk*PCBLCK + bin*PCIN + TM + outst*OUTSTATE
prob <- ( pnorm( (log(DOOP/DELTA) - rho*Ef - Xb) / sep) / sep)/DOOP*(
( pnorm ( ( DELTA-DBAR)/sed ) /sed) * gama/((1-(1-gama)*DOOV)^2)*(1-a)+a)
b <- prob<0.1^(10)
prob <- prob*(1-b) + b*0.1^(10)
# L gets a minus sign because S-PLUS minimizes a function
L <- -log(prob)
}
-----------------------------------------------------------------------
This message was distributed by s-news@wubios.wustl.edu. To unsubscribe
send e-mail to s-news-request@wubios.wustl.edu with the BODY of the
message: unsubscribe s-news
|