Create symmetrical matrix from list of values

Question

I need to make a matrix (in the form of a numpy array) by taking a list of parameters of length N and returning an array of dimensions N+1 x N+1 where the off-diagonals are symmetric and each triangle is made up of the values given. The diagonals are equal to 0 - the row sum of the off-diagonals.

I don't really like having to instantiate an empty array, fill the top triangle, then create a new matrix with the bottom triangle filled in. I'm wondering if I can do this in fewer lines (while still being readable) and also wondering if it can be made faster. vals is never going to be very large (between 2 and 100, most likely) but the operation is going to be repeated many times.

def make_sym_matrix(dim, vals):
 my_matrix = np.zeros([dim,dim], dtype=np.double)
 my_matrix[np.triu_indices(dim, k=1)] = vals
 my_matrix = my_matrix + my_matrix.T
 my_matrix[np.diag_indices(dim)] = 0-np.sum(my_matrix, 0)
 return my_matrix

• hpaulj: If you are doing this many times for the same dim, try to reuse xy=np.triu_indices(dim,1). For dim=4, that takes about half the time; so reuse may double the overall speed.

Answer 1

First of all I would use np.zeros() to initialize your matrix. Using np.empty() can create a matrix with large values relative to your values on the diagonal which will affect the computation of 0-np.sum(my_matrix, 0) due to numeric underflow.

For example, just run this loop and you'll see it happen:

for i in xrange(10):
 print make_sym_matrix(4, [1,2,3,4,5,6])

Secondly, you can avoid taking the transpose by re-using the triu indices. Here is the code combining these two ideas:

def make_sym_matrix(n,vals):
 m = np.zeros([n,n], dtype=np.double)
 xs,ys = np.triu_indices(n,k=1)
 m[xs,ys] = vals
 m[ys,xs] = vals
 m[ np.diag_indices(n) ] = 0 - np.sum(m, 0)
 return m

• C_Z_: Already switched from np.empty to np.zeros, thanks for the suggestions. This should be faster than the code I have, right? Because it doesn't involve creating another matrix (like I have in the line: my_matrix = my_matrix + my_matrix.T)
• ErikR: I would imagine so, but ultimately the only way to know for sure is to benchmark it.
• Janne Karila: 0 - np.sum(m, 0) could be just -np.sum(m, 0), and m[ np.diag_indices(n) ] -= np.sum(m, 0) seems faster still.
• hpaulj: At least on this small sample array I'm not seeing much difference in time. Timewise a transpose is a trivial operation.