matlab - 如何在MATLAB中模拟间歇性需求？

Question

上下文:间歇性需求是指(需求)时间序列数据，其特征是存在零需求时段。这些非零需求通常很小但变化很大。间歇性需求有很多应用，对许多传统的预测和控制系统提出了挑战。

如何在 MATLAB 中生成间歇需求的随机样本？

在 R 中 tsintermittent包中，有一个模拟间歇性需求的函数叫做simID() .
我不知道有任何内置的 MATLAB 函数可以实现这一点(而且我非常怀疑它们是否存在)。

相关参数包括可变性和间歇性。

Answer 1

可以使用 nbinrnd() 使用统计工具箱在 MATLAB 中复制该函数和 binornd() .

语法是[d] = simID(n, obs, adi, cv2, level) 和

n:要创建的需求流数量
obs:模拟的时间段数
adi:平均需求区间(时间段)
CV2:coefficient of variation平方
水平:非零需求的平均值
d: 'n' x 'obs' 矩阵与指定的需求流

这明确假设非零需求到达遵循 Bernoulli distribution非零要求 negative binomial distribution .^[1]

_{[1] Petropoulos F.、Makridakis S.、Assimakopoulos V. 和 Nikolopoulos K.(2014 年)“需求预测中的‘类(class)用马’”，欧洲运筹学杂志 ，卷。 237，第 1 期，第 152-163 页。}

% Test in MATLAB 
D = simID(1,50000,2.5,2.5,3);

mean(D)                          % 1.2089
std(D)                           % 3.3913

# Test in R
# install.packages('tsintermittent')
library('tsintermittent')

D = simID(1,50000,2.5,2.5,3)
mean(D$ts.1)                      # 1.19632
sd(D$ts.1)                        # 3.261225

由于没有相同的随机种子，值不完全一致。

---

MATLAB 代码:(使用 R2019a 测试)

function [d] = simID(n, obs, adi, cv2, level)
% n: number of demand streams to create
% obs: number of time periods to simulate
% adi: average demand interval
% cv2: coefficient of variation squared
% level: average of nonzero demand 
%
% d: 'n' x 'obs' matrix with the specified demand streams
%
% If 'level' empty, defaults to sample from Uniform(10,100). 
% 
% MATLAB implementation of simID() from 'tsintermittent' package for R
% https://www.rdocumentation.org/packages/tsintermittent/versions/1.9/topics/simID
% https://CRAN.R-project.org/package=tsintermittent
%
% This simulator assumes that non-zero demand arrivals follow a bernoulli distribution 
% and the non-zero demands a negative binomial distribution. 
% Petropoulos F., Makridakis S., Assimakopoulos V. & Nikolopoulos K. (2014) 
% "'Horses for Courses' in demand forecasting", 
% European Journal of Operational Research, Vol. 237, No. 1, pp. 152-163

% Input Error Checking ****************************************************
narginchk(3,5)
if isempty(level), level = 10 + 90*rand(); end
if cv2<=0, error('cv2 must be positive'), end
if n<0 || obs<0 || adi < 0, error('n, obs, and adi must be positive'), end
if adi < 1, error('adi must be greater than 1'), end
% End (Input Error Checking) **********************************************

m = level - 1;
if ( cv2 ~= 0)
    p = (m / ( cv2 * ( (m+1)^2) ) ); % calculates the p for negative binomial function
    r = m * p / (1 - p); % calculates the r for the negative binomial function
    d = binornd(1, 1 / adi, n, obs) .* (nbinrnd(r, p, n, obs) + 1);
else
    d = binornd(1, 1/adi, n, obs) .* repmat(m + 1, n, obs);
end

end