В этом примере показано, как выполнить аппроксимирование кривыми нелинейного метода наименьших квадратов с помощью Основанного на проблеме Рабочего процесса Оптимизации.
Уравнение модели для этой проблемы
где , , , и неизвестные параметры, ответ, и время. Проблема требует данных в течение многих времен tdata и (шумные) измерения ответа ydata. Цель состоит в том, чтобы найти лучшее и , значение тех значений, которые минимизируют
Как правило, у вас есть данные для проблемы. В этом случае сгенерируйте искусственные зашумленные данные для проблемы. Используйте A = [1,2] и r = [-1,-3] как базовые значения и использование 200 случайных значений от 0 к 3 как данные времени. Постройте получившиеся точки данных.
rng default % For reproducibility A = [1,2]; r = [-1,-3]; tdata = 3*rand(200,1); tdata = sort(tdata); % Increasing times for easier plotting noisedata = 0.05*randn(size(tdata)); % Artificial noise ydata = A(1)*exp(r(1)*tdata) + A(2)*exp(r(2)*tdata) + noisedata; plot(tdata,ydata,'r*') xlabel 't' ylabel 'Response'
Данные, кажется, являются шумными. Поэтому решение, вероятно, не будет совпадать с исходными параметрами A и r очень хорошо.
Найти параметры оптимальной подгонки A и r, сначала задайте переменные оптимизации с теми именами.
A = optimvar('A',2); r = optimvar('r',2);
Создайте выражение для целевой функции, которая является суммой квадратов, чтобы минимизировать. Поскольку ответ имеет экспоненциальный термин, запишите анонимную функцию для ответа и преобразуйте ответ на выражение оптимизации с помощью fcn2optimexpr. Смотрите преобразуют нелинейную функцию в выражение оптимизации.
fun = @(A,r) A(1)*exp(r(1)*tdata) + A(2)*exp(r(2)*tdata); response = fcn2optimexpr(fun,A,r); obj = sum((response - ydata).^2);
Создайте задачу оптимизации с целевой функцией obj.
lsqproblem = optimproblem("Objective",obj);Для подхода, основанного на проблеме задайте начальную точку как структуру с именами переменных как поля структуры. Задайте начальный A = [1/2,3/2] и начальный r = [-1/2,-3/2].
x0.A = [1/2,3/2]; x0.r = [-1/2,-3/2];
Рассмотрите формулировку задачи.
show(lsqproblem)
OptimizationProblem :
Solve for:
A, r
minimize :
sum((anonymousFunction1(A, r) - extraParams{1}).^2)
where:
anonymousFunction1 = @(A,r)A(1)*exp(r(1)*tdata)+A(2)*exp(r(2)*tdata);
extraParams{1}:
2.9278
2.7513
2.7272
2.4199
2.3172
2.3961
2.2522
2.1974
2.1666
2.0944
1.9566
1.7989
1.7984
1.7540
1.8318
1.6745
1.6874
1.5526
1.5229
1.5680
1.4784
1.5280
1.3727
1.2968
1.4012
1.3602
1.2714
1.1773
1.2119
1.2033
1.2037
1.1729
1.1829
1.1602
1.0448
1.0320
1.0397
1.0334
1.0233
1.0275
0.8173
0.9373
1.0202
0.8896
0.9791
0.9128
0.7763
0.7669
0.6579
0.7135
0.7978
0.7164
0.7071
0.6429
0.6676
0.6782
0.6802
0.6328
0.6301
0.7406
0.4908
0.7126
0.6136
0.6269
0.4668
0.4963
0.5007
0.5226
0.3764
0.4824
0.3930
0.4390
0.4665
0.4490
0.4841
0.4539
0.3698
0.3974
0.3356
0.3045
0.4131
0.3561
0.3506
0.3960
0.3625
0.3446
0.3778
0.3565
0.3187
0.2677
0.2664
0.3572
0.2129
0.2919
0.1764
0.3210
0.3016
0.2572
0.2514
0.1301
0.2825
0.1372
0.1243
0.2421
0.1888
0.2547
0.2559
0.2632
0.1801
0.2309
0.2134
0.2495
0.2332
0.2512
0.1875
0.1861
0.2397
0.0803
0.1579
0.1196
0.1541
0.1978
0.2034
0.1095
0.1332
0.1567
0.1345
0.1635
0.1661
0.0991
0.1366
0.0387
0.1922
0.1031
0.0714
0.1178
0.0568
0.1255
0.0957
0.2313
0.1710
-0.0148
0.1316
0.0385
0.0946
0.1147
0.1436
0.0917
0.1840
0.0786
0.1161
0.1327
0.1026
0.1421
0.1142
0.0553
0.0036
0.1866
0.0634
0.0974
0.1203
0.0939
0.0429
0.0640
0.0811
0.1603
0.0427
0.1244
0.0993
0.0696
0.0264
0.0641
0.0703
0.0010
0.0793
0.0267
0.0625
0.0834
0.0204
0.0507
0.0826
-0.0272
0.1161
0.1832
0.1100
0.0453
0.0826
0.0079
0.1531
0.1052
0.0965
0.0132
0.0729
0.0287
0.0410
0.0280
0.0049
0.0102
0.0442
-0.0343
Решите задачу.
[sol,fval] = solve(lsqproblem,x0)
Solving problem using lsqnonlin. Local minimum found. Optimization completed because the size of the gradient is less than the value of the optimality tolerance.
sol = struct with fields:
A: [2x1 double]
r: [2x1 double]
fval = 0.4724
Постройте получившееся решение и исходные данные.
figure responsedata = evaluate(response,sol); plot(tdata,ydata,'r*',tdata,responsedata,'b-') legend('Original Data','Fitted Curve') xlabel 't' ylabel 'Response' title("Fitted Response")
График показывает, что подходящие данные совпадают с исходными зашумленными данными довольно хорошо.
Смотрите, как тесно подходящие параметры совпадают с исходными параметрами A = [1,2] и r = [-1,-3].
disp(sol.A)
1.1615
1.8629
disp(sol.r)
-1.0882 -3.2256
Подходящие параметры выключены приблизительно на 15% в A и 8% в r.