Eigenvector calculation involving three matrices

Code performs well for small matrices, but gives memory error for large matrices of size 63000x93000, even when I am using server with memory of 64GB, space 100GB, swap 40GB.

# -*- coding: utf-8 -*-
from numpy import array, argmax, concatenate
import csv

UFILE = open(Filename, 'r')
with UFILE:
for i in READER:
yield i

U = []
for k in ufinal:
TEMP = []
for j in k:
TEMP.append(eval(j))
U.append(TEMP)
del TEMP
del k

S = []
for s in sfinal:
for t in s:
S.append(eval(t))
del t

EIGVAL = []
for e in EIGVALfinal:
for l in e:
EIGVAL.append(eval(l))
del l

EIGVEC = []
for e in EIGVECfinal:
for l in e:
EIGVEC.append(eval(l))
del l

NEW_U = []
NEW_EIGVAL = []
NEW_EIGVEC = []
for i in range(len(S)):
if S[i] > 0.001:
NEW_U.append(U[i])
NEW_EIGVAL.append(EIGVAL[i])
NEW_EIGVEC.append(EIGVEC[i])
del EIGVAL
del EIGVEC
# Select the first r columns of u corresponding to the r principle Eigenvector of
# MatrixForEigenvalues
TRANSPOSE_NEW_EIGVEC = array(NEW_EIGVEC).real.T   # Done for easy accesiblity of matrix
TRANSPOSE_MATRIX_U = array(NEW_U).real.T  # Done for easy accesiblity of matrix
FINAL_ARRAY_EIGVAL = []
FINAL_ARRAY_EIGVEC = []
FINAL_MATRIX_U = []
for i in range(len(array(NEW_EIGVAL).real)):
j = argmax(array(NEW_EIGVAL).real)
FINAL_ARRAY_EIGVAL.append((array(NEW_EIGVAL).real)[j])
FINAL_ARRAY_EIGVEC.append(TRANSPOSE_NEW_EIGVEC[j])
FINAL_MATRIX_U.append(TRANSPOSE_MATRIX_U[j])
TRANSPOSE_NEW_EIGVEC = concatenate((TRANSPOSE_NEW_EIGVEC[:j], TRANSPOSE_NEW_EIGVEC[j+1:]))
TRANSPOSE_MATRIX_U = concatenate((TRANSPOSE_MATRIX_U[:j], TRANSPOSE_MATRIX_U[j+1:]))
RESULT_MATRIX_U = array(FINAL_MATRIX_U).T     # This is the actual R
r_file = open('r.csv', 'w')
with r_file:
WRITER = csv.writer(r_file)
WRITER.writerows(RESULT_MATRIX_U)
print(RESULT_MATRIX_U)
del FINAL_MATRIX_U


Performance

Here's where all the action happens:

# Select the first r columns of u corresponding to the r principle Eigenvector of
# MatrixForEigenvalues
TRANSPOSE_NEW_EIGVEC = array(NEW_EIGVEC).real.T   # Done for easy accessibility of matrix
TRANSPOSE_MATRIX_U = array(NEW_U).real.T  # Done for easy accessibility of matrix
FINAL_ARRAY_EIGVAL = []
FINAL_ARRAY_EIGVEC = []
FINAL_MATRIX_U = []
for i in range(len(array(NEW_EIGVAL).real)):
j = argmax(array(NEW_EIGVAL).real)
FINAL_ARRAY_EIGVAL.append((array(NEW_EIGVAL).real)[j])
FINAL_ARRAY_EIGVEC.append(TRANSPOSE_NEW_EIGVEC[j])
FINAL_MATRIX_U.append(TRANSPOSE_MATRIX_U[j])
TRANSPOSE_NEW_EIGVEC = concatenate((TRANSPOSE_NEW_EIGVEC[:j], TRANSPOSE_NEW_EIGVEC[j+1:]))
TRANSPOSE_MATRIX_U = concatenate((TRANSPOSE_MATRIX_U[:j], TRANSPOSE_MATRIX_U[j+1:]))
RESULT_MATRIX_U = array(FINAL_MATRIX_U).T     # This is the actual R


The first thing I notice is multiple calls to array(), specifically for array(NEW_EIGVAL).real. This is inefficient; doing all the necessary conversions at the beginning and then referencing the variables would be much more efficient. Next, I get the strong suspician that j = argmax(array(NEW_EIGVAL).real) could be moved out of the loop, since array(NEW_EIGVAL).real isn't changing between iterations of the for loop. This in turn makes me strongly suspect that the line appending to FINAL_ARRAY_EIGVAL is buggy in some manner, since it's just getting the same element of NEW_EIGVALUE each time. Another thing I notice is that neither FINAL_ARRAY_EIGVAL nor FINAL_ARRAY_EIGVEC is referenced outside of the loop, making them functionally useless at the moment. If you plan to use them later, they can stay, but otherwise, they are just unnecessary cycles. Finally, concatenate is a very hacky and non idiomatic way to remove one element of an array; I would probably try to use delete() instead.

Since I'm unfamiliar with matrix eigenvector calculations, I unfortunately can't help you much more with optimization, unless I did some more research (and I probably would need some more information about your datasets). Anyway, I have some other suggestions to help you improve your Python coding style in general:

Descriptive variable names are useful

I recommend looking at PEP 8 variable naming conventions. It would improve the readability of your code.

