[R] constrOptim and further arguments

Uwe Ligges ligges at statistik.tu-dortmund.de
Thu Dec 22 09:40:09 CET 2011



On 21.12.2011 18:08, Michael Griffiths wrote:
> Dear List,
>
> I have the code below, where I am using the constrained optimisation
> package, 'constrOptim.nl' to find the values of two values, b0 and b1.
>
> I have no problems when I enter further variable information DIRECTLY into
> the functions, fn, and heq. In this instance I require fn to have -0.0075
> appended to it, and in the case of heq, h[1] has -0.2.
>
>
> library(alabama)
>
> fn<-function(x)
> (((1/(1+exp(-x[1]+x[2]))+(1/(1+exp(-x[1])))+(1/(1+exp(-x[1]-x[2])))))/3)-0.0075
>
> heq<- function(x) {
> h<- rep(NA,2)
> h[1]<- (1-(1/(1+exp(-x[1]))/(1/(1+exp(-x[1]-x[2])))))-0.2    #Constraint
> 1: diff = 20%
> h[2]<- -log((1-0.0075)/0.0075)-x[1]                #Constraint 2: set
> b0
> h
> }
>
> p0<-c(-4,0.3)
> ans<-constrOptim.nl(par=p0,fn=fn,heq=heq)
>
> This code does indeed work as expected and the results have been verified
> elsewhere.
>
> However, I do have problems when I wish to enter this information
> indirectly via the constrOptim function, e.g.
>
> fn<-function(x,RR)
> (((1/(1+exp(-x[1]+x[2]))+(1/(1+exp(-x[1])))+(1/(1+exp(-x[1]-x[2])))))/3)-RR
>
> heq<- function(x,Diff) {
> h<- rep(NA,2)
> h[1]<- (1-(1/(1+exp(-x[1]))/(1/(1+exp(-x[1]-x[2])))))-Diff    #Constraint
> 1: diff = 20%
> h[2]<- -log((1-RR)/RR)-x[1]
> #Constraint 2: set b0
> h
> }
>
> p0<-c(-4,0.3)
> ans<-constrOptim.nl(par=p0,fn=fn,heq=heq, RR=0.0075, Diff=0.2)
>
>
> In this instance I get the following error,
>
> 'Error in heq(theta, ...) : unused argument(s) (RR = 0.0075)'
>
> I do not appear to be able to input further arguments in this manner, any
> help from the list would be gratefully received.


Reason is that you need to specify fn and heq in a way that will accept 
these arguments!

So the defitionions should start with

fn <- function(x,RR,Diff)

and

heq <- function(x,RR,Diff)

respectively.

Uwe Ligges





>
> Regards
>
> Mike Griffiths
>
>
>
> On Tue, Dec 20, 2011 at 11:53 AM, Michael Griffiths<
> griffiths at upstreamsystems.com>  wrote:
>
>> Dear List,
>>
>> I am using constrOptim to solve the following
>>
>> fr1<- function(x) {
>>      b0<- x[1]
>>      b1<- x[2]
>>      ((1/(1+exp(-b0+b1))+(1/(1+exp(-b0)))+(1/(1+exp(-b0-b1)))))/3
>> }
>>
>> As you can see, my objective function is
>> ((1/(1+exp(-b0+b1))+(1/(1+exp(-b0)))+(1/(1+exp(-b0-b1)))))/3 and I would
>> like to solve for both b0 and b1.
>>
>> If I were to use optim then I would derive the gradient of the function
>> (grr) as follows:
>>
>> fr2<-
>> expression(((1/(1+exp(-b0+b1))+(1/(1+exp(-b0)))+(1/(1+exp(-b0-b1)))))/3)
>> grr<- deriv(fr2,c("b0","b1"), func=TRUE)
>>
>> and then simply use optim via
>>
>> optim(c(-5.2,0.22), fr1, grr)
>>
>> My problem is that I wish to place constraints (b0>=-0.2 and b1>= 0.1)
>> upon the values of b0 and b1. I can set the constraints matrix and boundary
>> values to
>>
>> ui=rbind(c(1,0),c(0,1)) and ci=c(-0.2,0.1), however, when I come to run
>> constrOptim function via
>>
>>
>> constrOptim(c(-0.1,0.2), fr1, grr, ui=rbind(c(1,0),c(0,1)), ci=c(-0.2,0.1))
>>
>> I get the following error message:
>>
>> "Error in .expr1 + b1 : 'b1' is missing"
>>
>> So, it seems to me that I am doing something incorrectly in my
>> specification of grr in constrOptim.
>>
>> If the list could help, I would be most grateful.
>>
>> Regards
>>
>> Mike Griffiths
>>
>>
>>
>> --
>>
>> *Michael Griffiths, Ph.D
>> *Statistician
>>
>> *Upstream Systems*
>>
>> 8th Floor
>> Portland House
>> Bressenden Place
>> SW1E 5BH
>>
>> <http://www.google.com/url?q=http%3A%2F%2Fwww.upstreamsystems.com%2F&sa=D&sntz=1&usg=AFrqEzfKYfaAalqvahwrpywpJDL9DxUmWw>
>>
>> Tel   +44 (0) 20 7869 5147
>> Fax  +44 207 290 1321
>> Mob +44 789 4944 145
>>
>> www.upstreamsystems.com<http://www.google.com/url?q=http%3A%2F%2Fwww.upstreamsystems.com%2F&sa=D&sntz=1&usg=AFrqEzfKYfaAalqvahwrpywpJDL9DxUmWw>
>>
>> *griffiths at upstreamsystems.com<einstein at upstreamsystems.com>*
>>
>> <http://www.upstreamsystems.com/>
>>
>
>
>



More information about the R-help mailing list