Vibrational energy flow in the villin headpiece subdomain: Master equation simulations

Temperature Dependence of Thermal Conductivity of Proteins: Contributions of Thermal Expansion and Grüneisen Parameter

The thermal conductivity of many materials depends on temperature due to several factors, including variation of heat capacity with temperature, changes in vibrational dynamics with temperature, and change in volume with temperature. For proteins some, but not all, of these influences on the variation of thermal conductivity with temperature have been investigated in the past. In this study, we examine the influence of change in volume, and corresponding changes in vibrational dynamics, on the temperature dependence of the thermal conductivity. Using a measured value for the coefficient of thermal expansion and recently computed values for the Grüneisen parameter of proteins we find that the thermal conductivity increases with increasing temperature due to change in volume with temperature. We compare the impact of thermal expansion on the variation of the thermal conductivity with temperature found in this study with contributions of heat capacity and anharmonic coupling examined previously. Using values of thermal transport coefficients computed for proteins we also model heating of water in a protein solution following photoexcitation.

Thermal Energy Transport through Nonbonded Native Contacts in Protein

Within the protein interior, where we observe various types of interactions, nonuniform flow of thermal energy occurs along the polypeptide chain and through nonbonded native contacts, leading to inhomogeneous transport efficiencies from one site to another. The folded native protein serves not merely as thermal transfer medium but, more importantly, as sophisticated molecular nanomachines in cells. Therefore, we are particularly interested in what sort of "communication" is mediated through native contacts in the folded proteins and how such features are quantitatively depicted in terms of local transport coefficients of heat currents. To address the issue, we introduced a concept of inter-residue thermal conductivity and characterized the nonuniform thermal transport properties of a small globular protein, HP36, using equilibrium molecular dynamics simulation and the Green-Kubo formula. We observed that the thermal transport of the protein was dominated by that along the polypeptide chain, while the local thermal conductivity of nonbonded native contacts decreased in the order of H-bonding > π-stacking > electrostatic > hydrophobic contacts. Furthermore, we applied machine learning techniques to analyze the molecular mechanism of protein thermal transport. As a result, the contact distance, variance in contact distance, and H-bonding occurrence probability during MD simulations are found to be the top three important determinants for predicting local thermal transport coefficients.

Vibrational Energy Landscapes and Energy Flow in GPCRs

We construct and analyze disconnectivity graphs to provide the first graphical representation of the vibrational energy landscape of a protein, in this study β2AR, a G-protein coupled receptor (GPCR), in active and inactive states. The graphs, which indicate the relative free energy of each residue and the minimum free energy barriers for energy transfer between them, reveal important composition, structural and dynamic properties that mediate the flow of energy. Prolines and glycines, which contribute to GPCR plasticity and function, are identified as bottlenecks to energy transport along the backbone from which alternative pathways for energy transport via nearby noncovalent contacts emerge, seen also in the analysis of first passage time (FPT) distributions presented here. Striking differences between the disconnectivity graphs and FPT distributions for the inactive and active states of β2AR are found where structural and dynamic changes occur upon activation, contributing to allosteric regulation.

Intra-protein interacting collective modes in the terahertz frequency region

Understanding how proteins work requires a thorough understanding of their internal dynamics. Proteins support a wide range of motions, from the femtoseconds to seconds time scale, relevant to crucial biological functions. In this context, the term "protein collective dynamics" refers to the complex patterns of coordinated motions of numerous atoms throughout the protein in the sub-picosecond time scale (terahertz frequency region). It is hypothesized that these dynamics have a substantial impact on the regulation of functional dynamical mechanisms, including ligand binding and allosteric signalling, charge transport direction, and the regulation of thermodynamic and thermal transport properties. Using the theoretical framework of hydrodynamics, the collective dynamics of proteins had previously been described in a manner akin to that of simple liquids, i.e. in terms of a single acoustic-like excitation, related to intra-protein vibrational motions. Here, we employ an interacting-mode model to analyse the results from molecular dynamics simulations and we unveil that the vibrational landscape of proteins is populated by multiple acoustic-like and low-frequency optic-like modes, with mixed symmetry and interfering with each other. We propose an interpretation at the molecular level of the observed scenario that we relate to the side-chains and the hydrogen-bonded networks dynamics. The present insights provide a perspective for understanding the molecular mechanisms underlying the energy redistribution processes in the interior of proteins.

