Сравните paretosearch и gamultiobj

Скрипт Open Live Script

Этот пример показывает, как создать набор точек на передней стороне Парето, использующей и paretosearch и gamultiobj. Целевая функция имеет две цели и двумерную контрольную переменную x. mymulti3 целевой функции доступен на вашем сеансе MATLAB®, когда вы нажимаете кнопку, чтобы отредактировать или попробовать этот пример. Также скопируйте код mymulti3 в свой сеанс. Для скорости вычисления векторизована функция.

type mymulti3
function f = mymulti3(x)
%
f(:,1) = x(:,1).^4 + x(:,2).^4 + x(:,1).*x(:,2) - (x(:,1).^2).*(x(:,2).^2) - 9*x(:,1).^2;
f(:,2) = x(:,2).^4 + x(:,1).^4 + x(:,1).*x(:,2) - (x(:,1).^2).*(x(:,2).^2) + 3*x(:,2).^3;

Основной пример и графики

Найдите Множества Парето для целевых функций с помощью paretosearch и gamultiobj. Установите опцию UseVectorized на true для добавленной скорости. Включайте функцию построения графика, чтобы визуализировать Множество Парето.

rng default
nvars = 2;
opts = optimoptions(@gamultiobj,'UseVectorized',true,'PlotFcn','gaplotpareto');
[xga,fvalga,~,gaoutput] = gamultiobj(@(x)mymulti3(x),nvars,[],[],[],[],[],[],[],opts);

Optimization terminated: average change in the spread of Pareto solutions less than options.FunctionTolerance.

optsp = optimoptions('paretosearch','UseVectorized',true,'PlotFcn',{'psplotparetof' 'psplotparetox'});
[xp,fvalp,~,psoutput] = paretosearch(@(x)mymulti3(x),nvars,[],[],[],[],[],[],[],optsp);

Pareto set found that satisfies the constraints. 

Optimization completed because the relative change in the volume of the Pareto set 
is less than 'options.ParetoSetChangeTolerance' and constraints are satisfied to within 
'options.ConstraintTolerance'.

Вычислите теоретически точные места на передней стороне Парето при помощи mymulti4. Функция mymulti4 доступна на вашем сеансе MATLAB®, когда вы нажимаете кнопку, чтобы отредактировать или попробовать этот пример.

type mymulti4
function mout = mymulti4(x)
%
gg = [4*x(1)^3+x(2)-2*x(1)*(x(2)^2) - 18*x(1);
    x(1)+4*x(2)^3-2*(x(1)^2)*x(2)];
gf = gg + [18*x(1);9*x(2)^2];

mout = gf(1)*gg(2) - gf(2)*gg(1);

Функция mymulti4 оценивает градиенты этих двух целевых функций. Затем, для области значений значений x(2), используйте fzero, чтобы определить местоположение точки, где градиенты точно параллельны, который является где вывод mout = 0.

a = [fzero(@(t)mymulti4([t,-3.15]),[2,3]),-3.15];
for jj = linspace(-3.125,-1.89,50)
    a = [a;[fzero(@(t)mymulti4([t,jj]),[2,3]),jj]];
end
figure
plot(fvalp(:,1),fvalp(:,2),'bo');
hold on
fs = mymulti3(a);
plot(fvalga(:,1),fvalga(:,2),'r*');
plot(fs(:,1),fs(:,2),'k.')
legend('Paretosearch','Gamultiobj','True')
xlabel('Objective 1')
ylabel('Objective 2')
hold off

gamultiobj находит точки с немного более широким распространением на пробеле целевой функции. Постройте решения на пробеле переменной решения, наряду с теоретической оптимальной кривой Парето и контурным графиком этих двух целевых функций.

[x,y] = meshgrid(1.9:.01:3.1,-3.2:.01:-1.8);
mydata = mymulti3([x(:),y(:)]);
myff = sqrt(mydata(:,1) + 39);% Spaces the contours better
mygg = sqrt(mydata(:,2) + 28);% Spaces the contours better
myff = reshape(myff,size(x));
mygg = reshape(mygg,size(x));

