[R] optim

[Remigijus Lapinskas]

Sunday, March 02, 2003, 9:48:26 PM, you wrote:

MK> ripley at stats.ox.ac.uk wrote:
>> You are not solving the same problem, which is more than a little `odd'.

MK> Why?

MK>  >>> want to choose the three parameters so that
MK> these three integrals are as close to, resp., 2300, 4600 and 5800 as
MK> possible. As I have three equations with three unknowns, I expect the
MK> exact fit, i.e., the SS (see below) to be zero.

Good evening to all!

In fact, I want to estimate three unknown parameters TETA[1],TETA[2]
(the scale parameter) and TETA[3]. What I know is 10 numbers
aa <- c(2300,4600,5800,
7100,...,37700) and these are the empirical counterparts of
the integrals over [0,0.1],...[0.9,1]. I don't know whether this is
correct but in order to find the parameters I minimize the
sum((integr-aa)^2) (kind of the "moment method"). As I had problems
with optim, I reduced the problem to three intervals.

I have to admit, I'm still fighting with optim. Following Prof
Ripley's suggestion, I replaced

res <- optim(init,LSS,aa=aa,method = "L-BFGS-B",lower=c(0,0,0.5))

by

res <- optim(init,LSS,aa=aa,method = "L-BFGS-B",lower=c(0,0,0.5),
control=list(parscale=c(1,1000,1)))

Now I have

> source("C:/Program Files/R/integral.R")
[1]     2.50      7000.00         0.84       # initial
[1]     2.052587 66734.476822    42.110597   # final
[1] 42820.43

instead of what I had earlier

>[1]    2.50      7000.00         0.84        # initial
>[1]    2.3487221 6999.9999823    0.5623628   # final
>[1] 75613.05

Now SS=42820, but it is still far from zero. I agree, I did not
implemented all what I was advised to, but it could take some time.

Many thanks for the help.
Remigijus

***************************************************
***************************************************


>Do read the help page.  It says:
>
> [...]
>
>`parscale' A vector of scaling values for the parameters.
>          Optimization is performed on `par/parscale' and these should
>          be comparable in the sense that a unit change in any element
>          produces about a unit change in the scaled value.

******************************************************
******************************************************

>>Dear all,
>>
>>I have a function MYFUN which depends on 3 positive parameters TETA[1],
>>TETA[2], and TETA[3]; x belongs to [0,1].
>>I integrate the function over [0,0.1], [0.1,0.2] and
>>[0.2,0.3] and want to choose the three parameters so that
>>these three integrals are as close to, resp., 2300, 4600 and 5800 as
>>possible. As I have three equations with three unknowns, I expect the
>>exact fit, i.e., the SS (see below) to be zero. However, the optim
>>function never gives me what I expect, the minimal SS value(=res$value)
>>never comes close to zero, the estimates of the parameters, res$par,
>>wildly depends on init etc.
>>I would be grateful for any comments on this miserable situation.
>>
>>aa <- c(2300,4600,5800)
>>init <- c(2.5,8000,0.84) # initial values of parameters
>>print(init)
>>###################
>>myfun <- function(x,TETA) TETA[2]*(((1-x)^(-1/TETA[3]))-
>>1)^(1/TETA[1])
>>###################
>>x <- seq(0,0.3,by=0.01)
>>plot(x,myfun(x,init),type="l")
>>###################
>>LSS <- function(teta,aa)
>>{
>>integr <- numeric(3)
>>   for(i in 1:3)
>>   {integr[i] <- 10*integrate(myfun,
>>   lower=(i-1)/10,upper=i/10,TETA=teta)$value
>>   }
>>SS <-  # SS=Sum of Squares
>>SS
>>}
>>####################
>>res <- optim(init,LSS,aa=aa,
>>method = "L-BFGS-B",lower=c(0,0,0.5))
>>print(res$par)
>>print(res$value)
>>
>>
>>
>>>source("C:/Program Files/R/integral.R")
>>
>>[1]    2.50      7000.00         0.84        # initial
>>[1]    2.3487221 6999.9999823    0.5623628   # final
>>[1] 75613.05                                 # minSS
>>
>>>source("C:/Program Files/R/integral.R")
>>
>>[1]     2.5      15000            0.84       # initial
>>[1]     2.125804 14999.999747     2.241179   # final
>>[1] 50066.35                                 # minSS
>>