Locating dynamic contributions to allostery via determining rates of vibrational energy transfer

Determining rates of energy transfer across non-covalent contacts for different states of a protein can provide information about dynamic and associated entropy changes during transitions between states. We investigate the relationship between rates of energy transfer across polar and nonpolar contacts and contact dynamics for the β₂-adrenergic receptor, a rhodopsin-like G-protein coupled receptor, in an antagonist-bound inactive state and agonist-bound active state. From structures sampled during molecular dynamics (MD) simulations, we find the active state to have, on average, a lower packing density, corresponding to generally more flexibility and greater entropy than the inactive state. Energy exchange networks (EENs) are computed for the inactive and active states from the results of the MD simulations. From the EENs, changes in the rates of energy transfer across polar and nonpolar contacts are found for contacts that remain largely intact during activation. Change in dynamics of the contact, and entropy associated with the dynamics, can be estimated from the change in rates of energy transfer across the contacts. Measurement of change in the rates of energy transfer before and after the transition between states thereby provides information about dynamic contributions to activation and allostery.

Time-resolved spectroscopic mapping of vibrational energy flow in proteins: Understanding thermal diffusion at the nanoscale

Vibrational energy exchange between various degrees of freedom is critical to barrier-crossing processes in proteins. Hemeproteins are well suited for studying vibrational energy exchange in proteins because the heme group is an efficient photothermal converter. The released energy by heme following photoexcitation shows migration in a protein moiety on a picosecond timescale, which is observed using time-resolved ultraviolet resonance Raman spectroscopy. The anti-Stokes ultraviolet resonance Raman intensity of a tryptophan residue is an excellent probe for the vibrational energy in proteins, allowing the mapping of energy flow with the spatial resolution of a single amino acid residue. This Perspective provides an overview of studies on vibrational energy flow in proteins, including future perspectives for both methodologies and applications.

Energy Transport and Its Function in Heptahelical Transmembrane Proteins

Photoproteins such as bacteriorhodopsin (bR) and rhodopsin (Rho) need to effectively dissipate photoinduced excess energy to prevent themselves from damage. Another well-studied seven transmembrane (TM) helices protein is the β₂ adrenergic receptor (β₂AR), a G protein-coupled receptor for which energy dissipation paths have been linked with allosteric communication. To study the vibrational energy transport in the active and inactive states of these proteins, a master equation approach [J. Chem. Phys.2020, 152, 045103] is employed, which uses scaling rules that allow us to calculate energy transport rates solely based on the protein structure. Despite their overall structural similarity, the three 7TM proteins reveal quite different strategies to redistribute excess energy. While bR quickly removes the energy using the TM7 helix as a "lightning rod", Rho exhibits a rather poor energy dissipation, which might eventually require the hydrolysis of the Schiff base between the protein and the retinal chromophore to prevent overheating. Heating the ligand adrenaline of β₂AR, the resulting energy transport network of the protein is found to change significantly upon switching from the active state to the inactive state. While the energy flow may highlight aspects of the inter-residue couplings of β₂AR, it seems not particularly suited to explain allosteric phenomena.

Energy Transport in Class B GPCRs: Role of Protein-Water Dynamics and Activation

We compute energy exchange networks (EENs) through glucagon-like peptide-1 receptor (GLP-1R), a class B G-protein-coupled receptor (GPCR), in inactive and two active states, one activated by a peptide ligand and the other by a small molecule agonist, from results of molecular dynamics simulations. The reorganized network upon activation contains contributions from structural as well as from dynamic changes and corresponding entropic contributions to the free energy of activation, which are estimated in terms of the change in rates of energy transfer across non-covalent contacts. The role of water in the EENs and in activation of GLP-1R is also investigated. The dynamics of water in contact with the central polar network of the transmembrane region is found to be significantly slower for both activated states compared to the inactive state. This result is consistent with the contribution of water molecules to activation of GLP-1R previously suggested and resembles water dynamics in parts of the transmembrane region found in earlier studies of rhodopsin-like GPCRs.

Computational Study on the Thermal Conductivity of a Protein

Protein molecules are thermally fluctuating and tightly packed amino acid residues strongly interact with each other. Such interactions are characterized in terms of heat current at the atomic level. We calculated the thermal conductivity of a small globular protein, villin headpiece subdomain, based on the linear response theory using equilibrium molecular dynamics simulation. The value of its thermal conductivity was 0.3 ± 0.01 [W m^-1 K^-1], which is in good agreement with experimental and computational studies on the other proteins in the literature. Heat current along the main chain was dominated by local vibrations in the polypeptide bonds, with amide I, II, III, and A bands on the Fourier transform of the heat current autocorrelation function.