figure;
hold on
contour(x,y,mygg,50)
contour(x,y,myff,50)
plot(xp(:,1),xp(:,2),'bo')
plot(xga(:,1),xga(:,2),'r*')
plot(a(:,1),a(:,2),'-k')
xlabel('x(1)')
ylabel('x(2)')
hold off

В отличие от решения paretosearch, решение gamultiobj имеет точки в экстремальных концах области значений на пробеле целевой функции. Однако решение paretosearch имеет больше точек, которые ближе к истинному решению и на пробеле целевой функции и на пробеле переменной решения. Число точек на передней стороне Парето отличается для каждого решателя, когда вы используете опции по умолчанию.

Переключенная проблема

Что происходит, если решение вашей проблемы имеет контрольные переменные, которые являются большими? Исследуйте этот случай путем сдвига проблемных переменных. Для неограниченной проблемы может перестать работать gamultiobj, в то время как paretosearch более устойчив к таким сдвигам.

Для более легкого сравнения задайте 35 точек на передней стороне Парето для каждого решателя.

shift = [20,-30];
fun = @(x)mymulti3(x+shift);
opts.PopulationSize = 100; % opts.ParetoFraction = 35
[xgash,fvalgash,~,gashoutput] = gamultiobj(fun,nvars,[],[],[],[],[],[],opts);

Optimization terminated: average change in the spread of Pareto solutions less than options.FunctionTolerance.

gamultiobj не удается найти полезное Множество Парето.

optsp.PlotFcn = [];
optsp.ParetoSetSize = 35;
[xpsh,fvalpsh,~,pshoutput] = paretosearch(fun,nvars,[],[],[],[],[],[],[],optsp);
Pareto set found that satisfies the constraints. 

Optimization completed because the relative change in the volume of the Pareto set 
is less than 'options.ParetoSetChangeTolerance' and constraints are satisfied to within 
'options.ConstraintTolerance'.

figure
plot(fvalpsh(:,1),fvalpsh(:,2),'bo');
hold on
plot(fs(:,1),fs(:,2),'k.')
legend('Paretosearch','True')
xlabel('Objective 1')
ylabel('Objective 2')
hold off

paretosearch находит распространение точек решения равномерно почти целой возможной областью значений.

Добавление границ, даже довольно свободных единиц, помогает и gamultiobj и paretosearch найти соответствующие решения. Установите нижние границы -50 в каждом компоненте и верхние границы 50.

opts.PlotFcn = [];
optsp.PlotFcn = [];
lb = [-50,-50];
ub = -lb;
[xgash,fvalgash,~,gashoutput] = gamultiobj(fun,nvars,[],[],[],[],lb,ub,opts);
Optimization terminated: average change in the spread of Pareto solutions less than options.FunctionTolerance.

[xpsh2,fvalpsh2,~,pshoutput2] = paretosearch(fun,nvars,[],[],[],[],lb,ub,[],optsp);
Pareto set found that satisfies the constraints. 

Optimization completed because the relative change in the volume of the Pareto set 
is less than 'options.ParetoSetChangeTolerance' and constraints are satisfied to within 
'options.ConstraintTolerance'.

figure
plot(fvalpsh2(:,1),fvalpsh2(:,2),'bo');
hold on
plot(fvalgash(:,1),fvalgash(:,2),'r*');
plot(fs(:,1),fs(:,2),'k.')
legend('Paretosearch','Gamultiobj','True')
xlabel('Objective 1')
ylabel('Objective 2')
hold off

В этом случае оба решателя находят хорошие решения.

Запустите paretosearch с решения gamultiobj

Получите подобную область значений решений от решателей путем запуска paretosearch с решения gamultiobj.

optsp.InitialPoints = xgash;
[xpsh3,fvalpsh3,~,pshoutput3] = paretosearch(fun,nvars,[],[],[],[],lb,ub,[],optsp);
Pareto set found that satisfies the constraints. 

Optimization completed because the relative change in the volume of the Pareto set 
is less than 'options.ParetoSetChangeTolerance' and constraints are satisfied to within 
'options.ConstraintTolerance'.

