2013-05-22 31 views
15

Çok boyutlu Gaussian'dan rastgele değişkenleri örneklemek ve rasgele değişkenlerin güç spektrumunu hesaplamak için Cholesky ayrışmasını kullanırım. numpy.linalg.cholesky'dan aldığım sonuç her zaman yüksek frekanslarda scipy.linalg.cholesky'dan daha yüksek güce sahiptir.Numune ve scipy'de koleski arasındaki fark nedir?

Bu iki işlev arasındaki bu sonuca neden olabilecek farklılıklar nelerdir? Hangisi daha sayısal olarak kararlı?

n = 2000 

m = 10000 

c0 = np.exp(-.05*np.arange(n)) 

C = linalg.toeplitz(c0) 

Xn = np.dot(np.random.randn(m,n),np.linalg.cholesky(C)) 

Xs = np.dot(np.random.randn(m,n),linalg.cholesky(C)) 

Xnf = np.fft.fft(Xn) 

Xsf = np.fft.fft(Xs) 

Xnp = np.mean(Xnf*Xnf.conj(),axis=0) 

Xsp = np.mean(Xsf*Xsf.conj(),axis=0) 
+0

Scipy faq'den [NumPy ve SciPy arasındaki fark nedir?] (Http://new.scipy.org/faq.html#what-is-the-difference-between-numpy-and-scipy) : "Her durumda, SciPy co lineer cebir modüllerinin daha fazla özellikli versiyonlarının yanı sıra diğer birçok sayısal algoritmayı da kullanır. " Ayrıca bkz. [Neden her ikisi de numpy.linalg' ve scipy.linalg'? Fark nedir?] (Http://new.scipy.org/faq.html#why-both-numpy-linalg-and-scipy-linalg-what-s-the-difference). –

cevap

19

scipy.linalg.choleskynp.linalg.cholesky oysa, size varsayılan olarak üst üçgen ayrışma veriyor size alt üçgen versiyonunu veriyor: İşte

kullandığım koddur. scipy.linalg.cholesky için Dokümanlar: Örneğin

cholesky(a, lower=False, overwrite_a=False) 
    Compute the Cholesky decomposition of a matrix. 

    Returns the Cholesky decomposition, :math:`A = L L^*` or 
    :math:`A = U^* U` of a Hermitian positive-definite matrix A. 

    Parameters 
    ---------- 
    a : ndarray, shape (M, M) 
     Matrix to be decomposed 
    lower : bool 
     Whether to compute the upper or lower triangular Cholesky 
     factorization. Default is upper-triangular. 
    overwrite_a : bool 
     Whether to overwrite data in `a` (may improve performance). 

:

:
>>> scipy.linalg.cholesky([[1,2], [1,9]]) 
array([[ 1.  , 2.  ], 
     [ 0.  , 2.23606798]]) 
>>> scipy.linalg.cholesky([[1,2], [1,9]], lower=True) 
array([[ 1.  , 0.  ], 
     [ 1.  , 2.82842712]]) 
>>> np.linalg.cholesky([[1,2], [1,9]]) 
array([[ 1.  , 0.  ], 
     [ 1.  , 2.82842712]]) 

iki kere ve bunun yerine linalg.cholesky(C,lower=True) kullanmak, sonra ben gibi cevaplar almak aynı rasgele matris uygulamak için kodunuzu değiştirirseniz

>>> Xnp 
array([ 79621.02629287+0.j, 78060.96077912+0.j, 77110.92428806+0.j, ..., 
     75526.55192199+0.j, 77110.92428806+0.j, 78060.96077912+0.j]) 
>>> Xsp 
array([ 79621.02629287+0.j, 78060.96077912+0.j, 77110.92428806+0.j, ..., 
     75526.55192199+0.j, 77110.92428806+0.j, 78060.96077912+0.j]) 
>>> np.allclose(Xnp, Xsp) 
True 
+0

, üst üçgeni numpy'nin cholesky işleviyle hesaplamak mümkün mü? Resmi belge bu soruya yardımcı görünmüyor (https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.linalg.cholesky.html). –

İlgili konular