Energy Bilocalization Effect and the Emergence of Molecular Functions in Proteins

Proteins are among the most complex molecular structures, which have evolved to develop broad functions, such as energy conversion and transport, information storage and processing, communication, and regulation of chemical reactions. However, the mechanisms by which these dynamical entities coordinate themselves to perform biological tasks remain hotly debated. Here, a physical theory is presented to explain how functional dynamical behavior possibly emerge in complex/macro molecules, thanks to the effect that we term bilocalization of thermal vibrations. More specifically, our approach allows us to understand how structural irregularities lead to a partitioning of the energy of the vibrations into two distinct sets of molecular domains, corresponding to slow and fast motions. This shape-encoded spectral allocation, associated to the genetic sequence, provides a close access to a wide reservoir of dynamical patterns, and eventually allows the emergence of biological functions by natural selection. To illustrate our approach, the SPIKE protein structure of SARS-COV2 is considered.

Activation-Induced Reorganization of Energy Transport Networks in the β2 Adrenergic Receptor

We compute energy exchange networks (EENs) through the β₂ adrenergic receptor (β₂AR), a G-protein coupled receptor (GPCR), in inactive and active states, based on the results of molecular dynamics simulations of this membrane bound protein. We introduce a new definition for the reorganization of EENs upon activation that depends on the relative change in rates of energy transfer across noncovalent contacts throughout the protein. On the basis of the reorganized network that we obtain for β₂AR upon activation, we identify a branched pathway between the agonist binding site and the cytoplasmic region, where a G-protein binds to the receptor when activated. The pathway includes all of the motifs containing molecular switches previously identified as contributing to the allosteric transition of β₂AR upon agonist binding. EENs and their reorganization upon activation are compared with structure-based contact networks computed for the inactive and active states of β₂AR.

Energy Transfer across Nonpolar and Polar Contacts in Proteins: Role of Contact Fluctuations

Molecular dynamics simulations of the villin headpiece subdomain HP36 have been carried out to examine relations between rates of vibrational energy transfer across non-covalently bonded contacts and equilibrium structural fluctuations, with focus on van der Waals contacts. Rates of energy transfer across van der Waals contacts vary inversely with the variance of the contact length, with the same constant of proportionality for all nonpolar contacts of HP36. A similar relation is observed for hydrogen bonds, but the proportionality depends on contact pairs, with hydrogen bonds stabilizing the α-helices all exhibiting the same constant of proportionality, one that is distinct from those computed for other polar contacts. Rates of energy transfer across van der Waals contacts are found to be up to 2 orders of magnitude smaller than rates of energy transfer across polar contacts.

A Time-Dependent Quantum Approach to Allostery and a Comparison With Light-Harvesting in Photosynthetic Phenomenon

The allosteric effect is one of the most important processes in regulating the function of proteins, and the elucidation of this phenomenon plays a significant role in understanding emergent behaviors in biological regulation. In this process, a perturbation, generated by a ligand in a part of the macromolecule (the allosteric site), moves along this system and reaches a specific (active) site, dozens of Ångströms away, with a great efficiency. The dynamics of this perturbation in the macromolecule can model precisely the allosteric process. In this article, we will be studying the general characteristics of allostery, using a time-dependent quantum approach to obtain rules that apply to this kind of process. Considering the perturbation as a wave that moves within the molecular system, we will characterize the allosteric process with three of the properties of this wave in the active site: (1) t_a, the characteristic time for reaching that site, (2) A_a, the amplitude of the wave in this site, and (3) B_a, its corresponding spectral broadening. These three parameters, together with the process mechanism and the perturbation efficiency in the process, can describe the phenomenon. One of the main purposes of this paper is to link the parameters t_a, A_a, and B_a and the perturbation efficiency to the characteristics of the system. There is another fundamental process for life that has some characteristics similar to allostery: the light-harvesting (LH) process in photosynthesis. Here, as in allostery, two distant macromolecular sites are involved—two sites dozens of Ångströms away. In both processes, it is particularly important that the perturbation is distributed efficiently without dissipating in the infinite degrees of freedom within the macromolecule. The importance of considering quantum effects in the LH process is well documented in literature, and the quantum coherences are experimentally proven by time-dependent spectroscopic techniques. Given the existing similarities between these two processes in macromolecules, in this work, we suggest using Quantum Mechanics (QM) to study allostery.

