将向量放入矩阵+变换中答案

【问题标题】：Put a vector inside a matrix + transformation将向量放入矩阵+变换中
【发布时间】：2018-11-06 09:35:29
【问题描述】：

所以我有一个用 Numpy 创建的向量，称为

V = [10 20 30 40  0  1]

我想要一个像这样的矩阵 M：

 [10.  0.  0.  0.  0.  0.]
 [20. 10.  0.  0.  0.  0.]
 [30. 20. 10.  0.  0.  0.]
 [40. 30. 20. 10.  0.  0.]
 [ 0. 40. 30. 20. 10.  0.]
 [ 1.  0. 40. 30. 20. 10.]
 [ 0.  1.  0. 40. 30. 20.]
 [ 0.  0.  1.  0. 40. 30.]
 [ 0.  0.  0.  1.  0. 40.]
 [ 0.  0.  0.  0.  1.  0.]
 [ 0.  0.  0.  0.  0.  1.]

为此，我使用 for 循环和向量提取，但它太长了，因为我的向量 V 有 500 列，而矩阵 M 有 500*2 -1 行和 500 列。

此外，对于不同的向量 V，我必须至少重复 100 000 次

是否有可能，使用矩阵计算并避免循环来得到这个结果？（尽可能快）

谢谢！

（我在 Spyder 上使用 Python 3.6）

编辑：我的解决方案：

t=480 n=10000

t1 = time.time()
for p in range(n):
    for j in range(M.shape[1]):
        M[j:j+t,j] = np.transpose(V[:])
print(time.time()-t1)

14 秒仅 10000 次...太长了

编辑 2：评论中的解决方案基准：

（这里的普拉特是 V）

t1 = time.time()
for p in range(n):
    for j in range(M.shape[1]):
        M[j:j+t,j] = np.transpose(Prate[:])
print(time.time()-t1) 


t1 = time.time()
for p in range(n):
    n = len(Prate)
    m = np.tile(np.concatenate((np.array(Prate), np.zeros(t))), t)[:2*t*t-t]
    result = m.reshape(t, -1).T
print(time.time()-t1)  

t1 = time.time()
for p in range(n):
    ind = np.arange(t)
    indices = ((ind[:,None] + ind).ravel() , np.repeat(ind, t))
    base = np.zeros((n1, t))
    base[indices] = np.tile(Prate, t)
print(time.time()-t1)

输出：

16.737313747406006
29.46031618118286 
3.6843104362487793

编辑 3：为了避免通过两倍于所需大小的数组，我提出了不同的问题：

我有一个向量（1x6）：

V = [1 20 5 0  0  9]

我想要一个像这样的矩阵 M (6x6)：

 [1.  20.  5.  0.  0.  9.]
 [0.  1.  20.  5.  0.  0.]
 [0.  0.  1.  20.  5.  0.]
 [0.  0.  0.  1.  20.  5.]
 [0.  0.  0.  0.  1.  20.]
 [0.  0.  0.  0.   0.  1.]

在每一行中，它是相同的向量 V（它的一部分），但有一个偏移量以获得三角矩阵。

如何在没有循环的情况下做到这一点？

（这只是一个简单的例子，但真正的向量 V 要大得多）

谢谢你：D

【问题讨论】：

标签： python python-3.x matrix

【解决方案1】：

如果您只想要一个快速的解决方案而不是避免循环的解决方案，您可以简单地使用Numba 或Cython。

示例

import numpy as np
import numba as nb

@nb.njit()
def create_mat(V):
  arr=np.zeros((V.shape[0],V.shape[0]),dtype=V.dtype)
  for i in range(V.shape[0]):
    for j in range(i,V.shape[0]):
      arr[i,j]=V[j-i]
  return arr

时间

V=np.random.rand(10000)
#The first call has an constant compilation overhead of about 0.2s, 
#which is neglected here.
create_mat: 0.35s

【讨论】：

我不知道与 Python for 循环相比，使用 njit 的 numba 这么快！谢谢

【解决方案2】：

这是一种矢量化方法：

In [74]: ind = np.arange(6)

In [75]: indices = ((ind[:,None] + ind).ravel() , np.repeat(ind, 6))

In [76]: base = np.zeros((11, 6))

In [77]: base[indices] = np.tile(V, 6)