figure
plot(fvalpsh3(:,1),fvalpsh3(:,2),'bo');
hold on
plot(fvalgash(:,1),fvalgash(:,2),'r*');
plot(fs(:,1),fs(:,2),'k.')
legend('Paretosearch','Gamultiobj','True')
xlabel('Objective 1')
ylabel('Objective 2')
hold off

Теперь решение paretosearch подобно решению gamultiobj, несмотря на то, что некоторые точки решения отличаются.

Запустите с Одно-объективных решений

Другой способ получить хорошее решение состоит в том, чтобы начать с точек, которые минимизируют каждую целевую функцию отдельно.

От многоцелевой функции создайте одно целевую функцию, которая выбирает каждую цель в свою очередь. Используйте переключенную функцию от предыдущего раздела. Поскольку вы даете точки хорошего начала решателям, вы не должны задавать границы.

nobj = 2; % Number of objectives

x0 = -shift; % Initial point for single-objective minimization
uncmin = cell(nobj,1); % Cell array to hold the single-objective minima
allfuns = zeros(nobj,2); % Hold the objective function values
eflag = zeros(nobj,1);
fopts = optimoptions('patternsearch','Display','off'); % Use an appropriate solver here
for i = 1:nobj
    indi = zeros(nobj,1); % Choose the objective to minimize
    indi(i) = 1;
    funi = @(x)dot(fun(x),indi);
    [uncmin{i},~,eflag(i)] = patternsearch(funi,x0,[],[],[],[],[],[],[],fopts); % Minimize objective i
    allfuns(i,:) = fun(uncmin{i});
end
uncmin = cell2mat(uncmin); % Matrix of start points

Запустите paretosearch с одно-объективных минимальных точек и обратите внимание, что он имеет полный спектр в своих решениях. paretosearch добавляет, что случайная начальная буква указывает на предоставленные единицы в порядке иметь генеральную совокупность, по крайней мере, людей options.ParetoSetSize. Точно так же gamultiobj добавляет случайные точки в предоставленные единицы, чтобы получить генеральную совокупность, по крайней мере, людей (options.PopulationSize)*(options.ParetoFraction).

optsp = optimoptions(optsp,'InitialPoints',uncmin);
[xpinit,fvalpinit,~,outputpinit] = paretosearch(fun,nvars,[],[],[],[],[],[],[],optsp);
Pareto set found that satisfies the constraints. 

Optimization completed because the relative change in the volume of the Pareto set 
is less than 'options.ParetoSetChangeTolerance' and constraints are satisfied to within 
'options.ConstraintTolerance'.

Теперь решите проблему с помощью gamultiobj, начинающего с начальных точек.

opts = optimoptions(opts,'InitialPopulationMatrix',uncmin);
[xgash2,fvalgash2,~,gashoutput2] = gamultiobj(fun,nvars,[],[],[],[],[],[],opts);
Optimization terminated: average change in the spread of Pareto solutions less than options.FunctionTolerance.

figure
plot(fvalpinit(:,1),fvalpinit(:,2),'bo');
hold on
plot(fvalgash2(:,1),fvalgash2(:,2),'r*');
plot(fs(:,1),fs(:,2),'k.')
plot(allfuns(:,1),allfuns(:,2),'gs','MarkerSize',12)
legend('Paretosearch','Gamultiobj','True','Start Points')
xlabel('Objective 1')
ylabel('Objective 2')
hold off

Оба решателя заполняют переднюю сторону Парето между экстремальными точками с довольно точными и хорошо распределенными решениями.

Просмотрите конечные пункты на пробеле переменной решения.

figure;
hold on
xx = x - shift(1);
yy = y - shift(2);
contour(xx,yy,mygg,50)
contour(xx,yy,myff,50)
plot(xpinit(:,1),xpinit(:,2),'bo')
plot(xgash2(:,1),xgash2(:,2),'r*')
ashift = a - shift;
plot(ashift(:,1),ashift(:,2),'-k')
plot(uncmin(:,1),uncmin(:,2),'gs','MarkerSize',12);
xlabel('x(1)')
ylabel('x(2)')
hold off

Смотрите также

|

Похожие темы

Документация Global Optimization Toolbox
Поддержка
Для просмотра документации необходимо авторизоваться на сайте