[R] Are the error messages of ConstrOptim() consisten with each other?
Dear Friends. I found something very puzzling with constOptim(). When I change the parameters for ConstrOptim, the error messages do not seem to be consistent with each other: constrOptim(c(0.5,0.3,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) Error in constrOptim(c(0.5, 0.3, 0.5), f = fit.error, gr = fit.error.grr, : initial value not feasible constrOptim(c(0.5,0.9,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) Error in constrOptim(c(0.5, 0.9, 0.5), f = fit.error, gr = fit.error.grr, : initial value not feasible constrOptim(c(0.3,0.5,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) Error in f(theta, ...) : argument lambda1 is missing, with no default I only changed the parameters, how come the lambda1 that is not missing in the first 2 cases suddently become missing? For your convenience, I put the complete code below: Best Wishes Yuchen Luo rm(list = ls()) mat=5 rint=c(4.33,4.22,4.27,4.43,4.43,4.44,4.45,4.65,4.77,4.77) tot=rep(13319.17,10) sh=rep(1553656,10) sigmae=c(0.172239074,0.188209271,0.193703774,0.172659891,0.164427247,0.24602361,0.173555309,0.186701165,0.193150456, 0.1857315601) ss=c(56.49,56.39,56.55,57.49,57.37,55.02,56.02,54.35,54.09, 54.67) orange=rep(21.25,10) apple2=expression(rint*(1.0-rec)*(1.0-(pnorm(-lambda/2.0+log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/lambda)-((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))*pnorm(-lambda/2.0-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/lambda))+(exp(rint*(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)*ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0+0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))-sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*! 1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(-sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0+0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))-(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0+0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt((lambda*lam! bda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh* 1000.0))-sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(-sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0+0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))/((pnorm(-lambda/2.0+log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/lambda)-((ss+(tot/sh*100! 0.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))*pnorm(-lambda/2.0-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/lambda))-(pnorm(-sqrt((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*mat+lambda*lambda)/2.0+log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda)))/sqrt((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*mat+lambda*lambda))-((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))*pnorm(-sqrt((sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*mat+lambda*lambda)/2.0-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0)/lbar*exp
Re: [R] Are the error messages of ConstrOptim() consisten with each other?
Dear Professor Murdoch. Thank you for your help! 1. I believe c(0.5,0.3,0.5) satisfies the constrain because I did the following experiment ui=-1*ui ci=-1*ci constrOptim(c(0.5,0.3,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) The same error message pops up. Any theta ( in this case, c(0.5,0.3,0.5)) cannot violate both ui%*%theta=ci and -ui%*%theta=-ci. 2. There is lambda1 available. The 0.3 in c(0.5,0.3,0.5) is lambda1. If you plug c(0.5,0.3,0.5) into fit.error and fit.error.grr by fit.error(0.5,0.3,0.5) fit.error.grr(0.5,0.3,0.5) It works. Best Wishes Yuchen Luo On 9/10/07, Duncan Murdoch [EMAIL PROTECTED] wrote: Yuchen Luo wrote: Dear Friends. I found something very puzzling with constOptim(). When I change the parameters for ConstrOptim, the error messages do not seem to be consistent with each other: constrOptim(c(0.5,0.3,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) Error in constrOptim(c(0.5, 0.3, 0.5), f = fit.error, gr = fit.error.grr , : initial value not feasible Not feasible means it doesn't satisfy the constraints. constrOptim(c(0.5,0.9,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) Error in constrOptim(c(0.5, 0.9, 0.5), f = fit.error, gr = fit.error.grr , : initial value not feasible constrOptim(c(0.3,0.5,0.5), f=fit.error, gr=fit.error.grr, ui=ui,ci=ci) Error in f(theta, ...) : argument lambda1 is missing, with no default This time your starting values satisfied the constraints, so your objective function was called, but you didn't pass it a value for lambda1. I only changed the parameters, how come the lambda1 that is not missing in the first 2 cases suddently become missing? For your convenience, I put the complete code below: Best Wishes Yuchen Luo rm(list = ls()) mat=5 rint=c(4.33,4.22,4.27,4.43,4.43,4.44,4.45,4.65,4.77,4.77) tot=rep(13319.17,10) sh=rep(1553656,10) sigmae=c(0.172239074,0.188209271,0.193703774,0.172659891,0.164427247, 0.24602361,0.173555309,0.186701165,0.193150456, 0.1857315601) ss=c(56.49,56.39,56.55,57.49,57.37,55.02,56.02,54.35,54.09, 54.67) orange=rep(21.25,10) apple2=expression(rint*(1.0-rec)*(1.0- (pnorm(-lambda/2.0+log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/lambda)-((ss+(tot/sh*1000.0)*lbar)/(tot/sh* 1000.0)/lbar*exp(lambda*lambda))*pnorm(-lambda/2.0-log(((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/lambda))+(exp(rint*(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh* 1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)*ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rint/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))-sqrt(0.25+2.0*rint/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 *(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*! 1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(-sqrt(0.25+2.0*rint/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+sqrt(0.25+2.0*rint/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 *(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))-(((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rint/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt((lambda*lam! bda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh* 1000.0))-sqrt(0.25+2.0*rint/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0*(sigmae*ss/(ss+lbar*(tot/sh* 1000.0)))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(-sqrt(0.25+2.0*rint/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+sqrt(0.25+2.0*rint
[R] What does it mean by initial value not available?
Dear friends. I use ConstrOptim( ) and got error message initial value not available. My understanding of initial value not available is that one of the following 3 cases happens: 1.The objective function is not well defined at the point of the initial value. 2. The differentiation of the objective function is not well defined at the point of the initial value. 3. The initial value violate the constrain of ui %*% theta - ci = 0 But my situation does not belong to any of the above cases. I have attached my code bellow and could you please help me take a look? Best Wishes Yuchen Luo rm(list = ls()) ui=matrix(c(1,-1,0,0,0,0,0,0,1,-1,0,0,0,0,0,0,1,-1),6,3) ci=c(0,-0.5,0,-2,0,-0.6) apple=expression(rr*(1.0-rec)*(1.0-(pnorm(-lambda/2.0+log(((ss+(tot/sh* 1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/lambda)-((ss+(tot/sh*1000.0)*lbar)/(tot/sh* 1000.0)/lbar*exp(lambda*lambda))*pnorm(-lambda/2.0-log(((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/lambda))+(exp(rr*(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh* 1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)*ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))-sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 *(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(-sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 *(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))-(((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))-sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 *(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(-sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0))+sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 *(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt((lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 ))/((pnorm(-lambda/2.0+log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/lambda)-((ss+(tot/sh*1000.0)*lbar)/(tot/sh* 1000.0)/lbar*exp(lambda*lambda))*pnorm(-lambda/2.0-log(((ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/lambda))-(pnorm(-sqrt((sigmae*ss/(ss+lbar*(tot/sh* 1000.0)))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*mat+lambda*lambda)/2.0+log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/sqrt((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0)))*mat+lambda*lambda))-((ss+(tot/sh* 1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda))*pnorm(-sqrt((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*mat+lambda*lambda)/2.0-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/sqrt((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*(sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*mat+lambda*lambda)))*exp(-rr*mat)-(exp(rr*(lambda*lambda/(sigmae*ss/(ss+lbar*(tot/sh* 1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0)*ss+(tot/sh*1000.0 )*lbar)/(tot/sh*1000.0)/lbar*exp(lambda*lambda))^(sqrt(0.25+2.0*rr/ (sigmae*ss/(ss+lbar*(tot/sh*1000.0)))/(sigmae*ss/(ss+lbar*(tot/sh*1000.0 +0.5)*pnorm(-log(((ss+(tot/sh*1000.0)*lbar)/(tot/sh*1000.0 )/lbar*exp(lambda*lambda)))/((sigmae*ss/(ss+lbar*(tot/sh*1000.0 )))*sqrt(mat+(lambda*lambda
Re: [R] What does it mean by initial value not available?
Dear Professor Mordoch. Thank you very much for your help! Your time is highly appreciated! I do intend to optimize over 3 parameters and the way I did it is constrOptim(c(0.5,0.3,0.5), fit.error, fit.error.grr, ui=ui,ci=ci). Also, There is a missing value for both sigmae and ss. (The last 3rd and 4th line) Please correct them as: sigmae=c(0.172239074,0.188209271,0.193703774,0.172659891,0.164427247,0.24602361,0.173555309,0.186701165,0.193150456, 0.1857315601) ss=c(56.49,56.39,56.55,57.49 ,57.37,55.02,56.02,54.35,54.09,54.67) Best Wishes Yuchen Luo On 9/9/07, Duncan Murdoch [EMAIL PROTECTED] wrote: On 09/09/2007 7:01 AM, Yuchen Luo wrote: Dear friends. I use ConstrOptim( ) and got error message initial value not available. My understanding of initial value not available is that one of the following 3 cases happens: 1.The objective function is not well defined at the point of the initial value. 2. The differentiation of the objective function is not well defined at the point of the initial value. 3. The initial value violate the constrain of ui %*% theta - ci = 0 But my situation does not belong to any of the above cases. I have attached my code bellow and could you please help me take a look? I haven't tried your code, but there's one obvious error: ... fit.error=function(rec,lambda,lbar) {sum((eval(apple)*1000-orange)^2/(orange^2)) } The function optimizes over the elements of the first parameter. You've got two other parameters there, and I think you're trying to optimize over them as well. Put them all into one vector. The documentation for constrOptim doesn't make this as clear as it should; I'll clarify (by copying the docs from ?optim). Duncan Murdoch fit.error.grr=function(rec,lambda, lbar) {drec=sum(eval(D(apple,'rec'))*(eval(apple)*1000-orange)/(orange^2)) dlambda=sum(eval(D(apple,'lambda'))*(eval(apple)*1000-orange)/(orange^2)) dlbar=sum(eval(D(apple,'lbar'))*(eval(apple)*1000-orange)/(orange^2)) c(drec,dlambda,dlbar) } rr=c(4.33,4.22,4.27,4.43,4.43,4.44,4.45,4.65,4.77,4.77) tot=rep(13319.17,10) sh=rep(1553656,10) sigmae=c(0.172239074,0.188209271,0.193703774,0.172659891,0.164427247, 0.24602361,0.173555309,0.186701165,0.193150456) ss=c(56.49,56.39,56.55,57.49,57.37,55.02,56.02,54.35,54.09) orange=rep(21.25,10) constrOptim(c(0.5,0.3,0.5), fit.error, fit.error.grr, ui=ui,ci=ci) [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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. __ R-help@stat.math.ethz.ch 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.
[R] 'initial value not feasible' in constrOptim
Dear friends. I am using function constrOptim(c(0.5,0.3,0.5), fit.error, fit.error.grr, ui=-1*ui,ci=-1*ci) and I am confronted with error message initial value not feasible I plug in the initial value of (0.5,0.3,0.5) to function fit.error and fit.error.grr and have pretty reasonable result. I inequality ui %*% theta - ci = 0 as suggested in the R manual and it is satisfied. In case that this is a typo of the manual, I let ui=-ui and ci=-ci and try constrOptim again but the same warning message pops up. Could you please point me a way out of this? I am actually trying to translate a fortran code to R and the function I want to replace is DBCPOL, which used 'the complex method described in Nelder and Mead (1965) and Gill et al. (1981). I believe contrOptim is better than it because it is newer, is it? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] 'singular gradient matrix’ when using nl s() and how to make the program skip nls( ) and run on
Dear friends. I use nls() and encounter the following puzzling problem: I have a function f(a,b,c,x), I have a data vector of x and a vectory y of realized value of f. Case1 I tried to estimate c with (a=0.3, b=0.5) fixed: nls(y~f(a,b,c,x), control=list(maxiter = 10, minFactor=0.5 ^2048),start=list(c=0.5)). The error message is: number of iterations exceeded maximum of 10 Case2 I then think maybe the value of a and be are not reasonable. So, I let nls() estimate (a,b,c) altogether: nls(y~f(a,b,c,x), control=list(maxiter = 10, minFactor=0.5 ^2048),start=list(a=0.3,b=0.5,c=0.5)). The error message is: singular gradient matrix at initial parameter estimates. This is what puzzles me, if the initial parameter of (a=0.3,b=0.5,c=0.5) can create 'singular gradient matrix', then why doesn't this 'singular gradient matrix' appear in Case1? I have tried to change the initial value of (a,b,c) around but the problem persists. I am wondering if there is a way out. My another question is, I need to run 220 of nls() in my program with different y and x. When one of the nls() encounter a problem, the whole program stops. In my case, the 3rd nls() runs into a problem. I would still need the program to run the remaining 217 nls( )! Is there a way to make the program skip the problematic nls() and complete the ramaining nls()'s? Your help will be highly appreciated! Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] How to signal the end of the table?
I am using a for loop to read a table row by row and I have to specify how many records are there in the table. I need to read row by row because the table is huge and the memory not large enough for the whole table.: number.of.records=100 fp=file(abc.csv,r) pos=seek(fp, rw=read) for (i in 1:number.of.record){ current.row=scan(file=fp, sep=',', what=list(count=1, cusip6=, idate=1, spread=1.1, vol252=1.1, vol1000=1.1, st_debt=1.1, lt_debt=1.1 , total_liab= 1.1, cr=1.1, shrout=1.1, prc=1.1, mkt_cap=1.1, rtng=1.1, sec=1.1, cr3m=1.1, cr5y=1.1, ust3m=1.1, ust5y=1.1), flush=TRUE, nlines=1,quiet=T) ... } I need to know the number of records in the table and put it in the variable named number.of.records. When I have a new table that I do not know how many records it has, I use excel to open the file to figure it out and put it in variable number.of.records. I often have many tables to try and every one of them has thousands of recordsit takes a lot of time and trouble to adjust the code every time I read a different table. I am wondering if I can change the for loop to a while loop: while (the end of the table has not been reached) { current.row=scan(file=fp, sep=',', what=list(count=1, cusip6=, idate=1, spread=1.1, vol252=1.1, vol1000=1.1, st_debt=1.1, lt_debt=1.1 , total_liab= 1.1, cr=1.1, shrout=1.1, prc=1.1, mkt_cap=1.1, rtng=1.1, sec=1.1, cr3m=1.1, cr5y=1.1, ust3m=1.1, ust5y=1.1), flush=TRUE, nlines=1,quiet=T) ... } The problem is how to articulate while (the end of the table has not been reached), or equivalently, how to signal the end of the table? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] How to let the program know the end of a file is reached?
Dear collegues. I am using scan( ) to read from a table (a csv file). I am wondering how to let the program know that the end of the file is reached? Your help will be highly appreciated! Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] How to store the parameter estimated by nls( ) to a variable?
Dear Colleagues. I believe this should be a problem encountered by many: nls( ) is a very useful and efficient function to use if we are just to display the estimated value on screen. What if we need R to store the estimated parameter in a variable? For example: x=rnorm(10, mean=1000, sd=10) y=x^2+100+rnorm(10) a=nls(y~(x^2+para),control=list(maxiter = 1000, minFactor=0.5 ^1024),start=list(para=0.0)) How to store the estimated value of para in this case, in a variable, say, b? It is easy to display a and find all the information. How ever, I need to fit a different set of x and y in every loop of my code and I need to store the estimated values for further use. I have checked both the online manual and several S-plus books but no example as such showed up. Your help will be highly appreciated! Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
Re: [R] How to store the parameter estimated by nls( ) to a variable?
Dear Professor Murdoch. Thank you so much Best Wishes Yuchen Luo On 8/12/07, Duncan Murdoch [EMAIL PROTECTED] wrote: Yuchen Luo wrote: Dear Colleagues. I believe this should be a problem encountered by many: nls( ) is a very useful and efficient function to use if we are just to display the estimated value on screen. What if we need R to store the estimated parameter in a variable? For example: x=rnorm(10, mean=1000, sd=10) y=x^2+100+rnorm(10) a=nls(y~(x^2+para),control=list(maxiter = 1000, minFactor=0.5 ^1024),start=list(para=0.0)) How to store the estimated value of para in this case, in a variable, say, b? It is easy to display a and find all the information. How ever, I need to fit a different set of x and y in every loop of my code and I need to store the estimated values for further use. I have checked both the online manual and several S-plus books but no example as such showed up. Your help will be highly appreciated! coef(a) will get what you want. coef() works for most modelling functions where it makes sense. Duncan Murdoch [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] How to write to a table column by column?
Dear friends. Every loop of my program will result in a list that is very long, with a structure similar to the one below: Lst - list(name=Fred, wife=Mary, daily.incomes=c(1:850)) Please notice the large size of daily.incomes. I need to store all such lists in a csv file so that I can easily view them in Excel. Excel cannot display a row of more than 300 elements, therefore, I have to store the lists as columns. It is not hard to store one list as a column in the csv file. The problem is how to store the second list as a second column, so that the two columns will lie side by side to each other and I can easily compare their elements. ( If I use 'appened=TRUE', the second time series will be stored in the same column. ) Thank you for your tine and your help will be highly appreciated!! Best Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] set up a blank csv file and write time series to it row by row
Dear Friends. Greetings! I have asked the question of how to set up a blank file and write a list to it as a row for many times, with the number of lists unknown. I have received many beautiful solutions. Thanks go to Professor *Murdoch, Professor *Menne, Professor Grothendieck and Dr. Olshansky. I have organized the solutions below: ## *Set up a blank table in harddrive and write to it row by row* *#Method 1* *blank = data.frame(name=character(0), wife=character(0), no.children=numeric(0))* write.csv(blank, 'file1.csv', row.names=FALSE) a1 = list(name=Tom, wife=Joy, no.children=9) a2 = list(name=Paul, wife=Alic, no.children=5) write.table(a1, file=file1.csv, sep=',', append=TRUE, row.names=FALSE, col.names=FALSE) write.table(a2, file=file1.csv, sep=',', append=TRUE, row.names=FALSE, col.names=FALSE)** * * *#Method 2* *blank = data.frame(name=character(0), wife=character(0), no.children=numeric(0))* write.csv(blank, 'file2.csv') a1 = list(name=Tom, wife=Joy, no.children=9) a2 = list(name=Paul, wife=Alic, no.children=5) write.table(a1, file=file2.csv, sep=',', append=TRUE, row.names=2, col.names=FALSE) write.table(a2, file=file2.csv, sep=',', append=TRUE, row.names=3, col.names=FALSE)** ### My problem now is, how to write a time series (instead of a list) to a csv file? Also, how to set up such a csv file to accept the time series? I know the length of the time series' but I do not know how many of them are going to come up. Examples are : bb1=c(1:10) bb2=c(101:110) How to write bb1 and bb2 to a csv file and how to set up blank csv file to accept such time series in the first place? Your help will be highly appreciated!!! Best Wishes! Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] Saving an expression to a file
Dear Friends. I have a very long expression and I use function D to find its derivative, which will be even longer. I save the resulting expression in a variavle, say bbb. But when I tried to display bbb on the screen, the R screen is not long enough for me to view it. Is there a way to save the expression to a file? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] Looping through a series of (csv) files
Dear Colleagues. This should be a very common operation and I believe there should be a nice way in R to handle it. I couldn't find it in the manual or by searching online. I am wondering if I could ask for some help in this community. I have 48 csv files; each stores the data for a specific month. The 48 corresponding months are consecutively from January, 2001 to December, 2004. I name the files A200101, A200102, .., A200112, A200201, ,etc. I want to process file A2000101 and store the result to a new file named B200101, process file A200102 and store the result to a new file named B200102 etc. I do not want to manually change a little bit of the code to read a different file and write to a different file every time. I want the program to be able to loop through all the files. My question is, how to loop through the 48 files? Your help will be highly appreciated! Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] Setting up a blank table with column names in the hard drive
Dear Friends. Greetings! This should be a very common operation and I believe there should be a nice way in R to handle it. I couldn't find it in the manual or by searching on line. I am wondering if I could ask for some help in this community. I am trying to record the results of my program to a csv file in the hard drive so as to save memory space and also to read the results in excel after running the program. Every loop of my program will result in a list with element belonging to different class. For example, things like a1 - list(name=Fred, wife=Mary, no.children=3) a2 - list(name=Tom, wife=Joy, no.children=9) a3 - list(name=Paul, wife=Alic, no.children=5) I want the columns to have titles, in the example above, I want to see the title name, wife and no.children in the excel file. To set up the table in the csv file, I need to add at least one row in the table up front. How ever, I do not have the first loop of the program completed yet at that time. If I add a row that is meaningless, how may I delete it after all the loops are completed and all the meaningful rows are added to the table? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] set up a talbe with column name and write to the table row by row
Dear Friends. I tried to set up a table with column name and write to it row by row. I have 3 rows to write to the table: a1 - list(name=Fred, wife=Mary, no.children=3) a2 - list(name=Tom, wife=Joy, no.children=9) a3 - list(name=Paul, wife=Alic, no.children=5) When I write to it, the first write command has to be different than others in order to put the column name in the table: write.table(a1, file = tt.csv, sep=',',row.names=FALSE,col.name=TRUE) write.table(a2, file=tt.csv, sep=',', append=TRUE, row.names=FALSE, col.names=FALSE) write.table(a3, file=tt.csv, sep=',', append=TRUE, row.names=FALSE, col.names=FALSE) Is there a way to seperate the operation of 1.putting the column name in the table and 2. writing the first row to the table ? Best Wishes Yuchen [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] if within a function
Dear Friends. I found a puzzling phenomenon in R when you use 'if' within a function: # defining a function aaa aaa=function(a) {if (a==1) {aaa=1}; if (a!=1) {aaa=2} } # using the function: b=20 bbb=aaa(b) bbb [1] 2 typeof(bbb) [1] double c=1 ccc=aaa(c) ccc NULL typeof(ccc) [1] NULL It seems that only the last 'if' phrase works. Is it an instrinsic weakness of R? Is there a way to get around it? ( I use 'elseif' to get around this when there are only two cases to choose from, but what if there are more than two cases to choose from?) Best Yuchen [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] How to lower the 'miniFactor' when using nls( ) function?
Dear Friends. I used nls( ) function and encountered error message step factor 0.000488281 reduced below 'minFactor' of 0.000976563. I then tried the following: 1) Put nls.control(minFactor = 1/(4096*128)) inside the brackets of nls, but the same error message shows up. 2) Put nls.control(minFactor = 1/(4096*128)) as a separate command before the command that use nls( ) function, again, the same thing happens, although the R responds to the nls.control( ) function immediately: - $maxiter [1] 50 $tol [1] 1e-05 $minFactor [1] 1.907349e-06 -- I am wondering how may I change the minFactor to a smaller value? The manual that comes with the R software about nls( ) is very sketchy --- the only relevant example I see is a separate command like 2). A more relevent question might be, is lower the 'minFactor' the only solution to the problem? What are the other options? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ [EMAIL PROTECTED] 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.
[R] nls.control( ) has no influence on nls( ) !
Dear Friends. I tried to use nls.control() to change the 'minFactor' in nls( ), but it does not seem to work. I used nls( ) function and encountered error message step factor 0.000488281 reduced below 'minFactor' of 0.000976563. I then tried the following: 1) Put nls.control(minFactor = 1/(4096*128)) inside the brackets of nls, but the same error message shows up. 2) Put nls.control(minFactor = 1/(4096*128)) as a separate command before the command that use nls( ) function, again, the same thing happens, although the R responds to the nls.control( ) function immediately: - $maxiter [1] 50 $tol [1] 1e-05 $minFactor [1] 1.907349e-06 -- I am wondering how may I change the minFactor to a smaller value? The manual that comes with the R software about nls( ) is very sketchy --- the only relevant example I see is a separate command like 2). A more relevent question might be, is lower the 'minFactor' the only solution to the problem? What are the other options? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
[R] How to manipulate the pointer of a file?
Dear friends. With file( ) to obtain a pointer of a file, every time we use scan ( ) to read one row of it, the pointer will point to the next row of the file. In the following example, d1 and d2 are obtained the same way but they correspond to different rows of the same file because the pointer of the file moves down a row when a row of the file is read. The following is an example: a1 - list(name=Fred, wife=Mary, no.children=3) a2 - list(name=Tom, wife=Joy, no.children=9) a3 - list(name=Paul, wife=Alic, no.children=5) write.table(a1, file = tt.csv, sep=',',row.names=FALSE,col.name=TRUE) write.table(a2, file=tt.csv, sep=',', append=TRUE, row.names=FALSE, col.names=FALSE) write.table(a3, file=tt.csv, sep=',', append=TRUE, row.names=FALSE, col.names=FALSE) fp=file(tt.csv,r) c=scan(file=fp, sep=',', what=list(c1=, c2=, c3=), flush=TRUE, nlines=1) d1=scan(file=fp, sep=',', what=list(name=, wife=, no.kids=0), flush=TRUE, nlines=1) d1 R output: --- $name [1] Fred $wife [1] Mary $no.kids [1] 3 - d2=scan(file=fp, sep=',', what=list(name=, wife=, no.kids=0), flush=TRUE, nlines=1) d2 R Output: $name [1] Tom $wife [1] Joy $no.kids [1] 9 -- My question is, how to manipulate the pointer of the file further? For example, what if I need the pointer to go back to the previous row? Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ [EMAIL PROTECTED] 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.
[R] Reading a large csv file row by row
Hi, my friends. When a data file is large, loading the whole file into the memory all together is not feasible. A feasible way is to read one row, process it, store the result, and read the next row. In Fortran, by default, the 'read' command reads one line of a file, which is convenient, and when the same 'read' command is executed the next time, the next row of the same file will be read. I tried to replicate such row-by-row reading in R.I use scan( ) to do so with the skip= xxx option. It takes only seconds when the number of the rows is within 1000. However, it takes hours to read 1 rows. I think it is because every time R reads, it needs to start from the first row of the file and count xxx rows to find the row it needs to read. Therefore, it takes more time for R to locate the row it needs to read. Is there a solution to this problem? Your help will be highly appreciated! Best Wishes Yuchen Luo [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.
Re: [R] Reading a csv file row by row
Hi, my friends. When a data file is large, loading the whole file into the memory all together is not feasible. A feasible way is to read one row, process it, store the result, and read the next row. In Fortran, by default, the 'read' command reads one line of a file, which is convenient, and when the same 'read' command is executed the next time, the next row of the same file will be read. I tried to replicate such row-by-row reading in R.I use scan( ) to do so with the skip= xxx option. It takes only seconds when the number of the rows is within 1000. However, it takes hours to read 1 rows. I think it is because every time R reads, it needs to start from the first row of the file and count xxx rows to find the row it needs to read. Therefore, it takes more time for R to locate the row it needs to read. Is there a solution to this problem? Your help will be highly appreciated! [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch 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.