Computing mit Sparse-Matrix in scipy Python

Question

Ich versuche, das Stoppkriterium described in this paper in einem Richardson-Lucy Dekonvolution Algorithmus zu implementieren. Mein Testbild hat 2.2 Mpx. Ich brauche zu berechnen:

estimator = numpy.dot(P_e_ortho, im_deconv.flatten()) 

wo

im_deconv = [[]] # Image L channel as a 2D array 
m = 2.2E6 # im_deconv.size 

P_e_ortho = scipy.sparse.identity(m, dtype='int8') - \ 
      1/m * np.ones((m, m), dtype='int8') 

Also im Grunde P_e_ortho hat 1 - 1/m auf den Diagonalen und - 1/m anderswo.

Jetzt gibt dieser Code einen Speicherfehler zurück (4.8 × 10² Zellen erforderlich). Wie kann ich vermeiden, mit der ganzen quadratischen Matrix in der Berechnung umzugehen?

scipy.sparse.identity(m, dtype='int8') * (1 - 1/m) 

funktioniert gut, um die Diagonale einzustellen, aber wie nicht diagonale Elemente zu ändern?

Answer 1

Lösung gefunden:

# Construct a Row-based linked list sparse matrix 
P_e_ortho = scipy.sparse.lil_matrix((m, m)) 

P_e_ortho[:, :] = -1/m # Sets all elements 
P_e_ortho.setdiag(1-1/m) # Sets the diagonal 
P_e_ortho = P_e_ortho.tocsr() # Convert to Compressed Sparse Row matrix 

estimator = P_e_ortho.dot(im_deconv.flatten()) # Compute the dot product