在不使用for循环的情况下生成迭代序列

Question

我想通过以下与非中心卡方分布相关的规则生成一个长度为N的随机序列x：

x_n+1~χ_Γ2(λx_n(

其中，Γ是表示自由度的给定常数，λ也是预先指定的，并且λ和x_n的乘积是非中心性参数，x₁应该是给定的。我写了以下代码来生成这样的序列，并在x_{1＝0.04、μ＝0.005、λ＝100和N＝1e5的情况下运行：}

tic;
N = 1e5;
x = zeros(1,N);
x(1) = 0.04;
nu = 0.005;
lambda = 100;
for i = 1:N-1
x(i+1) = ncx2rnd(nu,lambda*x(i));
end
toc;

为了说明我的问题，我测试了另一个与上面不同的例子。在这里，我考虑从分布χ_Γ²(λ(生成N=1e5个样本，其中Γ=0.005和λ=100:

tic;
N = 1e5;
x = zeros(1,N);
nu = 0.005;
lambda = 100;
for i = 1:N
x(i) = ncx2rnd(nu,lambda);
end
toc;
tic;
N = 1e5;
nu = 0.005;
lambda = 100;
x = ncx2rnd(nu,lambda*ones(1,N));
toc;

这两种方法等效。然而，事实证明，不用于循环的第二种方法比第一种方法快得多。两个例子之间的区别在于，在第二个例子中，生成某个样本的规则不需要先前样本的信息，而在第一个例子中情况并非如此，因此可以在不使用for循环的情况下同时生成所有样本。基于此，我想知道避免for循环是否会加速代码执行。那么，当对先前样本的依赖性规则明确时，是否有任何MATLAB内置函数可以生成第一个示例中显示的随机序列，而不使用for循环？如果规则是线性的，我知道函数filter是一个可能的选择，那么像第一个例子这样的情况呢？

Answer 1

从逻辑上讲，不进行迭代就不可能计算出迭代的东西。如果x(n+1)依赖于x(n)，则必须首先计算x(n)；聪明的把戏；在这里

这只会让我们在循环中优化计算，特别是ncx2rnd。与大多数MATLAB内置函数一样，它已经相当简洁和高性能，但仍有一些需要考虑的地方。注意，我要建议的是使用edit ncx2nrd来查看这个内置函数的内部，该函数包含MathWorks版权下的代码，我只是注意到它的工作原理。

1. 有一些输入检查可以处理大小不正确的输入和/或具有负值的输入。如果你可以自己承担验证的负担(即你知道你的输入是有效的(，那么你可以将函数简化为单一的数学运算：

   % function r = ncx2rnd(v,delta)
   r = 2.*randg(poissrnd(delta./2, sizeOut)) + 2.*randg(v./2,sizeOut);
   

   运行此独立程序可节省大约20%的处理时间，这些时间用于输入验证(使用标称N=1e5(。

2. 在MathWorks语法中，delta等于您的lambda*x(i)，包括v在内的另一个术语独立于x，因此您可以在循环之外计算它，即对randg的一个调用进行矢量化。再次使用N=1e5，这使总时间节省约25%。

结果将意味着对您的示例进行此更改：

% Common inputs
N = 1e5;
nu = 0.1;
lambda = 0.1;
% Baseline example
x = zeros(1,N);
x(1) = 0.04;
for i = 1:N-1
x(i+1) = ncx2rnd(nu,lambda*x(i));
end
% ~25% faster alternative, with no input validation and partially vectorised
x = zeros(1,N);
x(1) = 0.04;
vTerm = 2.*randg(nu./2, [1,N]);
for i = 1:N-1
x(i+1) = 2.*randg(poissrnd(lambda*x(i)./2, [1,1])) + vTerm(i);
end