Exploring data visualization techniques in molecular dynamics: Analysis of dynamic cross-correlation maps as dynamic networks

Dynamic cross correlation maps as dynamic networks

In MD simulations, DCCM analysis is typically restricted to the Alpha Carbons ( C α ) of the protein. Such restriction facilitates the investigation of large molecular systems through dimensionality reduction, while preserving essential structural information. This simplification is standard practice in the field.^6–12 Mathematically, the method is defined as the covariance of three-dimensional displacement vectors normalized by the product of their root-mean-square fluctuations. ⁹ In this context, atomic fluctuations in the same direction relative to the average structure contribute positively to the correlation coefficient, whereas fluctuations in opposite directions contribute negatively. Therefore, this method quantifies the synergistic dynamic behavior of residue pairs relative to their mean positions. ⁹

The spatial arrangement of these vectors complicates interpretation compared to two-dimensional environments, as the third dimension introduces additional complexity to the linear correlations.^7,13 However, normalizing the data by the product of the root-mean-square fluctuations ensures that the final method captures the essential covariances. ⁹ This normalization bounds the coefficients between − 1 and + 1 ,^3,8,9,11,12 highlighting the most significant dynamic contributions. A value of − 1 indicates perfect anti-correlation, where displacement vectors are oriented in opposite directions; this corresponds to concerted motions where residues either move away from each other or approach one another. ¹¹ In contrast, a positive correlation ( + 1 ) indicates that the displacement vectors share the same direction relative to the mean structure, indicating the particles move in a coordinated, similar fashion. ¹¹

The mathematical formulation of the method is derived from the variance (equation (2)) and covariance (equation (1)) of atomic displacements. The complete Dynamic Cross-Correlation Map function is presented in equation (3). In this context, covariance quantifies the joint variability between two three-dimensional vectors, while variance measures the fluctuation of a single vector. In these equations, the angle brackets 〈…〉 denote the ensemble average (or time average) over the simulation trajectory.

c ( i , j ) = 〈 Δ r i · Δ r j 〉

(1)

c ( i , i ) = 〈 Δ r i · Δ r i 〉 = 〈 Δ r i 2 〉

(2)

C ( i , j ) = c ( i , j ) [ c ( i , i ) c ( j , j ) ] 1 / 2 = 〈 Δ r i · Δ r j 〉 〈 Δ r i 2 〉 1 / 2 〈 Δ r j 2 〉 1 / 2

(3)

It is important to note that, although this statistical method reveals linear correlations and anti-correlations, the interpretation of this information in a biological context is the Structural Bioinformatics researcher’s responsibility. Here, we explore the possibilities of visualizing these correlations in terms of time. While Bach ¹⁴ employed correlation networks for visualization with Matrix Cubes in the context of different tasks, including correlation between signal connectivity related to brain functions, in the present study we chose Structural Bioinformatics as the applied knowledge field.

Graph Theory offers a robust framework to improve the understanding of complex biological phenomena, ³ particularly in the interpretation of three-dimensional conformational dynamics. This information can be represented as a dynamic network, where atoms are defined as nodes and their interactions constitute the edges. ³ Graph-theoretic metrics, such as connectivity, modularity, and centrality, are crucial for characterizing the behavior of these molecular complexes. Despite their potential, these visualization strategies remain underutilized in current structural bioinformatics research.

A fundamental representation of network topology is the adjacency matrix, ¹⁵ where rows and columns correspond to the network’s nodes. In this structure, diagonal elements represent self-loops. Each matrix entry indicates the existence of an edge between the corresponding nodes, with the numerical value representing the associated weight. Specifically, this study employs an undirected weighted graph model, which inherently results in a symmetric adjacency matrix.

The use of data visualization within dynamic networks

The taxonomy proposed by Bach et al. (2017) ¹⁶ outlines the main data visualization classes for space-time cubes: extraction, flattening, geometry transformation, and content transformation. The authors discuss the advantages and limitations of this technique, facilitating the description, critical assessment of visualization challenges, and comparison of different tasks applied to temporal data. Regarding extraction operations, point and plane extraction from the volume are performed via orthogonal slicing. Within geometry transformation operations, time translation is addressed. Finally, the content transformation class employs labeling and data filtering techniques. By employing these operations in a composite manner, ¹⁶ it is possible to construct an interactive environment that enables users to explore data efficiently. This approach aligns with the methodology developed by Hong et al., ¹⁷ which emphasizes user interaction within the application.

