2012-03-12 41 views
2

我正在嘗試設置一個優化腳本,該腳本將查看一組模型,將曲線擬合到模型,然後對它們進行優化,並受到一些參數的限制。限制R中的優化問題

從本質上講,我的收入是成本的函數,在一個遞減函數中,我有多個投資組合,比如說4或5.作爲一個輸入,我有成本和收入數字,以設定的增量。我想要做的就是對收益= A *成本^ B的投資組合曲線進行擬合,然後對不同投資組合進行優化,以找出每個投資組合之間針對集合預算的最優成本分配。下面的代碼(我爲它的不雅感到歉意,我確信有很多改進要做!)基本上讀取我的數據,在這種情況下,模擬,創建必要的數據幀(這是可能出現我的不雅行爲時),計算每個模擬曲線的必要變量,並生成圖形以檢查數據的擬合曲線。

我的問題是,現在我有如下形式的5條曲線:

收入= A *成本^ B(不同的A,B和每個功能的成本)

而且我想知道,給出的5個變量,應該怎麼拆我的成本兩者之間,所以我要到5條曲線主題的總和優化到

成本< =預算

我知道我需要使用constrOptim,但我花了幾個小時把我的頭撞到我的桌子上(實際上是幾個小時,而不是字面上敲我的頭......),我仍然無法弄清楚如何設置功能這樣它就可以最大限度地增加收入,但會受到成本的限制...

在這裏的任何幫助將不勝感激,這一直困擾着我幾個星期。

謝謝!

豐富

## clear all previous data 

rm(list=ls()) 
detach() 
objects() 

library(base) 
library(stats) 

## read in data 

sim<-read.table("input19072011.txt",header=TRUE) 
sim2<-data.frame(sim$Wrevenue,sim$Cost) 

## identify how many simulations there are - here you can change the 20 to the number of steps but all simulations must have the same number of steps 

portfolios<-(length(sim2$sim.Cost)/20) 

## create a matrix to input the variables into 

a<-rep(1,portfolios) 
b<-rep(2,portfolios) 
matrix<-data.frame(a,b) 

## create dummy vector to hold the revenue predictions 

k<-1 
j<-20 

for(i in 1:portfolios){ 

test<-sim2[k:j,] 

rev9<-test[,1] 
cost9<-test[,2] 

ds<-data.frame(rev9,cost9) 

rhs<-function(cost, b0, b1){ 
b0 * cost^b1 

m<- nls(rev9 ~ rhs(cost9, intercept, power), data = ds, start = list(intercept = 5,power = 1)) 

matrix[i,1]<-summary(m)$coefficients[1] 
matrix[i,2]<-summary(m)$coefficients[2] 

k<-k+20 
j<-j+20 

} 

## now there exists a matrix of all of the variables for the curves to optimise 

matrix 
multiples<-matrix[,1] 
powers<-matrix[,2] 
coststarts<-rep(0,portfolios) 

## check accuracy of curves 

k<-1 
j<-20 

for(i in 1:portfolios){ 

dev.new() 

plot(sim$Wrevenue[k:j]) 
lines(multiples[i]*(sim$Cost[k:j]^powers[i])) 

k<-k+20 
j<-j+20 

} 
+0

相關問題:http://stackoverflow.com/questions/9592369/in-r-how-do-i-find-the-optimal-variable-to-maximize-or-minimize-correlation-bet/9593809#9593809 – 2012-03-12 12:17:53

回答

5

如果你想找到 值cost[1],...,cost[5] 最大化revenue[1]+...+revenue[5] 受到限制cost[1]+...+cost[5]<=budget (和0 <= cost[i] <= budget), 可以參數化的一組可行解 如下

cost[1] = s(x[1]) * budget 
cost[2] = s(x[2]) * (budget - cost[1]) 
cost[3] = s(x[3]) * (budget - cost[1] - cost[2]) 
cost[4] = s(x[4]) * (budget - cost[1] - cost[2] - cost[3]) 
cost[5] = budget - cost[1] - cost[2] - cost[3] - cost[4] 

其中x[1],...,x[4]是查找 (對它們沒有約束)的參數 和s是實線R和段(0,1)之間的任何雙射。

# Sample data 
a <- rlnorm(5) 
b <- rlnorm(5) 
budget <- rlnorm(1) 

# Reparametrization 
s <- function(x) exp(x)/(1 + exp(x)) 
cost <- function(x) { 
    cost <- rep(NA,5) 
    cost[1] = s(x[1]) * budget 
    cost[2] = s(x[2]) * (budget - cost[1]) 
    cost[3] = s(x[3]) * (budget - cost[1] - cost[2]) 
    cost[4] = s(x[4]) * (budget - cost[1] - cost[2] - cost[3]) 
    cost[5] = budget - cost[1] - cost[2] - cost[3] - cost[4] 
    cost 
} 

# Function to maximize 
f <- function(x) { 
    result <- sum(a * cost(x)^b) 
    cat(result, "\n") 
    result 
} 

# Optimization 
r <- optim(c(0,0,0,0), f, control=list(fnscale=-1)) 
cost(r$par) 
+0

謝謝!我將在明天早上花費大量時間試着理解你寫的內容,如果我有任何問題(我不可避免)會回覆你。但非常感謝您的幫助! – 2012-03-12 17:38:20

+0

嗨,謝謝,這真的很有用,它現在可以工作。我現在唯一想做的工作是根據我有多少個投資組合(可以是2到100)進行優化。但我相信我會弄清楚:)。謝謝! – 2012-03-15 13:56:20

+0

你看過LIM包嗎? http://cran.r-project.org/package=LIM – user2030503 2014-03-11 17:43:04