Energy flow and intersubunit signalling in GSAM: A non-equilibrium molecular dynamics study

Non-equilibrium molecular dynamics simulations of vibrational energy flow induced by the imposition of a thermal gradient have been performed on the μ2-dimeric enzyme glutamate-1-semialdehyde aminomutase (GSAM), the key enzyme in the biosynthesis of chlorophyll, in order to identify energy transport pathways and to elucidate their role as potential allosteric communication networks for coordinating functional dynamics, specifically the negative cooperativity observed in the motion of the two active site gating loops. Fully atomistic MD simulations of thermal diffusion were executed with a GROMACS simulation package on a fully solvated GSAM enzyme by heating various active site target ligands (initially, catalytic intermediates and cofactors) to 300K while holding the remainder of the protein and the solvent bath at 10K and monitoring the temperature T(t) of all the enzyme residues as a function of time over a 1ns observation window. Energy is observed to be deposited in a relatively small number of discrete chains of residues most of which contribute to specific structural or biochemical functionality. Thermal linkages between all thermally active chains were established by isolating a specific pair of chains and performing a thermal diffusion simulation on the pair, one held at 300K and the other at 10K, with the rest of the protein frozen in its initial atomic configuration and thus thermally unresponsive. Proceeding in this way, it was possible to map out multiple pathways of vibrational energy flow leading from one of the active sites through a network of contiguous residues, many of which were evolutionarily conserved and linked by hydrogen bonds, into the other active site and ultimately to the other gating loop.

Role of atomic contacts in vibrational energy transfer in myoglobin

Heme proteins are ideal systems to investigate vibrational energy flow at the atomic level. Upon photoexcitation, a large amount of excess vibrational energy is selectively deposited on heme due to extremely fast internal conversion. This excess energy is redistributed to the surrounding protein moiety and then to water. Vibrational energy flow in myoglobin (Mb) was examined using picosecond time-resolved anti-Stokes ultraviolet resonance Raman (UVRR) spectroscopy. We used the Trp residue directly contacting the heme group as a selective probe for vibrationally excited populations. Trp residues were placed at different position close to the heme by site-directed mutagenesis. This technique allows us to monitor the excess energy on residue-to-residue basis. Anti-Stokes UVRR measurements for Mb mutants suggest that the dominant channel for energy transfer in Mb is the pathway through atomic contacts between heme and nearby amino acid residues as well as that between the protein and solvent water. It is found that energy flow through proteins is analogous to collisional exchange processes in solutions.

Recent developments in the computational study of protein structural and vibrational energy dynamics

Recent developments in the computational study of energy transport in proteins are reviewed, including advances in both methodology and applications. The concept of energy exchange network (EEN) is discussed, and a recent calculation of EENs for the allosteric protein FixL is reviewed, which illustrates how residues and protein regions involved in the allosteric transition can be identified. Recent work has examined relations between EENs and protein dynamics as well as structure. We review some of the computational studies carried out on several proteins that explore connections between energy conductivity across polar contacts in proteins and between proteins and water and equilibrium dynamics of the contacts, and we discuss some of the recent experimental work that addresses this topic.

Normal mode analysis and beyond

Normal mode analysis provides a powerful tool in biophysical computations. Particularly, we shed light on its application to protein properties because they directly lead to biological functions. As a result of normal mode analysis, the protein motion is represented as a linear combination of mutually independent normal mode vectors. It has been widely accepted that the large amplitude motions throughout the entire protein molecule can be well described with a few low-frequency normal modes. Furthermore, it is possible to represent the effect of external perturbations, e.g., ligand binding, hydrostatic pressure, as the shifts of normal mode variables. Making use of this advantage, we are able to explore mechanical properties of proteins such as Young's modulus and compressibility. Within thermally fluctuating protein molecules under physiological conditions, tightly packed amino acid residues interact with each other through heat and energy exchanges. Since the structure and dynamics of protein molecules are highly anisotropic, the flow of energy and heat should also be anisotropic. Based on the harmonic approximation of the heat current operator, it is possible to analyze the communication map of a protein molecule. By using this method, the energy transfer pathways of photoactive yellow protein were calculated. It turned out that these pathways are similar to those obtained via the Green-Kubo formalism with equilibrium molecular dynamics simulations, indicating that normal mode analysis captures the intrinsic nature of the transport properties of proteins.

