पायथन

Question

में घने चोलस्की अपडेट मुझे कोई लाइब्रेरी/कोड में इंगित कर सकता है जो मुझे पाइथन (numpy) में Cholesky अपघटन पर निम्न रैंक अपडेट करने की इजाजत देता है? मैटलैब इस कार्यक्षमता को 'चोलुपडेट' नामक फ़ंक्शन के रूप में प्रदान करता है। LINPACK में भी यह कार्यक्षमता है, लेकिन इसमें (मेरे ज्ञान के लिए) अभी तक LAPACK को पोर्ट नहीं किया गया है और इसलिए यह उपलब्ध नहीं है। scipy।

मुझे पता चला कि scikits.sparse CHOLMOD के आधार पर एक समान कार्य प्रदान करता है, लेकिन मेरे matrices घने हैं।

क्या कोई कोड 'cholupdate' की कार्यक्षमता के साथ पाइथन के लिए उपलब्ध है जो numpy के साथ संगत है?

धन्यवाद!

Answer 1

This guy विज्ञान और numpy/scipy का उपयोग कर कुछ ऐसा कर रहा है।

Answer 2

यहाँ एक अजगर पैकेज कि रैंक 1 अद्यतन और Cholesky कारकों पर downdates Cython का उपयोग कर करता है: https://github.com/jcrudy/choldate

उदाहरण:

from choldate import cholupdate, choldowndate 
import numpy 

#Create a random positive definite matrix, V 
numpy.random.seed(1) 
X = numpy.random.normal(size=(100,10)) 
V = numpy.dot(X.transpose(),X) 

#Calculate the upper Cholesky factor, R 
R = numpy.linalg.cholesky(V).transpose() 

#Create a random update vector, u 
u = numpy.random.normal(size=R.shape[0]) 

#Calculate the updated positive definite matrix, V1, and its Cholesky factor, R1 
V1 = V + numpy.outer(u,u) 
R1 = numpy.linalg.cholesky(V1).transpose() 

#The following is equivalent to the above 
R1_ = R.copy() 
cholupdate(R1_,u.copy()) 
assert(numpy.all((R1 - R1_)**2 < 1e-16)) 

#And downdating is the inverse of updating 
R_ = R1.copy() 
choldowndate(R_,u.copy()) 
assert(numpy.all((R - R_)**2 < 1e-16)) 

Answer 3

यह एक रैंक -1 अद्यतन या downdate NumPy सरणी पर क्या करना चाहिए आर और एक्स को '+' या '-' अपडेट या डॉन्डेंडेट के अनुरूप है। (विकिपीडिया पृष्ठ पर MATLAB cholupdate से पोर्ट किया गया: http://en.wikipedia.org/wiki/Cholesky_decomposition):

def cholupdate(R,x,sign): 
    import numpy as np 
     p = np.size(x) 
     x = x.T 
     for k in range(p): 
     if sign == '+': 
      r = np.sqrt(R[k,k]**2 + x[k]**2) 
     elif sign == '-': 
      r = np.sqrt(R[k,k]**2 - x[k]**2) 
     c = r/R[k,k] 
     s = x[k]/R[k,k] 
     R[k,k] = r 
     if sign == '+': 
      R[k,k+1:p] = (R[k,k+1:p] + s*x[k+1:p])/c 
     elif sign == '-': 
      R[k,k+1:p] = (R[k,k+1:p] - s*x[k+1:p])/c 
     x[k+1:p]= c*x[k+1:p] - s*R[k, k+1:p] 
     return R