Optimization of randomized neural networks for transfer operator approximation

Technical Summary: Optimization of Randomized Neural Networks for Transfer Operator Approximation

Problem Statement
Data-driven approximation of transfer operators (such as the Koopman and Perron–Frobenius operators) is essential for analyzing the global behavior of complex dynamical systems, including metastable states in molecular dynamics and coherent structures in fluid flows. Traditional methods like Extended Dynamic Mode Decomposition (EDMD) rely on fixed dictionaries of basis functions (e.g., polynomials, Gaussians), the optimal selection of which is system-dependent and often an open problem. While neural networks (NNs) can learn these dictionaries, fully trained networks require computationally expensive backpropagation, suffering from issues like slow convergence, local minima, and gradient instability.

Randomized neural networks, such as the previously proposed RaNNDy, address these computational costs by fixing the weights and biases of hidden layers and training only the output layer via a closed-form solution. However, this approach has a critical limitation: the quality of the approximation is entirely dependent on the initial random selection of weights and biases. If the randomly generated basis functions fail to represent the system's eigenfunctions accurately, the resulting operator approximation is poor, and there is no mechanism within the standard RaNNDy framework to correct this.

Methodology
The authors propose an algorithm to optimize the randomized basis of RaNNDy by tuning the hyperparameters of the activation function and the distributions from which weights and biases are sampled, while keeping the network's structural parameters (weights and biases) fixed once sampled.

1. Parametric Activation and Basis: The method introduces tunable parameters $\omega$ (e.g., the scale of the distribution for weights/biases, or parameters of a parametric activation function like $\tanh$) into the basis function definition: $\psi(x, \omega) = \sigma(W(\omega_W)x + b(\omega_b), \omega_a)$.
2. Variational Optimization: The algorithm utilizes the Rayleigh variational principle to define a loss function. For self-adjoint operators (reversible systems), the goal is to maximize the sum of the eigenvalues of the projected operator. For non-self-adjoint operators (non-reversible systems), the method maximizes the trace of the projected forward–backward operator to approximate singular values.
   • The loss function is formulated as $L(\omega) = \text{tr}(\hat{A}(\omega))$, where $\hat{A}(\omega)$ is constructed from data-driven covariance and cross-covariance matrices ($\hat{C}_{00}, \hat{C}_{01}$, etc.).
3. Optimization Loop: The hyperparameters $\omega$ are updated iteratively using gradient descent or Bayesian optimization to maximize the trace of the operator approximation. Once $\omega$ is optimized, the output layer weights are computed via a closed-form eigenvalue problem, preserving the computational efficiency of randomized networks.

Key Contributions

1. Novel Optimization Algorithm: The paper introduces a method to optimize the randomized basis of RaNNDy by tuning activation function hyperparameters and distribution scales, bridging the gap between fully randomized and fully trained networks.
2. Application to Diverse Systems: The approach is validated on three distinct benchmark problems:
   • Random walks on Graphons: Demonstrating the ability to reconstruct graphons and transition densities from data.
   • Bickley Jet: A fluid dynamics benchmark for detecting coherent sets.
   • Protein NuG2: High-dimensional molecular dynamics data to analyze folding–unfolding transitions.

Results
The numerical experiments demonstrate that optimizing the activation parameters significantly improves approximation accuracy compared to using a fixed, poorly initialized random basis:

• Graphons: The algorithm converged in fewer than 20 epochs, correctly identifying a spectral gap corresponding to three metastable sets. The optimized basis allowed for the accurate reconstruction of the graphon and transition density. In contrast, a grid search for fixed distributions was less effective, and fully trained networks (VAMPnets) with higher learning rates diverged or got stuck in local minima.
• Bickley Jet: The method successfully detected coherent sets by clustering the optimized singular functions. The algorithm converged in fewer than 10 iterations, revealing the system's coherent structures more effectively than the initial random basis.
• Protein NuG2: The optimized basis accurately captured the dominant eigenvalues and eigenfunctions, distinguishing between folded and unfolded states of the protein. The initial random basis yielded poor approximations, highlighting the necessity of the proposed optimization.

Significance and Claims
The authors position this work as a "compromise" between fully randomized and fully trained neural networks. The significance lies in retaining the computational advantages of RaNNDy (closed-form solution, no backpropagation for hidden layers, lower training costs) while mitigating its primary weakness: sensitivity to the initial random selection of basis functions.

The paper claims that tuning the distributions of hidden weights and biases is easier due to a simpler loss landscape compared to training full networks, potentially allowing for the use of algorithms like Bayesian optimization. The authors modestly suggest that while the current work focuses on specific activation parameters, future work could explore activation functions with more hyperparameters and conduct systematic runtime comparisons on more complex dynamical systems to further assess the trade-offs of randomized approaches.