In [78]: base
Out[78]: 
array([[10.,  0.,  0.,  0.,  0.,  0.],
       [20., 10.,  0.,  0.,  0.,  0.],
       [30., 20., 10.,  0.,  0.,  0.],
       [40., 30., 20., 10.,  0.,  0.],
       [ 0., 40., 30., 20., 10.,  0.],
       [ 1.,  0., 40., 30., 20., 10.],
       [ 0.,  1.,  0., 40., 30., 20.],
       [ 0.,  0.,  1.,  0., 40., 30.],
       [ 0.,  0.,  0.,  1.,  0., 40.],
       [ 0.,  0.,  0.,  0.,  1.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  1.]])

【讨论】：

您好，感谢您提供此解决方案。在我的电脑上，它比前两个慢……但学习其他方法很好！
@rniii 你如何对它们进行基准测试？请将您的代码添加到问题中。
错误的计时方法，如果我没有做错，你的方法会更快（见编辑2）！太棒了：D
但这很奇怪，因为在我的时间里 n=1000 比 n=10000 慢...你觉得呢？
@mii：我跑了banchmarks，实际上确实说这种方法（奇怪的是）慢了，但也许我用错了方法？

【解决方案3】：

是的，对于数组v，我们可以先用len(v)extra0`s构造一个numpy数组：

n = len(V)
m = np.tile(np.concatenate((V, np.zeros(n))), n)[:2*n*n-n]
result = m.reshape(n, 1).T

对于给定的V 列表，这给了我们：

>>> np.tile(np.concatenate((V, np.zeros(n))), n)[:2*n*n-n].reshape(n,-1).T
array([[10.,  0.,  0.,  0.,  0.,  0.],
       [20., 10.,  0.,  0.,  0.,  0.],
       [30., 20., 10.,  0.,  0.,  0.],
       [40., 30., 20., 10.,  0.,  0.],
       [ 0., 40., 30., 20., 10.,  0.],
       [ 1.,  0., 40., 30., 20., 10.],
       [ 0.,  1.,  0., 40., 30., 20.],
       [ 0.,  0.,  1.,  0., 40., 30.],
       [ 0.,  0.,  0.,  1.,  0., 40.],
       [ 0.,  0.,  0.,  0.,  1.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  1.]])

对于一个包含 500 个元素和 10'000 次运行的 numpy 数组，我们得到：

>>> timeit(f, number=10000)
5.285840999999891

因此，在我的机器上转换具有 500 个元素的单个数组需要 0.5 毫秒。因此，构建所有这些数组大约需要 52.86 秒。

编辑：我实现了三个尝试如下：

def wil():
    n = len(V)
    return np.tile(np.concatenate((V, np.zeros(n))), n)[:2*n*n-n].reshape(n,-1).T


def mii():
    n = len(V)
    M = np.zeros((2*n-1, n))
    for j in range(n):
        M[j:j+n,j] = np.transpose(V[:])
    return M


def kas():
    n = len(V)
    ind = np.arange(n)
    indices = ((ind[:,None] + ind).ravel() , np.repeat(ind, n))
    base = np.zeros((2*n-1, n))
    base[indices] = np.tile(V, n)
    return base

并生成一个包含 500 个元素的随机数组：

V = np.random.randn(500)

然后我们运行测试：

>>> timeit(wil, number=10000)
3.461620999999923
>>> timeit(mii, number=10000)
13.704932000000099
>>> timeit(kas, number=10000)
159.63497699999994

【讨论】：

您好，谢谢您的回答。但是，我的解决方案需要更多时间...
@mii：那么分享您的解决方案当然是有意义的。
不是更快，但我猜大致相同（我发布了我的解决方案）
@mii：如果我使用timeit 运行此程序，您将在 21.44 秒内获得本地结果。
我会在我的电脑上查看！我没有找到每种方法的那些时间。谢谢！

【解决方案4】：

使用以下内容：

from scipy.linalg import toeplitz
res=np.tril(toeplitz(V).T

输出：

res
>>array([[10, 20, 30, 40,  0,  1],
         [ 0, 10, 20, 30, 40,  0],
         [ 0,  0, 10, 20, 30, 40],
         [ 0,  0,  0, 10, 20, 30],
         [ 0,  0,  0,  0, 10, 20],
         [ 0,  0,  0,  0,  0, 10]])

【讨论】：