The work proposed by Hong et al. ¹⁷ seeks to clarify the integration of information across visualizations of varying dimensionality, for instance, enabling the integration of 3D and 2D environments. It maps the usage of distinct visual data representation methods that utilize techniques analogous to those developed by Bach. ⁴ This recent mapping guides the deployment of data visualization techniques, highlighting the design patterns applied in this research.

The literature presents robust techniques for visualizing dynamic networks through the aggregation of time slices in three-dimensional environments. Gohnert et al. ¹⁸ propose 3D DynNetVis, a web-based tool (built with Three.js) that stacks graph layouts along a third temporal dimension, facilitating the inspection of network evolution in accordance with Shneiderman’s “mantra”. ¹⁹ Similarly, Bach et al. ²⁰ introduce MultiPiles, which utilizes the stacking of adjacency matrices. While focusing on dense and weighted networks, such as those observed in brain activity, the authors group similar time slices to reduce redundancy and visual clutter, validating the use of dense time-series weighted networks in complex scientific domains.

Beyond the temporal structure, data comprehension depends on interaction techniques and the identification of topological patterns. MatrixExplorer ²¹ established the importance of coupling matrix and graph (node-link) visualizations in a synchronized manner, allowing for reordering to reveal clusters and connectivities. Recent advances focus on the exploration of changes and explainability: Wen et al. ²² developed DiffSeer, which visually highlights weight differences between adjacent time slices, while Shu et al. ²³ propose interactive approaches to assist users in identifying and interpreting complex visual patterns in matrices, demonstrating superiority over static learning methods. The development of such tools must consider usability barriers and the various stages of exploration. Alkadi et al. ²⁴ define an eight-step framework for Network Information Visualization, encompassing every step from data formatting and import to the final interpretation of network concepts. This methodological workflow facilitates the identification of limitations and informs design decisions, ensuring that the transition from data modeling to visual exploration, whether via matrices or diagrams, meets the researcher’s cognitive needs.

Application of dynamic cross correlation maps in structural bioinformatics

The literature widely employs DCCM as a comparative tool to elucidate functional alterations in biological systems by contrasting distinct molecular states. Wang et al. ¹¹ and Dash et al. ⁸ apply this method to compare native structures (wild-type) with mutant variants, demonstrating how point mutations propagate dynamic effects throughout the structure and may induce anti-correlated motions that are functionally relevant. Extending the analysis to molecular interactions, Parida et al. ¹⁰ employ the technique to assess the impact of ligands, enabling the identification of molecules that stabilize the protein or induce perturbations within its internal motion network.

Beyond static comparisons, more recent and foundational approaches emphasize the temporal dimension and methodological robustness. Okamoto and Ando ⁹ propose Time-Dependent Dynamic Cross-Correlation (TDDCC) to track signal propagation, addressing the need for tools that account for the temporal evolution of data. This perspective is complemented by Kasahara et al., ⁷ who discuss the limitations of standard DCCM in capturing transient states, as well as by the influential work of Hünenberger et al. ⁶ The latter highlights the slow convergence of correlation maps, underscoring the necessity for extended simulations and suggesting that visual tools can facilitate the monitoring of pattern stabilization.

The work by Cavatão et al. ¹² investigates the Melanocortin-1 Receptor (MC1R) through DCCM analysis across six distinct systems, comprising the native protein and mutants in both the presence and absence of ligands. The authors identify how specific mutations alter the transition dynamics between active and inactive states, thereby explaining the loss of molecular function.

Dataset

The dataset used in this research is derived from the study by Cavatao et al., ¹² which investigates the effects of distinct mutations on the mechanism of action of the Melanocortin-1 Receptor (MC1R). It is important to note that the detailed biological investigation of this mechanism is outside the scope of the present work.

The dataset comprises three molecular complexes: the wild-type (WT) and two mutant variants. Each system was simulated under two conditions: the protein structure in isolation (apo state) and the protein bound to a specific ligand (holo state). For each condition, five independent simulation replicas were performed, each spanning 1  μ s (10,000 frames). This yields a total of 30 simulations ( 2 conditions × 3 complexes × 5 replicas ), amounting to 300,000 frames. These data were obtained from the study by Cavatão et al. ¹²