In particular, instead of using names like U and K, I would recommend using more descriptive names. Also, instead of continually phasing names out and deleting them, it may be better to create temp variables to store the new value and then delete those, so the same variable has the same name throughout, but that's a matter of personal style (I'll address this more in the next section.)

Touch-ups

I notice a lot of del statements in your code. This can be avoided by taking advantage of the for-loop scope and not naming one use variables:

U = []
for k in uread('U.csv'):
# ...

S = []
for s in uread('s.csv'):
# ...

EIGVAL = []
for e in uread('EIGVAL.csv'):
# ...

EIGVEC = []
for e in uread('EIGVEC.csv'):
# ...


With statements can use an as clause to specify the name of the file so you don't have to define it on a separate line:

def uread(filename):
with open(filename, 'r') as UFILE:
for i in READER:
yield i


List comprehensions can significantly simplify your code, and eliminate unnecessary (explicit) variable assignment:

Original:

for k in ufinal:
TEMP = []
for j in k:
TEMP.append(eval(j))
U.append(TEMP)
del TEMP
del k


Improved:

for k in ufinal:
U.append([eval(i) for i in k])
del k


Better yet, you can just turn most of your initial loops into list comprehensions:

Original:

U = []
for k in uread('U.csv'):
U.append([eval(i) for i in k])

S = []
for s in uread('s.csv'):
for t in s:
S.append(eval(t))

EIGVAL = []
for e in uread('EIGVAL.csv'):
for l in e:
EIGVAL.append(eval(l))

EIGVEC = []
for e in uread('EIGVEC.csv'):
for l in e:
EIGVEC.append(eval(l))


Improved:

U = [[eval(j) for j in i] for i in uread('U.csv')]
S = [eval(j) for i in uread('s.csv') for j in i]
EIGVAL = [eval(j) for i in uread('EIGVAL.csv') for j in i]
EIGVEC = [eval(j) for i in uread('EIGVEC.csv') for j in i]


Final improved code

Warning: this is untested, so one or two things may be broken (but I've tried my very best to avoid that.) Since I don't know your datasets, I thought I'd leave testing up to you.

# -*- coding: utf-8 -*-
from numpy import array, argmax, delete
import csv

with open(filename, 'r') as UFILE:
for i in READER:
yield i

U = [[eval(j) for j in i] for i in uread('U.csv')]
S = [eval(j) for i in uread('s.csv') for j in i]
EIGVAL = [eval(j) for i in uread('EIGVAL.csv') for j in i]
EIGVEC = [eval(j) for i in uread('EIGVEC.csv') for j in i]

NEW_U = []
NEW_EIGVAL = []
NEW_EIGVEC = []
for i in range(len(S)):
if S[i] > 0.001:
NEW_U.append(U[i])
NEW_EIGVAL.append(EIGVAL[i])
NEW_EIGVEC.append(EIGVEC[i])
del U
del EIGVAL
del EIGVEC
# Select the first r columns of u corresponding to the r principle Eigenvector of
# MatrixForEigenvalues
NEW_EIGVAL = array(NEW_EIGVAL).real
NEW_EIGVEC = array(NEW_EIGVEC)
TRANSPOSE_NEW_EIGVEC = NEW_EIGVEC.real.T   # Done for easy accessibility of matrix
TRANSPOSE_MATRIX_U = array(NEW_U).real.T  # Done for easy accessibility of matrix
FINAL_ARRAY_EIGVAL = []
FINAL_ARRAY_EIGVEC = []
FINAL_MATRIX_U = []
j = argmax(NEW_EIGVAL)
for i in range(len(NEW_EIGVAL)):
FINAL_ARRAY_EIGVAL.append((NEW_EIGVAL)[j])
FINAL_ARRAY_EIGVEC.append(TRANSPOSE_NEW_EIGVEC[j])
FINAL_MATRIX_U.append(TRANSPOSE_MATRIX_U[j])
TRANSPOSE_NEW_EIGVEC = delete(TRANSPOSE_NEW_EIGVEC, j)
TRANSPOSE_MATRIX_U = delete(TRANSPOSE_MATRIX_U, j)
RESULT_MATRIX_U = array(FINAL_MATRIX_U).T     # This is the actual R
with open('r.csv', 'w') as r_file:
WRITER = csv.writer(r_file)
WRITER.writerows(RESULT_MATRIX_U)
print(RESULT_MATRIX_U)

• Thanks Graham. Using delete is not giving all the lists appended to FINAL_MATRIX_U. Instead, it is giving the first list and last item of other list. Also, we cannot use writerows from csv with RESULT_MATRIX_U. – An student Aug 23 '18 at 10:08

It looks like you are creating duplicate numpy.array objects.

for i in range(len(array(NEW_EIGVAL).real)):
j = argmax(array(NEW_EIGVAL).real)
FINAL_ARRAY_EIGVAL.append((array(NEW_EIGVAL).real)[j])
#...


Each iteration creates array(NEW_EIGVAL) twice, yet the value of NEW_EIG is not changing. The first is discarded, but the last must be kept around since a slice of it is being appended into another object. You should instead ...

NEW_EIGVAL_R = array(NEW_EIGVAL).real
for i in range(len(NEW_EIGVAL_R)):
j = argmax(NEW_EIGVAL_R)
FINAL_ARRAY_EIGVAL.append(NEW_EIGVAL_R[j])


Wait. It doesn’t look like j will be changing in the loop either. It could be moved out of the loop too.

Minor: you don’t delete U after creating NEW_U.