将创建连续幂矩阵的向量化

Question

将创建连续幂矩阵的向量化

matlabmatrixvectorization

12

令 x=1:100，N=1:10。我想要创建一个矩阵 x^N，使得第 i 列包含条目 [1 i i^2 ... i^N]。

我可以使用 for 循环轻松完成这个任务。但是是否有一种使用向量化代码的方法呢？

- alext87

4

非常感谢能够有机会看到这么多不同的方法，给个赞。 - zellus

5个回答

6

听起来像是一个范德蒙矩阵。因此使用 vander 函数：

A = vander(1:100);
A = A(1:10, :);

- Oliver Charlesworth

唯一的问题是VANDER创建的是方阵，因此你需要做很多额外的工作。 - gnovice

@gnovice：这不是一个很棘手的问题，可以参考上面的代码解决。 - Oliver Charlesworth

实际上，由于OP在评论中提到了使用方阵进行操作，因此这个函数看起来是一个不错的选择。你需要做的唯一更改就是使用ROT90旋转它，以便第i列包含i的幂：A = rot90(vander(1:N)); - gnovice

5

由于你的矩阵并不是那么大，最直接的方法是使用MESHGRID和逐元素幂运算符.^：

[x,N] = meshgrid(1:100,0:10);
x = x.^N;

这将创建一个11行100列的矩阵，其中每一列i包含[i^0; i^1; i^2; ... i^10]。

- gnovice

2

很好的meshgrid使用方法。我对这种方法和这里提出的类似方法唯一的问题是它们不以低成本使用第i行来计算第i+1行。cumprod变量确实可以做到，但有些美学缺失。 - Adrien

@Adrien：对于大矩阵，你说得没错，通过连续乘以基础值计算值可能会更快。但对于像这个例子中相对较小的矩阵，选择易于阅读和理解但略微低效的方法似乎是一个合理的权衡。 ;) - gnovice

我喜欢这种方法... 在我的实际编程中，我正在创建 10000x10000 的数据。所以也许我应该使用更有效的方法。不过知道这种方式还是很好的。谢谢！ - alext87

@alext87：这些矩阵非常大，你很可能会遇到内存问题。实际上，它们太大了，以至于当我创建一个那么大的矩阵时，甚至没有足够的内存来执行基本操作。例如，clear all; V=ones(10000); V=V.';会给我一个内存不足的错误。此外，将数字提高到如此大的幂几乎没有任何意义。即使只是2^1000，结果也是1.0715e+301，接近双精度浮点数的最大值。 - gnovice

@Adrien：实际上，Oli Charlesworth建议的函数VANDER利用了您提到的高效计算。您可以在MATLAB命令窗口中键入type vander来查看代码。 - gnovice

2

不确定它是否真正符合您的问题。

bsxfun(@power, cumsum(ones(100,10),2), cumsum(ones(100,10),1))

编辑: 正如Adrien所指出的那样，我的第一次尝试并不符合问题要求。

xn = 100;
N=10;
solution = [ones(1,xn); bsxfun(@power, cumsum(ones(N,xn),2), cumsum(ones(N,xn),1))];

- zellus

计算矩阵的第一行应该为1，在您的情况下它是1:100。 - Adrien

1

为什么不使用易于理解的for循环？

c = [1:10]'; %count to 100 for full scale problem
for i = 1:4; %loop to 10 for full scale problem
    M(:,i) = c.^(i-1)
end

要理解人们展示的聪明向量化版本的这段代码，需要更多的思考。我的方法比较野蛮，但是任何人读它都能理解。

我喜欢易于理解的代码。

（是的，我可以预先分配内存。但对于这种小型情况降低清晰度不值得。）

- MatlabDoug

2

如果您使用循环，则应逐行填充矩阵，而不是使用.^运算符。相反，用1填充第一行，然后通过将前一行乘以[1:100]来创建每个新行。这样既具有可读性又具有效率。 - Adrien

2

Adrien的建议提高了12%的速度。对于正在讨论的问题规模来说，这是无法察觉的。我必须运行这个程序100,000次才能使整个过程从2.5秒变为2.8秒。对我来说，我的版本更容易阅读，因为没有M的初始化步骤，但口味因人而异！ :) - MatlabDoug

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Adrien · Accepted Answer

我会选择：

x = 1:100;
N = 1:10;
Solution = repmat(x,[length(N)+1 1]).^repmat(([0 N])',[1 length(x)]);

另一种解决方案（可能更有效）：

Solution = [ones(size(x)); cumprod(repmat(x,[length(N) 1]),1)];

或者甚至这样：

 Solution = bsxfun(@power,x,[0 N]');

希望这可以帮到你。