In the work by Cavatão et al., ¹² DCCM analysis was performed using GROMACS ³² through a procedure consisting of the following steps: gmx trjcat was used to concatenate system replicas; gmx trjconv generated a reference file with initial coordinates; and gmx covar calculated the covariance matrix. Subsequently, MDAnalysis^33,34 was employed to normalize the matrix by the product of data variances. The final result is a DCCM derived from the concatenation of all replicas for a single molecular system.

Preprocessing

Data preprocessing is required to prepare the inputs for the visualization tool. This process involves removing solvent molecules and retaining only Alpha Carbon ( C α ) atoms. The GROMACS protocols used for this step are detailed in the accompanying README.md file. The required inputs are the trajectory (.trr) and index (.ndx) files obtained from Cavatão et al. ¹² The output consists of processed files containing only the C α data: a compressed trajectory file (.xtc) and a structure file (.gro), which serve as inputs for the MDTraj library. Finally, trajectory segmentation is performed to calculate the DCCM for each time interval, structuring the data for visual representation.

Computation of sliced dynamic cross correlation maps

The Python-based component processes the simulation data, utilizing the MDTraj library to parse the input files generated during preprocessing. These data comprise topological information regarding the macromolecular structure; specifically, one file stores trajectory coordinates, while the other contains essential system attributes.

The raw data utilized in this work are available at https://linkly.link/2Rkqo. The algorithm automatically detects the number of systems and replicas within a target directory, loading them into memory via the processTrajectory function. This function leverages the MDTraj package, ³⁰ a Python library optimized for the manipulation of MD simulation data. Additionally, the slice size parameter used for segmentation is defined within this routine.

Coordinates corresponding to each simulation frame are extracted using MDTraj and stored in a list. A global alignment of the entire trajectory to a reference structure (the initial frame) is performed using the library’s superpose function. This procedure eliminates global rotational and translational motions by minimizing the RMSD, ensuring that the calculated correlations reflect solely the protein’s internal dynamics. Subsequently, the algorithm partitions the aligned trajectory and computes the standard DCCM for each segment. This calculation is implemented in Python, leveraging the NumPy library ³⁵ for efficient matrix operations.

Time segmentation is performed to capture the temporal evolution inherent in MD trajectories; therefore, the system is modeled by partitioning the simulation into discrete time windows. This approach generates a sequence of graph states, effectively representing the DCCM as a dynamic graph. ¹⁸ Such dynamic network representations provide a robust foundation for advanced data visualization. ³⁶

For each selected trajectory segment, the standard DCCM method is applied. Initially, the mean position of each atom across the segment’s conformations is determined. Subsequently, fluctuation vectors for each atom in every conformation are calculated by subtracting the mean position from the instantaneous positions. The core step involves calculating the covariance matrix (using NumPy’s np.einsum) to measure the coordination of atomic movements. The notation ‘tij,tkj->ik’ defines the tensor contraction rule, instructing the algorithm to multiply the fluctuations of atom ( i ) by those of atom ( k ). The absence of indices t (time) and j (spatial coordinates) in the output implies a summation over these axes. This effectively calculates the dot product of the displacement vectors over the entire trajectory slice, yielding the raw covariance matrix. Finally, this matrix is normalized by dividing the covariance by the product of the standard deviations of the corresponding atoms (derived from the square root of the diagonal variances). This results in a final DCCM matrix with coefficients ranging from − 1 to 1 . An alternative implementation is provided using the CuPy library, ³⁷ allowing for GPU acceleration during correlation calculations.

The application of DCCM analysis to discrete time windows of the trajectory may lead to noisy results. Short windows yield calculations that are not truly significant correlations between atoms, impacting the correct interpretation of the results. In this work, we do not evaluate the significance of the correlation obtained using the selected time windows. Selecting the most accurate trajectory size is relative to each biological focus and is a responsibility of the researcher. Although we assume time-dependent noise, our aim in this work is to provide a method that separates an MD trajectory into discrete time windows and computes DCCM for each slice. More detailed information on Section Constraints on Temporal Granularity, Statistical Robustness, and Method Linearity.

Data storage

Subsequently, the generated data is stored in an optimized format to minimize storage overhead by exploiting the symmetry of the matrix to eliminate redundancy. This optimization facilitates efficient data transfer between the processing module and the visualization interface. The dataset generated for the visualization examples is available at https://linkly.link/2Rljq.

