罗兰关于MATLAB的艺术

将想法转化为MATLAB

求最优值

著者罗兰美国舒尔，

274年ビュー(過去30日間)| いいね| 13コメント
你曾经需要解决一个有局部极小值的优化问题吗?你用什么策略来解决它，试图找到“最佳”答案?今天我要讲一个简单的策略，在 全局优化工具箱 ．

解决一个简单的问题

或者至少让我们试试。我有一些数据，我想用一条特定形式的曲线来拟合它。首先让我们看看药代动力学数据。参考: 基于全局优化的非线性代数模型参数估计。 计算机与化学工程，第22卷，增编1,1998年3月15日，S213-S220页。
数据是时间与注意力的对比
T = [3.92, 7.93, 11.89, 23.90, 47.87, 71.91, 93.85, 117.84]
t = 1×8
3.9200 7.9300 11.8900 23.9000 47.8700 71.9100 93.8500 117.8400
C = [0.163, 0.679, 0.679, 0.388, 0.183, 0.125, 0.086, 0.0624]
c = 1×8
0.1630 0.6790 0.6790 0.3880 0.1830 0.1250 0.0860 0.0624
我喜欢查看数据，一方面是为了确保我没有输入错误，另一方面是为了了解整个系统。实际上，让我们把数据可视化。
情节(t、c、“o”）
包含(“时间”）
ylabel (“浓度”）

三室模型

在参考文献中，我们拟合一个3室模型，3个衰减指数的和。
c=b1 e(-b4 t)+b2 e(-b5 t)+b3 e(-b6 t)" style="vertical-align:-6px"> c ＝ b 1 e （ - b 4 t ） + b 2 e （ - b 5 t ） + b 3. e （ - b 6 t ）
我们可以把这个模型表示为 匿名函数 的 t (时间)及模型参数[ (1) b(2)……b (6)] ．
模型= @ (b, t) b (1) * exp (t - b (4) *) + b (2) * exp (t - b (5) *) + b (3) * exp (t - b (6) *)
模型=function_handle与价值:@ (b, b t) (1) * exp (t - b (4) *) + b (2) * exp (t - b (5) *) + b (3) * exp (t - b (6) *)

定义优化问题

接下来，我们定义优化问题来解决使用 具体问题具体分析公式 ．这允许我们选择我们想要的求解器，提供数据，并自然地表达约束和选项。
问题= createOptimProblem (“lsqcurvefit”，．..
“目标”、模型．..
“xdata”t“ydata”c．..
“x0”的(6),．..
“磅”， [-10 -10 -10 0 0 0]，．..
乌兰巴托的， [10 10 10 0.5 0.5 0.5]，．..
“选项”optimoptions (“lsqcurvefit”，．..
“OutputFcn”@curvefittingPlotIterates,．..
“显示”，“没有”））
问题=结构体字段:
摘要目的:@ (b, b t) (1) * exp (t - b (4) *) + b (2) * exp (t - b (5) *) + b (3) * exp (t - b (6) *) x0: [1 1 1 1 1 1] xdata: [3.9200 7.9300 11.8900 23.9000 47.8700 71.9100 93.8500 117.8400] ydata:[0.1630 0.6790 0.6790 0.3880 0.1830 0.1250 0.0860 0.0624]磅:[-10 -10 -10 0 0 0]乌兰巴托:[10 10 10 0.5000 0.5000 0.5000]解算器:“lsqcurvefit”选项:[1×1 optim.options.Lsqcurvefit]

解决这个问题

首先直接解决问题一次。
b = lsqcurvefit(问题)
b = 1×6
0.1842 0.1836 0.1841 0.0172 0.0171 0.0171
您会注意到，模型没有很好地拟合数据，甚至没有遵循数据的形状。

MultiStart

我们看看能不能从不同的点开始。
女士= MultiStart;
ms.Display =“通路”；
rng默认的
数字
抽搐
[~，fval,exitflag,output,金宝搏官方网站solutions] = run(ms, problem, 50)
运行当地当地当地一阶指数exitflag f (x) # iter F-count最优1 3 0.222 8 0.000154 21 154 0.001864 63 0.0006396 - 2 3 3 3 0.009442 44 315 0.01989 4 3 e-05 34 245 0.002586 1.462 1.454 e-05 19 140 1.079 e-05 6 3 1.475 e-05 24 175 0.009445 357 0.002266 0.006883 7 3 8 3 1.495 e-05 32 231 0.006853 9 3 1.466 e-05 35 252 0.004783 0.00944 80 567 0.01042 11 3 1.471 e-05 175 0.001576 287 1.566 0.005472 12 3 e-05 24 13 3 0.009439 24 175 0.0005121 14 3 0.009451 41 294 0.02935 189年15 3 0.009493 26日0.004837 16 3 1.476 e-05 287 0.008288 287 1.494 0.006352 17 3 e-05 40 18 62 441 0.009446 1.457 0.01296 19 3 e-05 22 161 1.488 0.001755 20 3 e-05 53 266 378 0.00944 0.007608 21 3 3722 3 1.464 e-05 24 175 0.006878 0.009449 0.003709 23 3 43 308 0.02515 24 336 0.007942 0.009447 47 25 3 1.455 e-05 23 168 0.009442 0.001621 26 3 32 231 0.01328 27日3 1.479 e-05 40 287 0.004821 28日3 1.479 e-05 18 133 0.0009721 0.006878 1.456 e-05 72 29 511 30 31 3 1.455 e-05 42 301 0.001122 336 0.01537 0.009441 47 32 3 0.009451 47 336 0.008942 33 30.0001729 315 0.01062 105 0.009442 0.0003276 34 3 44 35 3 0.0001751 21 154 1.509 7.71 e-05 36 3 e-05 26 189 0.009458 0.009896 37 3 39 280 0.02208 1.454 e-05 24 38 175 0.009441 427 0.0107 7.815 e-08 39 3 40 245 0.002981 1.472 e-05 34 41 3 1.503 e-05 22 161 0.00585 42 3 0.00952 112 0.009439 0.008492 43 3 21 154 0.000769 1.462 e-05 44 64 4550.0005576 45 3 0.009439 17 126 0.0001567 46 3 1.471e-05 30 217 0.001973 47 3 0.009444 38 273 0.02022 48 3 1.474e-05 24 175 0.004799 49 3 1.522e-05 42 301 0.008228
MultiStart从所有起始点完成运行。所有50个本地求解器都以一个正的本地求解器退出标志运行。
fval = 1.4540 e-05
exitflag = 1
输出=结构体字段:
funcCount: 12726 localSolverTotal: 50 localSolverSuccess: 50 localSolverIncomplete: 0 localSolverNoSolution: 0 message: 'MultiStart completed the runs from all start points。↵↵所有50个本地求解器的运行都以一个正的本地求解器退出标志汇合。
金宝搏官方网站解决方案= 1×50对象
1 2 3. 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20. 21 22 23 24 25 26 27 28 29 30.
1 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution 1×1 GlobalOptimSolution
serialTime = toc;

设想最佳解决方案

第50个解，也就是上面画的，不一定是最好的。幸运的是, MultiStart 按最好到最坏的顺序排列金宝搏官方网站解决方案。所以我们只需要看第一个。
curvefittingPlotIterates(金宝搏官方网站解决方案)
你现在可以看到，第50个不是最好的解决方案，因为最后一个显示的均方误差要高出10倍。

MultiStart与并行计算

现在，我将看看是否可以在本地使用所有4个内核作为并行工作器来提高性能。
ms.UseParallel = true;
gcp;
抽搐;
rng默认的
运行(女士,问题,50);
并行运行本地求解器。运行当地当地当地一阶指数exitflag f (x) # iter F-count最优1 3 8 0.222 63 0.0006396 567 0.00944 80 0.01042 252 0.00478 1.466 e-05 35 2 3 154 21 0.000154 1.476 0.001864 16 3 e-05 189 0.004837 287 0.009493 0.006352 15 3 26 14 3 0.009451 41 294 0.02935 315 0.01989 0.009442 44 22 175 0.003709 1.464 e-05 24 21 30．009443.72660．006878 20 3 1.488e-05 53 378 0.007608 4 3 1.462e-05 34 245 0.002586 8 3 1.495e-05 32 231 0.006853 7 3 0.009445 50 357 0.002266 6 3 1.475e-05 24 175 0.006883 5 3 1.454e-05 19 140 1.079e-05 28 3 1.479e-05 18 133 0.006878 27 3 1.479e-05 40 287 0.004821 26 3 0.009442 32 231 0.01328 25 3 1.455e-05 23 168 0.001621 13 3 0.009439 24 175 0.0005121 12 3 1.566e-05 24 175 0.001576 11 3 1.471e-05 40 287 0.005472 34 3 0.009442 44 315 0.01062 33 3 0.0001729 14 105 0.0003276 32 3 0.009451 47 336 0.008942 19 3 1.457e-05 22 161 0.001755 18 3 0.009446 62 441 0.01296 17 3 1.494e-05 40 287 0.008288 24 3 0.009447 47 336 0.007942 23 3 0.009449 43 308 0.02515 31 3 0.009441 47 336 0.01537 30 3 1.455e-05 42 301 0.001122 29 3 1.456e-05 72 511 0.0009721 42 3 0.00952 15 112 0.008492 45 3 0.009439 17 126 0.0001567 37 3 0.009458 39 280 0.02208 36 3 1.509e-05 26 189 0.009896 35 3 0.0001751 21 154 7.71e-05 40 3 1.472e-05 34 245 0.002981 39 3 0.009441 60 427 0.0107 38 1 1.454e-05 24 175 7.815e-08 41 3 1.503e-05 22 161 0.00585 49 3 1.522e-05 42 301 0.008228 47 3 0.009444 38 273 0.02022 44 3 1.462e-05 64 455 0.0005576 43 3 0.009439 21 154 0.000769 48 3 1.474e-05 24 175 0.004799 50 3 0.009445 40 287 0.02166 46 3 1.471e-05 30 217 0.001973 MultiStart completed the runs from all start points. All 50 local solver runs converged with a positive local solver exit flag.
parallelTime = toc;

计算加速

由于水池启动时间的关系，在第二轮之前速度可能不会明显加快。自从我更早开始，我就看到了不错的速度。
加速= serialTime / parallelTime
加速= 2.5014

您是否存在解决方案对起点敏感的问题?

告诉我们 关于你对问题解决方案空间的探索。如果解决方案对起始位置敏感，则可以考虑使用 MultiStart 以及全局优化工具箱中的其他技术。
版权所有:The MathWorks, Inc.

附录

下面是绘制迭代的代码。
dbtypecurvefittingPlotIterates
1 function stop = curvefittingPlotIterates(x,optimValues,state) 2%用于绘制优化算法迭代图的输出函数。3 4%版权所有2010 The MathWorks, Inc. 5 6 persistent x0 r;7 if nargin == 1 8 showPlot(x(1).X,x(1).X0{:}，x(1).Fval) 9 else 10 switch state 11 case 'init' % store initial point for later use 12 x0 = x;13 case 'done'优化后，显示解决方案在plot title 16 r = optimValues.resnorm;7 showPlot(x,x0,r) 18 end 19 end 20 end 21 if nargout > 0 22 stop = false;28 if = @(b,x) b(1)*exp(-b(4).*x) + b(2).*exp(-b(5).*x) +…29 b(3)。* exp (- b (6) . * x);30 31 persistent h ha 32 if is空(h) || ~isvalid(h) 33 x = [3.92, 7.93, 11.89, 23.90, 47.87, 71.91, 93.85, 117.84];34 y = [0.163, 0.679, 0.679, 0.388, 0.183, 0.125, 0.086, 0.0624]; 35 plot(x,y,'o'); 36 xlabel('t') 37 ylabel('c') 38 title('c=b_1e^{-b_4t}+b_2e^{-b_5t}+b_3e^{-b_6t}') 39 axis([0 120 0 0.8]); 40 h = line(3:120,f(b,3:120),'Color','r','Tag','PlotIterates'); 41 42 else 43 set(h,'YData',f(b,get(h,'XData'))); 44 end 45 s = sprintf('Starting Value Fitted Value\n\n'); 46 47 for i = 1:length(b) 48 s = [s, sprintf('b(%d): % 2.4f b(%d): % 2.4f\n',i,b0(i),i,b(i))]; 49 end 50 s = [s,sprintf('\nMSE = %2.4e',r)]; 51 52 if isempty(ha) || ~isvalid(ha) 53 % Create textbox 54 ha = annotation(gcf,'textbox',... 55 [0.5 0.5 0.31 0.32],... 56 'String',s,... 57 'FitBoxToText','on',... 58 'Tag','CoeffDisplay'); 59 end 60 ha.String = s; 61 drawnow 62 63 end
|
カテゴリ:
优化
  • 打印
  • 发送电子邮件

另请参阅

コメント

コメントを残すには,ここをクリックしてMathWorksアカウントにサインインするか新しいMathWorksアカウントを作成します。