Mise en œuvre du prélèvement multivariée gaussienne avec la matrice de covariance qui est circulante sur Python

Question

Je veux mettre en œuvre la solution indiquée ici: https://math.stackexchange.com/questions/2234019/fastest-way-of-sampling-multivariable-guassian-with-covariance-matrix-that-is-ci/

Lequel cela devrait-je utiliser?
https://docs.scipy.org/doc/scipy/reference/linalg.html

Ai-je besoin d'écrire mon propre code? Quelle méthode donne le temps de calcul le plus efficace?

Answer 1

I traduit une de C++ pour python http://people.sc.fsu.edu/~jburkardt/cpp_src/toeplitz_cholesky/toeplitz_cholesky.html

Ceci calcule la factorisation de Cholesky d'une matrice de Toeplitz définie non négative symétrique.

import numpy as np 


def toeplitz_cholesky_lower(n, g): 
    # G[0][0:N] contains the first row. 
    # G[1][1:N] contains the first column. 

    arr = [[0] * n for _ in range(n)] 
    for entry in range(n): 
     arr[entry][0] = g[0][entry] 
    g[0][1:] = g[0][:-1] 
    g[0][0] = 0.0 

    for i in range(1,n): 
     rho = - g[1][i]/g[0][i] 
     div = np.sqrt((1.0 - rho)*(1.0 + rho)) 
     for j in range(i,n): 
      g1j = g[0][j] 
      g2j = g[1][j] 
      g[0][j] = (g1j + rho*g2j)/div 
      g[1][j] = (rho*g1j + g2j)/div 
     for j in range(i,n): 
      arr[j][i] = g[0][j] 
     for j in range(n-1, i, -1): 
      g[0][j] = g[0][j-1] 
     g[0][i] = 0.0 
    return arr