global search with fmincon returns [] when a particular constraint is added

3 次查看（过去 30 天）
Mic
Mic 2016-1-18
编辑： Matt J 2016-1-19
Hello, I'm running a global search algorithm to determine the constraint optimum of a loglikelihood function. I'm repeating the same optimization for all k < R, in order to populate the matrix xfinals
% Create the control matrix X with all the required control variables:
jj = ones((C-2),1); % vector of constants
Xscale = [ jj, matrix(:,4),matrix(:,5),matrix(:,6),matrix(:,7),matrix(:,8),matrix(:,9),matrix(:,11),matrix(:,12),matrix(:,13),matrix(:,14),matrix(:,15),matrix(:,16),matrix(:,18), matrix(:,19),matrix(:,20), matrix(:,23), matrix(:,24), matrix(:,25), matrix(:,26) ]; % Fixed effect for scale parameter lambda
Xshape = [ jj, matrix(:,4),matrix(:,15),matrix(:,16),matrix(:,23), matrix(:,24), matrix(:,25), matrix(:,27) ]; % Fixed effect for shape parameter k
% Fixed effect for scale parameter lambda
%Xscale = [ jj, matrix(:,4), matrix(:,5),matrix(:,7), matrix(:,9), matrix(:,11), matrix(:,12), matrix(:,13), matrix(:,14), matrix(:,15),matrix(:,16),matrix(:,18), matrix(:,19), matrix(:,20), matrix(:,21),matrix(:,24), matrix(:,23), matrix(:,26),matrix(:,27) ];
% Fixed effect for shape parameter k
%Xshape = [ jj, matrix(:,14), matrix(:,19), matrix(:,20), matrix(:,21),matrix(:,24), matrix(:,23), matrix(:,26),matrix(:,27) ];
[RXscale, CXscale] = size(Xscale);
[RXshape, CXshape] = size(Xshape);
% Total number of unknown per alpha, beta combination:
unknowns = CXscale + CXshape + 1;
% The initial points (theta0) is a matrix with 3 rows (one per parameter)
% and a number of columns equal to the number the number of variables in X.
theta0 = rand(unknowns,1) * .00005;%*.000001; % initial points (column vector) * 0.001
Ae = [];
be = [];
Aeq = [];
beq = [];
lb = - ones(1, unknowns) * 100;
ub = ones(1, unknowns) * 100;
for k=1:R
if count_min_v(k,1) > 80
xfinals(k,:) = zeros(1,unknowns);
loglikelihood(k,1) = 0;
exitreason(k,1) = 6;
else
% make sure to transpose the vector of bids to be sent to the MLE
% function
% Set MaxIter and MaxfunEvals to 1000
values = (v_pos(k,1:(C-2)))';
f = @(theta) fun_MLE_controls( values, theta, Xscale, CXscale, Xshape, CXshape );
opts = optimoptions(@fmincon,'Display','final-detailed','MaxIter', 20000,...
'MaxfunEvals',20000,'TolX',1e-18,'TolFun',1e-18,'TolCon',1e-18,'Algorithm','interior-point','Hessian','bfgs');
mycon = @(theta) constraint( values, theta, Xscale, CXscale, Xshape, CXshape ) ;
problem = createOptimProblem('fmincon','objective',f,'x0',theta0,'lb',lb,'ub',ub,'nonlcon',mycon,'options',opts) ;
gs = GlobalSearch('Display','iter','TolFun',1e-18,'TolX',1e-18,'StartPointsToRun','bounds-ineqs');
[xfinal, fval, exitflag, output] = run(gs,problem);
%[xfinal, fval, exitflag] = fmincon(f,theta0,Ae,be,Aeq,beq,lb,ub,mycon,options); % no nonlinear constraint []
xfinals(k,:) = xfinal;
exitreason(k,1) = exitflag;
loglikelihood(k,1) = fval; %fun_MLE_controls( values, xfinals(k,:), Xscale, CXscale, Xshape, CXshape );
end
end
Since I want the expression in c(2) and c(3) to be larger than zero for all values, I tell the code to make the smallest value of the vector operation larger than zero. If the smallest value is larger than zero, then all the values in the vector are lareger than zero.
function [ c,ceq ] = constraint( v, parameters, Xscale, CXscale, Xshape, CXshape )
% LAMBDA and K as in fun_MLE_controls must be larger than zero
% Define the Gamma-Weibull parameters
lambda = parameters(1:CXscale,1); % scale
k = parameters((CXscale+1):(CXshape+CXscale),1); % shape
zi = parameters((CXscale+CXshape+1),1); % unobserved heterogeneity
LAMBDA = exp(sum(Xscale * lambda));
K = exp(sum(Xshape * k));
lk = v.^K; % lk = v_i^k
sum_k = sum(lk);
THETA = LAMBDA^(-K);
c(1) = - ((K + zi)/K);
pdf = K * THETA^(1+(zi/K)) * (v.^(zi+K-1)) .* exp(-THETA * v.^K) / gamma(1+(zi/K));
c(2) = -min(pdf);
c(3) = -min((v.^(-1)).*((zi + K - 1) - K .* THETA .* (v.^K)) + ...
(1./(gammainc(THETA .* v.^K,(1 + (zi/K)),'upper'))));
ceq = [];
I do not report the MLE function, but there is no problem with that because the code works fine without constraint c(3). When I include c(3) xfinal = [] and the program stops with the following message:
Subscripted assignment dimension mismatch.
Error in fioretti_auction_Charity_controls_GlobalSearch (line 270)
xfinals(k,:) = xfinal;
Could anyone help me please?

请先登录，再回答此问题。

回答（1 个）

Matt J
Matt J 2016-1-18
You need to decide what goes in xfinals when no feasible solution is found. Something like this would make sense:
if ~isempty(xfinal)
xfinals(k,:) = xfinal;
else
xfinals(k,:)=nan;
end
  2 个评论
Mic
Mic 2016-1-18
Thank you Matt for your input, but I do not understand why it returns []. It should at least return the starting points right? Also, the problem is in the "gammainc" function, because if I remove that function the code goes through without problem.
Matt J
Matt J 2016-1-18
编辑：Matt J 2016-1-19
Thank you Matt for your input, but I do not understand why it returns []. It should at least return the starting points right?
What exitflag value is returned in the cases where xfinal=[]? That will help us interpret what is being returned for xfinal.
From the documentation on GlobalSearch, it doesn't look to me that a solution will be accepted if fmincon does not return a positive exitflag, indicating convergence. Thus, it stands to reason that if no point satisfying the constraints exists, the set of solutions returned by run() could be empty.
Also, the problem is in the "gammainc" function, because if I remove that function the code goes through without problem.
That suggests that you are over-constraining your problem, i.e., the addition of the gammainc constraint makes the whole problem infeasible. But, we need the exitflag to know what's really going on.
Another thing I now notice is that your tolerances TolCon, TolFun, etc... have been set very tight (1e-18). That could make it very difficult to achieve the convergence criteria, which has the same effect as if the problem is infeasible.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Solver-Based Nonlinear Optimization 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by