Mögliche Optimierung in Cython: numpy Array

Question

Das folgende ist mein Cython-Code zum Zeichnen von multivariaten Normalverteilung. Ich benutze Schleife, weil jedes Mal ich unterschiedliche Dichte habe. (conLSigma ist der Cholesky-Faktor)

Dies dauert eine Menge Zeit, weil ich für jede Schleife inverse und Cholesky-Zerlegung nimmt. Es ist schneller als reiner Python-Code, aber ich habe mich gefragt, ob es irgendeinen Weg gibt, wie ich die Geschwindigkeit steigern kann.

from __future__ import division 

import numpy as np 

cimport numpy as np 

ctypedef np.float64_t dtype_t 

cimport cython 
@cython.boundscheck(False) 
@cython.wraparound(False) 

def drawMetro(np.ndarray[dtype_t, ndim = 2] beta, 
       np.ndarray[dtype_t, ndim = 3] H, 
       np.ndarray[dtype_t, ndim = 2] Sigma, 
       float s): 

    cdef int ncons = betas.shape[0] 
    cdef int nX = betas.shape[1] 
    cdef int con 

    cdef np.ndarray betas_cand = np.zeros([ncons, nX], dtype = np.float64) 
    cdef np.ndarray conLSigma = np.zeros([nX, nX], dtype = np.float64) 

    for con in xrange(ncons): 
     conLSigma = np.linalg.cholesky(np.linalg.inv(H[con] + Sigma)) 
     betas_cand[con] = betas[con] + s * np.dot(conLSigma, np.random.standard_normal(size = nX)) 

    return(betas_cand) 

Answer 1

Die Cholesky-Zerlegung erzeugt eine untere dreieckige Matrix. Dies bedeutet, dass fast die Hälfte der Multiplikationen, die in np.dot erledigt werden, nicht ausgeführt werden müssen. Wenn Sie die Zeile

betas_cand[con] = betas[con] + s * np.dot(conLSigma, np.random.standard_normal(size = nX)) 

in

tmp = np.random.standard_normal(size = nX) 
for i in xrange(nX): 
    for j in xrange(i+1): 
     betas_cand[con,i] += s * conLSigma[i,j] * tmp[j] 

jedoch ändern, werden Sie auch ändern müssen

cdef np.ndarray betas_cand = np.zeros([ncons, nX], dtype = np.float64) 

in

cdef np.ndarray betas_cand = np.array(betas) 

konnte Sie von Scheiben natürlich Nutzung für die Multiplikationen, aber ich bin nicht Sicher, ob es schneller sein wird als ich es vorgeschlagen habe. Wie auch immer, hoffentlich kommst du auf die Idee. Ich glaube nicht, dass es viel mehr gibt, was Sie tun können, um dies zu beschleunigen.

Answer 2

Was ist mit der Berechnung der Cholesky-Zerlegung zuerst und Invertierung der unteren Dreiecksmatrix nach der Rücksubstitution. Dies sollte schneller sein als linalg.cholesky (linalg.inv (S)).