Learning stochastic dynamics and predicting emergent behavior using transformers

We show that a neural network originally designed for language processing can learn the dynamical rules of a stochastic system by observation of a single dynamical trajectory of the system, and can accurately predict its emergent behavior under conditions not observed during training. We consider a lattice model of active matter undergoing continuous-time Monte Carlo dynamics, simulated at a density at which its steady state comprises small, dispersed clusters. We train a neural network called a transformer on a single trajectory of the model. The transformer, which we show has the capacity to represent dynamical rules that are numerous and nonlocal, learns that the dynamics of this model consists of a small number of processes. Forward-propagated trajectories of the trained transformer, at densities not encountered during training, exhibit motility-induced phase separation and so predict the existence of a nonequilibrium phase transition. Transformers have the flexibility to learn dynamical rules from observation without explicit enumeration of rates or coarse-graining of configuration space, and so the procedure used here can be applied to a wide range of physical systems, including those with large and complex dynamical generators.

Introduction-Learning the dynamics governing a simulation or experiment is a demanding task because the number of possible dynamical transitions increases exponentially with the physical size of the system.For large systems these transitions are too numerous be enumerated explicitly.Here we show that it is possible to circumvent this restriction by using a neural network to parameterize a stochastic dynamics.In particular, we show that a recently-introduced neural network called a transformer [1], popular in the fields of natural-language processing and computer vision [2][3][4][5][6], can express a general dynamics in an efficient way.It can be trained offline, i.e. by observation only [7], to learn the dynamical rules of a model, even when those rules are numerous and nonlocal.Forward-propagated trajectories of the trained transformer can then be used to reproduce the behavior of the observed model, and to predict its behavior when applied to conditions not seen during training.
Previous work has shown that it is possible to learn the rules of deterministic dynamics, such as deterministic cellular automata [8][9][10], or of stochastic dynamics for small state spaces, using maximum-likelihood estimation on the rates of the generator [11].Similar methods have been used to learn the form of intermolecular potentials that influence the dynamical trajectories of particles undergoing Brownian motion [12][13][14].Several approaches exist in which high-dimensional dynamical systems are approximated by lower-dimensional ones, such as Markov-state models [15,16].In some cases the coarse-graining procedures used to produce such models involve variational methods [17] and neural networks [18].Coarse-graining methods have also been used to obtain deterministic hydrodynamic equations from stochastic trajectories of active matter, allowing for the extraction of hydrodynamic transport coefficients [19,20].Our work complements these approaches by showing that it is possible to learn the dynamical rules of stochastic systems without explicit enumeration of rates or coarse-graining of configuration space, thereby allowing treatment of large and complex systems.From the observation of a single dynamical trajectory a transformer can identify how many classes of process exist and what are their rates, providing physical insight into the dynamics and allowing it to be simulated in new settings, where new phenomena can be discovered.
We focus on the case of a lattice model of active matter, simulated using continuous-time Monte Carlo dynamics [21] (in the Supplemental Information we show that the transformer can be used to treat a second class of model, one realization of which has nonlocal dynamical rules).We allow the transformer to know that the rates for this dynamics are independent of time, and that possible moves consist of single particles rotating in place or translating one lattice site at a time (both restrictions can be relaxed within our framework).However, we do not allow the transformer to know the rates for each move, and, because each rate could in principle depend on the state of the entire system, explicit enumeration of rates would require a generator with many more than 10 100 entries for the system size considered.From observation of a single trajectory of the model, carried out at a density at which its steady state comprises small, dispersed clusters, the transformer learns that particle moves fall into a small number of classes, and accurately determines the associated rates.Forward-propagated trajectories of the trained transformer at the training density reproduce the model's behavior.Moreover, forward-propagated trajectories of the transformer carried out at densities higher than that used in training exhibit motility-induced phase separation [22][23][24][25][26].The details of this phase separation match those of the original model, although that information was not available to the transformer during training.
The trained transformer is therefore able to accurately extrapolate a learned dynamics to predict the existence and details of an emergent phenomenon that it had not previously observed.Given that the transformer is expressive enough to represent a nonlocal dynamics, these results indicate the potential of such devices to learn dynamical rules and study emergent phenomena from observations of dynamical trajectories in a wide variety of settings.
Learning dynamics by observation -Imagine that we are given a dynamical trajectory ω of total time T .The trajectory starts in configuration (microstate) C 0 , and visits K additional configurations C k (Fig. 1(a)).In configuration C k it is resident for time ∆t C k .Schematically, We are told that ω was generated by a dynamics whose rates W C→C for passing between configurations C and C we do not know.We will call this unknown dynamics the original dynamics.
Here we show it is possible to efficiently learn the original dynamics offline, i.e. solely by observation of ω.We start by constructing a synthetic dynamics, which consists of a set of allowed configuration changes {C → C }, which must include those observed in ω, and associated rates W (θ) C→C .Without prior knowledge of the system we should allow the rates for these moves to depend, in principle, on the entire configuration of the system.The number of possible rates grows exponentially with system size, and so treating a system of appreciable size requires the use of an expressive parameterization of the synthetic dynamics.Here we parameterize the rates W (θ) C→C of the synthetic dynamics using the weights θ of a neural network.
One way to learn the original dynamics is to propagate the synthetic dynamics and alter its parameters θ until the dynamical trajectories it generates resemble ω.One drawback of this approach is that original and synthetic dynamics are stochastic, and so comparison of trajectories can be made only in a statistical sense, potentially requiring the generation of many synthetic trajectories at each stage of training.In addition, a comparison of this nature would require the introduction of additional order parameters, different combinations of which may result in different outcomes of training.Instead, we train the syn-  with which it would have generated ω [11].We consider continuous-time Monte Carlo dynamics, in which case see Section S1.Training proceeds by adjusting the parameters θ of the neural network until U (θ) ω no longer increases; see Fig. 1(b).The synthetic dynamics obtained in this way -the learned dynamics -is then the best approximation to the original dynamics that our choice of allowed configuration changes and method of training allows: for a sufficiently long trajectory ω, the dynamics that maximizes U (θ) ω is the original dynamics, W C→C .Other types of dynamics, such as Langevin dynamics, can be studied analogously, by choosing an appropriate replacement for the trajectory log-likelihood.
Original dynamics-The original dynamics we consider is a lattice model of active matter simulated using continuous-time Monte Carlo [21].It consists of a two-dimensional periodic square lattice of size L 2 , occupied by n volume-excluding particles.Each particle α ∈ {1, . . ., n} possesses a unit orientation vector e α that points toward one of the four neighboring sites.The orientation vector of each particle rotates π/2 clockwise or counter-clockwise with rate D. A particle moves to a vacant adjacent lattice site with rate v + if it points toward that lattice site, and with rate v 0 otherwise.The steady state of this model depends on the particle density φ = n/L 2 .At small values of φ, typical configurations consist of small clusters of particles.Upon increasing φ, for sufficiently large v + , the system undergoes the nonequilibrium phase transition called motility-induced phase separation (MIPS).We shall show that the existence of this phase transition can be deduced by observation of a single trajectory obtained at a value of φ at which MIPS is not present.
Training synthetic dynamics -We introduce a general synthetic dynamics using a neural-network architecture called a transformer [1].We allow the transformer to know only that the dynamics is time-independent and consists of single-particle translations and rotations, although these restrictions can be lifted within this framework [27].In microstate C, the transformer represents the transition rates W (θ) C→C to each of the microstates C connected to C through translation or rotation of a single particle (Fig. 1(c)).The transformer learns which particle interactions are relevant to each of these moves, and what their rates are.To train the transformer we perform gradient ascent on its weights using backpropagation in order to maximize the log-likelihood U (θ) ω , Eq. ( 1), with which it would have generated ω.This trajectory is of length T = 5 × 10 3 , on a 30 × 30 lattice, with parameters φ = 0.124, v + = 10, v 0 = 1, and D = 0.1.MIPS is not present at these parameter values; see Fig. 3.
During training we operate the transformer in one of two modes.In Mode 1, the transformer freely predicts ln W W is a hyperparameter that constrains the complexity of the learned dynamics, and provides a measure of the number of dis- W up to a value N W .The value N W provides insight into the structure of the generator of the original dynamics, signaling, for instance, the presence of translational invariance.Section S2 provides additional details of the architectures of both types of neural-network dynamics and their optimization.We have used lattice models in this paper, but the transformer architecture can be directly applied to off-lattice models in any dimension.
In Fig. 2(a During training we did not assume that the dynamical rules are local, nor that some processes (those that violate volume exclusion) are suppressed.The transformer was able to learn both things.If we know that interactions are of finite range then such knowledge can be used to reduce the number of transformer parameters required to learn dynamics (Section S4).Transformers can also learn long-ranged interactions if they are present (Section S5).We also note that learned rates for forbidden processes (inset Fig. 2(b)) are small and decrease with training time, but are not exactly zero: the result is that in forward-propagated trajectories a small fraction of particles can experience overlaps.If volume exclusion is suspected then it can be imposed directly.In addition, with Monte Carlo methods it is possible to determine that the rate of a forbidden process is exactly zero, even given a finite-length training trajectory; see Section S3.
Discovering a phase transition using learned dynamics -In Fig. 3 we show that trajectories generated by the trained transformer can be used to determine the existence of a nonequilibrium phase transition not seen during training.We randomly initialize a configuration at a chosen density φ and propagate the transformer dynamics for fixed time T (Fig. 1(d)).At the training density φ = 0.124, the model's steady state consists of small clusters, but trajectories generated by the transformer at larger values of φ show motility-induced phase separation: the transformer has therefore "discovered" this emergent phenomenon.In Fig. 4 we quantify the details of this phase separation.We measure the fraction of particles with four neighboring occupied sites f 4 , and the variance of that quantity, as well as the number of clusters n c and the average cluster size s c .The time averages of these observables are shown as a function of φ for trajectories obtained with the transformer, both in Mode 1 and Mode 2. For comparison, we show the same quantities from trajectories generated using the original dynamics.The agreement between original and learned dynamics is good, and slightly better using Mode 2, indicating that the transformer, trained under conditions for which no phase separation is observed (see the vertical line in the figure), has predicted the existence and details of a non-equilibrium phase transition.
Conclusions -We have shown that the stochastic dynamics of a many-body system can be efficiently determined using machine-learning tools developed for language processing.A neural network called a transformer can function as an expressive ansatz for the generator of a many-body dynamics, for systems large enough that its possible rates are too numerous to represent explicitly.For instance, for the lattice model of active matter considered here, a 30 × 30 lattice at density φ = 0.1 admits 900 90 ∼ 10 125 arrangements of particles.Each particle takes 1 of 4 rotational states, can move in 4 directions and undergo 2 types of rotation, meaning that there are in principle ∼ 10 182 possible rates.Trained on this model, the transformer learns its dynamics, correctly identifying its local and translationally-invariant nature, and the numerical values of the associated rates.Forward-propagated trajectories of the transformer, carried out at higher densities than that observed during training, show motility-induced phase separation.The details of this nonequilibrium phase transition "discovered" by the transformer agree with those of the original model.Our work shows that it is possible to learn the dynamical rules of stochastic systems without explicit enumeration of rates or coarse-graining of configuration space, complementing existing papers on learning dynamics and pointing the way to the treatment of large and complex systems.

S1. PATH WEIGHT OF A CONTINUOUS-TIME MONTE CARLO DYNAMICS
Consider a dynamical trajectory ω of total time T .The trajectory starts in configuration (microstate) C 0 , and visits K additional configurations C k .In configuration We are told that ω was generated by a continuoustime Monte Carlo dynamics [29,30] whose rates W C→C for passing between configurations C and C we do not know.We will call this unknown dynamics the original dynamics.
In order to learn the original dynamics we introduce a new continuous-time Monte Carlo model called the synthetic dynamics.The synthetic dynamics consists of a set of allowed configuration changes {C → C }, which must include those observed in ω, and associated rates Rates are parameterized by a vector θ = {θ 1 , . . ., θ N } of N numbers (in the main text these numbers corresponds to the weights of the transformer).
We train the synthetic dynamics by maximizing the log-likelihood U where C k →C , the sum running over all transitions allowed from C k .The probability density [31] with which the synthetic dynamics would have chosen the associated residence time ∆t The product of transition-and residence-time factors is Noting that the probability of the final portion of the trajectory, C K the log-likelihood with which the synthetic dynamics would have generated ω is The sum in (S6) is taken over the trajectory ω, i.e. over all configuration changes and corresponding residence times.To train the synthetic dynamics we adjust its parameters θ until (S6) no longer increases.The synthetic dynamics obtained in this way -the learned dynamics -is then the best approximation to the original dynamics that our choice of allowed configuration changes and method of training allows: for a sufficiently long trajectory ω, the dynamics that maximizes (S6) is the original dynamics W C→C .
In sections S3, S4 and S5, we illustrate the training procedure in situations of increasing complexity.

S2. TRANSFORMER AND TRAINING DETAILS
The neural network used to treat the active-matter model described in the main text (and the FA models described in Section S4 and Section S5) is a transformer [1].A transformer relies on an attention mechanism that allows it to learn which parts of a configuration are relevant for a particular process.This generality ensures that it is not biased toward learning local interactions, as is the case for e.g. a convolutional neural network.
The first step in calculating the transition rates is a learned representation of the current state of the system.We first embed particle positions and orientations (or spin states for the FA model) as d h -dimensional vectors using trainable weight matrices; d h is a hyperparameter controlling the expressivity of our neural-net model.We then sum the representations of the position and spin for each particle, which serve as the input to the transformer.We do not impose the boundary conditions of our lattice models; the transformer has to learn these through its positional embedding.For the lattice active matter model, we do not use the empty sites for computational efficiency.Instead, the transformer must learn which neighboring sites are occupied.
Next, we calculate the attention matrix for the entire configuration using multi-head attention [1], and for each particle obtain a d h -dimensional vector containing a weighted sum of the properties of all other particle (the weighting being a measure of the attention paid to each particle).These vectors are then processed using fullyconnected neural networks.We apply this alternating process of attention and application of fully-connected neural networks n l times.The final output of these blocks ✓ T (5) 2 ⇥ 10 6 (11) x/L is used to calculate the transition rate for each possible particle update (spin flip for the FA model or particle rotation or translation for the active-matter model).
Training in Mode 1, the rates are obtained by applying a fully-connected neural network to the output vectors of the transformer.We apply the same network for each particle.This fully-connected neural network has one output node for each possible particle update, the value of ln W assigned to the corresponding transition.Training in Mode 2, we first classify the transformer's output vectors using a fully-connected neural network with N W output nodes and a softmax activation function, again for each possible particle update.The class with the highest probability is sent, as a one-hot vector, to another fully-connected neural network with one output node, which calculates the value of ln W for each of the N W classes. Picking the highest-probability class is not a differentiable operation, and so we use a straightthrough estimator to obtain the gradients to optimize these neural networks [32] .
The results in this paper were obtained with the hyperparameters d h = 64 and n l = 2.We used the AdaBelief optimizer [33] with a learning rate of 10 −4 to optimize the transformer's weights.To obtain a baseline for the trajectory log-likelihood U (θ) ω , we first train a Mode 1 neural-network dynamics on the provided trajectory.For efficiency we train for several epochs on smaller sections of the trajectory; during the final stages of training we use the entire trajectory to obtain more accurate gradients of the trajectory log-likelihood.Next, we train a Mode 2 neural-network dynamics to gain insight into the model's generator.We initialize the first layers of the neural network (the embedding and transformer layers) with the weights obtained with the Mode 1 dynamics, which leads to much faster convergence.

S3. LEARNING A LOCAL DYNAMICS WITH KNOWLEDGE OF ITS LOCALITY
In this section we consider an original dynamics whose rules are local, and assume that we know that its rules are local.The learning procedure then amounts to identifying the correct numerical values for the rates of each local process.This is a conceptually simple case, but worth considering because it illustrates the precision with which rates can be learned, the fact that it is possible to learn the existence of forbidden processes, and also demonstrates some of the convenient features of learning dynamics offline, without propagating new trajectories.
We consider the original dynamics to be the onedimensional Fredrickson-Andersen (FA) model with periodic boundary conditions [34].The FA model is a lattice model of a supercooled liquid whose dynamical rules give rise to slow relaxation and complex space-time behavior [35,36].On each site of a lattice lives a binary spin that can be down (0) or up (1), intended to model immobile and mobile regions of a supercooled liquid.The dynamical rules for the model are nearest-neighbor ones: spins can only flip if at least one of their nearest neighbors is up; if so, down spins flip up with rate c, and up spins flip down with rate 1 − c.We choose c = 0.3, giving the rates shown in Table I.We use the shorthand 001, 101, etc. to denote the 8 possible configuration changes J = 1 (8) We start by generating a single FA model trajectory ω of length T = 10 6 , using a model with 15 lattice sites.This is the trajectory from which we attempt to learn the rates of the model that generated it.A segment of this trajectory of length T /500 is shown in Fig. S1(a).
To construct the synthetic dynamics we make the assumption that only local processes are allowed, i.e. that only one spin at a time can flip.We also assume that the dynamical rules of the model are independent of time.In this section we further assume that the dynamical rules are local, and that these rules are translationally invariant.The synthetic dynamics constrained in this way contains 8 parameters θ = {θ 1 , . . ., θ i , . . ., θ 8 } that correspond to the rates of the 8 processes shown in Table I.
To find the rates that maximize (S6) we proceed as follows.We initialize the rates θ i by choosing them to be random numbers uniformly distributed on (0, 1], and then apply the following Monte Carlo algorithm.At each 0 0.25 S3.As Fig. S1(d), but now the red curve results from calculations in which the value 10 −2 was added to the learned rates 000 and 010 in order to introduce an intentional error.The result -substantially different to the exact answer -shows that comparing the fluctuations of learned and true dynamics is a discriminating test of the learning process, and highlights the precision of that process.
step of the learning procedure we propose new parameters and accept the proposal if (S6) increases or remains the same.This Monte Carlo algorithm is equivalent, for small values of the proposal-size parameter σ, to noisy clipped gradient ascent on the function U (θ) ω , described (for θ i > 0) by the Langevin equation [37] Here n is the number of steps of the learning algorithm, ∇ is the N -dimensional gradient in the coordinates θ, and η is a Gaussian white noise with zero mean and covariance We set σ = 0.05.We carried out 100 independent learning simulations, each begun from different random initial rates θ i , and each trained on the single trajectory ω.Each simulation converged, within about 10 3 steps, to the same value of (S6); we show 5 examples in Fig. S1(b) (colored lines).This value is equal to (S6) evaluated using the true rates of the FA model (black dashed line), although that information was not available to the algorithm during training.The rates produced by one learning simulation are shown in Fig. S1(c) and Table I: the learned rates are numerically close to those of the FA model.
Notably, the learning process has correctly identified that the rates 000 and 010 of the model used to produce the trajectory of Fig. S1(a) are exactly zero.Observing that neither transition occurs in a trajectory of finite length allows us only to bound the rates with which those (12) x/L processes occur.The learning procedure has done better, determining that vanishing rates lead to the largest attainable value of (S6), and has therefore identified the existence of forbidden processes in the FA model.
The non-vanishing rates shown in Fig. S1(c) and Table I are close to those of the FA model, but not identical.A natural question to ask is how well the learned models can approximate the dynamics generated by the original one.Doing so requires the introduction of order parameters.A popular order parameter for the FA model is the activity a = A/T , the number of configuration changes A per unit time T [38,39].Forward-propagated trajectories of duration T = 10 6 of the learned models indeed have the same typical activity as the FA model trajectory, a 0 ≈ 3.2.
However, a more discriminating measure of similarity is to compare the fluctuations of the learned and original dynamics.Fluctuations of the activity can be characterized by the dynamical large-deviation rate function J(a) = lim T →∞ −T −1 ln ρ T (A), a measure of the logarithmic probability ln ρ T (A) of observing, for a trajectory of length T , a particular value of A [40,41].We used the VARD method described in Refs.[42,43] to calculate J(a) for the learned dynamics (we used the neuralnetwork ansatz described in [42]), with each of the 100 learned models constrained to produce a different value of a. Values of J(a) calculated in this way are shown as green circles in Fig. S1(d); they closely approximate the exact rate function of the FA model (black dashed line), which we calculated by numerically diagonalizing the model's rate matrix [41].
The numerical similarity between J(a) of the original and learned models indicates that the precision indicated in Table I and in panels (a) and (b) of Fig. S1 is sufficient to produce an essentially exact description of the behav-ior of the original model, at least as far as the activity is concerned.(In Fig. S3 we show results in which we have made an intentional error by adding 10 −2 to the learned rates 000 and 010; in that case there is a clear discrepancy between the learned and true rate functions.)This comparison also indicates that there is sufficient information in ω to calculate the likelihood with which the original dynamics would produce trajectories never seen in ω.The Conway-Maxwell-Poisson (CMP) function shown in gray in Fig. S1(d) is an upper bound on J(a) inferred by sampling trajectories containing only typical configurations [42,44], such as those seen in Fig. S1(a).The large discrepancy between the CMP bound and the true rate function indicates that rare trajectories of the FA model are dominated by configurations very different to those seen in ω, and therefore never observed during training.Nonetheless, rates inferred from observation of ω can be used to calculate the probability with which long trajectories containing these previously unseen, rare configurations will be observed.
The precision of learning increases with the length of the trajectory ω.In Fig. S2 we show the mean-squared difference between the N = 8 rates of the FA model θ i and the rates θ i learned from an FA model trajectory ω of length T (a single trajectory is used for each value of T ).The precision of the learning process increases with T over the range shown.We end this section by highlighting a notable feature of the "offline" learning process.During training, no trajectories of the synthetic dynamics are propagated.Given a synthetic dynamics we propose a change of its parameters, calculate the likelihood with which that dynamics would have generated the original trajectory ω, and accept the proposal if (S6) does not decrease.Some of the proposals accepted in this way would have been rejected had we trained by comparing ω with forward-propagated trajectories of the synthetic dynamics.In Fig. S4 we propagate some of the synthetic models produced after n steps of training (here we consider an FA model of 100 sites).Next to each trajectory we display a panel comparing the learned rates (blue squares) and true rates (green circles).Here, synthetic rates were initially set to unity.As with the smaller FA model, the training procedure converges to a good approximation of the true rates, and again the rates of the forbidden transitions 000 and 010 are correctly identified to be zero.In general terms the synthetic trajectories look increasingly like those of the FA model as n increases, but there are some notable exceptions.At step n = 80 we encounter an absorbing state of all 0s.This configuration is not accessible to the FA model from any configuration with at least one up spin, and if we were to train by comparison of trajectories generated by original-and synthetic dynamics then this synthetic dynamics would be rejected.However, such a comparison is never made during offline training, and after additional training steps the synthetic dynamics converges to the true dynamics.

S4. LEARNING A LOCAL DYNAMICS WITHOUT KNOWLEDGE OF ITS LOCALITY
In Section S3 we assumed that ω was generated by a dynamics whose rules are local and translationally invariant, allowing us to define a synthetic dynamics with few parameters.If we relax these restrictions then the number of possible rates increases exponentially with system size, and so direct representation of the rates of the synthetic dynamics becomes impractical for large systems.Instead, we can express a general synthetic dynamics using a neural network, which here we take to be a transformer.
We now relax the assumptions of locality and translational invariance in its interaction rules, and so the transformer must determine which features of the configuration are relevant to each process (here we use the version of the FA model whose rates are proportional to the number of nearest neighbors in the up state).We generate a training trajectory ω of length T = 10 6 using a model with N = 100 lattice sites and rate parameter c = 0.3.For each configuration of the trajectory the transformer must calculate the rate of flipping each spin, and therefore has to represent N × 2 N possible transition rates.
We first trained a transformer in Mode 1 (making no restrictions on the number of distinct rate values), and observed that the transformer has the capacity to represent the transition rates for this large state space: the trajectory log-likelihood obtained with the synthetic dynamics rapidly converges to that obtained with the orig- W classes, and a second neural network determines the rate for each of these classes.
In Fig. S5, we show how the optimized trajectory loglikelihood depends on the value of N (θ) W .As the number of distinct rates is increased, the trajectory log-likelihood increases until N (θ) W ≥ 5, at which point it remains constant.The small number of rates required to model these dynamics tells us that the original dynamics is translationally invariant.We also gain insight into the number of distinct processes present in the original dynamics (because multiple processes can have the same transition rate, this value is a lower bound on the number of distinct processes).In this case, the value N W = 5 is equal to the number of rates in the original dynamics; we compare the exact and learned rates in the inset of Fig. S5.
While restricting the interactions to each spin's nearest neighbors limited the number of rates in Section S3 such that they can be modeled explicitly, prior knowledge of the interactions can also be built into our neural-network framework.For instance, if a maximum interaction range is assumed, local attention [45] can be used which limits the attention calculation to a fixed number of neighboring sites, reducing the computational cost and accelerating convergence.

S5. LEARNING A NONLOCAL DYNAMICS
The dynamical rules in the lattice active matter model described in the main text and the FA model described in the previous section are local, depending only on nearestneighbor particles.In this section we demonstrate that the transformer can likewise learn nonlocal dynamical rules, where rates depend on the states of distant particles and the number of distinct rates is large.We consider a conditioned dynamics of the FA model, whose trajectories are biased towards having atypical values of the activity.The generator of this model is given by where W s is obtained by multiplying the off-diagonal elements of the FA-model generator W by e −s , θ(s) is the largest eigenvalue of W s , and L is the corresponding left eigenvector of W s as a diagonal matrix [46,47].Although the sets of forbidden and allowed transitions are the same for both W and W doob s , the transition rates of W doob s depend on the entire lattice configuration.Recently, the conditioned dynamics of several prototypical dynamical systems have been uncovered using neural-network methods [48,49] and reinforcement learning [50,51].
In Fig. S6, we show the result of training a transformer in Mode 1 on a trajectories of conditioned FA models of length T = 10 7 with N = 14 lattice sites.As values for the conditioning field, we use s = s ≈ 0.017 for which the susceptibility χ(s) = θ (s) is maximal for this lattice size, and s = −0.1.Typical trajectories of the conditioned dynamics for s = −0.1 have larger activity than typical trajectories of the unconditioned FA model.Trajectories generated by the trained transformers are shown in Figure S6

FIG. 1 .
FIG. 1. Schematic of our dynamics-learning procedure.(a)We are provided with a trajectory ω, a time series of configurations, and wish to learn the dynamics that created it.For the lattice-based active-matter model studied here, red or blue indicates a particles whose orientation vector points toward an occupied or empty site, respectively.(b) We parameterize a general dynamics using a neural network called a transformer.Rates connecting configurations depend on the weights of the transformer, which are adjusted during training in order to maximize the log-likelihood with which it would have generated ω.(c) The transformer receives the position and orientation of all particles, and must calculate the transition rates to translate or rotate each particle.To do so, it must learn which interactions affect these rates (line thickness denotes attention given to each particle), and their numerical values.(d) Once trained, the neural-network dynamics can be forward-propagated to generate new trajectories, even under conditions not observed in ω.The transformer calculates the rates for all possible transitions C k → {C k }, represented by the blue blobs, at each step.

FIG. 2 .
FIG. 2. (a) Training of a transformer in Mode 1 (unrestricted rates) to maximize the log-likelihood U (θ) ω , Eq. (1), of the training trajectory ω.The horizontal black line denotes the value of the path weight associated with the original model.(b) Dependence of U (θ) ω for a transformer trained in Mode 2 on N (θ) W , the number of distinct classes of move it was asked to identify.This procedure allows us to identify the existence of N W = 4 distinct rates.Inset: Evolution of the rates during training in Mode 2, with N (θ) W = 4.The horizontal black lines denote the values of the rates in the original dynamics.

φ 3 φ =0. 3 φ =0. 5 φ =0. 5 FIG. 3 .
FIG. 3. Trajectories of the lattice active-matter model generated using the dynamics learned by the transformer.The top row shows time-ordered snapshots of a trajectory generated at density φ = 0.124, the value used during training.The two middle rows use densities φ = 0.3 and φ = 0.5; here, motility-induced phase separation can be seen.For comparison, the bottom row shows a trajectory generated with the original dynamics at φ = 0.5.
C→C for each possible transition.In Mode 2, the transformer assigns each transition to one of an integer number N (θ) W of classes, and a second neural network assigns a value ln W (θ) C→C to each class.N (θ)

FIG. 4 .
FIG. 4. Comparison of learned and original dynamics.(a) Time-averaged fraction of particles with four neighboring occupied sites, f4, as a function of density φ, measured over trajectories generated using a transformer trained in Mode 1 (crosses) and Mode 2 (plusses).Angle brackets denote time averages.For each density, 10 trajectories of length T = 10 4 were created using each of the dynamics.Error bars are smaller than the symbol size.Training was done only at φ = 0.124 (vertical dashed line).Squares denote results obtained using the original dynamics.The remaining panels have the same format and show (b) the variance of f4, (c) the number of clusters nc, and (d) the averaged cluster size sc.
) we show the results of training in Mode 1.The trajectory log-likelihood U (θ) ωincreases with the number of observations (epochs) of the trajectory ω, and converges to the value U ω that is obtained using the original dynamics.This value, not available to the transformer during training, indicates that the learned transition rates W (θ) are numerically very close to those of the original dynamics, W .In Fig.2(b) we show the results of training in Mode 2, for several values of N (θ) W .These results show that N W = 4, indicating that the transformer has correctly learned the degree of complexity of the original model, whose dynamical rules are translationally invariant and consist of 4 distinct rates.The inset to Fig. 2(b) shows the evolution with training time of the values of the 4 rates, compared with their values in the original model.
would have generated ω.To calculate U (θ) ω we start by considering the portion C k ∆t C k − −− → C k+1 (S1) of ω, which involves a transition C k → C k+1 and a residence time ∆t C k .The probability with which the synthetic dynamics would have generated the transition