Topology, landscapes, and biomolecular energy transport

While ubiquitous, energy redistribution remains a poorly understood facet of the nonequilibrium thermodynamics of biomolecules. At the molecular level, finite-size effects, pronounced nonlinearities, and ballistic processes produce behavior that diverges from the macroscale. Here, we show that transient thermal transport reflects macromolecular energy landscape architecture through the topological characteristics of molecular contacts and the nonlinear processes that mediate dynamics. While the former determines transport pathways via pairwise interactions, the latter reflects frustration within the landscape for local conformational rearrangements. Unlike transport through small-molecule systems, such as alkanes, nonlinearity dominates over coherent processes at even quite short time- and length-scales. Our exhaustive all-atom simulations and novel local-in-time and space analysis, applicable to both theory and experiment, permit dissection of energy migration in biomolecules. The approach demonstrates that vibrational energy transport can probe otherwise inaccessible aspects of macromolecular dynamics and interactions that underly biological function.

Understanding vibrational energy transfer in macromolecules has been challenging to both theory and experiment. Here the authors use non-equilibrium molecular dynamics to reveal the relationship between heat transport in a model peptide, emergent nonlinearity, and the underlying free energy landscape.

Energy Transport Pathways in Proteins: A Non-equilibrium Molecular Dynamics Simulation Study

To facilitate the observation of biomolecular energy transport in real time and with single-residue resolution, recent experiments by Baumann et al. ( Angew. Chem. Int. Ed. 2019 , 58 , 2899 , DOI: 10.1002/anie.201812995 ) have used unnatural amino acids β-(1-azulenyl)alanine (Azu) and azidohomoalanine (Aha) to site-specifically inject and probe vibrational energy in proteins. To aid the interpretation of such experiments, non-equilibrium molecular dynamics simulations of the anisotropic energy flow in proteins TrpZip2 and PDZ3 domains are presented. On this account, an efficient simulation protocol is established that accurately mimics the excitation and probing steps of Azu and Aha. The simulations quantitatively reproduce the experimentally found cooling times of the solvated proteins at room temperature and predict that the cooling slows by a factor 2 below the glass temperature of water. In PDZ3, vibrational energy is shown to travel from the initially excited peptide ligand via a complex network of inter-residue contacts and backbone transport to distal regions of the protein. The supposed connection of these energy transport pathways with pathways of allosteric communication is discussed.

Vibrational Energy in Proteins Correlates with Topology

The exchange of vibrational energy in proteins is crucial for their function. Here, we establish a connection between quantities related to it with geometry-based properties such as the proteins' residues coordination number. This relation is proven by molecular simulation in a neuro-pharmacologically relevant transmembrane receptor. The connection demonstrated here paves the way to studies of protein allostery and conformational changes based solely on protein structure.

Computational study to understand the energy transfer pathways within amicyanin

Vibrational energy diffusivities between the residues present in Amicyanin copper protein are calculated and presented in form of communication map. From those results energy flow pathways from the copper metal ion to the inter protein residue Glu31 are identified. Our finding suggests many different pathways are possible and copper metal ion in oxidized and reduced state switches the pathways. Our finding also suggests the cooperative nature of surrounding residues and water molecules towards selecting the pathways. The major transport channels in the oxidised state are, Cu²⁺---> MET28---> LYS29---> TYR30---> GLU31 and Cu²⁺---> MET98---> TYR30--- > GLU31. And in the reduced state Cu⁺---> CYS9---> TYR30---> GLU31 and Cu⁺---> MET28---> LYS2---> TYR30---> GLU31. We studied further the interaction energies between the copper ion and neighbouring residues using B3LYP/QZVP method. Both the methods complement each other in predicting the energy flow pathways and the cooperative nature of residues.

Protein structural robustness to mutations: an in silico investigation

