450 likes | 669 Views
Towards molecular simulations of biological cells . time scale maximal system size Molecular Dynamics 1ns - 1 s 100.000 atoms = (10 nm) 3 Brownian Dynamics 1 s – 1ms 1 – 100 rigid proteins = (100 nm) 3 Random Walk 1 s – 1ms
E N D
Towards molecular simulations of biological cells time scale maximal system size Molecular Dynamics 1ns - 1s 100.000 atoms = (10 nm)3 Brownian Dynamics 1s – 1ms 1 – 100 rigid proteins = (100 nm)3 Random Walk 1s – 1ms Diffusion equation 1s – 1ms cell subcompartments (1 - 10 m)3 (e.g. Virtual Cell) Dynamic Monte Carlo 1 ns – 1 s 10.000 reactions Network models no time scale no length scale Bioinformatics III
Molecular dynamics: internal dynamics of proteins Bioinformatics III
Molecular dynamics: folding of alpha-helix (ca. 100 ns) http://www.stanford.edu/group/pandegroup/folding/results.html Bioinformatics III
Molecular dynamics: protein unfolding http://www.stanford.edu/group/pandegroup/folding/results.html Bioinformatics III
Molecular dynamics: membrane dynamics http://www.mpibpc.mpg.de/abteilungen/071/bgroot/gallery.html Bioinformatics III
Molecular dynamics: membrane dynamics http://www.mpibpc.mpg.de/abteilungen/071/bgroot/gallery.html Bioinformatics III
Protein dynamics time scale maximal system size protein diffusion 10 fs = 10-14 s fastest bond vibrations, 10-5 cm2 s-1 duration the catalytic step of a chemical reaction 1 ps = 10-12 s rotational correlation time of a water molecule frequency of ring flips of Tyr and Phe rings < 1 ns = 10-9 s < life-times of hydrogen bonds 1 ns - 1 dynamics of protein loops, protein-protein association 1s – 1ms crossing of membrane? 1 ms – 1 s protein folding/unfolding Bioinformatics III
Supramolecular systems The changes in the electrostatic potential of the 70S ribosome occurring during the ratchet-like motion. Potential was obtained by solving the Poisson-Boltzmann equation with APBS for various backbone structures of the ribosome derived from normal mode analysis. Red color denotes potential isosurface of -5kT/e and blue +1kT/e. http://chemcca51.ucsd.edu/gallery/018.html Bioinformatics III
Supramolecular systems Cowpea chlorotic mottle virus (CCMV) is used as a model system to investigate the organization of the RNA genome inside the virus capsid. A coarse grain model, in which each nucleic acid is treated as a sphere with negative charge, is necessary to simulate a system of the size of 34,200 amino acid residues and about 3,000 RNA nucleic acids. The electrostatic potential of the capsid is calculated using the Adaptive Poisson-Boltzmann Solver (APBS), and the potential grid is used to evaluate the electrostatic interaction energy of the RNA spheres inside the capsid. In the crystal structure, there are some short strands of RNAs present. These RNAs are kept intact in the simulation. In other words, they are included in the electrostatic potential grid, and they do not move in the Monte Carlo simulation. The following conditions are used in the simulation: RNA sphere size 7.5 A, charge -0.25e, dielectric constant 4.0, Electrostatic potential grid size 2.5 A. Total energy of the system is the electrostatic interations between the RNAs and their electrostic energy within the potential grid of the capsid. The image above is a stereo view of the conformation with the lowest energy, as viewed from outside the virus. The virus capsid is shown together with pre-existing RNA (red spheres = low-energy conformations during MC simulation, orange spheres = ribonucleotides that appear in the crystal structure). http://chemcca51.ucsd.edu/gallery/018.html Bioinformatics III
Integrated Model of Epidermal Growth Factor Receptor Trafficking and Signal Transduction Diagram showing the compartments involved in receptor trafficking and the receptor movement pathways within the cell. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Integrated Model of Epidermal Growth Factor Receptor Trafficking and Signal Transduction Cells of living organism sense their environment and respond to environmental stimuli. Cellular signaling mechanisms govern how information from the environment is decoded, processed and transferred to the appropriate locations within the cell. Signaling through the receptor tyrosine kinase (RTK) family of receptors regulates a wide range of biological phenomena, including cell proliferation and differentiation. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Integrated Model of Epidermal Growth Factor Receptor Trafficking and Signal Transduction Signaling pathways of various RTKs are reasonably well characterized. Common features: - receptor self-phosphorylation on tyrosine residues - subsequent interaction with molecules containing SH2 and phospho-Tyr residues. The signal from the receptor is transmitted to downstream effector molecules through a series of protein-protein interactions, such as the MAP kinase cascade. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Integrated Model of Epidermal Growth Factor Receptor Trafficking and Signal Transduction The EGF receptor can be activated by the binding of any one of a number of different ligands. Each ligand stimulates a somewhat different spectrum of biological responses. The effect of different ligands on EGFR activity is quite similar at a biochemical level the mechanisms responsible for their differential effect on cellular responses are unkown. After binding of any of its ligands, EGFR is rapidly internalized by endocytosis. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Integrated Model of Epidermal Growth Factor Receptor Trafficking and Signal Transduction Diagram showing the compartments involved in receptor trafficking and the receptor movement pathways within the cell. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Integrated Model of Epidermal Growth Factor Receptor Trafficking and Signal Transduction Different EGFR ligands vary in their ability to bind to EGFR as a function of receptor microenvironment such as intravesicular pH. After endocytosis, receptor-ligand complexes pass through several different compartments that vary in their intravesicular milieu. Receptor movement among cellular compartments („receptor trafficking“) can exert a significant effect on the activity of the complexes. The different intracellular compartments also vary in their access to some of the substrates of the EGFR kinase. This coupled relationship between substrate access and ligand-dependent activity in different endocytic compartments suggests that trafficking could function to „decode“ the information unique to each ligand. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
3 functions of trafficking (1) controlling the magnitude of the signal (2) controlling the specificity of the response (3) controlling the duration of the response. Understanding the relative contribution of these 3 aspects for any given combination of cells, conditions, and ligands is very difficult use computational models, hooray! Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Computational modelling of EGF receptor system (1) trafficking and ligand-induced endocytosis (2) signaling through Ras or MAP kinases This work combines both aspects for the first time into a single model. Most approaches to building computational kinetic models have severe drawbacks when representing spatially heterogenous processes on a cellular scale. Traditional approach: - formulate set of coupled ODEs (reaction rate equations) for the time-dependent concentration of chemical species - use integrator to propagate the concentrations as a function of time given the rate constants and a set of initial concentrations. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Dynamic Monte Carlo (DMC) Gillespie showed 1977 that this formal deterministic approach can be translated to a stochastic scheme, termed Dynamic Monte Carlo. Because the molecules forming the physical system are chemical entities, they must participate in the reactions as integer species. Traditional approach based on continuum treatment of chemical kinetics ignores the discrete nature of the problem. This is most problematic when the number of reacting molecules is small use of a discrete representation is more appropriate in kinetic simulation studies in cellular systems. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Multiple time scale problem In DMC, reactions are considered events that occur with certain probabilities over set intervals of time. The event probabilities depend on the rate constant of the reaction and on the number of molecules participating in the reaction. In many interesting natural problems, the time scale of the events are spread over a large spectrum. Then it is very inefficient to treat all processes at the time scale of the fastest individual reaction. In the EGFR signaling network, - receptor phosphorylation after ligand binding occurs almost instantaneously - vesicle formation or sorting to lysosomes require many minutes. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Solution to multiple time scale problem Resat et al. (2001) introduced Probability-Weighted DMC to speed-up the simulation by factor 20 – 100. Different processes are only tested at variant times depending on their probabilities = very unlikely processes compute MC decision very infrequently. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Signal transduction model of EGF receptor signaling pathway Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Species in the EGF receptor signaling model Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Receptor and ligand group definitions Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Rate constants of the ligand:receptor interactions Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Early endosome inclusion coefficients These are adjusted to yield the experimentally determined rates of ligand-free and ligand-bound receptor internalization. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Time course of phosphorylated EGF receptors (a) Total number of phosphorylated EGF receptors in the cell. Curves represent the number of activated receptors when the cell is stimulated with different ligand doses at the beginning. The y axis represents the number of receptors in thousands. (b ) Ratio of the number of phosphorylated receptors that are internalized to that of the phosphorylated surface receptors (in the text this ratio is referred as the In/Sur ratio of the phosphorylated receptors). (c) Ratio of the number of internalized receptors to the number of surface receptors. Curves are colored as: [L] = 0.2 (magenta), 1 (blue), 2 (green), and 20 (red) nM. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Distribution of the receptors among cellular compartments Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Stimulation of EGFR signaling pathway by different ligands Comparison of the results when the EGFR signaling pathway is stimulated with its ligands EGF (red) and TGF- (green). (a ) Total number of receptors in the cell as a function of time after 20 nM ligand is added to the system. Red diamond (EGF) and green square (TGF-) points show the experimental results. The experimental results at short times were normalized to overlap them with the computational results. (b) Distribution of the receptors between intravesicular compartments and the cell membrane. (c) Distribution of the phosphorylated receptors between intravesicular compartments and the cell membrane. In the figures, y axes represent the number of receptors in thousands. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Ratio of internal/surface receptors The ratio of the In/Sur ratios when the EGFR signaling pathway is stimulated with its ligands EGF and TGF- at 20 nM ligand concentration. Comparison of computational (solid lines) and experimental (points) results. Ratio of the ratios for the phosphorylated (i.e., activated) (blue), and total (phosphorylated + unphosphorylated) number (magenta) of receptors. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Assume that there is a ligand of EGFR which has exactly the same association to the receptor properties of EGF. However, it dissociates from the receptor at a higher (e.g. 2.06 times faster) rate when the pH is at the pH level of the endosomes. It was assumed that the only rates that change are the rates of ligand dissociation from the receptor (reaction 1 in Fig.2) and of the reaction that separates a receptor dimer into monomers (reaction 2). Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Comparison of the computed results for the hypothetical ligand (blue) that was investigated in the third set of simulations with the results for the ligands EGF (red) and TGF- (green). Ligand concentration is 20 nM. excessive ligand dissociation in the intravesicular compartments can only partially explain the signaling differences between EGF and TGF-. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Number of SOS molecules recruited to the receptor (a) Dependence on the ligand type, EGF (red) and TGF- (green) at 20-nM ligand concentration. (b) Dependence on the EGF ligand concentration, [L] = 0.2 (magenta), 1 (blue), 2 (green), and 20 (red) nM. In the figures, y axes represent the number of receptor:Sos complexes in thousands. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Number of SOS molecules recruited to the receptor Number of Sos molecules recruited to the receptor complexes through (a) the Grb2 and (b) the Shc/Grb2 branch of the EGFR signal transduction network. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Findings (1) The model predictions for the receptor downregulation and receptor distribution among cellular compartments are in very good agreement with experiment. (2) The distribution of the total (phosphorylated plus unphosphorylated) number of receptors among cellular compartments has a monotonic dependence on the stimulating ligand concentrations. Due to the recycling and degradation patterns of the receptors, the distribution of the total number of receptors among cellular compartments does not saturate until large ligand doses are provided. ... Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Summary Large-scale simulations of the kinetics of biological signaling networks are becoming feasible. Here, the model consists of hundreds of distinct compartments and ca. 13.000 reactions/events that occur on a wide spatial-temporal range. Methods like the probability-weighted DMC are promising tools for studying complex cellular systems using molecular quanta. Resat et al. Biophys Journal 85, 730 (2003) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Components of the autocrine EGFR signaling system and their interactions in the model. Transport, covalent modifications, and protein binding processes are shown by solid arrows. Kinetic acceleration of a process by a species is depicted by a dashed arrow. EGFR is receptor of ligand L; Grb2 is adaptor protein; Sos is activator protein; Raf, MEK, and ERK are protein kinases; p denotes phosphorylated form of a protein, and pp is a double-phosphorylated form. Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Qualitatively different, stable spatial distributions of autocrine signaling. (a) Central section through spherically symmetric distributions of active ERK inside the cell (color coded) and of EGFR ligand outside the cell (grayscale). (b) Central section of axially symmetric, spatially self-organized distributions of the same signaling proteins as in a. (c ) Color coding used. ERK activity is measured by the fraction of the double-phosphorylated form of ERK. The boundary between the colored and the grayscale zones is the cell boundary; the nonshaded circle in the center is the cell nucleus. (Other parameter values: Rtotal = 5 x 103µm-2 in a and 2 x 104µm-2 in b; G1 = 180 µm-2 min-1; DL = 10-8 cm2 s-1; V9 = V10 = 900 nM min-1; kps = 0; [Grb]total = 30 nM; [Sos]total = 34 nM; Linf = 0; k1 = 1.8 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Abrupt transitions between qualitatively different stable distributions of autocrine signaling caused by slow variation of a control parameter, here, the surface density of EGF receptors. The solid and dashed curves represent the steady-state fractions of activated, ligand-bound forms of EGFR at the opposite poles of the cell. (Other parameter values: G1 = 1.25 x 105µm-2 min-1, DL = 10-8 cm2 s-1, V9 = V10 = 900 nM min-1, kps = 0, [Grb]total = 3 nM, [Sos]total = 3.4 nM, Linf = 0; k1 = 36 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Sections through the multidimensional parameter space of the model that show the domains of the parameter values under which the stable states of the system differ qualitatively. In the blue domain, the stable state is no signaling. In the red domain, the stable state is uniform signaling with a spherically symmetric distribution of protein species. In the green domain, the stable state is spatially self-organized signaling that possesses only axial symmetry. Numbers in the graphs mark locations in the parameter space that are discussed in the text. (a) Receptor density is the total density of receptors and all their complexes at the membrane (Rtotal). Ligand release rate is the maximum rate (G1). (b) The phosphatase activity is the maximum rate of ERK dephosphorylation (V9 = V10). (c) The strength of the negative feedback is measured by the rate constant of Sos phosphorylation, kps. The strength of the positive feedback is measured by the rate constant of EGF release, G1. (The abscissae in a and c refer to the same numerical parameter of the model in its different mechanistic roles.) (d ) On the axes, [Grb] and [Sos] are the total concentrations of Grb2 and Sos in the cell, including the Grb2 or Sos components of the different protein complexes. (Other parameter values: Linf = 0; k1 = 36 min-1, a; k1 = 1.8 min-1, b–d.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Collapse of the spatially self-organized state of the autocrine system of the cell subjected to paracrine stimulation. On the abscissa axis is the concentration of EGFR ligand in the medium far from the autocrine cell. The solid and broken curves are the steady-state fractions of activated, ligand-bound forms of EGFR at the opposite poles of the cell. (Other parameter values: Rtotal = 1.75 x 103µm-2, G1 = 3.6 x 103µm-2 min-1, DL = 10-8 cm2 s-1, V9 = V10 = 900 nM min-1, kps = 0, [Grb]total = 30 nM, [Sos]total = 30 nM, k1 = 1.8 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Collapse of the spatially self-organized state of the autocrine system of the cell subjected to paracrine stimulation. On the abscissa axis is the concentration of EGFR ligand in the medium far from the autocrine cell. The solid and broken curves are the steady-state fractions of activated, ligand-bound forms of EGFR at the opposite poles of the cell. (Other parameter values: Rtotal = 1.75 x 103µm-2, G1 = 3.6 x 103µm-2 min-1, DL = 10-8 cm2 s-1, V9 = V10 = 900 nM min-1, kps = 0, [Grb]total = 30 nM, [Sos]total = 30 nM, k1 = 1.8 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Collapse of the spatially self-organized state of the autocrine system of the cell subjected to paracrine stimulation. On the abscissa axis is the concentration of EGFR ligand in the medium far from the autocrine cell. The solid and broken curves are the steady-state fractions of activated, ligand-bound forms of EGFR at the opposite poles of the cell. (Other parameter values: Rtotal = 1.75 x 103µm-2, G1 = 3.6 x 103µm-2 min-1, DL = 10-8 cm2 s-1, V9 = V10 = 900 nM min-1, kps = 0, [Grb]total = 30 nM, [Sos]total = 30 nM, k1 = 1.8 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Collapse of the spatially self-organized state of the autocrine system of the cell subjected to paracrine stimulation. On the abscissa axis is the concentration of EGFR ligand in the medium far from the autocrine cell. The solid and broken curves are the steady-state fractions of activated, ligand-bound forms of EGFR at the opposite poles of the cell. (Other parameter values: Rtotal = 1.75 x 103µm-2, G1 = 3.6 x 103µm-2 min-1, DL = 10-8 cm2 s-1, V9 = V10 = 900 nM min-1, kps = 0, [Grb]total = 30 nM, [Sos]total = 30 nM, k1 = 1.8 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III
Application of elementary modesMetabolic network structure of E.coli determineskey aspects of functionality and regulation Collapse of the spatially self-organized state of the autocrine system of the cell subjected to paracrine stimulation. On the abscissa axis is the concentration of EGFR ligand in the medium far from the autocrine cell. The solid and broken curves are the steady-state fractions of activated, ligand-bound forms of EGFR at the opposite poles of the cell. (Other parameter values: Rtotal = 1.75 x 103µm-2, G1 = 3.6 x 103µm-2 min-1, DL = 10-8 cm2 s-1, V9 = V10 = 900 nM min-1, kps = 0, [Grb]total = 30 nM, [Sos]total = 30 nM, k1 = 1.8 min-1.) Maly et al. Biophys Journal 86, 10 (2004) Bioinformatics III