FIG
FIG. S1.(a) Portion of the FA model trajectory ω from which we attempt to learn the FA model dynamical rules.Space is vertical, time is horizontal, and 1s and 0s blue and white, respectively.(b) Log-likelihood (Eq.(S6) divided by trajectory length T ) with which 5 distinct trained synthetic models realize the trajectory shown in panel (a), as a function of the number of training steps n.(c) Rates for one of the trained synthetic dynamics as a function of n.Rates are color-coded according to their values in the original model: blue, green, and red rates have values in the original model of 0.7, 0.3, and 0, respectively.(d) The large-deviation rate function calculated from forward propagation of 100 distinct trained synthetic models (green) compared with the rate function of the true FA model (black).
FIG. S4.Comparison of the true FA model dynamics (top) with the forward-propagated dynamics of the synthetic models produced after n steps of training.Small panels show the 8 learned rates (blue squares) and true rates (green circles).
FIG. S5.Training a transformer in Mode 2 to maximize the log-likelihood U (θ) ω with which it would have generated a trajectory of an FA model with 100 lattice sites.The value ofU (θ) ω is maximized if N θ W ≥ 5,which is equal to the number of distinct rates of the original model.The horizontal black line denotes the exact value.Inset: Evolution of the rates during training, for the case N (θ) W = 5.The horizontal black lines denote the values of the rates in the original dynamics.

-
FIG. S6.Learning the long-ranged dynamics of the conditioned FA model.The dynamical rules of the conditioned generator are nonlocal and consist of a large number of distinct rates.(a) Trajectory of length T = 2000 generated by a transformer trained on a trajectory of the conditioned FA model with N = 14 lattice sites, c = 0.3, and s = s ≈ 0.017.Space is vertical, time is horizontal, and 1s and 0s are blue and white, respectively.(b) Optimization of a transformer to maximize the log-likelihood of the training trajectory.The horizontal black line denotes the value of the path weight associated with the original model.(c) Probability of the lattice containing n spins in state 1 during a trajectory of length T = 10 6 , both for the original (exact conditioned) dynamics and the transformer dynamics.(d) Histogram of the activity a of 1000 trajectories of length T = 10 4 .(e-g) As (a-d), but now for s = −0.1.

TABLE I .
a nearest-neighbor dynamics: each triplet indicates the process in which the central spin changes state.Thus 011 indicates the process 011 → 001, while 000 indicates the process 000 → 010, etc.The rates for the processes 000 and 010 in the FA model are zero, a feature responsible for the model's complex dynamical behavior.Comparison of the true FA model rates and those learned from the trajectory shown in Fig.S1(a). of