The present invention relates to Method for simulation of biological pathways, chemical reaction pathways, biomolecules, and nano-molecular-machines by using electronic circuits.
This general field is known as “Molecular Biology” (MB), “Cell Biology” (CB), “Nano-technology” (NT), “Chemistry” (CH). When used for analysis of disease processes, this field is referred to as “Pathology”. When used for the examination of the effect of pharmaceutical agents on biological systems, this field is referred to as “Drug evaluation” (DE) or “Drug testing” (DT).
The Need for Simulation of Biological Pathways and Chemical Reaction pathways by a new approach
Simulation of the kinetics of biological pathways and chemical reaction pathways facilitates the understanding of biological processes, complex chemical reactions, disease processes, and the effect of mutations on biological systems. It can also be used for testing the effect of drugs/molecules on biological systems and disease processes. These are described in “Simulation of prokaryotic genetic circuits”, H. H. McAddams, and A. Arkin. Annu. Rev. Biophys. Biomol. Struct. 27, 199-224 (1998).
So far, computational methods appear to be the only practically used approach for simulation of biological pathways as described in “A mathematical model of caspase function in apoptosis”, M. Fussenegger, J. E. Bailey, and J. Varner. Nature Biotech. 18, 768-774 (2000); “Mathematical modeling of epidermal growth factor receptor signaling through the phospholipase C pathway: Mechanistic insights and predictions for molecular interventions”, J. M. Haugh, A. Wells, D. A. Lauffenburger. Biotech Bioeng 70, 225-238 (2000); “Computational analysis of biochemical systems, a practical guide for biochemists and molecular biologists”, E. O. Voit. Cambridge Press (2000). In these works, the kinetics of binding and reaction events in a pathway, under certain approximations, have also been modeled by a set of first-order linear differential equations.
The kinetics of binding and reaction events in a pathway have also been described by a set of first-order nonlinear ordinary differential equations as proposed in “Systemic properties of ensembles of metabolic networks: application of graphical and statistical methods to simple unbranched pathways”, R. Alves, and M. A. Savageau. Bioinformatics 16, 534-547 (2000); and in “Inferring qualitative relations in genetic networks and metabolic pathways”, T. Akusu, S. Miyano, and S. Kuhara. Bioinformatics 16, 727-734 (2000).
While they have been successfully used in studying individual signaling and metabolic pathways involving dozens of proteins and other molecules, in the foreseeable future, computational methods are not expected to have the capacity for simulation of large and complex pathways or multiple pathways involving hundreds or more number of proteins as described in “Metabolic networks: a signal-oriented approach to cellular models” J. W. Lengeler. Biol. Chem. 381, 991-920 (2000).
Cellular events frequently involve “crosstalk” interactions among multiple pathways that may include thousands or more proteins and other molecules, and the analysis of such complex biochemical networks is further complicated by the great number of feedback and feedforward loops and of regulatory networks involved in cellular control, as described in “Metabolic networks: a signal-oriented approach to cellular models” J. W. Lengeler. Biol. Chem. 381, 991-920 (2000) and in “Mathematical modeling of epidermal growth factor receptor signaling through the phospholipase C pathway: Mechanistic insights and predictions for molecular interventions”, J. M. Haugh, A. Wells, D. A. Lauffenburger, Biotech Bioeng 70, 225-238 (2000).
The kinetics of such a large network is described by up to thousands of first order nonlinear differential equations. Significantly more numbers of equations may be needed for simulation of cellular and multi-cellular systems.
Given that it is unlikely that the computer resource in the near future is capable of solving the “long-time” kinetics of such a large number of nonlinear equations, a new approach needs to be introduced.
The Need for Structural Optimization and Molecular Dynamics Simulation of Biomolecules and Nano-molecular-systems by a New Approach
Molecular dynamics simulation is a widely used method for simulation and modeling of structures, motions, binding, and thermodynamics of biomolecules and nano-molecular-systems. These are described in “Computer simulation of biomolecular systems: Theretical and experimental applications”, Eds. W. F. van Gunsteren, P. K. Weiner, and A. J. Wilkinson. ESCOM Leiden (1993); “Computational nanotechnology”, R. C. Merkle, Nanotechnology, 2, 134-141 (1991).; nanostructure engineering as described in “Simulated engineering of nanostructures”, D. W. Brenner, S. B. Sinnott, J. A. Harrison, and O. A. Shenderova. Nanotechnology 7 161-167 (1996). It has been used in facilitating the study of a variety of important problems such as protein folding as described in “Pathways to a protein folding intermediate observed in a 1-microsecond simulation in aqueous solution”, Y. Duan and P. A. Kollman. Science 282, 740-744 (1998); drug design as described in “Binding pathway of retinal to bacterio-opsin: A prediction by molecular dynamics simulations”, B. Isralewitz, S. Izrailev, and K. Schulten. Biophys. J. 73, 2972-2979; nanostructure engineering as described in “Simulated engineering of nanostructures”, D. W. Brenner, S. B. Sinnott, J. A. Harrison, and O. A. Shenderova. Nanotechnology 7 161-167 (1996).
Nano-molecular-systems refer to: nano-sized or nano-structured materials (such as molecular films, nanotubes, nanoscopic particles with specific electronic, optical or magnetic properties), nano-scale molecular mechanical or manufacturing systems (capable of guiding reactive molecules to 0.1 nm precision), other nano-scale devices (molecular motors, carriers, containers, pumps, circuits, tools etc.). These are described in “Nanosystems: molecular machinery, manufacturing, and computation” by K. Eric Drexler (Wiley 1992).
So far, computational methods are the only practically used approach for molecular dynamics simulation as described in “Computer simulation of biomolecular systems: Theretical and experimental applications”, Eds. W. F. van Gunsteren, P. K. Weiner, and A. J. Wilkinson. ESCOM Leiden (1993). The molecular dynamics simulation of biomolecules can be described by a set of second-order nonlinear ordinary differential equations as indicated in “A second generation force field for simulation of proteins, nucleic acids, and organic molecules”, W. D. Cornell et. al., J. Am. Chem. Soc. 117, 5179-5197 (1995); and in “CHARMM: A program for macromolecular energy, minimization, and dynamics calculations”, B. R. Brooks, et. al. J. Comp. Chem. 4, 187-217 (1983). While they have been successfully used in studying some problems of proteins and other molecules, in the foreseeable future, computational methods are not expected to have the capacity for simulation of proteins and nucleic acids for orders of magnitude longer than the microsecond range.
The estimated time scale for folding/unfolding events in proteins, ligand-receptor binding/dissociation, and interaction of nano-molecular-machine with its substrates ranges from hundreds of milliseconds to tens of seconds as described in “Kinetic analysis of the unfolding and refolding of ribonuclease T1 by a stopped-flow double-mixing technique”, L. M. Mayr et. al., Biochemistry 35, 5550-5561 (1996); and in “Studies of the binding of actinomycin and related compounds to DNA”, W. Muller, and D. M. Crothers, J. Mol. Biol. 35, 251-290 (1968). In contrast, the longest achieved simulation time is 1 microsecond for a small protein, consisting of only 36 amino acids, on a 256-node CRAY T3E supercomputer, as described in “Pathways to a protein folding intermediate observed in a 1-microsecond simulation in aqueous solution”, Y. Duan and P. A. Kollman. Science 282, 740-744 (1998). Hence there is a big gap (5 to 7 orders of magnitude) between simulation time scale and that of certain molecular events. Given that it is unlikely the computer approach in the near future is capable of conducting “long” enough time-scale molecular dynamics simulation, a new approach needs to be introduced.
The present invention seeks to address the problems above, and in particular to provide a new approach for the simulation of biological and/or chemical reaction pathway, and for the molecular dynamic simulation of biomolecules and nano-molecular systems by comprising constructing set of binding, reaction, concentration and/or molecular dynamics equations and electronic circuits corresponding to such equations.
According to a first embodiment, the invention provides a method for simulation of at least one biological and/or chemical reaction pathway comprising:
In particular, the binding, reaction and/or concentration equation may be a linear or non-linear first- or second-order ordinary differential equation (ODE).
The electronic circuit may comprise at least one of the following circuit units: linear protein/molecule concentration ±kxi unit; protein/molecule concentration power ±kxja unit; ligand-protein concentration product ±Lxi unit; ligand-protein concentration power product ±kLxia unit; protein-protein concentration product ±kxixj unit; protein-protein concentration power product ±kxiaxjb unit; stochastic rate constant generator unit that replaces a forward/reverse binding/reaction rate constant k by a stochastic value; and
According to an aspect, the method of the invention further comprises maintaining the voltage level of the circuit between two fixed voltage values. For example, by: multipling scaling factors to the concentration equation; applying at least one resistor and/or amplifier at one or more connection point of the circuit, thereby scaling-down or scaling-up the voltage of one or more segment of the circuit; and/or applying automatic gain control circuits.
According to another aspect, the method of the invention further comprises adding at least one unit circuit comprising an electronic random noise generator and/or a multiplicator amplifier at one or more connection point of the circuit.
The method of the invention may further comprise determining the effect of at least one drug comprising adding at least one circuit unit associated with a receptor protein at one or more connection point of the circuit.
The method may further comprise determining deficiency, mutation and/or deletion of at least one protein of the biological pathway, comprising setting a voltage ceiling, a voltage range and/or a fixed voltage at one or more connection point of the circuit.
The invention also provides an electronic circuit system for the simulation of at least one biological and/or chemical reaction pathway, comprising at least one electronic circuit representing a set of binding, reaction and/or concentration equation. In particular, an electronic circuit system, wherein the binding, reaction and/or concentration equation is a linear or non-linear first- or second-order ordinary differential equation (ODE).
The electronic circuit may comprise at least one the following circuit units: linear protein/molecule concentration ±kxi unit; protein/molecule concentration power ±kxja unit; ligand-protein concentration product ±Lxi unit; ligand-protein concentration power product ±kLxja unit; protein-protein concentration product ±kxixj unit; protein-protein concentration power product ±kxiaxjb unit; stochastic rate constant generator unit that replaces a forward/reverse binding/reaction rate constant k by a stochastic value; and
wherein, k is a forward/reverse binding/reaction rate constant, L the concentration of a ligand, xi and xj are the concentration of protein/molecule xi and xj respectively, a and b are order of power of xi and xj respectively.
According to another embodiment, the invention provides a method for molecular dynamics simulation of biomolecules and/or nano-molecular systems comprising:
The equation may be a linear or non-linear second order ordinary differential equation (ODE).
According to an aspect, in the method of the invention, the electronic circuit comprises at least one atom-position circuit unit, wherein the atom-position circuit unit represents the position of an atom of a molecule or a molecular system.
In particular, the atom-position circuit unit comprises at least one atom-atom interaction circuit subunit, the atom-atom interaction circuit subunit representing a sub-unit of atom-atom interactions within a molecule or a molecular system and comprising internal bond stretch, angle bending, torsion, non-bonded units; and/or bond stretch, angle bending, and/or torsion units; between at least two nearest sub-units of a molecule.
Further, the atom-atom interaction circuit subunit may represent a term in the molecular dynamics equation, and wherein the atom-atom interaction circuit subunit comprises at least one of the following: bond stretch x unit, bond stretch y unit, bond stretch z unit, angle bending x type-A unit, angle bending x type-B unit, angle bending y type-A unit, angle bending y type-B unit, angle bending z type-A unit, angle bending z type-B unit, torsion x type-A unit, torsion x type-B unit, torsion y type-A unit, torsion y type-B unit, torsion z type-A unit, torsion z type-B unit, non-bonded x unit, non-bonded y unit, non-bonded z unit, hydrogen-bond x unit, hydrogen-bond y unit, and hydrogen-bond z unit; and wherein x, y, and z represent the coordinates of each atom of the molecule, and type-A represents the case of the atom being in the middle-position of an angle bending or torsion connection with other atoms, and type-B represents the case of the atom being in the end-position of an angle bending or torsion connection with other atoms.
The method may further comprise maintaining the voltage level in the circuit between two fixed voltage values. For example, by wherein x, y and z represent the coordinate of the molecule, and the voltage level of the circuit is maintained between two fixed voltage values by:
The invention also provides a circuit group representing the interaction pattern in the chemical structure of a molecule or a sub-unit of interaction pattern in the chemical structure of a molecule comprising:
The electronic circuit may comprise at least one circuit unit, the circuit unit comprising at least one circuit group, the circuit group representing a sub-unit of interaction pattern in the chemical structure of a molecule and comprising internal bond stretch, angle bending, torsion, non-bonded units; and/or bond stretch, angle bending, and/or and torsion units; between at least two nearest sub-unit of a molecule.
The circuit unit may represent a term in the molecular dynamic equation, and wherein the circuit unit comprises at least one of the following: bond stretch x unit, bond stretch y unit, bond stretch z unit, angle bending x type-A unit, angle bending x type-B unit, angle bending y type-A unit, angle bending y type-B unit, angle bending z type-A unit, angle bending z type-B unit, torsion x type-A unit, torsion x type-B unit, torsion y type-A unit, torsion y type-B unit, torsion z type-A unit, torsion z type-B unit, non-bonded x unit, non-bonded y unit, non-bonded z unit, hydrogen-bond x unit, hydrogen-bond y unit, and hydrogen-bond z unit; and wherein x, y, and z represent the coordinates of each atom of the molecule, and type-A unit represents the case of the atom being in the middle-position of an angle bending or torsion connection with other atoms, and type-B represents the case of the atom being in the end-position of an angle bending or torsion connection with other atoms.
The invention further relates to a method for the manufacture of an electronic circuit representing at least one biomolecule and/or nano-molecular system, the electronic circuit comprising at least a unit circuit comprising at least a circuit group, wherein the circuit group represents a sub-unit of interaction pattern in the chemical structure of a molecule or a molecular system comprising:
and
Method for Simulation of Biological and/or Chemical Reaction Pathway
According to a first embodiment, the present invention provides a method for simulation of a single or multiple biological and/or chemical reaction pathways.
It has particular application in simulation of large networks of pathways that are beyond the scope of computer methods.
It provides a potentially useful means for facilitating the study of cellular and multi-cellular processes, disease processes, effect of gene or protein mutations on biological systems, and effect of drugs on biological systems.
A single or a set of biological pathways and/or chemical reaction pathways can be represented by a map composed of a network of boxes or labels (reaction objects) connected by arrows (binding/reaction paths). This kind of representation is widely used in popular biological pathway databases such as SPAD (http://www.qrt.kyushu-u.ac.jp/spad/) and KEGG (http://www.genome.ad.jp/kegg/kegg2.html) and as described in “Mathematical modeling of epidermal growth factor receptor signaling through the phospholipase C pathway: Mechanistic insights and predictions for molecular interventions”, J. M. Haugh, A. Wells, D. A. Lauffenburger. Biotech Bioeng 70, 225-238 (2000). Each binding/reaction event in these pathways can be described by a binding/reaction equation, whose form depends on whether the event is reversible or irreversible.
The set of binding and/or reaction (binding/reaction) equations for pathways can be modeled by a set of linear and/or non-linear first-order and/or second-order ordinary differential equations. For example, non-linear first-order ordinary differential equations (ODEs) derived from the balance equation and generalized mass-action kinetics as described in “Understanding the control of metabolism” D. Fell, Portland Press, London, UK (1996) and in “The regulation of cellular systems” R. Heinrich, and S Schuster, Chapman & Hall, New York (1996). Each ODE describes the time-dependent behaviour of the concentration of a protein/molecule in the pathways. Further, equations representing the concentration of the molecules of the pathway map of the biological and/or chemical reaction pathway are also provided. The concentration equations can also be at least a set of linear and/or non-linear first-order and/or second-order ordinary differential equations. For example, non-linear first-order ordinary differential equations (ODEs).
According to this embodiment of the invention, at least one set of electronic circuits represents the at least one set of the above equations. In particular, binding, reaction and/or concentration non-linear ODEs.
Electronic circuits have been proposed to solve certain first-order and second-order nonlinear ordinary differential equation as described in “The transition from solitons to chaos in the solution of the logistic equation” M. I. Sobhy, and A. S. Burman. Int. J. Bifurcation Chaos 10, 2823-2829 (2000) as well as linear ordinary differential equation (up to second-order) as described in “The analog computer as an aid to the teaching elementary quantum mechanics”, M. K. Summers. Phys. Educ 13, 22-27 (1996).
However, electronic circuits have never been disclosed nor suggested for the purpose of simulation of biological and/or chemical reaction pathway.
Operational amplifiers, function circuits and other electronic components can be used to construct operations of differentiation, integration, sum, subtraction, multiplication, inversion, exponential, logarithm, power, and others as described in “Electronic circuits and design” D. A. Neamen. McGraw Hill (2001), “Function circuits: Design and applications” Y. J. Wong and W. E. Ott, McGraw-Hill (1976), “Modern operational circuit design” J. L. Smith, Wiley-Interseicne (1971). Hence, a set of nonlinear or linear ordinary differential equations composed of these terms can be represented by a set of electronic circuits composed of these operational amplifiers and other electronic components. The time-dependent behaviour of the concentration of a protein/molecule can be measured by the voltage at a specific point in the circuits.
Accordingly, for simulation of biological and/or chemical reaction pathways, the time-dependent behavior of the concentration of one or more reaction objects, defined as proteins/nucleic acids/ligands/substrates/reactants/products, in the pathways can be derived by measuring the voltage at specific points in the circuits.
In order to accurately describe the microscopic kinetics of binding or reaction events in pathways, all possible states for each protein or molecule need to be considered. Examples of these states include active, inactive, ligand/substrate-bound, dimmer, and any other state. The corresponding circuit section for the concentration in a particular state is typically composed of up to a few dozen operational amplifiers, a similar number of connection points (two for measuring its voltage, a few for connection to external voltage sources, and several for connection to and from other circuit units), and a number of small components such as switches etc. The circuits for a large pathway or multiple pathways may thus contain up to thousands or more such sections connected to each other and to the external voltage sources, which is feasible to construct.
Accordingly, the first embodiment provides a method for simulation of at least one biological and/or chemical reaction pathway comprising:
The binding, reaction and/or concentration equation may be a linear or non-linear first- or second-order ordinary differential equation (ODE). Preferably, the binding, reaction and/or concentration equation is a non-linear first-order ordinary differential equation (ODE).
In a further aspect of this invention, methods are disclosed for maintaining the voltages in the circuits (concentrations of molecules and functions of concentrations) within the allowed range. Scaling factors (which are non-zero real numbers) are applied to the concentrations and to the kinetic equations so as to keep these concentrations and each term in the equations within the allowed range. In addition, resistors and amplifiers can be used to scale-down and scale-up the voltages at certain segment of the circuits. Moreover automatic gain control circuits can also be used to ensure the voltages are within the required range. The automatic gain control circuits are described in “Function circuits: Design and applications” Y. J. Wong and W. E. Ott, McGram-Hill, 1976.
Accordingly, the invention provides a method further comprising maintaining the voltage level of the circuit between two fixed voltage values.
In particular, the voltage level of the circuit is maintained between two fixed voltage values by:
In a further aspect of this invention, stochastic effect on the binding/reaction events, which is important in cases involving low concentration or short-time event, can be taken into account by replacing the relevant binding/reaction rate constants by stochastic values (random variations around average rate constants determined by experiments or by theoretical methods). This can be accomplished in the circuits by adding at specific sites special circuit units composed of an electronic random noise generator device, such as Wayne's random noise generators, and a multiplication amplifier.
Accordingly, the method of the invention further comprises adding at least one unit circuit comprising an electronic random noise generator and/or a multiplicator amplifier at one or more connection point of the circuit.
In a further aspect of the invention, design plans for several circuit units are introduced to represent one or more, or even all possible terms in a circuit equation. Accordingly, in the method of the invention, the electronic circuit may comprise at least one of the following circuit units: linear protein/molecule concentration ±kxi unit; protein/molecule concentration power ±kxia unit; ligand-protein concentration product ±Lxi unit; ligand-protein concentration power product ±kLxia unit; protein-protein concentration product ±kxixj unit; protein-protein concentration power product ±kxiaxjb unit; stochastic rate constant generator unit that replaces a forward/reverse binding/reaction rate constant k by a stochastic value; and
By introducing additional circuit sections or by modifying the circuits, it is possible to simulate events such as the effect of drugs, protein or gene mutations, and diseases on a pathway or a group of pathways.
Therefore, in a further aspect of this invention, effect of drugs can be simulated by adding one or more circuit units (the design plan for each of these units is disclosed) at specific sites in the circuits that are associated with the respective receptor proteins. These added units represent the forward and reverse binding/reaction process between these drugs and their respective receptor proteins.
Accordingly, the method of the invention, further comprises determining the effect of at least one drug comprising adding at least one circuit unit associated with a receptor protein at one or more connection point of the circuit.
In a further aspect of this invention, deficiency or mutation of one or more proteins can be simulated by setting a voltage ceiling, or a voltage range, or a fixed voltage value at specific points in the circuits.
Accordingly, the method of the invention further comprises determining deficiency, mutation and/or deletion of at least one protein of the biological pathway, comprising setting a voltage ceiling, a voltage range and/or a fixed voltage at one or more connection point of the circuit.
The method of claim 1, wherein the biological pathway comprises at least one object, and wherein the object is protein, nucleic acid, ligand, substrate, inhibitor or antagonist, activator or agonist, reactant and/or reaction product.
According to the first embodiment of the invention, it is also provided an electronic circuit system for the simulation of at least one biological and/or chemical reaction pathway, comprising at least one electronic circuit representing a set of binding, reaction and/or concentration equation.
In the electronic circuit system, the binding, reaction and/or concentration equation may be a linear or non-linear first- or second-order ordinary differential equation (ODE). Preferably, a non-linear first-order ODE.
According to a particular aspect, the invention relates to the disclosure of design plans for the construction of circuit units that represent one or more or even all possible terms in the ODEs describing a biological or a chemical reaction pathway.
Accordingly, in the electronic circuit system of the invention, the electronic circuit may comprise at least one, or even all, the following circuit units: linear protein/molecule concentration ±kxi unit; protein/molecule concentration power ±kxia unit; ligand-protein concentration product ±Lxi unit; ligand-protein concentration power product ±kLxia unit; protein-protein concentration product ±kxixj unit; protein-protein concentration power product ±kxjaxjb unit; stochastic rate constant generator unit that replaces a forward/reverse binding/reaction rate constant k by a stochastic value; and
Referring now to the drawings and, in particular, to
The process 100 begins at a state 110, wherein the pathways are presented as a map of a network of boxes or labels (proteins or molecules) connected by arrows (binding/reaction paths).
The process 100 then moves to a sate 120 wherein the binding or reaction equation for every process in the pathway is constructed from the pathway map. Based on whether it is reversible or irreversible, the equation of an individual binding/reaction step takes the form:
where R1, R2, . . . are binding/reaction objects (proteins/molecules etc) and P1, P2, . . . are product objects (proteins/molecules etc), K+and K−are forward and reverse binding/reaction rate constant respectively.
The first few reaction equations for the illustrative caspase pathway in
where L=Fas ligand, x1=ligand-free Fas surface receptor, x2=ligand-bound Fas surface receptor, x3=clustered ligand-bound Fas surface receptor, x4=FADD protein, x5=complex of FADD-Fas receptor, . . . , and ki+and ki are the forward and reverse reaction rate constant in the i-th path respectively. These rate constants can either be obtained from experimental data, or from theoretical computations, or estimated by empirical means based on statistical analysis of available experimental data of available of relevant reactions.
In general, the reaction equation for the i-th path is:
where aij and bik are the stoichiometric coefficient of the j-th and k-th reaction objects respectively; j=1, . . . , Nj; k=1, . . . , Nk; and Nj and Nk are number of reaction objects and number of product objects respectively.
The process 100 then moves to a state 130 where in the corresponding kinetic equations for the concentration of each protein or molecule are derived from the balance equation and generalized mass-action kinetics as described in “Understanding the control of metabolism” D. Fell, Portland Press, London, UK (1996) and in “The regulation of cellular systems” R. Heinrich, and S Schuster, Chapman & Hall, New York (1996). Each equation is a non-linear first order ordinary differential equation (ODE) of the concentration of a protein/molecule in the pathways. In general, for a set of equations involving the n-th reaction object:
the corresponding ODE is:
dXn/dt=ki−Πk(Xk)blk+ki′+Πj′(Xj′)aij′+ . . . ki+(xn)ainΠj(Xj)aij−ki′(xn)ai′nΠjk′(Xk′)bi′k′. . . +ωn
where xn is the concentration of object xn, ωn is its rate of expression/production, aij and bij are the power of the concentration in cases where the reaction rate is determined by power law representation.
The first few ODEs for the illustrative caspase pathway are:
dx1/dt=k1−X2−k1+L x1+ω1 dx3/dt=k2 X2+k3−X5−k3+x3x4
dx2/dt=k1+L X1−(k1−+k2)x2 dx4/dt=k3−X5+k5X6−k3+X3X4−k5+X4 X5+ω4
where the ligand concentration L and rate of expression ωn n of protein xn are assumed to be time-dependent function.
Referring back to
The process 100 then moves to a simulation state 150 wherein the time-dependent kinetics of the concentration of one or more proteins or molecules in a pathway or a group of pathways can be determined by measurement of voltages at specific sites of the circuits. For the circuit illustrated in
Construction Plan for the Specially Designed Circuit Units Representing Different Terms in an ODE For Pathway Simulation
Referring now to
Methods for Maintaining the Voltages in the Circuits Within Allowed Range
Normally the concentration of a molecule can be scaled into the allowed voltage range of a circuit. In case that this may not be so, or in case that the derivative of concentration exceeds the voltage range, or in case that some terms in the kinetic equation may exceed the voltage range, the following methods can be used to maintain the voltage level of the circuit between two voltage values:
The effect of drugs on a pathway or a group of pathways can be simulated by adding ligand-protein concentration product units to the circuit sections of the drug receptor proteins. This corresponds to case 3 in
Illustrative Example: Electronic Pathway Emulator of a Caspase Pathway
As mentioned before in the description of the first embodiment of the invention, a pathway can be described by a set of binding/reaction equations. This set of binding and/or reaction equations can be further modeled by a set of ODEs. An illustrative biological pathway is part of caspase pathway in
where the ligand concentration L and rate of expression Con can change with time (a time-dependent function). Variable x1 to x20 represent the protein concentration at the respective paths. This set of nonlinear ODEs can be represented by a corresponding set of electronic circuits composed of different function blocks. Thus, an electronic biological pathways emulator can be built to mimic and study the biological process.
Among the given 20 simultaneous ODEs, Eqn. (I) to (IV) are of typical types. These typical equations can be emulated using the circuits shown in
The emulator circuit has been constructed and the measured voltages, representing variations of protein concentrations in the illustrative biological pathway, are shown in
Molecular Dynamics Simulations of Biomolecules and/or Nano-Molecular-System
According to a second embodiment, the invention relates to molecular dynamics simulations of biomolecules and/or nano-molecular-systems.
Nano-molecular-systems refer to: nano-sized or nano-structured materials (such as molecular films, nanotubes, nanoscopic particles with specific electronic, optical or magnetic properties), nano-scale molecular mechanical or manufacturing systems (capable of guiding reactive molecules to 0.1 nm precision), other nano-scale devices (molecular motors, carriers, containers, pumps, circuits, tools etc.). These are described in “Nanosystems: molecular machinery, manufacturing, and computation” by K. Eric Drexler (Wiley 1992).
The molecular dynamics of a protein/nucleic acid/organic molecule can be represented by a set of linear and/or non-linear second-order ordinary differential equations (ODEs) derived from the Newton's second law of motion given the energy functions, force fields and a starting 3D structure of the system, as described in “A second generation force field for simulation of proteins, nucleic acids, and organic molecules”, W. D. Cornell et. al., J. Am. Chem. Soc. 117, 5179-5197 (1995); and in “CHARMM: A program for macromolecular energy, minimization, and dynamics calculations”, B. R. Brooks, et. al. J. Comp. Chem. 4, 187-217 (1983). Each ODE describes the time-dependent change of the x- or y- or z- coordinate of each atom in a studied molecule.
According to the second embodiment of the invention, a corresponding set of electronic circuits is proposed to represent this set of linear or non-linear ODEs. This set of electronic circuits is composed of operational amplifiers, function circuits and other electronic circuit elements. A combination of one or more specially designed circuit units (composed of operational amplifiers, function circuits and other electronic components together with differentiator and summing amplifier) is used to represent each term of the ODEs. The time-dependent position of each atom in a protein/nucleic acid/organic molecule can be measured by the voltage at a specific point in the circuits.
Electronic circuits have been proposed to solve certain first-order and second-order nonlinear ordinary differential equation as described in “The transition from solitons to chaos in the solution of the logistic equation” M. I. Sobhy, and A. S. Burman. Int. J. Bifurcation Chaos 10, 2823-2829 (2000) as well as linear ordinary differential equation (up to second-order) as described in “The analog computer as an aid to the teaching elementary quantum mechanics”, M. K. Summers. Phys. Educ 13, 22-27 (1996). However, electronic circuits have never been disclosed nor suggested for the purpose of molecular dynamics simulations of biomolecules and/or nano-molecular systems.
If a 3D structure is not available, a starting structure of a molecule can be generated from its molecular bond connectivity profile and from the standard average bond length, angle and torsion. Solution of these molecular dynamics simulation equations gives the time evolution of the position of each atom of a studied molecule. All the terms in this equation can be represented by a combination of existing operational amplifiers, function circuits and other electronic components. Hence it is feasible to use electronic circuits to solve these equations and thus conduct molecular dynamics simulations.
A protein normally consists of thousands of atoms. Each atom has three equations describing its x-, y-, or z- coordinate respectively. The corresponding circuit section for each coordinate is typically composed of up to a few dozen bond-stretch/angle/torsion units plus up to thousands of non-bonded units. If a force field with separate hydrogen bond terms is used, then this circuit section also contains a few additional units for hydrogen bonds. If explicit solvent is used, then this circuit section also contains additional non-bonded units related to the non-bonded interactions with water molecules. Each unit consists of approximately a few dozen operational amplifiers and function amplifiers, a similar number of connection points (two for measuring its voltage, a few for connection to external voltage sources, and the rest for connection to and from other circuit units), and a number of small components such as switches etc. The circuits for a typical protein may thus contain up to thousands or more such sections connected to each other and to the external voltage sources, which is feasible to construct.
Operational amplifiers, function circuits and other electronic components can be used to construct operations of differentiation, integration, sum, subtraction, multiplication, inversion, exponential, logarithm, power, and others as described in “Electronic circuits and design” D. A. Neamen. McGraw Hill (2001), “Function circuits: Design and applications” Y. J. Wong and W. E. Ott, McGraw-Hill (1976), “Modern operational circuit design” J. L. Smith, Wiley-Interseicne (1971). Hence, a set of nonlinear or linear ordinary differential equations composed of these terms can be represented by a set of electronic circuits composed of these operational amplifiers and other electronic components. The time-dependent behavior of the concentration of a protein/molecule can be measured by the voltage at a specific point in the circuits.
The present invention provides a method for molecular dynamics simulation of biomolecules and/or nano-molecular systems comprising:
The equations may be linear and/or non-linear second order ordinary differential equations (ODEs). Preferably, non-linear second order ordinary differential equations (ODEs).
The biomolecule and/or nano-molecular system of the invention comprise in general, but not exclusively, amino acids, nucleotides, organic molecules, and/or inorganic molecules.
According to one aspect of the invention, construction plan is disclosed for circuit group of each unit in a biomolecule (amino acid or nucleotide) or a nano-molecular-system to represent all its internal interaction bonded and non-bonded interactions and bonded interactions to its nearest neighbouring amino acid or nucleotide. Each circuit group is composed of at least one of: (1) internal bond stretch, angle bending, torsion, non-bonded units; and (2) bond stretch, angle bending, and torsion units to its two nearest neighbouring amino acid or nucleotide. In a molecular dynamics simulation of a protein or nucleic acid, these circuit groups can be combined together according to the sequence of the protein or nucleic acid. Additional non-bonded units representing inter-residue non-bonded interactions can be added to each circuit group accordingly.
Accordingly, the electronic circuit of the invention may comprise at least one atom-position circuit unit, wherein the atom-position circuit unit represents the position of an atom of a molecule or a molecular system. The atom-position circuit unit comprises at least one atom-atom interaction circuit subunit, the atom-atom interaction circuit subunit representing a sub-unit of atom-atom interactions within a molecule or a molecular system and comprising at least one of: internal bond stretch, angle bending, torsion, non-bonded units; bond stretch, angle bending, and torsion units; between at least two nearest sub-unit of a molecule.
Each atom-atom interaction circuit subunit represents a term in the molecular dynamics equation and comprises at least one of the following: bond stretch x unit, bond stretch y unit, bond stretch z unit, angle bending x type-A unit, angle bending x type-B unit, angle bending y type-A unit, angle bending y type-B unit, angle bending z type-A unit, angle bending z type-B unit, torsion x type-A unit, torsion x type-B unit, torsion y type-A unit, torsion y type-B unit, torsion z type-A unit, torsion z type-B unit, non-bonded x unit, non-bonded y unit, non-bonded z unit, hydrogen-bond x unit, hydrogen-bond y unit, and hydrogen-bond z unit; and
According to a particular aspect, all the twenty subunit elements listed above may be used for the construction of the circuit units and electronic circuits of the invention.
In a further aspect of this invention, methods are disclosed for maintaining the voltages in the circuits (xyz coordinates and functions of xyz coordinates) within the allowed range. Scaling factors are applied to x, y, z coordinates and to the molecular dynamics equations so as to keep the x, y, z coordinates and each term in the equations within the allowed range. In addition, resistors and amplifiers can be used to scale-down and scale-up the voltages at certain segment of the circuits. Moreover automatic gain control circuits can also be used to ensure the voltages are within the required range. The automatic gain control circuits are described in “Function circuits: Design and applications” Y. J. Wong and W. E. Ott, McGram-Hill, 1976.
Accordingly, the method further comprises maintaining the voltage level in the circuit between two fixed voltage values. In particular, x, y and z represent the coordinates of the molecule, and the voltage level of the circuit is maintained between two fixed voltage values by:
According to another aspect of the second embodiment, it is provided a circuit group representing the interaction pattern in the chemical structure of a molecule or a sub-unit of interaction pattern in the chemical structure of a molecule comprising at least one of the following:
A circuit unit according to the invention may comprise at least one circuit group as disclosed.
The invention also provides an electronic circuit comprising at least one circuit unit, the circuit unit comprising at least one circuit group, the circuit group representing a sub-unit of interaction pattern in the chemical structure of a molecule and comprising internal bond stretch, angle bending, torsion, non-bonded units; and/or bond stretch, angle bending, and/or and torsion units; between at least two nearest sub-unit of a molecule.
The electronic circuit represents a term in the molecular dynamic equation, and wherein the circuit unit comprises at least one of the following: bond stretch x unit, bond stretch y unit, bond stretch z unit, angle bending x type-A unit, angle bending x type-B unit, angle bending y type-A unit, angle bending y type-B unit, angle bending z type-A unit, angle bending z type-B unit, torsion x type-A unit, torsion x type-B unit, torsion y type-A unit, torsion y type-B unit, torsion z type-A unit, torsion z type-B unit, non-bonded x unit, non-bonded y unit, non-bonded z unit, hydrogen-bond x unit, hydrogen-bond y unit, and hydrogen-bond z unit; and wherein x, y, and z represent the coordinates of each atom of the molecule, and type-A unit represents the case of the atom being in the middle-position of an angle bending or torsion connection with other atoms, and type-B represents the case of the atom being in the end-position of an angle bending or torsion connection with other atoms.
According to a particular aspect, the invention relates to design plans for the construction of circuit units representing all possible terms in the ODEs for molecular simulation of proteins, nucleic acids and organic molecules. The units are those represented above.
The invention also provides a method for the manufacture of an electronic circuit representing at least one biomolecule and/or nano-molecular system, the electronic circuit comprising at least a unit circuit comprising at least a circuit group, wherein the circuit group represents a sub-unit of interaction pattern in the chemical structure of a molecule or a molecular system comprising:
The equations of motion of molecular dynamics simulation of a molecule can be derived from Newton's second law of motion given the potential energy functions, force fields and a starting 3D structure. The potential energy functional for proteins, nucleic acids and organic molecules is given by:
V=½ΣbondsKr(R−Req)2+½ΣanglesKθ(θ−θeq)2+½ΣtorsionsVn[1+cos(n−γ)]+
Σnon bonded[Aij/rij12−Bij/rij6+qiqj/εrrij]+ΣH bondsVH(r) (1)
as described in “A second generation force field for simulation of proteins, nucleic acids, and organic molecules”, W. D. Cornell et. al., J. Am. Chem. Soc. 117, 5179-5197 (1995); and in “CHARMM: A program for macromolecular energy, minimization, and dynamics calculations”, B. R. Brooks, et. al. J. Comp. Chem. 4, 187-217 (1983).
The hydrogen bond potential VH(r) is different for different force fields. The following is a list of some of the VH(r) appeared in the literatures:
All force field parameters for each atom or atom-atom pair in a biomolecule or nano-molecular-machine are given in the literature. For example, AMBER force fields are described in “A second generation force field for simulation of proteins, nucleic acids, and organic molecules”, W. D. Cornell et. al., J. Am.
Chem. Soc. 117, 5179-5197 (1995). As for organic molecules, all parameters other than partial charges are available. Partial charges of an organic molecule can be computed from quantum chemistry software such as Gaussian.
From Newton's second law, the atomic position of the i-th atom, in terms of its xyz coordinates (xi, yi, zi) is given by:
where mi is the mass of the i-th atom; Brxi, Bryi, BrZi are the x, y z components of the B-matrix for bond stretching; Bθxi, Bθyi, Bθzi are the x, y, z components for the B-matrix for bond angle bending; Bφxi, Bφyi, Bφzi are the x, y, z components for the B-matrix for torsion. As described in “Vibrational states”, S. Califano. John Wiley & Sons, New York (1976), these B-matrix elements are given by:
where rij, Aij, Bij, Cij, qa, qc, pa, pc, Eb, Ec, Fb, Fc, Gb, Gc, s1, S2, S3, S4, S5, S6 are given by:
rij=[(X1−xj)2+(yi−yj)2+(zi−zj)2])1/2 (15)
Aij=(xj−xi)/rij (16)
Bij=(yj−yi)/rij (17)
Cij=(zj−zi)/rij (18)
qa=cOsθa=AijAik+BijBik+CijCik (19)
qC=cos θc=−(AikAkl+BikBkl+CikCkl) (20)
pa=sinθa=(1−qa)1/2 (21)
pc=sinθc=(1−qc)1/2 (22)
Eb=(Aijqa−Aik)/rijpa (23)
Ec=(Aikqa−Aij)/rikpa (24)
Fb=(Bijqa−Bik)/rijpa (25)
Fc=(Bikqa−Bij)/rikPa (26)
Gb=(Cijqa−Cik)/rijpa (27)
Gc=(Cikqa−Cij)/rikPa (28)
S1=AikBkl−BikAkl (29)
S2=AikBij−BikAij (30)
S3=AklCik−CklAik (31)
S4=AijCik−CijAik (32)
S5=BikCkl−CikBik (33)
S6=BikCij−CijBik (34)
This set of nonlinear ODEs contain terms of rij, (xi-xj)×(xk−xj), (xi−xj)/rij sin(nφ−γ) etc., which can be represented by a corresponding set of electronic analog circuits composed of known operational amplifiers, function circuits and other circuit components.
Referring now to
Construction Plan for the Specially Designed Circuit Groups Representing all the Internal Terms and Bonded terms with nearest neighbors in the ODEs of a sub-unit of a biomolecule or in a nano-molecular-machine
Referring now to
Referring now to
In these figures, each box represents the x, y and z circuit sections of each atom of an amino acid or nucleotide. An arrow represents a connection between two atoms such that the output voltage xj, yj, or Zj from the x, y or z circuit section of an atom is directed to a bond-stretch, angle bending or torsion unit in the corresponding section of the atom to which the arrow is pointed. For clarity purpose, the connection profile for bond stretch, angle bending, and torsion is displayed in separate Figures. In reality, these should be put together in the same circuit group configuration.
A circuit group is in general constructed by the following rules:
From the chemical structure of an amino acid or nucleotide:
Referring now to
Some circuit units in
Methods for Maintaining the Voltages in the Circuits Within Allowed Range
In large proteins, the distance between two atoms can be larger than a few hundred angstroms. Thus the range of the xyz coordinates in these proteins may exceed the allowed voltages in a circuit. In addition, some terms in the molecular dynamics simulation equations can become very large if two atoms get too close or too far apart from each other, thereby causing the voltages at certain parts of the circuits to exceed the allowed range.
These problems can be solved by maintaining the voltage level in the circuit between two fixed voltage values, by using at least one of the following methods:
Structural optimization and molecular dynamics of a protein can be described by a set of ODEs. An illustrative amino acid is governed by the collection of ODEs of each constituent atom (Equation 3-5). The circuit for each atom is illustrated in
Illustrative Example: Electronic Nano-Molecular-System Emulator
There have been efforts for developing nano-molecular-systems. One example is a mimic of an enzyme ribonuclease shown in
Structural optimization and molecular dynamics of such a nano-molecular-machine can be described by a set of ODEs. An illustrative sub-unit of cyclodextrin shown in