Function gpfa_log_likelihood 

fn gpfa_log_likelihood(
    y: &[f64],
    c: &[f64],
    d: &[f64],
    r_diag: &[f64],
    k_all: &[Vec<f64>],
    n_neurons: usize,
    n_bins: usize,
    n_latents: usize,
) -> f64

Exact marginal Gaussian log-likelihood via the Woodbury identity.

The marginal covariance Σ = A K Aᵀ + (I_T ⊗ R) is never formed densely. The Woodbury identity and the matrix-determinant lemma route both the quadratic form and the log-determinant through the n_state × n_state posterior precision M = K⁻¹ + AᵀR⁻¹A (Cholesky-factored):

yᵀ Σ⁻¹ y = yᵀ R⁻¹ y − (AᵀR⁻¹y)ᵀ M⁻¹ (AᵀR⁻¹y)
log|Σ|   = log|M| + log|K| + log|R_big|

This is the structured estimator of Yu et al. (2009): the exact marginal likelihood of the regularised model (the GP kernels carry the same 1e-6 jitter as the E-step), not an approximation.