MATLAB中含求和项的时变函数建模方法咨询

我尝试在MATLAB中建模如下时变函数：

$y(t) = \sum_{i=1}^{m} \gamma y(t-i) + \sum_{k=1}^{n} \beta u(t-k) + \delta$

其中$\gamma$、$\beta$、$\delta$为常数。我不确定处理求和项的最佳方式，虽然可以将$u(t)$和$y(t)$设为函数句柄，但symsum仅适用于符号运算。难道只能用for循环手动累加$y(t-i)$和$u(t-k)$？

我目前生成$y(t)$的代码如下：

N = 100;
y_out = zeros(N,1);
u = zeros(N,1);

for t = 1:1:N
u(t) = cos(beta*t);
end

y_out(1) = NARMA(beta,delta,gamma,true);
for t=2:N
    y_out(t) = NARMA(beta,delta,gamma);
end


function y_out = NARMA(beta,delta,gamma,first_call)
    % Set first_call to true only the first time you call the function
    persistent y
    if nargin>3 && first_call
        y = 0;
    end

    for i=1:1:m
    y = y + gamma*u + delta;
    end

    y_out = y;

end

解决方案 #

1. 手动循环累加（直接高效） #

对于这类时域递推模型，手动循环累加是最直接且高效的方式，完全不需要依赖符号运算。你的现有代码存在m/u未定义、persistent变量易出问题等缺陷，以下是修正后的可运行版本：

% 定义模型参数
N = 100;
m = 2; % y项的求和阶数，可按需修改
n = 3; % u项的求和阶数，可按需修改
gamma = 0.5;
beta = 0.3;
delta = 0.1;

% 向量化生成输入u(t)，替代循环更高效
u = cos(beta*(1:N)');

% 初始化y序列，处理初始条件（前m个值无历史项，仅保留delta）
y_out = zeros(N,1);
y_out(1:m) = delta;

% 递推计算后续y(t)
for t = m+1:N
    % 计算y的滑动求和项
    sum_y = gamma * sum(y_out(t-m : t-1));
    % 计算u的滑动求和项，确保索引不越界
    sum_u = beta * sum(u(max(1, t-n) : t-1));
    % 生成当前y(t)
    y_out(t) = sum_y + sum_u + delta;
end

2. 利用卷积简化求和（代码更简洁） #

求和项本质是滑动窗口求和，可借助conv函数实现，代码更简洁：

% 定义模型参数
N = 100;
m = 2;
n = 3;
gamma = 0.5;
beta = 0.3;
delta = 0.1;

u = cos(beta*(1:N)');
y_out = zeros(N,1);

% 初始条件处理
y_out(1:m) = delta;

% 生成求和对应的卷积核
kernel_y = gamma * ones(m,1);
kernel_u = beta * ones(n,1);

% 递推计算
for t = m+1:N
    % 卷积计算y的求和，取有效结果的最后一个值
    sum_y = conv(y_out(1:t-1), kernel_y, 'valid');
    sum_y = sum_y(end);
    % 卷积计算u的求和
    sum_u = conv(u(1:t-1), kernel_u, 'valid');
    sum_u = sum_u(end);
    y_out(t) = sum_y + sum_u + delta;
end

关键说明 #

• 符号运算工具（如symsum）更适合公式推导，数值计算场景下循环或卷积的效率更高
• 初始条件需根据实际场景调整：若t <= m时无历史y(t-i)值，可设为0或仅保留常数项$\delta$
• MATLAB对向量化代码优化更好，优先用向量化操作（如直接生成u序列）替代循环

内容的提问来源于stack exchange，提问作者Jonathan Frutschy