Hi May you help me correct my loop function. I want optim to estimates al_j; au_j; sigma_j; b_j by looking at 0 to 20, 21 to 40, 41 to 60 data points. The final result should have 4 columns of each of the estimates AND 4 rows of each of 0 to 20, 21 to 40, 41 to 60. ###MY code is n=20 runs=4 out=matrix(0,nrow=runs) llik = function(x) { al_j=x[1]; au_j=x[2]; sigma_j=x[3]; b_j=x[4] sum(na.rm=T, ifelse(a$R_j< 0, -log(1/(2*pi*(sigma_j^2)))- (1/(2*(sigma_j^2))*(a$R_j+al_j-b_j*a$R_m))^2, ifelse(a$R_j>0 , -log(1/(2*pi*(sigma_j^2)))- (1/(2*(sigma_j^2))*(a$R_j+au_j-b_j*a$R_m))^2, -log(pnorm(au_j,mean=b_j*a$R_m,sd=sqrt(sigma_j^2))- pnorm(au_j,mean=b_j*a$R_m,sd=sqrt(sigma_j^2))))) ) } start.par = c(0, 0, 0.01, 1) out1 = optim(llik, par=start.par, method="Nelder-Mead") for (i in 1: runs) { index_start=20*(i-1)+1 index_end= 20*i out[i]=out1[index_start:index_end] } out Thank you in advance Edward UCT ####My data R_j R_m -0.0625 0.002320654 0 -0.004642807 0.033333333 0.005936332 0.032258065 0.001060848 0 0.007114057 0.015625 0.005581558 0 0.002974794 0.015384615 0.004215271 0.060606061 0.005073116 0.028571429 -0.006001279 0 -0.002789594 0.013888889 0.00770633 0 0.000371663 0.02739726 -0.004224228 -0.04 0.008362539 0 -0.010951605 0 0.004682924 0.013888889 0.011839993 -0.01369863 0.004210383 -0.027777778 -0.04658949 0 0.00987272 -0.057142857 -0.062203157 -0.03030303 -0.119177639 0.09375 0.077054642 0 -0.022763619 -0.057142857 0.050408775 0 0.024706076 -0.03030303 0.004043701 0.0625 0.004951088 0 -0.005968731 0 -0.038292548 0 0.013381097 0.014705882 0.006424728 -0.014492754 -0.020115626 0 -0.004837891 -0.029411765 -0.022054654 0.03030303 0.008936428 0.044117647 8.16925E-05 0 -0.004827246 -0.042253521 0.004653096 -0.014705882 -0.004222151 0.029850746 0.000107267 -0.028985507 -0.001783206 0.029850746 -0.006372981 0.014492754 0.005492374 -0.028571429 -0.009005846 0 0.001031683 0.044117647 0.002800551 -- View this message in context: http://r.789695.n4.nabble.com/loop-in-optim-tp3643230p3643230.html Sent from the R help mailing list archive at Nabble.com.
Hi Edward, At least for me, your llik() function returns Inf for the starting values specified, so optim() never gets to estimate anything. You need to alter llik() or find starting parameters that work before worrying about getting the for loop working. Cheers, Josh On Mon, Jul 4, 2011 at 2:34 AM, EdBo <n.bowora at gmail.com> wrote:> Hi > > May you help me correct my loop function. > > I want optim to estimates al_j; au_j; sigma_j; ?b_j by looking at 0 to 20, > 21 to 40, 41 to 60 data points. > > The final result should have 4 columns of each of the estimates AND 4 rows > of each of 0 to 20, 21 to 40, 41 to 60. > > ###MY code is > > n=20 > runs=4 > out=matrix(0,nrow=runs) > > llik = function(x) > ? { > ? ?al_j=x[1]; au_j=x[2]; sigma_j=x[3]; ?b_j=x[4] > ? ?sum(na.rm=T, > ? ? ? ?ifelse(a$R_j< 0, -log(1/(2*pi*(sigma_j^2)))- > ? ? ? ? ? ? ? ? ? ? ? ? ? (1/(2*(sigma_j^2))*(a$R_j+al_j-b_j*a$R_m))^2, > ? ? ? ? ifelse(a$R_j>0 , -log(1/(2*pi*(sigma_j^2)))- > ? ? ? ? ? ? ? ? ? ? ? ? ? (1/(2*(sigma_j^2))*(a$R_j+au_j-b_j*a$R_m))^2, > > -log(pnorm(au_j,mean=b_j*a$R_m,sd=sqrt(sigma_j^2))- > ? ? ? ? ? ? ? ? ? ? ? ? ? pnorm(au_j,mean=b_j*a$R_m,sd=sqrt(sigma_j^2))))) > > ? ? ? ) > > ? } > > start.par = c(0, 0, 0.01, 1) > out1 = optim(llik, par=start.par, method="Nelder-Mead") > > > for (i in 1: runs) > { > ?index_start=20*(i-1)+1 > ?index_end= 20*i > ?out[i]=out1[index_start:index_end] > } > out > > > Thank you in advance > > Edward > UCT > ####My data > > R_j ? ? ? ? ? ? R_m > -0.0625 ? ? ? ? 0.002320654 > 0 ? ? ? ? ? ? ? -0.004642807 > 0.033333333 ? ? 0.005936332 > 0.032258065 ? ? 0.001060848 > 0 ? ? ? ? ? ? ? 0.007114057 > 0.015625 ? ? ? ?0.005581558 > 0 ? ? ? ? ? ? ? 0.002974794 > 0.015384615 ? ? 0.004215271 > 0.060606061 ? ? 0.005073116 > 0.028571429 ? ? -0.006001279 > 0 ? ? ? ? ? ? ? -0.002789594 > 0.013888889 ? ? 0.00770633 > 0 ? ? ? ? ? ? ? 0.000371663 > 0.02739726 ? ? ?-0.004224228 > -0.04 ? ? ? ? ? 0.008362539 > 0 ? ? ? ? ? ? ? -0.010951605 > 0 ? ? ? ? ? ? ? 0.004682924 > 0.013888889 ? ? 0.011839993 > -0.01369863 ? ? 0.004210383 > -0.027777778 ? ?-0.04658949 > 0 ? ? ? ? ? ? ? 0.00987272 > -0.057142857 ? ?-0.062203157 > -0.03030303 ? ? -0.119177639 > 0.09375 ? ? ? ? 0.077054642 > 0 ? ? ? ? ? ? ? -0.022763619 > -0.057142857 ? ?0.050408775 > 0 ? ? ? ? ? ? ? 0.024706076 > -0.03030303 ? ? 0.004043701 > 0.0625 ? ? ? ? ?0.004951088 > 0 ? ? ? ? ? ? ? -0.005968731 > 0 ? ? ? ? ? ? ? -0.038292548 > 0 ? ? ? ? ? ? ? 0.013381097 > 0.014705882 ? ? 0.006424728 > -0.014492754 ? ?-0.020115626 > 0 ? ? ? ? ? ? ? -0.004837891 > -0.029411765 ? ?-0.022054654 > 0.03030303 ? ? ?0.008936428 > 0.044117647 ? ? 8.16925E-05 > 0 ? ? ? ? ? ? ? -0.004827246 > -0.042253521 ? ?0.004653096 > -0.014705882 ? ?-0.004222151 > 0.029850746 ? ? 0.000107267 > -0.028985507 ? ?-0.001783206 > 0.029850746 ? ? -0.006372981 > 0.014492754 ? ? 0.005492374 > -0.028571429 ? ?-0.009005846 > 0 ? ? ? ? ? ? ? 0.001031683 > 0.044117647 ? ? 0.002800551 > > > > > > > > > > > > > -- > View this message in context: http://r.789695.n4.nabble.com/loop-in-optim-tp3643230p3643230.html > Sent from the R help mailing list archive at Nabble.com. > > ______________________________________________ > R-help at 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. >-- Joshua Wiley Ph.D. Student, Health Psychology University of California, Los Angeles http://www.joshuawiley.com/
Hi I have re-worked on my likelihood function and it is now working(#the code is below#). May you help me correct my loop function. I want optim to estimates al_j; au_j; sigma_j; b_j by looking at 0 to 20, 21 to 40, 41 to 60 data points. The final result should have 4 columns of each of the estimates AND 4 rows of each of 0 to 20, 21 to 40, 41 to 60. #likelihood function a=read.table("D:/hope.txt",header=T) attach(a) a llik = function(x) { al_j=x[1]; au_j=x[2]; sigma_j=x[3]; b_j=x[4] sum(na.rm=T, ifelse(a$R_j< 0, log(1/(2*pi*(sigma_j^2)))- (1/(2*(sigma_j^2))*(a$R_j+al_j-b_j*a$R_m))^2, ifelse(a$R_j>0 , log(1/(2*pi*(sigma_j^2)))- (1/(2*(sigma_j^2))*(a$R_j+au_j-b_j*a$R_m))^2, log(ifelse (( pnorm (au_j, mean=b_j * a$R_m, sd= sqrt(sigma_j^2))- pnorm(al_j, mean=b_j * a$R_m, sd=sqrt (sigma_j^2) )) > 0, (pnorm (au_j,mean=b_j * a$R_m, sd= sqrt(sigma_j^2))- pnorm(al_j, mean=b_j * a$R_m, sd= sqrt(sigma_j^2) )), 1)) )) ) } start.par = c(-0.01,0.01,0.1,1) out1 = optim(llik, par=start.par, method="Nelder-Mead") out1 -- View this message in context: http://r.789695.n4.nabble.com/loop-in-optim-tp3643230p3645031.html Sent from the R help mailing list archive at Nabble.com.
2 comments below. On 07/07/2011 06:00 AM, r-help-request at r-project.org wrote:> Date: Wed, 6 Jul 2011 20:39:19 -0700 (PDT) > From: EdBo <n.bowora at gmail.com> > To: r-help at r-project.org > Subject: Re: [R] loop in optim > Message-ID: <1310009959045-3650592.post at n4.nabble.com> > Content-Type: text/plain; charset=us-ascii > > I have one last theoretical question, I did not adjust my code prior so that > it maximise the likehood function. I googled that to make optim maximise you > multiply fn by -1. > > In my code, would that be the same as saying "-sum" on the "sum" part of my > code (see below)? > > llik = function(x) > { > al_j=x[1]; au_j=x[2]; sigma_j=x[3]; b_j=x[4] > sum(na.rm=T, > > ifelse(a$R_j< 0, log(1 / ( sqrt(2*pi) * sigma_j) )- >The optimx package has a "maximize" control because we felt that the fnscale approach, while perfectly correct, is not comfortable for users and is not standard across other optimization tools. Note that this package is undergoing a fairly extensive overhaul at the moment (the development version is on R-forge in the project 'optimizer') to include some safeguards on functions that return NaN etc. as well as a number of other changes -- hopefully improvements. A second comment on this looping: Why do you not use the parameters from the last estimation as the starting parameters for the next? Unless you are expecting very extreme changes over the moving window of data, this should appreciably speed up the optimization. John Nash