Given the symmetric nature of correlation matrices (where ( i , j ) equals ( j , i ) ), the raw data is highly redundant. To minimize storage requirements, the method extracts only the upper triangular elements of each matrix using the np.triu_indices function from the NumPy/CuPy library. This operation linearizes the matrix, flattening the data into a one-dimensional array. These unique values are subsequently serialized into a compact binary format, significantly reducing the storage footprint. The resulting file consists of a header, which specifies the number of slices and atoms, and a payload, which stores the residue names and the sliced correlation values.

Visualization

To incorporate and visualize DCCM information, the framework adopted in this study adheres to the elementary taxonomy of space-time cube operations proposed by Bach et al. ⁵ We implemented this technique as a web application, leveraging graph theory and data visualization concepts to explore data derived from MD simulations. Figure 1 illustrates this approach utilizing space-time cubes. Specifically, the operations selected from Bach et al., ¹⁶ as described earlier, guided the functional development of the tool.

OPEN IN VIEWER

Figure 1.

Matrix cube visualization.

The visualization interface is implemented in JavaScript utilizing the Three.js library. ³¹ This library facilitates the construction of a 3D environment for rendering scenes derived from the processed data. The software architecture employs a modular design (using functions and classes) to isolate the visual component logic. This approach not only streamlines development and dynamic scene reconstruction but also ensures efficient memory usage on the client device. Functionally, the application retrieves binary data from the server, while the computational burden of rendering is offloaded to the client-side hardware. This architecture effectively mitigates processing overhead on the hosting server.

Binary file reading

The DccmFunctions class, defined in DccmFunctions.js, includes methods responsible for data processing and binary file parsing. Specifically, the asynchronous function loadBinaryDCCM handles the extraction of all data required for the projection construction. Data is read sequentially from the binary file in the following order: the number of segments, the number of protein atoms, the residue names corresponding to the C α atoms, and the segmented correlations. The underlying algorithm calculates the byte offsets to correctly locate these data blocks. Notably, the full correlation matrix is reconstructed from the stored upper triangle during this step, facilitating data interpretation and subsequent visualization.

After parsing the binary file, the asynchronous function returns the number of generated slices, the total number of atoms, an ordered list of residue names, and the correlation matrices for each residue pair. In addition to this data, the function returns two closures, getDCCMValue and getSliceAsMatrix, which facilitate data access and are utilized throughout the visualization code. The getDCCMValue closure retrieves the correlation value between two atoms in a specific slice based on their indices. Meanwhile, getSliceAsMatrix exports a slice in a matrix format optimized for the rendering process.

Tool development

The software is built upon a modular architecture, where a main function defined in SceneManager.js controls the entire scene. Scene subjects are organized in a separate directory and are instantiated on demand by the main function. The SceneManager acts as an application controller, applying configuration changes and updates to the system being visualized.

The application design is inspired by the “Visual Information Seeking Mantra” proposed by Shneiderman. ¹⁹ This principle advocates a specific workflow for constructing data visualizations: overview first, zoom and filter, and then details-on-demand. This paradigm serves as a reference for various data visualization tools^4,18,25,38 as it establishes a foundational taxonomy for the field of Information Visualization.

The update function, defined in main.js, is primarily responsible for executing the rendering cycle. At the conclusion of each cycle, it invokes the renderer.render(scene, camera) method. This method instructs the Three.js WebGLRenderer to draw the current state of the scene to the screen from the camera’s perspective. This operation effectively presents the updated frame to the user.

The scene is initialized with a uniform gray background via the buildScene function. Subsequently, the remaining components, which are managed by the SceneManager, are added. These components include ambient lighting and a color-coded axis indicator. This indicator consists of line segments representing the x, y, and z axes, which intersect at the origin (0,0,0).

The user can manipulate the camera’s rotation and translation via cursor interaction. This functionality is enabled by the Three.js OrbitControls class, instantiated in the buildControls function, which facilitates scene navigation. The camera itself, initialized in the buildCamera function, is an instance of the THREE.PerspectiveCamera class. This class renders the visualization using a perspective projection that updates dynamically in response to user interaction.

A contextual background is added to the visualization. This structure takes the form of a parallelepiped that encapsulates the Matrix Cube, ⁵ maintaining a buffer distance to prevent overlap with the data. This geometric shape is dynamic; its dimensions scale to match the visualization, specifically adjusting its depth based on the number of slices selected for rendering. To ensure the observer’s view is not obstructed during navigation, the planes composing this background are rendered such that only their interior-facing surfaces are visible. Additionally, a subtle black border, defined via custom shader programming, outlines these planes to enhance the spatial perception of the visualization volume.

