[R] optim seems to be finding a local minimum

[Dimitri Liakhovitski]

Hello!

I am trying to create an R optimization routine for a task that's
currently being done using Excel (lots of tables, formulas, and
Solver).
However, otpim seems to be finding a local minimum.
Example data, functions, and comparison with the solution found in
Excel are below.
I am not experienced in optimizations so thanks a lot for your advice!
Dimitri

### 2 Inputs:
IV<-data.frame(IV=c(0,6672895.687,13345791.37,20018687.06,26691582.75,33364478.44,40037374.12,46710269.81,53383165.5,60056061.18,66728956.87,73401852.56,80074748.24,86747643.93,93420539.62,100093435.3,106766331,113439226.7,120112122.4,126785018.1,133457913.7,140130809.4,146803705.1,153476600.8,160149496.5,166822392.2,173495287.9,180168183.5,186841079.2,193513974.9,200186870.6))
DV<-data.frame(DV=c(0,439.8839775,829.7360945,1176.968757,1487.732038,1767.147276,2019.49499,2248.366401,2456.78592,2647.310413,2822.109854,2983.033036,3131.661246,3269.352233,3397.276321,3516.446162,3627.741311,3731.928591,3829.679009,3921.581866,4008.156537,4089.862363,4167.106955,4240.253215,4309.625263,4375.513474,4438.178766,4497.856259,4554.75841,4609.077705,4660.988983))

## Function "transformIV" transforms a data frame column "IV" using
parameters .alpha & .beta:
## It returns a data frame column IV_transf:
transformIV = function(.alpha,.beta) {
 IV_transf <- as.data.frame(1 - (1/exp((IV/.beta)^.alpha)))
 return(IV_transf)
}

### Function "mysum" calculates the sum of absolute residuals after a
regression with a single predictor:
mysum<- function(myIV,myDV){
  regr<-lm(myDV[[1]] ~ 0 + myIV[[1]])
  mysum<-sum(abs(regr$resid))
  return(mysum)
}

### Function to be optimized;
### param is a vector of 2 values (.alpha and .beta)
myfunc <- function(param){
  myalpha<-param[1]
  mybeta<-param[2]
  IVtransf<-transformIV(myalpha, mybeta)
  sumofdevs<-mysum(myIV=IVtransf,myDV=DV)
  return(sumofdevs)
}

# Optimizing using optim:
myopt <- optim(fn=myfunc, par=c(0.1,max(IV)), method="L-BFGS-B", lower=0)
(myopt)
myfunc(myopt$par)
## Comparing this solution to Excel Solver solution:
myfunc(c(0.888452533990788,94812732.0897449))

-- 
Dimitri Liakhovitski
marketfusionanalytics.com