Your function will not evaluate as coded. i.e., llfn(start.par) doesn't "work",
as there
are unequal numbers of arguments. Also, while R allows you to use variables
that are not
explicitly defined for a function, I've almost always got into trouble if I
don't pass
them VERY carefully.
Finally, as author of the code used to put CG in optim, I'll advise against its
use. One
of my least successful pieces of code. Rcgmin is better, but you really do need
analytic
derivatives to make it sing. For 5 parameters, use NM, or better the nmk from
dfoptim package.
Best, JN
On 05/15/2012 06:00 AM, r-help-requ...@r-project.org wrote:
> Message: 13
> Date: Mon, 14 May 2012 04:21:57 -0700 (PDT)
> From: infinitehorizon <barisvar...@hotmail.com>
> To: r-help@r-project.org
> Subject: Re: [R] Discrete choice model maximum likelihood estimation
> Message-ID: <1336994517063-4629910.p...@n4.nabble.com>
> Content-Type: text/plain; charset=us-ascii
>
> Hello again,
>
> I changed the name to tt.
> and for a and tt actually I was getting them from data, I didn't put them
> here in the question. Now I restructured my code and below I copy the full
> code, I tried many things but still getting the same error, I don't
> understand where is the mistake.
>
> I also defined one more variable to increase comprehension of the problem.
> Instead of data, I define three representative vectors in the beginning. If
> you run this code, you will see the error message.
>
> # Variables, a: discrete choice variable-dependent, x and tt are independent
> variables, tt is binary
> a <- c(1,1,2,3,2,3,1,2,2,3,1,1)
> x <- c(23,26,12,27,10,30,13,20,23,44,17,15)
> tt <- c(1,0,0,1,1,1,0,0,1,0,1,1)
>
> # First, just to see, the linear model
>
> LM <-lm(a ~ x+tt)
> coefLM <- coefficients(LM)
> summary(LM)
>
> # Probabilities for discrete choices for a=3, a=2 and a=1 respectively
> P3 <- function(bx,b3,b,tt) {
> P <- exp(bx*x+b3+b*(tt==1))/(1-exp(bx*x+b3+b*(tt==1)))
> return(P)
> }
> P2 <- function(bx,b2,b,tt) {
> P <- exp(bx*x+b2+b*(tt==1))/(1-exp(bx*x+b2+b*(tt==1)))
> return(P)
> }
> P1 <- function(bx,b1,b,tt) {
> P <- exp(bx*x+b1+b*(tt==1))/(1-exp(bx*x+b1+b*(tt==1)))
> return(P)
> }
>
> # Likelihood functions for discrete choices for a=3, a=2 and a=1
> respectively
>
> L3 <- function(bx,b1,b2,b3,b,tt) {
> P11 <- P1(bx,b1,b,tt)
> P22 <- P2(bx,b2,b,tt)
> P33 <- P3(bx,b3,b,tt)
>
> L3l <- P11*P22*P33
> return(L3l)
> }
>
> L2 <- function(bx,b1,b2,b,tt) {
> P11 <- P1(bx,b1,b,tt)
> P22 <- P2(bx,b2,b,tt)
>
> L2l <- P11*P22
> return(L2l)
> }
>
> L1 <- function(bx,b1,b,tt) {
> P11 <- P1(bx,b1,b,tt)
>
> L1l <- P11
> return(L1l)
> }
>
> # Log-likelihood function
>
> llfn <- function(param) {
>
> bx <- param[1]
> b1 <- param[2]
> b2 <- param[3]
> b3 <- param[4]
> b <- param[5]
>
> lL1 <- log(L1(bx,b1,b2,b,tt))
> lL2 <- log(L2(bx,b1,b2,b3,b,tt))
> lL3 <- log(L3(bx,b1,b2,b3,b,tt))
>
> llfn <- (a==1)*lL1+(a==2)*lL2+(a==3)*lL3
> }
> start.par <- c(0,0,0,0,0)
> est <- optim(param=start.par,llfn,x=x, a=a, tt=tt, method =
> c("CG"),control=list(trace=2,maxit=2000), hessian=TRUE)
>
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