A configuration panel is created to handle parameter manipulation. This interface is managed by the GUI class and is instantiated within the createPanel function. The panel provides users with control over the visualization through the following ordered elements:

Following the binary file parsing, the data flow is managed by the SceneManager. This manager utilizes the closures defined in loadBinaryDCCM alongside the DccmSlice class. The DccmSlice constructor is critical for the scene’s display; once the raw data is loaded, it is responsible for generating the point cloud. To achieve this, the getSliceAsMatrix closure is invoked to convert the correlation data into attribute buffers formatted for WebGL. These buffers are encapsulated within a Three.js THREE.BufferGeometry object and subsequently instantiated in the scene as a THREE.Points object, establishing the initial visualization.

The updateFromSettings function, triggered by GUI events, orchestrates the dynamic reconfiguration of the visualization. When the user modifies parameters, such as correlation thresholds, the underlying filtering logic is re-executed to recalculate point attributes, including color and visibility. Subsequently, the BufferGeometry attributes are updated in place and flagged for the renderer (setting needsUpdate to true). This efficient approach ensures the next render call displays the filtered points without the performance overhead of destroying and re-instantiating the THREE.Points object.

These rendered points represent the values derived from the symmetric correlation matrices. Their color is assigned by the function gradientColorForCorrelationForParticles, which implements a scalar-to-RGB mapping using piecewise linear interpolation. The visualization employs a diverging color scheme defined by three control points: blue, white, and red. The interpolation is applied in two distinct segments: negative values are interpolated between blue and white, while positive values are interpolated between white and red. This function is executed for each point in the visualization, determining its color based on the correlation value of the specific C α residue pair at the given simulation slice.

Reproducibility and availability

To ensure the software can be reproduced, modified, and audited across different contexts, it is developed as an open-source tool. Consequently, the source code is published on GitHub and released under the GNU General Public License (GPL). This licensing framework promotes broad distribution and guarantees the continued availability and development of the tool as free software.

To further ensure the tool’s reproducibility, the Frontend component employs Docker containerization. This approach encapsulates the application within a self-contained environment, automating the build process and strictly managing the dependency versions required for consistent rendering across different systems.

In contrast, the Backend is designed to operate on pre-processed data, specifically supporting trajectories that contain exclusively the C α atoms of the protein within the analyzed complex. A crucial factor in achieving the research objectives is the web-based availability of this case study. Providing online access allows future investigators to examine both the source code—ensuring reproducibility—and the interactive tool itself. To facilitate this, the Vite library ³⁹ is utilized to build the production environment required for publishing the application’s frontend. Consequently, researchers can grasp the described concepts by accessing the hosted application at https://linkly.link/2RfcK and the source code on GitHub at https://anonymous.4open.science/r/DCCM-Visualization-5B89. This transparency ensures the tool remains open and auditable, serving as a scientific resource for the community.

Overview of the data

To establish the initial data overview, the camera is automatically positioned to frame the visualization centrally within the viewport. This arrangement provides a comprehensive global view of the dataset, as illustrated in Figure 2. Furthermore, this stage supports immediate interaction, allowing the user to explore the visualization using mouse-based navigation controls.

OPEN IN VIEWER

Figure 2.

Overview of the matrix cube visualization applied to sliced DCCM.

Interactive operations

To facilitate data exploration, the tool provides interactive operations aligned with the principles proposed by Shneiderman. ¹⁹ These operations enable the user to manipulate the camera attributes described previously, specifically by allowing rotation, translation, and zooming around a point initially located at the center of the visualization volume.

Data filtering is implemented through two distinct mechanisms. The first is temporal filtering, ¹⁶ which allows the user to isolate and view specific time slices. By utilizing a slider interface, the user can seamlessly traverse the timeline, selecting distinct slices for detailed analysis as exemplified in Figure 3. Crucially, this interaction is optimized to ensure high-performance rendering. Even for systems with a relative large number of residues and slices, the transition occurs without perceptible lag, creating a practical environment for real-time data exploration.

OPEN IN VIEWER

Figure 3.

Visualization of a specific time slice selected from the composite representation, with additional point filtering applied based on normalized correlation thresholds.

