-2

I would like to find a maximum economic stress scenario restricted by a limit of the mahalanobis distance of this scenario. For this, I have to consider two functions in the optimization.

To make it easier, we can work with a simplifying problem: We have a simple linear model: y=a+bx. For this I want to minimize: sum(a+bx-y)^2. But also, I have for example the restriction that: (ab*5)/2<30.

To calculate this problem with the excel solver is not a problem. But, how I get this in r?

Klaus
  • 23
  • 3
  • 1
    what have you tried so far? SO is not a "please write code for me" site. You should look at making a good [MRE](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) – tbradley Oct 11 '17 at 12:41
  • Maybe post the steps of the Excel Solver and what have you *tried so far* and then and only then people might start contributing – amonk Oct 11 '17 at 12:42
  • 2
    StackOverflow is not a homework website please read https://stackoverflow.com/help/how-to-ask – F0XS Oct 11 '17 at 12:42
  • In all optimize functions in r like maxLik it is just possible to restrict parameters to an upper or an lower bound. For instance, like x1+x2+x3>=2. But unfortunatly, I do not find a function where I can restrict the parameters in a more sophisticated way like in the simple example in my question. It is enough, if somebody ever seen such a function and can give me a hint. – Klaus Oct 11 '17 at 13:10

1 Answers1

0

You could try to incorporate the constraint into the objective function, like this

# example data whose exact solution lies outside the constraint
x <- runif(100, 1, 10)
y <- 3 + 5*x + rnorm(100, mean=0, sd=.5) 

# big but not too big
bigConst <- sum(y^2) * 100

# if the variables lie outside the feasible region, add bigConst
f <- function(par, x, y) 
    sum((par["a"]+par["b"]*x-y)^2) + 
    if(par["a"]*par["b"]>12) bigConst else 0

# simulated annealing can deal with non-continous objective functions
sol <- optim(par=c(a=1, b=1), fn=f, method="SANN", x=x, y=y)

# this is how it looks like
plot(x,y)
abline(a=sol$par["a"], b=sol$par["b"])
Karsten W.
  • 17,826
  • 11
  • 69
  • 103