曲线拟合问题 MATLAB

计算科学 matlab 曲线拟合
2021-12-04 11:31:29

我有一个带有一些未知参数的微分方程系统，我需要找到最适合我拥有的数据的最佳参数集。

对于每个参数集选择（使用 for 循环），我的代码以数值方式求解 DE 系统，然后从时间向量中提取与数据时间值一致的时间子向量。下一步需要为每个与时间子向量匹配的变量提取一个新的输出子向量（以计算误差）。

这是我的代码的外观：

OO=[];
y2=(y(:,2))';
t1=t';
FF=[];
y1=(y(:,1))';

global "parameters and other global variables"

for parameter_1
    for parameter_2
        [t,y]=ode45('system',[time interval], [initial conditions]);
        time; % "time" is a code that extracts the time subvector

        for i=1:length(time)
            for j=1:length(t1)
                if time(i)==t1(j);
                    FF(end+1) = y1(j);
                end
            end
        end

        for k=1:length(time)
            for m=1:length(t1)
                if time(k)==t1(m);
                    OO(end+1) = y2(m);
                end
            end
        end

    end
end

此时代码还没有结束，但我希望它返回OO向量和FF用于循环的最后一个参数集选择的向量。但是，当代码运行（没有任何错误）时，它会产生向量OO并且FF比向量时间长：length(time) < {length(OO),length(FF)}.

1个回答

这是一个比通常简单地描述为“曲线拟合”的更一般的非线性优化问题。我建议你这样处理它。你需要定义一些东西：

  1. 您正在优化的参数。
  2. 要最小化的数量。
  3. 帮助函数来简化事情。

Matlab 有可以提供帮助的工具，它将帮助您将问题分解为如下几个部分：

主脚本

clear

% load the curve and time values you want the ODE to fit
[ tReference, yReference ] = loadKnownValues();

% load the initial conditions for the ODE:
initialConditions = loadInitialConditions();

% initial guess for parameter values (a column vector)
parametersGuess = [ 1; 6.83; ... ];

% define an objective function for the optimization
% I will assume you want to minimize the sum of the errors squared
% lsqnonlin is in the Optimization toolbox
%
% You could also use fminsearch. In that case you would want to use the
% objective function F = @(param) sum( objectiveFunction(param).^2 );
objectiveFunction = @(parameters) yReference ...
    - solveODE(tReference,initialConditions,parameters);

optimizedParameters = lsqnonlin(objectiveFunction,parametersGuess);

这个主脚本依赖于其他一些执行特定操作的函数。显然，您需要加载您尝试拟合的曲线和 ODE（loadKnownValues()和loadInitialConditions()）的初始条件。您还需要一个函数，在给定一组初始条件和参数的情况下求解 ODE，并在所需时间返回值。前两个应该很容易，最后一个可以如下完成：

求解ODE.m

function yVals = solveODE(tReference,initialConditions,parameters)
% yVals = solveODE(tReference,initialConditions,parameters)
% Solves the ODE for your problem and with the given initialConditions
% and parameters and returns the values at the specified times in
% tReference. The first time in tReference must correspond to the
% initial time.

% function that returns the right hand side of the ODE,
% use the parameters in the input vector "parameters"
rightHandSide = @(t,y) ... ;

[~, yVals] = ode45(rightHandSide, tReference, initialConditions);

显然，如果不了解更多有关系统的信息，我无法评论优化的收敛性或性能，但这应该可以帮助您入门。

其它你可能感兴趣的问题
上一篇测量图中顶点密度的差异？ 下一篇使用 P2/P1 元素时如何更新速度以包括压力