The second filtering mechanism is based on correlation magnitude. Implemented within the control panel, this function allows the user to set a threshold using a slider ranging from 0 to 1. This threshold applies to the absolute value of the correlation; any data point with a correlation magnitude below this limit is excluded from the rendering pipeline. This effectively declutters the visualization, isolating and highlighting only those residue pairs with strong correlations. This filtering capability is demonstrated in Figure 3.

Detailed information

To provide details-on-demand, the system employs tooltips, contextual labels that appear when the cursor hovers over a specific data point. The selection mechanism utilizes the THREE.Raycaster class, which projects a ray from the camera through the mouse position to detect intersections with the 3D point cloud. Upon selection, the tooltip displays precise metadata: the correlation coefficient, the identifiers of the interacting residue pair, and the corresponding time slice index, as shown in Figure 4.

OPEN IN VIEWER

Figure 4.

Detail view showing the interactive tooltip mechanism.

To enhance interpretability, residue identifiers are arranged along the base axes of the visualization, aligned with their corresponding columns as shown in Figure 4. Additionally, the display includes a label indicating the specific time window (segmentation) currently being viewed. This label provides essential temporal context, identifying the precise simulation interval during which the observed correlations occurred.

Case study

In the following examples, we explore the application of the proposed method, using traditional DCCM analysis as a baseline for comparison. To enhance visual clarity, a correlation magnitude threshold of 0.5 was applied. This filtering strategy reduces visual clutter, ensuring that only the strongest and most significant correlations are highlighted.

Regarding temporal segmentation, window sizes of 800 and 200 frames were adopted for Figures 6 and 7, respectively. These window sizes directly dictate the visualization’s depth along the time axis; specifically, a smaller frame count per segment results in a larger number of slices representing the data. These values were selected to demonstrate the tool’s flexibility in handling varying levels of temporal granularity and to illustrate the inherent trade-off regarding statistical noise ⁶ :

The standard DCCM analysis employed in the MC1R study by Cavatão et al. ¹² is presented in Figure 5. Crucially, the three-dimensional visualizations developed in the present research are generated from this identical dataset. This shared data foundation allows for a direct comparison between the two methods. Consequently, it becomes possible to identify distinct characteristics that emerge in the proposed visualization, thereby reinforcing the hypothesis that time-resolved analysis reveals correlations typically obscured by standard averaging methods.

OPEN IN VIEWER

Figure 5.

DCCM analysis of the wild-type (WT) protein-ligand complex.