Proteins possess qualities of robustness and adaptability to perturbations such as mutations, but occasionally fail to withstand them, resulting in loss of function. Herein, the structural impact of mutations is investigated independently of the functional impact. Primarily, we aim at understanding the mechanisms of structural robustness pre-requisite for functional integrity. The structural changes due to mutations propagate from the site of mutation to residues much more distant than typical scales of chemical interactions, following a cascade mechanism. This can trigger dramatic changes or subtle ones, consistent with a loss of function and disease or the emergence of new functions. Robustness is enhanced by changes producing alternative structures, in good agreement with the view that proteins are dynamic objects fulfilling their functions from a set of conformations. This result, robust alternative structures, is also coherent with epistasis or rescue mutations, or more generally, with non-additive mutational effects and compensatory mutations. To achieve this study, we developed the first algorithm, referred to as Amino Acid Rank (AAR), which follows the structural changes associated with mutations from the site of the mutation to the entire protein structure and quantifies the changes so that mutations can be ranked accordingly. Assessing the paths of changes opens the possibility of assuming secondary mutations for compensatory mechanisms.

Water-Mediated Energy Dynamics in a Homodimeric Hemoglobin

We examine energy dynamics in the unliganded and liganded states of the homodimeric hemoglobin from Scapharca inaequivalvis (HbI), which exhibits cooperativity mediated by the cluster of water molecules at the interface upon ligand binding and dissociation. We construct and analyze a dynamic network in which nodes representing the residues, hemes, and water cluster are connected by edges that represent energy transport times, as well as a nonbonded network (NBN) indicating regions that respond rapidly to local strain within the protein via nonbonded interactions. One of the two largest NBNs includes the Lys30-Asp89 salt bridge critical for stabilizing the dimer. The other includes the hemes and surrounding residues, as well as, in the unliganded state, the cluster of water molecules between the globules. Energy transport in the protein appears to be controlled by the Lys30-Asp89 salt bridge critical for stabilizing the dimer, as well as the interface water cluster in the unliganded state. Possible connections between energy transport dynamics in response to local strain identified here and allosteric transitions in HbI are discussed.

Scaling Rules for Vibrational Energy Transport in Globular Proteins

Computational studies of vibrational energy flow in biomolecules have to date mapped out transport pathways on a case-by-case basis. To provide a more general approach, we derive scaling rules for vibrational energy transport in a globular protein, which are identified from extensive nonequilibrium molecular dynamics simulations of vibrational energy flow in the villin headpiece subdomain HP36. We parametrize a master equation based on inter-residue, residue-solvent, and heater-residue energy-transfer rates, which closely reproduces the results of the all-atom simulations. From that fit, two scaling rules emerge, one for energy transport along the protein backbone which relies on a diffusion model and another for energy transport between tertiary contacts, which is based on a harmonic model. Requiring only the calculation of mean and variance of relatively few atomic distances, the approach holds the potential to predict the pathways and time scales of vibrational energy flow in large proteins.

Energy exchange network of inter-residue interactions within a thermally fluctuating protein molecule: A computational study

Protein function is regulated not only by the structure but also by physical dynamics and thermal fluctuations. We have developed the computer program, CURrent calculation for proteins (CURP), for the flow analysis of physical quantities within thermally fluctuating protein media. The CURP program was used to calculate the energy flow within the third PDZ domain of the neuronal protein PSD-95, and the results were used to illustrate the energy exchange network of inter-residue interactions based on atomistic molecular dynamics simulations. The removal of the α3 helix is known to decrease ligand affinity by 21-fold without changing the overall protein structure; nevertheless, we demonstrated that the helix constitutes an essential part of the network graph.

Allosteric communication pathways and thermal rectification in PDZ-2 protein: a computational study

Allosteric communication in proteins is a fundamental and yet unresolved problem of structural biochemistry. Previous findings, from computational biology ( Ota, N.; Agard, D. A. J. Mol. Biol. 2005 , 351 , 345 - 354 ), have proposed that heat diffuses in a protein through cognate protein allosteric pathways. This work studied heat diffusion in the well-known PDZ-2 protein, and confirmed that this protein has two cognate allosteric pathways and that heat flows preferentially through these. Also, a new property was also observed for protein structures: heat diffuses asymmetrically through the structures. The underling structure of this asymmetrical heat flow was a normal length hydrogen bond (∼2.85 Å) that acted as a thermal rectifier. In contrast, thermal rectification was compromised in short hydrogen bonds (∼2.60 Å), giving rise to symmetrical thermal diffusion. Asymmetrical heat diffusion was due, on a higher scale, to the local, structural organization of residues that, in turn, was also mediated by hydrogen bonds. This asymmetrical/symmetrical energy flow may be relevant for allosteric signal communication directionality in proteins and for the control of heat flow in materials science.