为什么非线性规划 (NLP) 求解器 Rsolnp 中的目标函数不符合要求？

Question

案例：

我有 3 个区域（巴西、新西兰、美国）（我的实际问题要大得多 - 31 个区域）。这三个地区通过迁移联系在一起。例如，如果有 10 个人从巴西移居美国（BRZ-USA），我们有移民到美国（人员流入）和从巴西移民（人员流出）。我有给定宇宙中所有可能的迁移流的迁移率数据集（3 * 2 = 6）。此外，我还有每个地区的人口数据集。当我将迁移率乘以人口时，我得到了迁移人数。然后我可以计算每个地区的移民人数和移民人数。从移民中减去移民得出净移民数（可以是正数或负数）。然而，由于我们有一个平衡的系统（每个地区的流入量等于流出量），所有地区的净移民总和应该为零。除了净移民率和人口外，我还有来自每个地区假设的未来情景的净移民人数。但是场景中的净移民人数与我可以从我的数据中计算出的人数不同。因此，我想向上和向下缩放 6 个迁移率（通过添加或减去一个固定数字），以使最终的净迁移计数符合场景值。我使用非线性规划 (NLP) 求解器 Rsolnp 来完成此任务（请参阅下面的示例脚本）。

问题：

我以最小二乘方程的形式指定了目标函数，因为它的目标是迫使 6 个标量尽可能接近零。此外，我正在使用等式约束函数来满足场景值。这一切都很好，求解器提供了标量，我可以将这些标量添加到导致迁移计数的迁移率中，以完美匹配场景值（请参阅脚本部分“测试是否达到目标”）。但是，我还想对目标函数应用权重（变量：w），以便某些标量上的更高值受到更强的惩罚。但是，无论我如何指定权重，我总是获得相同的解决方案（请参阅“不同权重的示例结果”）。因此，求解器似乎不遵守目标函数。有谁知道为什么会这样以及我如何更改目标函数以便可以使用权重？非常感谢您的帮助！

library(Rsolnp)

# Regions
regUAll=c("BRZ","NZL","USA") # "BRZ"=Brazil; "NZL"=New Zealand; "USA"=United States

#* Generate unique combinations of regions
uCombi=expand.grid(regUAll,regUAll,stringsAsFactors=F) 
uCombi=uCombi[uCombi$Var1!=uCombi$Var2,] # remove same region combination (e.g., BRZ-BRZ)
uCombi=paste(uCombi$Var2,uCombi$Var1,sep="-")

#* Generate data frames
# Migration rates - rows represent major age groups (row1=0-25 years, row2=26-50 years, row3=51-75 years)
dfnm=data.frame(matrix(rep(c(0.01,0.04,0.02),length(uCombi)),ncol=length(uCombi),nrow=3),stringsAsFactors=F) # generate empty df
names(dfnm)=uCombi # assign variable names

# Population (number of people) in region of origin
pop=c(rep(c(20,40,10),2),rep(c(4,7,2),2),rep(c(30,70,50),2))
dfpop=data.frame(matrix(pop,ncol=length(uCombi),nrow=3),stringsAsFactors=F) # generate empty df
names(dfpop)=uCombi # assign variable names

#* Objective function for optimization
# Note: Least squares method to keep the additive scalers as close to 0 as possible
#       The sum expression allows for flexible numbers of scalars to be included but is identical to: w[1](scal[1]-0)^2+w[2](scal[2]-0)^2+w[3](scal[3]-0)^2+w[4](scal[4]-0)^2+w[5](scal[5]-0)^2+w[6](scal[6]-0)^2
f.main=function(scal,nScal,w,dfnm,dfpop,regUAll){
  sum(w*(scal[1:nScal]-0)^2)
}

#* Equality contraint function
f.equal=function(scal,nScal,w,dfnm,dfpop,regUAll){

  #* Adjust net migration rates by scalar
  for(s in 1:nScal){
    dfnm[,s]=dfnm[,s]+scal[s]
  }

  #* Compute migration population from data
  nmp=sapply(dfpop*dfnm,sum) # sums migration population across age groups

  nmd=numeric(length(regUAll)); names(nmd)=regUAll                        # generate named vector to be filled with values
  for(i in 1:length(regUAll)){
    colnEm=names(nmp)[grep(paste0("^",regUAll[i],"-.*"),names(nmp))]      # emigration columns
    colnIm=names(nmp)[grep(paste0("^.*","-",regUAll[i],"$"),names(nmp))]  # immigration columns
    nmd[regUAll[i]]=sum(nmp[colnIm])-sum(nmp[colnEm])                     # compute net migration population = immigration - emigration
  }
  nmd=nmd[1:(length(nmd)-1)] # remove the last equality constraint value - not needed because we have a closed system in which global net migration=0

  return(nmd)
}

#* Set optimization parameters
cpar2=list(delta=1,tol=1,outer.iter=10,trace=1) # optimizer settings
nScal=ncol(dfnm)                  # number of scalars to be used
initScal=rep(0,nScal)             # initial values of additive scalars
lowScal=rep(-1,nScal)             # lower bounds on scalars
highScal=rep(1,nScal)             # upper bounds on scalars
nms=c(-50,10)                     # target values: BRZ=-50, NZL=10, USA=40; last target value does not need to be included since we deal with a closed system in which global net migration sums to 0
w=c(1,1,1,1,1,1)    # unity weights
#w=c(1,1,2,2,1,1)    # double weight on NZL
#w=c(5,1,2,7,1,0.5)  # mixed weights


#* Perform optimization using solnp
solRes=solnp(initScal,fun=f.main,eqfun=f.equal,eqB=nms,LB=lowScal,UB=highScal,control=cpar2,
             nScal=nScal,w=w,dfnm=dfnm,dfpop=dfpop,regUAll=regUAll)

scalSol=solRes$pars # return optimized values of scalars

# Example results for different weights
#[1]  0.101645349  0.110108019 -0.018876993  0.001571639 -0.235945755 -0.018134294 # w=c(1,1,1,1,1,1)  
#[1]  0.101645349  0.110108019 -0.018876993  0.001571639 -0.235945755 -0.018134294 # w=c(1,1,2,2,1,1)
#[1]  0.101645349  0.110108019 -0.018876993  0.001571639 -0.235945755 -0.018134294 # w=c(5,1,2,7,1,0.5)

#*** Test if target was reached
# Adjust net migration rates using the optimized scalars
for(s in 1:nScal){  
  dfnm[,s]=dfnm[,s]+scalSol[s]
}

# Compute new migration population
nmp=sapply(dfpop*dfnm,sum) # sums migration population across age groups

nmd=numeric(length(regUAll)); names(nmd)=regUAll                        # generate named vector to be filled with values
for(i in 1:length(regUAll)){
  colnEm=names(nmp)[grep(paste0("^",regUAll[i],"-.*"),names(nmp))]      # emigration columns
  colnIm=names(nmp)[grep(paste0("^.*","-",regUAll[i],"$"),names(nmp))]  # immigration columns
  nmd[regUAll[i]]=sum(nmp[colnIm])-sum(nmp[colnEm])                     # compute net migration population = immigration - emigration
}

nmd # should be -50,10,40 if scalars work correctly

Answer 1

尝试与零不同的起始值 (initScal)； sum( w * 0^2) = 0 代表所有 w。