仅在非零条目不同结果时评估 arrayfun

Question

给定两个向量 x,y 并创建距离矩阵

B=pdist2(x,y);

现在我评估某个函数 f(x)，

s=5;
if s-x > 0
    y=(1-x/s)^4*(1+4*x/s)/20;
else
    y=0;
end

其中 s 是一些参数。通过预定义 s 并使用

A=arrayfun(@f,B);

矩阵 A 是正定的。

由于矩阵是稀疏的，我只想在将返回非零的值处评估 f(x)。我的尝试是

B=pdist2(x,y);
B(B>s)=0;
B=B/s;
indexB=find(B);
Atmp=arrayfun(@f,B(indexB));
A(indexB)=Atmp;

但是现在 A 不是正定的，并且比以前具有更少的非零项。有人能解释一下我在做什么/正在做什么吗？

谢谢大家，我发现了错误。在 B 的对角线上有一些零条目，我必须单独计算。

Answer 1

我想知道您为什么要尝试索引B，然后将一些函数应用于每个单独的元素。相反，让我们矢量化您的代码

% Same setup
B = pdist2(x,y); s = 5;
B(B >= s) = 0;  B = B/s;
% Don't use arrayfun, instead use vectorised code
% The previous line already dealt with the 'else' case from f(x) 
B = ((1-B/s).^4).*(1+4*B/s)/20;

请注意，每次我们在两个非标量之间进行操作时，我们必须使用 element-wise 操作 .^ 或 .*。其他除法和乘法都是标量，但我们可以安全并养成使用它们的习惯！

B = ((1-B./s).^4).*(1+4.*B./s)./20;

---

关于 MATLAB 术语的注释：小心地说，您在这里有一个“稀疏矩阵”，即 MATLAB 中的 special variable type。你有一个数学上稀疏的标准矩阵！对稀疏矩阵进行操作可能会略有不同。

Answer 2

您应该预先分配与B 相同大小的矩阵A，否则MATLAB 不知道A 的正确大小：

A=zeros(size(B));

注意find 返回linear indices，所以A 变成了向量而不是矩阵。

Answer 3

您可以使用spfun 代替arrayfun 将函数应用于非零稀疏矩阵元素：

A = spfun(@f, B);

虽然矢量化可能更有效。

• spfun 也适用于密集矩阵和稀疏矩阵，但它的输出是稀疏矩阵

感谢@Wolfie 提供了一些测试和澄清。