[R] problem with optim()

Jonathan Daily biomathjdaily at gmail.com
Thu May 19 15:54:58 CEST 2011


I cc'd r-help to get this message back on the board.

Again, you haven't really answered my questions. Comments inline.

On Thu, May 19, 2011 at 9:45 AM, chirine wolley
<wolley.chirine at hotmail.com> wrote:
> First, thank you for ur response...
> Actually I didn't write the entire code ...X, Y and p_tilde are of
> course all defined at the beginning of my code.

And without them, your code example is not very enlightening, since it
doesn't run.

> Normally the optimal values that I get should let me create a model for
> prediction. However the model currently obtained is not good at all....

You still have provided neither a) what you expect nor b) what optim()
is returning. What context is normal?

> And I think that the fact  it gives me the exactly same values if I
> change the function to maximize shows that there's a problem in my code....

While I can't say if this makes a difference or if you accounted for
it, optim minimizes by default as outlined in the ?optim page.

>
>> Date: Thu, 19 May 2011 09:29:02 -0400
>> Subject: Re: [R] problem with optim()
>> From: biomathjdaily at gmail.com
>> To: wolley.chirine at hotmail.com
>> CC: r-help at r-project.org
>>
>> What do you mean when you say "wrong results"? What do you expect for
>> the output? Your code doesn't work for me because it references X in
>> places and X is not defined.
>>
>> Have you tested your functions to make sure they return reasonable values?
>>
>> On Thu, May 19, 2011 at 9:17 AM, chirine wolley
>> <wolley.chirine at hotmail.com> wrote:
>> >
>> > Dear R-users,
>> >
>> > I would like to maximize the function g above which depends on 4
>> > parameters (2 vectors, 1 real number, and 1 matrix)  using optim() and BFGS
>> > method. Here is my code:
>> >
>> > # fonction to maximize
>> >
>> > g=function(x)
>> > {
>> > x1 = x[1:ncol(X)]
>> > x2 = x[(ncol(X)+1)]
>> > x3 =
>> > matrix(x[(ncol(X)+2):(ncol(X)+1+ncol(X)*ncol(Y))],nrow=ncol(X),ncol=ncol(Y))
>> > x4 = x[(ncol(X)+1+ncol(X)*ncol(Y)+1):length(x)]
>> > res1=rep(0,nrow(X))
>> > res2=matrix(0,nrow=nrow(X),ncol=ncol(Y))
>> > log.res2=matrix(0,nrow=nrow(X),ncol=ncol(Y))
>> > res2.b=rep(0,nrow(X))
>> > res3 = rep(0,nrow(X))
>> > res3.b = rep(0,nrow(X))
>> > for (i in 1:nrow(X))
>> > {
>> > res1[i]=1/(1+exp(-t(x1)%*%X[i,]-x2))
>> > for (t in 1:ncol(Y))
>> > {
>> > res2[i,t] =
>> > ((1-(1+exp(-t(x3[,t])%*%X[i,]-x4[t]))^(-1))^(abs(Y[i,t]-Yb[i])))*(((1+exp(-t(x3[,t])%*%X[i,]-x4[t]))^(-1))^(1-abs(Y[i,t]-Yb[i])))
>> > log.res2[i,t]=log(res2[i,t])
>> > res2.b[i]=res2.b[i]+log.res2[i,t]
>> > }
>> > res3[i] = p_tilde[i]*log(res1[i])
>> > res3.b[i] = p_tilde[i]*(res2.b[i])
>> > }
>> > -(ncol(Y)*sum(res3)+sum(res3.b))
>> >
>> > }
>> >
>> > ##### Gradiants:
>> >
>> > gr=function(x)
>> > {
>> > x1 = x[1:ncol(X)]
>> > x2 = x[(ncol(X)+1)]
>> > x3 =
>> > matrix(x[(ncol(X)+2):(ncol(X)+1+ncol(X)*ncol(Y))],nrow=ncol(X),ncol=ncol(Y))
>> > x4 = x[(ncol(X)+1+ncol(X)*ncol(Y)+1):length(x)]
>> > gr1 = rep(0,ncol(X))
>> > gr4 = rep(0,ncol(Y))
>> > gr3 = matrix(0,nrow=ncol(X),ncol=ncol(Y))
>> > gr1.b = matrix(0,nrow=nrow(X),ncol=ncol(X))
>> > gr2.b = rep(0,nrow(X))
>> > eta = matrix(0,nrow=nrow(X),ncol=ncol(Y))
>> > d.eta.3 = array(0,dim=c(nrow(X),ncol(X),ncol(Y)))
>> > d.eta.4 = matrix(0,nrow=nrow(X),ncol=ncol(Y))
>> > gr3.b1 = array(0,dim=c(nrow(X),ncol(X),ncol(Y)))
>> > gr4.b1 = matrix(0,nrow=nrow(X),ncol=ncol(Y))
>> >
>> > #Gradiant of alpha and beta
>> >
>> > for (i in 1:nrow(X))
>> > {
>> > gr1.b[i,] =
>> > (2*p_tilde[i]-1)*((exp(-t(x1)%*%X[i,]-x2)*X[i,])/(1+exp(-t(x1)%*%X[i,]-x2))^2)
>> > gr2.b[i] =
>> > (2*p_tilde[i]-1)*((exp(-t(x1)%*%X[i,]-x2))/(1+exp(-t(x1)%*%X[i,]-x2))^2)
>> > }
>> > for (j in 1:ncol(X))
>> > {
>> > gr1[j] = sum(gr1.b[,j])
>> > }
>> > gr2 = sum(gr2.b)
>> >
>> >
>> > #Gradiant de w et gamma
>> > for (i in 1:nrow(X))
>> > {
>> > for (t in 1:ncol(Y))
>> > {
>> > eta[i,t] = 1/(1+exp(-t(x3[,t])%*%X[i,]-x4[t]))
>> > d.eta.3[i,,t] = eta[i,t]*(1-eta[i,t])*X[i,]
>> > d.eta.4[i,t] = eta[i,t]*(1-eta[i,t])
>> > gr3.b1[i,,t] =
>> > p_tilde[i]*((-abs(Y[i,t]-Yb[i]))*(1-eta[i,t])^(-1)+(1-abs(Y[i,t]-Yb[i]))*
>> > (eta[i,t])^(-1))*d.eta.3[i,,t]
>> > gr4.b1[i,t] =
>> > p_tilde[i]*((-abs(Y[i,t]-Yb[i]))*(1-eta[i,t])^(-1)+(1-abs(Y[i,t]-Yb[i]))*
>> > (eta[i,t])^(-1))*d.eta.4[i,t]
>> > }
>> > }
>> > for (t in 1:ncol(Y))
>> > {
>> > for (j in 1:ncol(X))
>> > {
>> > gr3[j,t] = sum(gr3.b1[,j,t])
>> > }
>> > gr4[t] = sum(gr4.b1[,t])
>> > }
>> > c(-gr1,-gr2,-gr3,-gr4)
>> > }
>> >
>> > opt = optim(c(alpha[,c+1],beta[c+1],w,gamma),g,gr,method="BFGS")
>> >
>> > The problem is that it gives me wrong results, and I have noticed that
>> > if I change my function to maximize (for example if, instead of
>> > -(ncol(Y)*sum(res3)+sum(res3.b)), I try to maximise -(ncol(Y)*sum(res3)), it
>> > gives me the exactly same results...which is not possible!
>> > So maybe I am using optim() in a wrong way...Does someone have an idea
>> > what could be wrong in my code ?
>> >
>> > Thank you very much in advance
>> >        [[alternative HTML version deleted]]
>> >
>> > ______________________________________________
>> > 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.
>> >
>>
>>
>>
>> --
>> ===============================================
>> Jon Daily
>> Technician
>> ===============================================
>> #!/usr/bin/env outside
>> # It's great, trust me.
>



-- 
===============================================
Jon Daily
Technician
===============================================
#!/usr/bin/env outside
# It's great, trust me.



More information about the R-help mailing list