MC1R is a transmembrane protein, and the coordinated movement of its transmembrane (Ts) components and N-terminal (N-Ter) region is associated with the protein’s activation and deactivation. Cavatão et al. ¹² describe this coordination of movements by employing a standard DCCM analysis (Figure 5). Standard DCCM analysis could identify which regions of MC1R are associated with an active state (red and blue points in Figure 5). Comparing the standard DCCM with the sliced approach proposed in this work reveals consistent data patterns, corroborating our results (red and blue points in the same positions (Figure 2). However, while the visualization by Cavatão et al. ¹² presents the data as a static 2D plane, effectively flattening the temporal dimension, the proposed method expands this view into three dimensions. This volumetric approach facilitates the identification of gradients where correlations transition from low to high magnitude, represented as solid volumes formed by temporally contiguous points (Figure 6). Consequently, it becomes possible to observe correlations and anti-correlations coexisting within specific time windows, thereby characterizing complex dynamic behaviors that are often obscured in planar representations.

OPEN IN VIEWER

Figure 6.

Details of the time-resolved DCCM visualization, contrasting with the standard static analysis shown in Figure 5.

Furthermore, a time-based description of the trajectory enables a more specific characterization of protein movements. For instance, Figures 6 and 7 show that T1 and T5 have a gain in correlated movements in the middle of the trajectory. It enables a time-detailed description of the MC1R activation process, which is not accessible in a standard DCCM analysis. Moreover, this approach allows the investigator to characterize related motions not only structurally distant but also temporally distant. As we observed in Figures 6 and 7, anti-correlation movements in the N-Ter region (blue region on base of both figures) start on the first slices, whereas the T5-T1 correlation (red region at center of both figures) starts only on the middle slices. It aggregates insights into MC1R’s time-dependent motion functionality, dependent on further analysis.

OPEN IN VIEWER

Figure 7.

Transition from an uncorrelated state to a period of significant correlation, highlighted by the pink bounding box.

Figure 6 highlights the volumetric clusters formed by the spatial neighborhood of contiguous points—a structural detail revealed only after applying the developed filters. To interpret these structures, the tooltip mechanism is essential, aiding in the precise identification of the constituent residues within a region of interest. Furthermore, the system’s interactive nature provides unconstrained navigational freedom, allowing the investigator to rigorously inspect these regions from multiple angles according to specific research needs.

By employing a higher temporal resolution (more slices) and filtering out weak signals, it is possible to observe distinct transitions between uncorrelated and correlated regions, as shown in Figure 7. This observation holds true for anti-correlations as well. In this specific case study, a 10,000-frame trajectory was segmented into 200-frame windows, yielding 50 discrete slices. Notably, a growing correlation emerges in the region spanning slices 25–36, highlighted by the pink bounding box. This distinct increase in correlation intensity may indicate a significant structural coupling event, such as a conformational change in a specific domain or the propagation of a dynamic effect following ligand accommodation. It is important to note that the same correlation and anti-correlation regions appear in both window-sizes, 800-frames and 200-frames (Figures 6 and 7). It suggests a stability correlation on those regions, once that the bigger statistical convergence in 800-frames found the same pattern in 200-frames window-size.

It is important to highlight that this visualization also aids in the identification of simulation artifacts-errors that may arise during the preparation of the molecular system ¹ or at other stages of the Structural Bioinformatics workflow. These artifacts can compromise the integrity of the research; therefore, their early detection is crucial for producing robust scientific work. In this context, an artifact often manifests as an abrupt, non-fluid transition in the observed correlations within the molecular complex. Ultimately, it is the investigator’s role to interpret these anomalies, determining their biological, physical, or chemical significance versus their nature as methodological errors.

These examples demonstrate that the proposed method enables a more detailed characterization of system dynamics. This capability allows the investigator to pinpoint the specific moments when conformational state transitions occur. Crucially, connecting these dynamic events to protein movements provides deeper insights into the system’s specific biological function.

More dataset exploration

We also applied our approach to a different dataset produced by Pinto et al. ⁴⁰ This dataset is an MD of PET30, an enzyme capable of plastic degradation. The original paper reported five trajectories of PET30 bound to a PET ligand, each describing a different interaction pattern between PET30 and the ligand. These different replica characteristics impact correlation patterns accessed by DCCM. The standard DCCM analysis yields mask correlations among most residues due to replica divergence in the trajectories (Supplemental Figure 1), as described in the original paper. ⁴⁰ However, when we split the first replica (10,000 frames) into 6 slices, we observe high correlation values forming in the last slice (Figures 8 and 9). These correlations inform the movement of residues not necessarily structure-related to the binding pocket. This information may indicate to investigators which interactions are involved in the protein-binding process and why they occur primarily in the final frames. As described by Pinto et al., ⁴⁰ this first replica demonstrates the PET ligand moving outside of the interaction pocket and then returning. Together with our results, this may indicate that the final slice pattern is related to PET30 adjustments for ligand interaction, although further analysis is required. Furthermore, 3D visualization using slice selection helps investigators overcome occlusion and enables a clear view of correlation transitions.

OPEN IN VIEWER

Figure 8.

General sliced PET30 DCCM overview.

OPEN IN VIEWER

Figure 9.

Time-sliced PET30 DCCM view.

Performance and scalability

A critical factor in the development of web-based visualization tools for bioinformatics is scalability with respect to the molecular system size. In the proposed approach, each point rendered in the three-dimensional scene corresponds to a single cell of the correlation matrix. Consequently, the total number of points ( P ) to be processed and rendered grows quadratically relative to the number of protein residues ( N ), obeying the relationship P ∝ N 2 .

For the MC1R system analyzed in this work Case Study section, it consists of 317 residues, each time slice generates a correlation matrix of 317 × 317 , resulting in 100,489 data points (accounting for the redundancy inherent in visualizing the full symmetric matrix). When visualizing multiple slices simultaneously, for instance, 50 slices, the application manages the rendering of approximately 5 million points. Under these conditions, and utilizing the implemented filters, the tool maintains a fluid frame rate, enabling real-time interaction on standard consumer hardware.

Simple performance metrics are shown in the Supplemental material, these were extracted from visualizations that ran in the following hardware specifications: Intel(R) Core(™) i5-9400F CPU @ 2.90 GHz, 24.0 GB of RAM and an NVIDIA GeForce GTX 1050 TI 4 GB . Supplemental Figure 2 shows the impact of threshold interaction on the number of points that are shown in the visualization disaggregated by the number of slices. For each system, a significant number of points are visible even when using a high threshold for correlation between atoms, due to the normalization of those values. The lower number of visible points at higher threshold reduces the computational cost, compared to rendering a larger number of points in the scene.

The augmentation of the time needed to render the scene regarding the augmentation of the number of visible points is given by Supplemental Figure 3, where the points form groups of the same slice size denoted by color, demonstrating that a high number of points takes a considerable time to render. So, for significantly larger molecular complexes, the O ( N 2 ) complexity imposes inherent technical limitations. Proteins exceeding 1000 residues would generate millions of points per individual slice, potentially leading to memory exhaustion (RAM) and GPU- and CPU-bound bottlenecks during the filtering and geometry update stages (BufferGeometry). Although the Three.js library leverages the GPU for efficient rasterization, the client-side pre-processing of color and position attributes remains computationally expensive. Consequently, the current architecture is validated for medium-sized proteins, with visualization of massive complexes requiring additional optimization strategies.

Frames Per Second (FPS) performance is affected by the threshold in this specific MD simulation. The boxplots (Supplemental Figure 4) show that a reasonable FPS value is observed in the visualizations with fewer than 100 slices, across different thresholds. While visualizations with more slices and smaller thresholds compromise FPS, this is important to grasp the idea of how many slices would be the limit for this specific simulation, giving a simple demonstration of when the application might fail. Disaggregating the data by threshold also shows a similar pattern in the boxplots, indicating that as the threshold decreases, the average FPS per slice decreases (i.e. more visible points).

Constraints on temporal granularity, statistical robustness, and method linearity

While the proposed tool facilitates the visual identification of dynamic patterns, it is critical to acknowledge the statistical limitations of calculating correlation matrices over short time windows. Specifically, reduced sampling periods may impede the statistical convergence of the covariance data. Consequently, what appears to be a dynamic correlation event could, in some cases, be an artifact of insufficient sampling rather than a genuine physical coupling.

The validity of the Pearson correlation, which forms the basis of the DCCM method, relies on adequate conformational sampling to ensure that the covariance between atomic displacements is statistically significant, a constraint highlighted by Hünenberger et al. ⁶ By segmenting the trajectory into slices of high granularity (such as 25-frame windows), the analysis incurs the risk that the recorded atomic motions are not representative of the global ensemble. Consequently, this insufficient sampling can introduce noise and generate spurious correlations.

The present work did not apply quantitative validation methods to independently verify the statistical significance of each individual slice. Consequently, the interpretation of the observed “transition correlations” requires caution, and the visualization should be regarded primarily as an exploratory instrument. In this context, the window size serves as a critical hyperparameter; the investigator must carefully adjust this value to strike an optimal balance between the desired temporal resolution and the statistical reliability of the data.

Crucially, while the proposed visualization enhances the temporal perception of the data, it remains subject to the inherent mathematical constraints of the standard DCCM method. Since the underlying algorithm relies on the Pearson correlation coefficient, the analysis is strictly limited to capturing linear relationships between atomic fluctuation vectors. Consequently, non-linear correlated motions remain undetectable by this approach, regardless of the quality of the three-dimensional representation or the granularity of the temporal slicing. Thus, the tool serves to optimize the visual interpretability of linear data.

2D and 3D view

Our 3D visualization, achieved by superimposing slices and changing the camera angle, provides an overall view of temporal slices and reveals general changes in correlation patterns. Furthermore, zoom and slice selection devices enable high-resolution imaging while maintaining detail and overcoming slice occlusion. All these features contribute to a clear understanding of the analysis and to user-directed manipulation. 2D visualizations can also overcome temporal occlusion. However, a side-by-side view of each segment introduces a limitation in the resolution and detail of the data. For instance, Hünenberger et al. ⁶ performed DCCM, separating the MD trajectory into four 200-ps slices and placing them side by side. However, longer MD trajectories would make it difficult to visualize detailed information in several 2D plots or to characterize collective patterns. SmallMultipile overcomes this limitation by stacking slices and enables user manipulation for slice visualization. ²⁰ However, this approach still makes it difficult for users to view individual slices’ contributions to the overall pattern. Although they have limitations, 2D or 3D methods positively contribute to data visualization. The application of those methods will depend on the area and the investigator’s interest.