在Numpy中生成对称矩阵

40

我正在尝试在numpy中生成对称矩阵。具体而言,这些矩阵应该有随机的条目和随机的内容。我们不关心主对角线上的条目是什么,因此我已经将其随机化了。

我采取的方法是首先生成一个nxn的全零矩阵,然后简单地循环遍历矩阵的索引。如何使用numpy更高效地完成这个过程?

import numpy as np
import random

def empty(x, y):
    return x*0

b = np.fromfunction(empty, (n, n), dtype = int)

for i in range(0, n):
    for j in range(0, n):
        if i == j:
            b[i][j] = random.randrange(-2000, 2000)
        else:
            switch = random.random()
            random.seed()
            if switch > random.random():
                a = random.randrange(-2000, 2000)
                b[i][j] = a
                b[j][i] = a
            else:
                b[i][j] = 0
                b[j][i] = 0
7个回答

65

您可以这样做:

import numpy as np

N = 100
b = np.random.random_integers(-2000,2000,size=(N,N))
b_symm = (b + b.T)/2

您可以在np.random或等效的scipy模块中选择任何分布。

更新:如果您正在尝试构建类似图形的结构,请务必检查networkx包:

http://networkx.lanl.gov

这个程序内置了多个例程来构建图表:

http://networkx.lanl.gov/reference/generators.html

如果你想添加一些随机放置的零,你可以生成一个随机索引集并将值替换为零。


1
谢谢!那是一个高效的解决方案。不过,有没有办法让它在随机位置放置零?这个矩阵应该代表图的某种邻接矩阵,因此具有随机分布的零的矩阵更可取。 - Ryan
8
@Ryan:你关心随机数的分布类型吗?如果你加上 b + b.T,你会得到一个集中在0附近的非均匀分布。 - unutbu
我正在验证一些矩阵的属性。这更多是为了提供一些数学性质的有力证据,因此这里的分布并不那么重要。谢谢! - Ryan
3
@unutbu,没错,请使用np.tril(a) + np.tril(a, -1).T - Ben Usman
np.random.randint(-5,5,size=(N,N)) seems to have been deprecated need to use randint - sushmit
它不会生成对称矩阵。 - Oalvinegro

44

我最好做:

a = np.random.rand(N, N)
m = np.tril(a) + np.tril(a, -1).T

因为在这种情况下,矩阵的所有元素都来自相同的分布(在此情况下是均匀分布)。


5
这是一种非常优雅的保持相同分布的方式! - Arash

4
import numpy as np

n = 5
M = np.random.randint(-2000,2000,(n,n))
symm = M@M.T
# test for symmetry
print(symm == symm.T)

这对我有用


3

矩阵中存在一种数学属性,允许轻松创建这样的结构:A.T * A,其中A是行向量,A.t是转置(列向量)。这总是返回一个正定对称方阵,它始终可逆,因此您不必担心空主元;)

# any matrix algebra will do it, numpy is simpler
import numpy.matlib as mt

# create a row vector of given size
size  = 3
A = mt.rand(1,size)

# create a symmetric matrix size * size
symmA = A.T * A


PS:对于矩阵而言,“always invertible”是正确的,但对于向量而言则不然。因此,使用mt.rand(size, size)会在数值上更加稳健。虽然这样会增加计算成本,但它仍会返回一个大小为size x size的矩阵。 - Bruno

2
如果您不介意对角线上有零,您可以使用以下代码片段:
def random_symmetric_matrix(n):
    _R = np.random.uniform(-1,1,n*(n-1)/2)
    P = np.zeros((n,n))
    P[np.triu_indices(n, 1)] = _R
    P[np.tril_indices(n, -1)] = P.T[np.tril_indices(n, -1)]
    return P

请注意,由于对称性,您只需要生成n*(n-1)/2个随机变量。

0

这里有一个优雅的答案,可以生成一个矩阵,其中所有条目都遵循相同的分布。然而,该答案会丢弃 (n-1)*n/2 个随机数而不使用它们。

如果您想让所有值都遵循相同的分布,一次性生成所有值并仅生成您要使用的值,那么可以运行以下命令:

>>> import numpy as np
>>> n = 5
>>> r = np.random.rand(n*(n+1)//2)
>>> sym = np.zeros((n,n))
>>> for i in range(n):
...     t = i*(i+1)//2
...     sym[i,0:i+1] = r[t:t+i+1]
...     sym[0:i,i] = r[t:t+i]
... 
>>> print(sym)
[[0.03019945 0.30679756 0.85722724 0.78498237 0.56146757]
 [0.30679756 0.46276869 0.45104513 0.28677046 0.10779794]
 [0.85722724 0.45104513 0.62193894 0.86898652 0.11543257]
 [0.78498237 0.28677046 0.86898652 0.13929717 0.45309959]
 [0.56146757 0.10779794 0.11543257 0.45309959 0.5671571 ]]

这里的想法是按照三角形数来确定随机向量中有多少个元素已被先前使用。给定此“t”值,将当前行填充到对角线位置(包括对角线),并将当前列填充到(但不包括)对角线位置。

0
我正在使用以下函数使矩阵在垂直和水平方向上对称:
def make_sym(a):
    w, h = a.shape
    a[w - w // 2 :, :] = np.flipud(a[:w // 2, :])
    a[:, h - h // 2:] = np.fliplr(a[:, :h // 2])

让我们来看看它是如何工作的:

>>> m = (np.random.rand(10, 10) * 10).astype(np.int)
>>> make_sym(m)
>>> m
array([[2, 7, 5, 7, 7, 7, 7, 5, 7, 2],
       [6, 3, 9, 3, 6, 6, 3, 9, 3, 6],
       [1, 4, 6, 7, 2, 2, 7, 6, 4, 1],
       [9, 2, 7, 0, 8, 8, 0, 7, 2, 9],
       [5, 5, 6, 1, 9, 9, 1, 6, 5, 5],
       [5, 5, 6, 1, 9, 9, 1, 6, 5, 5],
       [9, 2, 7, 0, 8, 8, 0, 7, 2, 9],
       [1, 4, 6, 7, 2, 2, 7, 6, 4, 1],
       [6, 3, 9, 3, 6, 6, 3, 9, 3, 6],
       [2, 7, 5, 7, 7, 7, 7, 5, 7, 2]])

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接