Atractores caóticos en rutas metabólicas

Atractores ca´oticos en rutas metab´olicas Chaotic attractors in metabolic pathways Alejandro Pliego Prenda Grado en F´ısica Trabajo fin de grado Tutora: Mar´ıa del Carmen Lemos Fern´andez Departamento de F´ısica de la Materia Condensada Junio del 2020 ´ Indice 1. Introducci´on 2 2. Metabolismo celular 4 2.1. Enzimas.................................... 4 2.2. Rutasmetab´olicas............................... 5 3. Modelizaci´on de reacciones enzim´aticas 8 3.1. Ley de acci´on de masas . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 3.2. Cin´etica de Michaelis-Menten . . . . . . . . . . . . . . . . . . . . . . . . 9 3.2.1. Reacciones reversibles . . . . . . . . . . . . . . . . . . . . . . . . 12 3.2.2. Inhibici´on ............................... 15 3.2.3. Reacciones bisustrato . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.3. Funci´ondeHill ................................ 19 3.4. Modelo de Monod-Wyman-Changeux . . . . . . . . . . . . . . . . . . . . 22 4. Modelizaci´on de rutas metab´olicas 34 4.1. Modelo para la parte alta de la gluc´olisis . . . . . . . . . . . . . . . . . . 34 4.2. Modelo de Decroly y Goldbeter . . . . . . . . . . . . . . . . . . . . . . . 41 5. Conclusiones 47 Referencias 48 Anexos: 51 A. Reacciones con inhibici´on 51 A.1. Inhibici´on competitiva . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 A.2. Inhibici´on acompetitiva . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 A.3. Inhibici´on no competitiva . . . . . . . . . . . . . . . . . . . . . . . . . . 57 B. C´odigos de MatLab utilizados para realizar algunas simulaciones 62 Resumen En este trabajo se aborda el estudio de las rutas metab´olicas mediante el an´alisis de las enzimas que las componen, poniendo especial ´enfasis en la fase inicial de la gluc´olisis. Primero, se introducir´a el concepto de enzima y se presentar´an las rutas como un conjunto de reacciones qu´ımicas cuya cin´etica se caracteriza por los tipos de enzimas que act´uen en ella como catalizador. Bajo esta premisa, se estudiar´an distintos modelos capaces de explicar el comportamiento de los diferentes tipos de enzimas. Por ´ultimo, con las herramientas anteriormente mencionadas, se construir´a un modelo para la simulaci´on de la parte alta de la gluc´olisis y se buscar´an los diferentes comportamientos que puedan darse en tal sistema. Abstract In this bachelor thesis, we will tackle the study of metabolic pathways through analysis of the enzymes that make them up, focusing mainly in the preparatory phase of glycolysis. First, we will introduce the concept of enzyme, and present pathways as chemical reactions whose kinetics are characterized by the types of enzymes that act as catalyst. Under this premise, we will study several models capable of explaining the behaviour of different types of enzymes. Finally, using the aforementioned tools, we will construct a model to simulate the upper part of glycolysis, and find the different processes that can take place within the system. 1 1. Introducci´on La biolog´ıa es la rama de la ciencia que estudia aquellos fen´omenos naturales en los que interviene alg´un organismo vivo. Tanto en la definici´on cl´asica de ser vivo, en la que se considera como tal a aquellos organismos capaces de realizar las funciones de nutrici´on, relaci´on y reproducci´on, como en la definici´on m´as moderna de Matuna y Valera [1], quienes entienden que un ser vivo es cualquier sistema de mol´eculas que posee una estructura organizada capaz de autorreplicarse y de cambiar su morfolog´ıa y composici´on para adaptarse al medio externo, el concepto de ser vivo implica un continuo intercambio y modificaci´on de la materia y la energ´ıa que lo compone. La suma de los cambios de materia y energ´ıa que realiza un ser vivo con el medio externo o entre distintas partes del mismo organismo se denomina metabolismo. Los estudios m´as tempranos sobre el metabolismo de los que se tienen constancia, fueron realizados en el siglo XIII por el m´edico ´arabe Ibn al-Nafis, quien una vez dijo que tanto el cuerpo humano como sus partes se encuentran en un estado continuo de disoluci´on y nutrici´on, por lo que inevitablemente est´an sometidos a un cambio permanente [2]. Sin embargo, no fue hasta el Renacimiento con la publicaci´on de Ars de statica medicina, de la mano del m´edico italiano Santorio Santorio, cuando se empez´o a dotar a la disciplina m´edica de cierto procedimiento experimental controlado y preciso [3]. Ya en el siglo XVIII, el qu´ımico franc´es Antoine-Laurent de Lavoisier, inspirado por los experimentos con plantas de Joseph Priestley, investig´o sobre el intercambio de gases entre los seres vivos y el medio, llegando a la conclusi´on de que la respiraci´on es una forma de combusti´on [4]. Posteriormente, unos de los descubrimientos m´as relevantes en el ´area de la biolog´ıa fueron los realizados sobre la fermentaci´on por Louis Pasteur, quien descubri´o que este proceso es causado por microorganismos e ide´o un m´etodo para destruir los microbios presentes en alimentos sin alterar su composici´on (pasteurizaci´on) [5]. Tambi´en destacan los estudios llevados a cabo por Friedrich W¨ohler, quien logr´o la s´ıntesis de la urea a partir de mol´eculas inorg´anicas [6], y los elaborados por Edward Buchner, quien consigui´o realizar el proceso de fermentaci´on alcoh´olica prescindiendo de c´elulas vivas [7], trabajos que supusieron un importante cambio de paradigma en el campo de la qu´ımi- 2 ca org´anica al invalidar la teor´ıa vitalista [8]. No fue hasta el a˜no 1940 cuando Gustav Embden, Otto Meyerhof y Jacob Parnas lograron la caracterizaci´on completa de la ruta glucol´ıtica [9], lo que supuso un avance significativo para la investigaci´on del metabolismo. Estos y otros muchos descubrimientos posteriores nos han dado la posibilidad, por un lado, de mejorar nuestra comprensi´on del mundo que nos rodea, y por otro, de lograr importantes progresos tanto en nutrici´on, permiti´endonos crear productos alimenticios cada vez m´as variados y espec´ıficos, como en medicina, donde ha servido para detectar m´as enfermedades as´ı como para aumentar y mejorar los tratamientos contra ellas [10], y en la industria farmac´eutica, permiti´endonos crear productos cada vez m´as eficaces, lo que nos permite disfrutar actualmente de una vida m´as duradera y de mejor calidad. En este trabajo se presenta un estudio cin´etico de las reacciones metab´olicas. Tras definir brevemente los conceptos de enzima y ruta metab´olica (cap´ıtulo 2), se introducir´a la modelizaci´on de la cin´etica enzim´atica en el cap´ıtulo 3, donde se obtendr´an los modelos para distintos tipos de enzimas desde las m´as simples, modelo de Michaelis-Menten, hasta las enzimas alost´ericas, modelo de Monod-Wyman-Changeux (MWC). Como aplicaci´on de los conceptos introducidos, en el cap´ıtulo 4 se analizar´an dos modelos propuestos para explicar la aparici´on de comportamientos oscilatorios complejos en las rutas metab´olicas: el modelo de la parte alta de la gluc´olisis y el modelo de Decroly y Goldbeter. Finalmente, tras la secci´on de conclusiones, se incluye un anexo en el que se deduce la cin´etica para distintos tipos de enzimas susceptibles a procesos de inhibici´on, procesos que si bien no se dan en las rutas estudiadas en el presente trabajo, son muy comunes en la naturaleza. 3 2. Metabolismo celular Se conoce por metabolismo al conjunto de procesos de car´acter fisicoqu´ımico que tiene lugar en los seres vivos gracias a los cuales, estos organismos son capaces de realizar determinadas funciones vitales [11]. En general, dichas funciones requieren, por un lado, de la s´ıntesis de biomol´eculas que realiza el organismo para poder renovar, cambiar o producir partes del mismo, y por otro lado, de la obtenci´on de materia y energ´ıa del medio externo. Tanto la formaci´on de biomol´eculas complejas como la producci´on de energ´ıa, son procesos que habitualmente dependen de una serie de reacciones qu´ımicas concretas, que ocurren de forma eficiente en los organismos vivos, y que se han ido optimizando por selecci´on natural durante miles de millones de a˜nos de evoluci´on [12]. 2.1. Enzimas Se denominan enzimas a aquellas prote´ınas que act´uan como catalizadores en ciertas reacciones qu´ımicas que tienen lugar en las c´elulas [13]. Estas mol´eculas est´an constituidas por una o varias cadenas polipept´ıdicas, es decir, secuencias de amino´acidos unidos mediante enlaces pept´ıdicos (apoenzima), aunque a veces contienen adem´as un grupo no prot´eico (grupo prost´etico o cofactor) que puede ser un cati´on met´alico o una mol´ecula org´anica (coenzima) [14]. Para que una reacci´on qu´ımica tenga lugar, primero los sustratos deben alcanzar el estado de activaci´on, un estado de transici´on energ´etica en el cual algunos de sus enlaces qu´ımicos son lo suficientemente d´ebiles como para que, posteriormente, se puedan formar unos enlaces distintos cambiando as´ı la naturaleza de los reactivos y dando lugar a los productos. La energ´ıa necesaria para que los sustratos alcancen el estado de transici´on se denomina energ´ıa de activaci´on. 4 El papel de las enzimas consiste en disminuir la energ´ıa de activaci´on de determinadas reacciones (figura 1) aumentando as´ı la velocidad a la que estas ocurren y, en ocasiones, logrando que dichas reacciones puedan ocurrir finalmente. Aunque muchas de las reacciones qu´ımicas que tienen lugar en las c´elulas pueden darse tambi´en en ambientes extracelulares con id´enticas condiciones de presi´on y temperatura, la acci´on enzim´atica es capaz de aumentar la velocidad de estos procesos en hasta 17 ´ordenes de magnitud [15]. Figura 1: Representaci´on del cambio de energ´ıa libre (∆G) en el curso de una reacci´on, en la que se establece la diferencia entre el balance energ´etico de una reacci´on qu´ımica con (curva roja) y sin (curva negra) la intervenci´on de catalizadores. Como puede observarse, la actuaci´on de la enzima consigue disminuir la energ´ıa de activaci´on, favoreciendo as´ı la reacci´on qu´ımica. Imagen tomada de [16]. Cada reacci´on catalizada es regulada por una enzima espec´ıfica. La especificidad de las enzimas se debe a su compleja estructura proteica, por la cual las enzimas presentan al menos una zona (centro activo o de uni´on) con la estructura concreta para acoplarse a un sustrato determinado y as´ı alcanzar el estado de transici´on. 2.2. Rutas metab´olicas Se puede definir una ruta metab´olica como una secuencia de reacciones enzim´aticas encadenadas en las que el producto de una reacci´on es, a su vez, el sustrato de la siguiente. Gracias a las rutas metab´olicas, unas mol´eculas iniciales pueden dar lugar a otras mol´eculas m´as complejas y ´utiles para el organismo (anabolismo), o bien, pueden producir otras mol´eculas m´as simples liberando energ´ıa en el proceso (catabolismo). Como ejemplo de una de estas rutas, en el transcurso del presente trabajo se estudiar´a la ruta glucol´ıtica, m´as concretamente aquella que se da en la levadura y el m´usculo esquel´etico. Esta ruta es una de las m´as estudiadas debido a su universalidad (ya que 5 es realizada por casi todos los seres vivos) y por ser precursora de diferentes procesos fisiol´ogicos imprescindibles para la vida. Se cree que la ruta glicol´ıtica es una de las rutas metab´olicas m´as antiguas que se conocen, puesto que esta puede darse perfectamente en un ambiente anaerobio (sin ox´ıgeno) como el que exist´ıa en la atm´osfera primitiva [17]. La gluc´olisis tiene lugar en el interior celular, concretamente en el citoplasma, y consta de diez reacciones qu´ımicas en las que la glucosa Glu es el sustrato inicial que se descompone finalmente en dos mol´eculas de piruvato Pyr, consumiendo en el transcurso dos mol´eculas de adenos´ın trifosfato o ATP en la parte alta del proceso, pero produciendo otras cuatro mol´eculas de ATP en su parte baja. En suma, por cada mol´ecula de glucosa que se consume en la c´elula se generan dos nuevas mol´eculas de ATP, las cuales almacenan en sus enlaces gran cantidad de energ´ıa que puede ser consumida para llevar a cabo otros muchos procesos metab´olicos. La parte alta de la gluc´olisis (1) consiste en la activaci´on de la glucosa para obtener energ´ıa de ella posteriormente. Esta parte consta de cuatro reacciones qu´ımicas catalizadas por diferentes enzimas. Glu ATP ADP G6P ATP ADP v5 v2F6P ATP ADP FBP v4 v1v3(1) En primer lugar, la enzima hexoquinasa fosforila a la glucosa para producir glucosa- 6-fosfato (G6P), a˜nadiendo un grupo fosfato al sexto carbono de la glucosa proveniente de la descomposici´on de una mol´ecula de ATP en una mol´ecula de adenos´ın difosfato o ADP. Luego, parte de la glucosa-6-fosfato (G6P) se puede consumir en otras rutas, pero otra parte continuar´a la ruta glicol´ıtica sufriendo un proceso de isomerizaci´on regulado por la enzima fosfoglucoisomerasa, que convertir´a la glucosa-6-fosfato en fructosa-6-fosfato (F6P), mol´ecula que contiene los mismos elementos que su predecesora pero formando distintos enlaces. La tercera reacci´on es una nueva fosforilaci´on, esta vez llevada a cabo por la enzima 6 De acuerdo con el esquema (27) la concentraci´on de complejo enzima-sustrato crece con velocidad v1y disminuye con velocidad v2. d[ES] dt =v1−v2=k1[E][S]−k−1[ES]−k2[ES] + k−2[E][P] (30) De nuevo, es posible aplicar la condici´on de estado cuasiestacionario (18). (k1[S] + k−2[P])[E] = (k−1+k2)[ES] (31) A continuaci´on, sustituyendo la concentraci´on de enzima libre por su equivalente mostrado en (21) en las ecuaciones (29) y (31), se obtienen (32) y (33) respectivamente. v= (k2+k−2[P])[ES]−k−2[Etotal][P] (32) (k1[S] + k−2[P])[Etotal]−(k1[S] + k−2[P])[ES]=(k−1+k2)[ES] (33) De esta ´ultima expresi´on, se puede despejar la concentraci´on de complejo enzima sustrato. [ES] = (k1[S] + k−2[P])[Etotal] k−1+k2+k1[S] + k−2[P](34) Ahora ya puede sustituirse [ES] en (32). v=(k2+k−2[P])(k1[S] + k−2[P])[Etotal] k−1+k2+k1[S] + k−2[P]−k−2[P](k−1+k2+k1[S] + k−2[P])[Etotal] k−1+k2+k1[S] + k−2[P] (35) Lo que es equivalente a: v= [Etotal]k2k1[S]−k−2k−1[P] k−1+k2+k1[S] + k−2[P](36) Esta es la ecuaci´on cin´etica para el caso reversible seg´un Michaelis-Menten, pero su expresi´on puede reducirse m´as dividiendo tanto el numerador como el denominador por la cantidad k−1+k2. v= [Etotal] k2k1 k−1+k2[S]−k−1k−2 k−1+k2[P] 1 + k1 k−1+k2[S] + k−2 k−1+k2[P] (37) 13 Recurriendo a la constante de Michaelis-Menten (26), y definiendo una nueva constante K−Man´aloga a la anterior (38), finalmente se alcanza la f´ormula (39). K−M=k−1+k2 k−2 (38) v= [Etotal] k2 KM[S]−k−1 K−M[P] 1 + [S] KM+[P] K−M (39) Para simplificar a´un m´as esta ecuaci´on, se puede recurrir a la velocidad m´axima Vmax, que se alcanzar´a cuando la concentraci´on de producto sea nula o la constante k−2 sea despreciable respecto a las dem´as, con lo que K−Mtender´ıa a infinito y la ecuaci´on de velocidad ser´ıa equivalente a la ecuaci´on cin´etica del caso Michaelis-Menten irreversible (25), es decir, Vmax =k2[Etotal]. De igual modo, existir´a una velocidad m´ınima Vmin que se alcanzar´a cuando el sustrato sea nulo o si la constante k1fuera despreciable, de este modo, tendr´ıamos Vmin =k−1[Etotal], lo cual es el equivalente a la velocidad m´axima en una ecuaci´on cin´etica tipo Michaelis-Menten que ocurre en sentido opuesto. v= Vmax KM[S]−Vmin K−M[P] 1 + [S] KM+[P] K−M (40) Figura 3: En esta gr´afica se muestra la velocidad de reacci´on (en sentido de la producci´on de producto) frente a la concentraci´on de sustrato. Las l´ıneas m´as azules indican una mayor concentraci´on de producto respecto a la constante K−M, en concreto, las cantidades de [P]/K−Maqu´ı representadas son 0, 0,1, 0,5, 1,5 y 3. Puede observarse que cuando la cantidad de producto es mayor, para bajas concentraciones de sustrato, la velocidad es negativa, esto es, que la reacci´on ocurre en sentido opuesto consumiendo producto para producir sustrato. En esta representaci´on, la constante KMtoma el valor 0,2M. 14 3.2.2. Inhibici´on Existen unas mol´eculas llamadas efectores, las cuales pueden influir en la acci´on enzim´atica alterando la cin´etica del proceso ya sea activando o, m´as com´unmente, inhibiendo a la enzima. Los inhibidores enzim´aticos reprimen la producci´on de una enzima cambiando sustancialmente su cin´etica. En el anexo A de este documento se estudian en detalle los tres casos m´as comunes de procesos de inhibici´on (inhibici´on competitiva, acompetitiva y no competitiva), deduciendo la ecuaci´on cin´etica para cada caso y representando algunos ejemplos de dichos comportamientos en funci´on de las concentraciones de sustrato e inhibidor presentes en el sistema. Un ejemplo de inhibici´on se encuentra en algunos tipos de hexoquinasa, como la hay en el tejido de la corteza cerebral donde existe una hexoquinasa capaz de ser inhibida por su propio producto, la glucosa-6-fosfato (G6P) [23]. 3.2.3. Reacciones bisustrato Existen enzimas con dos centros catal´ıticos en las que cada uno de estos centros es espec´ıfico para una especie qu´ımica diferente. Las reacciones reguladas por este tipo de enzimas se denominan reacciones bisustrato, y suelen escribirse de la siguiente forma: S1 S2P2 P1(41) Al tratarse de una enzima con dos centros de uni´on diferentes, existen dos caminos que la enzima puede tomar para unirse a los sustratos y dar lugar a los productos: unirse en primer lugar al sustrato S1para formar el complejo enzim´atico ES1y luego unirse al otro sustrato S2o, por el contrario, formar el complejo enzim´atico ES2en primer lugar y luego de eso unirse al sustrato S1[24]. 15 E+S1+S2 k1 k−1 ES1+S2 k2 k−2 ES1S2 k0 1 k0 −1 ES2+S1 k0 2 k0 −2 k3E+P1+P2(42) Al igual que en los casos anteriores, en primer lugar se calculan las velocidades para cada uno de los tramos. v1=k1[E][S1]−k−1[ES1] (43) v2=k2[ES1][S2]−k−2[ES1S2] (44) v0 1=k0 1[E][S2]−k0 −1[ES2] (45) v0 2=k0 2[ES2][S1]−k0 −2[ES1S2] (46) v=v3=k3[ES1S2] (47) En esta ocasi´on, v3se considera la velocidad de la reacci´on, pues es la que marca la tasa de producci´on. Aplicando la condici´on cuasiestacionaria a los distintos complejos enzim´aticos, se tienen las siguientes relaciones entre velocidades: d[ES1] dt =v1[S2]−v2= 0 →k1[E][S1][S2] + k−2[ES1S2]=(k−1+k2)[ES1][S2] (48) d[ES2] dt =v0 1[S1]−v0 2= 0 →k0 1[E][S1][S2] + k0 −2[ES1S2]=(k0 −1+k0 2)[ES2][S1] (49) d[ES1S2] dt =v2+v0 2−v3→k2[ES1][S2] + k0 2[ES2][S1] = (k−2+k0 −2+k3)[ES1S2] (50) 16 A continuaci´on, puede despejarse la concentraci´on del complejo enzima-sustrato ES1 en la ecuaci´on (48) y la concentraci´on del complejo ES2en la ecuaci´on (49): [ES1] = k1[S1][S2][E] + k−2[ES1S2] (k−1+k2)[S2](51) [ES2] = k0 1[S1][S2][E] + k0 −2[ES1S2] (k0 −1+k0 2)[S1](52) Sustituyendo [ES1] (51) y [ES2] (52) en (50), se tiene una expresi´on de [E] y [ES1S2]. k2(k0 −1+k0 2)(k1[S1][S2][E] + k−2[ES1S2]) + k0 2(k−1+k2)(k0 1[S1][S2][E] + k0 −2[ES1S2]) = = (k−1+k2)(k0 −1+k0 2)(k−2+k0 −2+k3)[ES1S2] (53) Despejando la cantidad [E], se encuentra la concentraci´on de enzima libre en funci´on de la concentraci´on del complejo enzim´atico ES1S2. [E] = k−1k−2(k0 −1+k0 2) + k0 −1k0 −2(k−1+k2)+(k−1+k2)(k0 −1+k0 2)k3 k1k2(k0 −1+k0 2)[S1][S2] + k0 1k0 2(k−1+k2)[S1][S2][ES1S2] (54) Asumiendo que la cantidad total de enzima involucrada en la reacci´on permanece constante, la concentraci´on de enzima total [Etotal] puede expresarse como la suma de la concentraci´on de enzima en estado libre [E] m´as la concentraci´on de enzima que forma parte de alg´un complejo enzim´atico [ES1], [ES2] y [ES1S2]: [Etotal] = [E]+[ES1]+[ES2]+[ES1S2] (55) En esta f´ormula se pueden sustituir ahora los valores de [ES1] y [ES2] obtenidos en (51) y (52) respectivamente. (k−1+k2)(k0 −1+k0 2)[S1][S2][Etotal]=(k−1+k2)(k0 −1+k0 2)[S1][S2][E]+ +k1(k0 −1+k0 2)[S1]2[S2][E] + k−2(k0 −1+k0 2)[S1][ES1S2] + k0 1(k−1+k2)[S1][S2]2[E]+ +k0 −2(k−1+k2)[S2][ES1S2]+(k−1+k2)(k0 −1+k0 2)[S1][S2][ES1S2] (56) 17 Y para eliminar la dependencia de [E], puede sustituirse esta por la ecuaci´on (54). (k−1+k2)(k0 −1+k0 2)(k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2))[S1]2[S2]2[Etotal] = = (k−1k−2(k0 −1+k0 2)+k0 −1k0 −2(k−1+k2)+(k−1+k2)(k0 −1+k0 2)k3)(k−1+k2)(k0 −1+k0 2)[S1][S2][ES1S2]+ +(k−1k−2(k0 −1+k0 2)+k0 −1k0 −2(k−1+k2)+(k−1+k2)(k0 −1+k0 2)k3)k1(k0 −1+k0 2)[S1]2[S2][ES1S2]+ +(k−1k−2(k0 −1+k0 2)+k0 −1k0 −2(k−1+k2)+(k−1+k2)(k0 −1+k0 2)k3)k0 1(k−1+k2)[S1][S2]2[ES1S2]+ + (k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2))k−2(k0 −1+k0 2)[S1]2[S2][ES1S2]+ + (k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2))k0 −2(k−1+k2)[S1][S2]2[ES1S2]+ + (k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2))(k−1+k2)(k0 −1+k0 2)[S1]2[S2]2[ES1S2] (57) Donde reordenando t´erminos y dividiendo ambos lados de la ecuaci´on por la cantidad (k−1+k2)(k0 −1+k0 2)[S1][S2], se tiene: (k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2))[S1][S2][Etotal] = = (k−1k−2(k0 −1+k0 2) + k0 −1k0 −2(k−1+k2))[ES1S2]+(k−1+k2)(k0 −1+k0 2)k3[ES1S2]+ + (k1k0 −1(k−2+k0 −2) + k−2k0 2(k1+k0 1))[S1][ES1S2]+(k0 −1+k0 2)k1k3[S1][ES1S2] + (k−1k0 1(k−2+k0 −2) + k2k0 −2(k1+k0 1))[S2][ES1S2]+(k−1+k2)k0 1k3[S2][ES1S2] + (k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2))[S1][S2][ES1S2] (58) De esta ecuaci´on puede despejarse [ES1S2] en funci´on de la concentraci´on total de enzima [Etotal]. Se llega as´ı a la expresi´on: [ES1S2] = [Etotal][S1][S2] KII +KS1[S1] + KS2[S2]+[S1][S2](59) Donde las constantes KS1,KS2yKII se definen como: KS1=(k1k0 −1(k−2+k0 −2) + k−2k0 2(k1+k0 1)) + (k0 −1+k0 2)k1k3 k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2)(60) KS2=(k−1k0 1(k−2+k0 −2) + k2k0 −2(k1+k0 1)) + (k−1+k2)k0 1k3 k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2)(61) KII =k−1k−2(k0 −1+k0 2) + k0 −1k0 −2(k−1+k2)+(k−1+k2)(k0 −1+k0 2)k3 k1k2(k0 −1+k0 2) + k0 1k0 2(k−1+k2)(62) 18 Volviendo ahora a la expresi´on para la velocidad (47), la f´ormula cin´etica puede expresarse como: v=Vmax[S1][S2] KII +KS1[S1] + KS2[S2]+[S1][S2](63) Donde la velocidad m´axima en este caso es Vmax =k3[Etotal]. 3.3. Funci´on de Hill Existen enzimas, con m´as de un centro activo, cuya cin´etica no puede explicarse utilizando los m´etodos hasta ahora descritos; estas enzimas se denominan enzimas alost´ericas. Uno de los primeros intentos de entender el comportamiento de este tipo de prote´ınas fue llevado a cabo por Archibald Vivian Hill a principios del siglo XX [25]. Se dice ligando a cualquier mol´ecula que pueda unirse temporalmente a una prote´ına. Se distinguen dos grupos de prote´ınas en funci´on de cu´antos ligandos pueda albergar: los mon´omeros, que pueden unirse a un s´olo ligando, y los olig´omeros (d´ımeros, tr´ımeros, tetr´ameros...), que poseen varios prot´omeros (o sitios de uni´on en el caso de las enzimas). Un ejemplo de tetr´amero es la hemoglobina, la cual no es una enzima. Esta prote´ına se une con varias mol´eculas de ox´ıgeno o de mon´oxido de carbono para transportarlas a trav´es del flujo sangu´ıneo [26]. En la secci´on dedicada a la cin´etica de Michaelis-Menten se estudi´o el caso de un mon´omero que admit´ıa unirse a un s´olo sustrato por vez (11). E+Sk1 k−1 ES (64) De acuerdo con la ley de acci´on de masas, la constante de uni´on vendr´a dada por: KB=k1 k−1 =[ES] [E][S]eq =K−1 D(65) 19 Se define, adem´as, la fracci´on de saturaci´on Yde una enzima como la cantidad de centros ocupados dividida entre la cantidad total de centros de todas las enzimas de la misma especie. En el caso de la reacci´on catalizada de Michaelis-Menten, ser´ıa: Y=[ES] [Etotal]=[ES] [ES]+[E]=KB[S] KB[S]+1 =[S] [S] + KD (66) Cuando se estudi´o el modelo cin´etico de Michaelis-Menten, se asumi´o que k1, k−1>> k2, por lo tanto, es posible aproximar la constante de disociaci´on a la constante de Michaelis-Menten KD≈KM. De esta forma, la f´ormula para la fracci´on de saturaci´on queda como: Y=[S] [S] + KM (67) Y comparando esta expresi´on para Ycon la ecuaci´on cin´etica (25), es posible encontrar una f´ormula que relacione la velocidad con la fracci´on de saturaci´on. Y=v Vmax (68) Si la enzima tiene varios centros de uni´on, entonces ocurre un fen´omeno de cooperaci´on. La cooperaci´on se da cuando un ligando se une a la enzima y eso influye (positiva o negativamente) en la futura uni´on de la enzima a otros ligandos. Para una enzima En, con ncentros de uni´on, se pueden dar las siguientes reacciones qu´ımicas: EnS+SEnS2 En+SEnS . . . EnSn−1+SEnSn (69) En este caso, la enzima tiene un ´unico ligando unido a ella en el complejo EnS, pero dos ligandos en el complejo EnS2, hasta llegar a los nligandos en el complejo enzim´atico EnSn, por lo que la fracci´on de saturaci´on para esta reacci´on ser´a: Y=[EnS] + 2[EnS2] + · · · +n[EnSn] n[Etotal]=[EnS] + 2[EnS2] + · · · +n[EnSn] n([En]+[EnS] + · · · + [EnSn]) (70) 20 Si la afinidad del ligando con una enzima es mayor cuanto m´as sustratos est´en unidos a esta (cooperaci´on positiva), entonces se puede asumir que cada complejo enzima-sustrato reaccionar´a con un nuevo sustrato en cuanto este se forme y que, por tanto, la concentraci´on de los complejos enzim´aticos intermedios puede suponerse nula ([EnSi]≈0 ∀i= 1, ..., n −1) [27]. Teniendo esto en cuenta, el conjunto de reacciones qu´ımicas (69) se reduce a s´olo una. En+nS EnSn(71) Donde la constante de uni´on es: KB=[EnSn] [En][S]neq (72) Y la fracci´on de saturaci´on, teniendo en cuenta la aproximaci´on queda como: Y=[EnSn] [En]+[EnSn]=KB[S]n 1 + KB[S]n(73) Recurriendo a la relaci´on entre la velocidad y la fracci´on de saturaci´on (68), puede obtenerse la ecuaci´on cin´etica. v=Vmax[S]n KD+ [S]n(74) Donde nes denominado coeficiente de Hill. Figura 4: La gr´afica de la izquierda muestra la velocidad de reacci´on frente a la concentraci´on de sustrato para coeficientes de Hill con valores 1 (curva negra), 2, 3 y 4 (curva m´as azul). La segunda gr´afica muestra el diagrama de Lineweaver-Burk correspondiente. Puede observarse c´omo para coeficientes de Hill mayores que n= 1 la curva toma forma sigmoidea. Para esta representaci´on se ha tomado KD= 0,5M. 21 3.4. Modelo de Monod-Wyman-Changeux El modelo de Hill supuso un importante paso para explicar los mecanismos que rigen el comportamiento de las enzimas alost´ericas. No obstante, con el tiempo ha sido necesario desarrollar otros modelos m´as precisos. Uno de estos, es el modelo MWC o modelo de Monod-Wyman-Changeux [28]. Este modelo supone que las enzimas alost´ericas est´an compuestas por varias subunidades en las que existen un sitio activo, al que puede unirse el sustrato, y un centro de regulaci´on, al que puede unirse una mol´ecula llamada efector o modulador enzim´atico. Cuando un efector se une a la enzima, esta sufre cambios conformacionales capaces de alterar significativamente la actividad enzim´atica [29]. La uni´on de un determinado efector a un centro de regulaci´on puede facilitar la uni´on de un sustrato con el centro activo de la subunidad correspondiente. El efector puede ser un activador A, si cuando este se une a la enzima se produce un fen´omeno de cooperaci´on positiva, o un inhibidor I, si la uni´on del efector al centro de regulaci´on dificulta que otros ligandos puedan acoplarse al sitio activo. Adem´as de lo ya mencionado, el modelo MWC se basa en las siguientes suposiciones: Cada enzima posee nsubunidades id´enticas entre s´ı. Cada subunidad puede presentar una forma activa Ro una forma inactiva T, que poseen distinta afinidad con el sustrato. En general, cada una de estas formas (activas e inactivas) comprenden un centro activo, al que puede acoplarse un sustrato, y uno o varios centros alost´ericos o reguladores, a los que pueden unirse un efector. Se asume que todas las subunidades de una misma enzima cambian de forma simult´aneamente y s´olo en ausencia de efectores ligados. La relaci´on de equilibrio entre las formas activas e inactivas viene dada por la constante alost´erica L. 22 En este punto, de acuerdo con (76), la concentraci´on de enzima libre en forma inactiva [T00] se puede relacionar con la concentraci´on de enzima libre en forma activa [R00] mediante las constantes cin´eticas k1yk2. [T00] = k1 k2 [R00] (106) Se tiene entonces que: [Tij] = n! (n−i)! i!a0 1[S] d0 1+k0in! (n−j)! j!a0 2[I] d0 2jk1 k2 [R00]∀j= 1, ..., n (107) Recordando ahora las definiciones (88) y teniendo en cuenta tanto (104) como (107), las expresiones que toman [Λi] y [Γi] son: [Λi] = n! (n−i)! i!a1[S] d1+kin X j=0 n! (n−j)! j!a2[A] d2j [R00] (108) [Γi] = n! (n−i)! i!a0 1[S] d0 1+k0in X j=0 n! (n−j)! j!a0 2[I] d0 2jk1 k2 [R00] (109) Para simplificar estas ecuaciones, resulta ´util introducir el teorema del binomio: (y+x)n= n X i=0 n! (n−i)! i!yn−ixi(110) Usando este teorema, cuando y= 1 se obtiene una relaci´on entre los factoriales similar a las ya obtenidas en (108) y en (109). n X i=0 n! (n−i)! i!xi= (1 + x)n(111) Luego: [Λi] = n! (n−i)! i!a1[S] d1+ki1 + a2[A] d2n [R00]∀i= 0, ..., n (112) [Γi] = n! (n−i)! i!a0 1[S] d0 1+k0i1 + a0 2[I] d0 2nk1 k2 [R00]∀i= 0, ..., n (113) 29 Para simplificar estas ecuaciones, conviene realizar algunos cambios de variables: KR≡d1 a1 α≡[S] KR ≡k d1 e≡1 1 + (114) A partir de estas definiciones, se puede demostrar la siguiente igualdad a la que se recurrir´a m´as adelante. a1[S] d1+k=a1[S] d1 d1 d1+k=α1 1 + =αe (115) Otros cambios de variables que se tomar´an son: KT≡d0 1 a0 1 c≡KR KT 0≡k0 d0 1 e0≡1 1 + 0(116) Con los que se puede demostrar otra expresi´on que se usar´a a continuaci´on. a0 1[S] d0 1+k0=a0 1[S] d0 1 d0 1 d0 1+k0=[S] KT 1 1 + 0=c[S] KR e0=cαe0(117) Hay dos constantes m´as que tambi´en ayudar´an a simplificar estas ecuaciones: la constante de disociaci´on del activador KAy del inhibidor KI. KA≡d2 a2 KI≡d0 2 a0 2 (118) Sustituyendo por las constantes que se acaban de definir en (112) y (113), y teniendo en cuenta las equivalencias (115) y (117): [Λi] = n! (n−i)! i!(αe)i1 + [A] KAn [R00]∀i= 0, ..., n (119) [Γi] = n! (n−i)! i!(cαe0)i1 + [I] KIn L[R00]∀i= 0, ..., n (120) Donde L≡k1/k2es la ya nombrada constante alost´erica, e indica la relaci´on de proporcionalidad que se da entre la forma activa e inactiva de la enzima libre. 30 Se busca ahora expresar estas cantidades en funci´on de la concentraci´on total de enzima [Etotal] disponible en el sistema. Para ello, en primer lugar, se sustituyen las expresiones obtenidas en (89). [Etotal] = n X i=0 n!(αe)i (n−i)! i!1 + [A] KAn [R00] + n X i=0 n!(cαe0)i (n−i)! i!1 + [I] KIn L[R00] (121) Y usando de nuevo la ecuaci´on hallada a partir del teorema del binomio (111): [Etotal] = (1 + αe)n1 + [A] KAn [R00] + (1 + cαe0)n1 + [I] KIn L[R00] (122) Ahora puede despejarse [R00] en funci´on de [Etotal]. [R00] = [Etotal] L(1 + cαe0)n1 + [I] KIn+ (1 + αe)n1 + [A] KAn(123) Volviendo a la ecuaci´on de velocidad (87): v=k n X i=0 i[Λi] + k0 kk n X i=0 i[Γi] (124) Es posible sustituir los valores de [Γi] y [Λi] obtenidos en (120) y (119) respectivamente. v=k n X i=0 in!(αe)i (n−i)! i!1 + [A] KAn [R00] + θk n X i=0 in!(cαe0)i (n−i)! i!1 + [I] KIn L[R00] (125) Donde θ≡k0/k se define como una proporci´on entre las constantes cin´eticas para la generaci´on de producto con la subunidad en forma activa e inactiva. Para eliminar los sumatorios de estas ecuaciones, se puede recurrir una vez m´as al teorema del binomio, derivando respecto a xla f´ormula (111). n(1 + x)n−1= n X i=0 in! (n−i)! i!xi−1(126) Y posteriormente, multiplicando por xambos t´erminos de la igualdad, obteniendo as´ı la siguiente f´ormula: n X i=0 in! (n−i)! i!xi=nx(1 + x)n−1(127) 31 Sustituyendo la expresi´on del sumatorio (127) en (125), puede expresarse la velocidad de producci´on del sistema sin necesidad de sumatorios. v=knαe(1 + αe)n−11 + [A] KAn [R00] + θkncαe0(1 + cαe0)n−11 + [I] KIn L[R00] (128) As´ı, se halla finalmente la ecuaci´on cin´etica, que es posible expresar en funci´on de la cantidad total de enzima gracias a (123). v=nk[Etotal] αe(1 + αe)n−11 + [A] KAn+θLcαe0(1 + cαe0)n−11 + [I] KIn L(1 + cαe0)n1 + [I] KIn+ (1 + αe)n1 + [A] KAn(129) Para obtenerse la expresi´on simplificada con la que suele expresarse esta ecuaci´on, es necesario definir la velocidad normalizada ν≡v KR. As´ı: v=KRσφ (130) Donde σes una cantidad constante distinta para cada tipo de enzima. σ=nk[Etotal] KR (131) Yφdepende de las concentraciones de las distintas sustancias implicadas en el proceso. φ= αe(1 + αe)n−11 + [A] KAn+θLcαe0(1 + cαe0)n−11 + [I] KIn L(1 + cαe0)n1 + [I] KIn+ (1 + αe)n1 + [A] KAn(132) Las ecuaciones que explican la cin´etica alost´erica de acuerdo al modelo MWC, suelen expresarse mediante el siguiente sistema normalizado de ecuaciones:    ˙α=ν−σφ ˙γ=qσφ −ksγ (133) Donde αes la concentraci´on de sustrato [S] dividida entre la constante de disociaci´on de los complejos en forma activa KR,γes la concentraci´on de producto [P] dividida entre la constante de disociaci´on de los complejos en forma inactiva KP,ν=v/KRes la velocidad normalizada y q=KR/KPes una constante de normalizaci´on. Lo que es igual al siguiente sistema de ecuaciones: 32    ˙ [S] KR=v KR −σφ ˙ [P] KP=KR KPσφ −ks[P] KP (134) La constante eindica la preferencia de una subunidad activa por desprenderse del sustrato (e≈1) frente a transformarlo en producto (e≈0), y la constante e0indica lo mismo para la forma inactiva. Aparece adem´as la constante θque determina la raz´on entre la constante cin´etica para la forma inactiva y activa, y que normalmente adquirir´a un valor peque˜no (θ << 1). El modelo MWC ofrece respuestas para el comportamiento de enzimas alost´ericas susceptibles al efecto de distintas mol´eculas activadoras e inhibidoras de la acci´on enzim´atica, por lo que suele ser muy usado para la modelizaci´on de rutas metab´olicas. Un ejemplo de esto puede encontrarse en los experimentos de Markus y Hess [30], quienes construyen un modelo simplificado para la ruta glucol´ıtica considerando ´unicamente dos reacciones enzim´aticas: la fosforilaci´on de la fructosa-6-fosfato (F6P) realizada por la enzima fosfofructoquinasa y la transferencia de un grupo fosfato del fosfoenolpiruvato (PEP) al ADP realizada por la enzima piruvato quinasa. En este modelo c´ıclico se suponen constantes las cantidades de ADP y ATP del sistema y se controla la entrada tanto de la fructosa-6-fosfato (F6P) como del fosfoenolpiruvato (PEP), siendo los productos finales de la reacci´on la fructosa-1,6-bisfosfato (FBP) y piruvato (Pyr). ATP ADP ATP + F6P FBP Pyr PEP Figura 5: Modelo simplificado para la gluc´olisis propuesto por Markus y Hess. Markus y Hess proponen el modelo mostrado en la figura 5 para explicar el comportamiento ca´otico observado en sus experimentos con extractos libres de Sacaromices cerevisae (levadura de la cerveza), bajo un suministro peri´odico de glucosa [31]. 33 4. Modelizaci´on de rutas metab´olicas En la presente secci´on se construir´an y analizar´an distintos modelos matem´aticos que describen la cin´etica enzim´atica de acuerdo a las ecuaciones obtenidas previamente. 4.1. Modelo para la parte alta de la gluc´olisis Para formar el modelo cin´etico de la parte alta de la gluc´olisis es necesario, en primer lugar, describir las reacciones principales de la ruta glicol´ıtica, mediante las cuales la glucosa se convierte finalmente en dos triosas (135), y en segundo lugar, a˜nadir las ecuaciones correspondientes a aquellas reacciones qu´ımicas que ocurren en el mismo entorno y que poseen la capacidad de alterar las concentraciones de los reactivos que intervienen en la ruta principal. En este caso, estas reacciones adicionales son la fosforilaci´on de las mol´eculas de ADP (136), el intercambio de grupos fosfato entre nucle´otidos (138) y el consumo de energ´ıa expresado como la p´erdida de grupos fosfato en las mol´eculas de ATP, que puede expresarse como una sola reacci´on (137). Glu ATP ADP G6P ATP ADP v5 v2F6P ATP ADP FBP v4 v1v3(135) ADP v6ATP (136) ATP v7ADP (137) ATP+AMP v82ADP (138) El sistema de ecuaciones diferenciales de primer orden que constituyen el sistema enzim´atico, y que describen la evoluci´on temporal de las concentraciones de las distintas sustancias en funci´on de las velocidades a las que se producen cada una de las reacciones que tienen lugar en el medio es: d[G6P] dt =v1−v2−v5(139) 34 d[F6P] dt =v2−v3(140) d[FBP] dt =v3−v4(141) d[ATP] dt =−v1−v3−v5+v6−v7−v8(142) d[ADP] dt =v1+v3+v5−v6+v7+ 2v8(143) d[AMP] dt =−v8(144) Las velocidades de reacci´on v4,v5,v6,v7yv8que aparecen en el sistema, se obtienen mediante la aplicaci´on directa de la ley de acci´on de masas para reacciones irreversibles en los primeros casos (8) y para reacciones reversibles en el ´ultimo de los casos (7). v4=k4[FBP] (145) v5=k5[G6P][ATP] (146) v6=k6[ADP] (147) v7=k7[ATP] (148) v8=k8[ATP][AMP]−k−8[ADP]2(149) Sin embargo, el resto de velocidades adquieren expresiones m´as complicadas al tratarse de reacciones catalizadas por enzimas espec´ıficas. La enzima que regula la primera reacci´on de la ruta glucol´ıtica es la hexoquinasa. Esta prote´ına posee dos centros activos diferentes, uno espec´ıfico para la glucosa y otro para el ATP. Cuando una mol´ecula de glucosa y otra de ATP se unen a esta enzima, inducen en ella un cambio conformacional que propicia el proceso de fosforilaci´on, dando como 35 resultado una mol´ecula de glucosa-6-fosfato (G6P) y otra de ADP. Se deduce entonces que la reacci´on que regula la hexoquinasa es una reacci´on bisustrato cuya velocidad vendr´ıa dada por la ecuaci´on (63). v1=V1max[Glu][ATP] KII +KGlu[Glu] + KATP [ATP]+[Glu][ATP](150) La segunda reacci´on de la ruta es regulada por la fosfoglucoisomerasa, una enzima que es capaz de alterar la estructura molecular de la glucosa-6-fosfato (G6P) para transformarla en fructosa-6-fosfato (F6P), pero tambi´en puede de realizar el procedimiento en sentido opuesto, partiendo de la F6Py obteniendo G6Pcomo producto. Se trata, entonces, de un ejemplo de cin´etica de tipo Michaelis-Menten reversible (40). v2= V2max KG6P[G6P]−V2min KF6P[F6P] 1 + [G6P] KG6P+[F6P] KF6P (151) Por ´ultimo, el paso de fructosa-6-fosfato (F6P) a fructosa-1,6-bisfosfato (FBP) es regulado por la fosfofructoquinasa, una enzima alost´erica que posee dos subunidades id´enticas. Cada una de estas subunidades cuenta con un centro activo al que se puede unir el sustrato y un centro alost´erico al que se puede ligar un efector. Si la enzima se encuentra en su forma activa, el AMP, que act´ua como activador, puede unirse al sitio alost´erico propiciando la formaci´on de producto. Si por el contrario, la enzima se encuentra en su forma inactiva, es el ATP lo que puede acoplarse a su centro alost´erico, pero este dificultar´a la formaci´on de producto, pues se trata de un inhibidor enzim´atico. Para encontrar la ecuaci´on de velocidad v3que describe la cin´etica de esta enzima se puede recurrir al modelo de Monod-Wyman-Changeux, teniendo en cuenta que la enzima tiene dos subunidades y que, en este caso, el centro catal´ıtico s´olo est´a disponible en la forma activa de la enzima, lo que en t´erminos del modelo MWC se traduce como n= 2 yc= 0. La funci´on φde la enzima ser´a: φ=αe(1 + αe) L1+ [I] KI 1+ [A] KA2 + (1 + αe)2 (152) 36 Entonces, la velocidad v3podr´ıa expresarse de acuerdo a la ecuaci´on cin´etica (130). No obstante, cuando se tratan sistemas enzim´aticos en los que el sustrato es capaz de unirse ´unicamente a la forma activa de la enzima, la curva experimental para la velocidad de reacci´on respecto a la concentraci´on de sustrato suele adquirir una forma sigmoidea. Esta curva se puede ajustar, en la mayor´ıa de los casos, a lo que se denomina la forma de Hill de la ecuaci´on cin´etica [32]. Para encontrar la forma de Hill de la ecuaci´on v3hay que suponer que la enzima es capaz de unirse a sus ligandos (efectores y sustratos) en un solo paso, como se mostr´o en el esquema (71) del cap´ıtulo 3. De este modo, las ´unicas reacciones qu´ımicas que deben considerarse son: el cambio reversible de la enzima de su forma activa R02 a su forma inactiva T02 (153), la captaci´on y liberaci´on de sustrato por la forma activa de la enzima (154) y la formaci´on de producto (155). 2ATP +R02 K2 AMP K2 ATP 2AMP +T02 (153) 2F6P+R02 a2 1 d2 1 R22 (154) R22 k2FBP +R02 (155) Aplicando la ley de acci´on de masas (9) a las reacciones que se acaban de mostrar y suponiendo que las dos primeras se encuentran en equilibrio, puede obtenerse el siguiente sistema de ecuaciones:          K2 AMP [ATP]2[R02] = K2 ATP [AMP]2[T02] K2 R[R22] = [F6P]2[R02] v3=k[R22] (156) Donde KR≡d1 a1se defini´o en (114). 37 Entones, la enzima fosfofructoquinasa (PFK) s´olo puede darse en tres formas diferentes: en su forma activa libre de sustrato R02, en su forma activa colmada de sustrato R22 y en su forma inactiva T02. Por tanto, la concentraci´on total de encima es: [PFK] = [R02]+[R22]+[T02] (157) Despejando [T02]y[R22] de la primera y segunda ecuaci´on del sistema (156) respectivamente y luego sustituyendo estas expresiones en (157), se tiene: [PFK]=[R02] + [F6P]2 K2 R [R02] + κ[ATP]2 [AMP]2[R02] (158) Donde κ≡K2 AMP K2 AT P . Puede despejarse la cantidad [R02] de la ecuaci´on (158). [R02] = [PFK] 1 + κ[ATP ]2 [AMP ]2+[F6P]2 K2 R (159) As´ı, volviendo al sistema de ecuaciones (156) es posible obtener una expresi´on para [R22] de la segunda ecuaci´on y sustituir su valor en la tercera ecuaci´on del sistema para obtener la velocidad v3en funci´on de [R02]. v3=k[F6P]2 K2 R [R02] (160) Por ´ultimo, sustituyendo la f´ormula (159) en (160), se obtiene finalmente la ecuaci´on cin´etica para la reacci´on alost´erica regulada por la fosfofructoquinasa (161). v3=V3max[F6P]2 K2 R1 + κ[ATP ] [AMP ]2+ [F6P]2 (161) Donde la velocidad m´axima V3max =k[PFK] viene dada como la constante cin´etica multiplicada por la cantidad total de enzima. 38 Figura 10: La gr´afica de la izquierda muestra la evoluci´on temporal de las concentraciones normalizadas α(en rojo), γ1(en verde) y γ2(en azul). La gr´afica de la derecha es una representaci´on en el espacio de fases de las concentraciones de las distintas sustancias que componen el sistema: En este caso, la constante cin´etica k3adquiere el valor 2,5s−1. Figura 11: La gr´afica de la izquierda muestra la evoluci´on temporal de las concentraciones normalizadas α(en rojo), γ1(en verde) y γ2(en azul). La gr´afica de la derecha es una representaci´on en el espacio de fases de las concentraciones de las distintas sustancias que componen el sistema: En este caso, la constante cin´etica k3adquiere el valor 3 s−1. Como puede observarse, el sistema es muy sensible al cambio de k3, constante asociada a la demanda de producto P2por parte de otras rutas, mostrando patrones de comportamiento totalmente diferentes. En concreto, para una baja demanda de producto, el sistema tiende a un ciclo l´ımite donde la concentraci´on de reactivos var´ıa con oscilaciones peri´odicas simples (figura 9), mientras que si la demanda de producto es alta, el sistema seguir´a realizando oscilaciones peri´odicas en las que los estados que este alcanza se repiten en el tiempo, pero siguiendo un patr´on mucho m´as complejo que en el caso anterior (figura 11). Sin embargo, existen ciertos valores intermedios de k3para los cuales el sistema muestra un comportamiento ca´otico, pudiendo observarse en algunos casos la formaci´on de un atractor extra˜no (figura 10). 45 En la siguiente imagen se muestran con m´as detalle el ciclo l´ımite y el atractor extra˜no que se obtienen para los valores de k3= 2 s−1yk3= 2,5s−1respectivamente. Figura 12: En la imagen de la izquierda se muestra un ciclo l´ımite en el espacio de fases (ampliaci´on de la figura 9), fruto de realizar la simulaci´on con el valor de la constante cin´etica k3= 2 s−1. En la imagen de la derecha se muestra el atractor ca´otico en el espacio de fases (ampliaci´on de la figura 10), que se obtiene con el valor k3= 2,5s−1. 46 5. Conclusiones El presente trabajo constituye una introducci´on al desarrollo de modelos matem´aticos capaces de explicar la cin´etica enzim´atica de las rutas metab´olicas. Para ello, se ha escogido la parte alta de la ruta glicol´ıtica como hilo conductor de este estudio, analizando la misma con herramientas de la f´ısica pero sin olvidar los aspectos relativos a la bioqu´ımica. En primer lugar, se han expuesto los conceptos te´oricos bioqu´ımicos necesarios para este estudio, explicando con especial inter´es la funci´on y el comportamiento de aquellas enzimas implicadas en la ruta metab´olica que es objeto de an´alisis. Posteriormente, se han deducido las ecuaciones que dictan la cin´etica de distintos tipos de enzimas partiendo como base de la ley de acci´on de masas. De este modo, se han hallado las ecuaciones para reacciones catalizadas por enzimas con un s´olo centro de uni´on, estudiadas por primera vez por Menten y Michaelis, y para las reacciones catalizadas por enzimas que admiten dos sustratos diferentes. Adem´as, el estudio de reacciones enzim´aticas se ha completado comentando la funci´on de Hill como propuesta para explicar los mecanismos alost´ericos de algunas enzimas y, seguidamente, obteniendo una ecuaci´on general para aquellas enzimas alost´ericas que son susceptibles a la acci´on de inhibidores y activadores, de acuerdo con el modelo expuesto por Monod, Wyman y Changeux. Por ´ultimo, se ha construido y simulado un modelo para la parte alta de la ruta glucol´ıtica, en el que se ha hallado, bajo un flujo de suministro de sustrato peri´odico, un comportamiento asociado al caos. Tambi´en se ha presentado el modelo de Decroly y Goldbeter, en el que se pueden observar atractores ca´oticos bajo un flujo de entrada constante de sustrato. 47 Referencias [1] Humberto Maturana y Francisco Varela. De Maquinas y Seres Vivos. Editorial Universitaria, 1 edition, 2003. [2] Richard Johnston and Max E. Valentinuzzi. Metabolism: The physiological powergenerating process: A history of methods to test human beings’ ”vital capacity”. IEEE pulse, 7(3):50–57, 2016. [3] Garabed Eknoyan. Santorio Sanctorius (1561–1636)–founding father of metabolic balance studies. American journal of nephrology, 19(2):226–233, 1999. [4] Charles A. Culotta. Respiration and the lavoisier tradition: theory and modification, 1777-1850. Transactions of the American Philosophical Society, 62(3):3–41, 1972. [5] Keith L. Manchester. Louis Pasteur (1822–1895)—chance and the prepared mind. Trends in biotechnology, 13(12):511–515, 1995. [6] E. Kinne-Saffran and R. K. H. Kinne. Vitalism and synthesis of urea. American journal of nephrology, 19(2):290–294, 1999. [7] Edward Buchner. Alcoholic fermentation without yeast cells. Cornish-Bowden (1997), pages 25–31, 1997. [8] Donald Voet and Judith G. Voet. Biochemistry. Wiley, 4th edition, 2010. [9] Byung Hong Kim and Geoffrey Michael Gadd. Bacterial Physiology and Metabolism. Cambridge University Press, 1 edition, 2008. [10] Paul N. Lithaw. Glycolysis: Regulation, Processes and Diseases. 2009. [11] P. Farreras Valent´ı y Ciril Rozman. Medicina Interna. Metabolismo y Nutrici´on Endocrinolog´ıa. Elsevier, 17 edition, 2014. [12] Nicolas Glansdorff, Ying Xu and Bernard Labedan. The last universal common ancestor: emergence, constitution and genetic legacy of an elusive forerunner. Biology direct, 3(1):29, 2008. [13] Alejandro G. Marangoni. Enzyme Kinetics: A Modern Approach. Wiley-Interscience, 1 edition, 2002. 48 [14] David L. Nelson, Michael M. Cox y Albert L. Lehninger. Lehninger: Principios de Bioqu´ımica. Omega, 5a edition, 2009. [15] Jeremy M. Berg, John L. Tymoczko y Lubert Stryer. Bioqu´ımica. Revert´e, 6a edition, 2008. [16] Edda Klipp, Axel Kowald, Wolfram Liebermeister and Christoph Wierling. Systems biology : a textbook. Wiley-VCH, 2nd edition, 2016. [17] Edda Klipp, Ralf Herwig, Axel Kowald, Christoph Wierling and Hans Lehrach. Systems Biology in Practice: Concepts, Implementation and Application. Wiley-VCH, 1 edition, 2005. [18] E. W. Lund. Guldberg and Waage and the Law of Mass Action. Journal of Chemical Education, 42(10):548, 1965. [19] Franklin C. McLean. Application of the law of chemical equilibrium (law of mass action) to biological problems. Physiological Reviews, 18(4):495–523, 1938. [20] Reinhart Heinrich and Stefan Schuster. The Regulation of Cellular Systems. Springer US, 1 edition, 1996. [21] S. C. Kou, Binny J. Cherayil, Wei Min, Brian P. English and X. Sunney Xie. Singlemolecule Michaelis-Menten equations. 2005. [22] Hans Lineweavera and Dean Burk. The determination of enzyme dissociation constants. Journal of the American chemical society, 56(3):658–666, 1934. [23] E. A. Newsholme, F. S. Rolleston and Kay Taylor. Factors affecting the glucose 6-phosphate inhibition of hexokinase from cerebral cortex tissue of the guinea pig. Biochemical Journal, 106(1):193–201, 1968. [24] Ross L. Stein. Kinetics of Enzyme Action: Essential Principles for Drug Hunters. Wiley, 1 edition, 2011. [25] Rudolf Gesztelyi, Judit Zsuga, Adam Kemeny-Beke, Balazs Varga, Bela Juhasz and Arpad Tosaki. The hill equation and the origin of quantitative pharmacology. Archive for history of exact sciences, 66(4):427–438, 2012. 49 [26] Archibald Vivian Hill. The combinations of haemoglobin with oxygen and carbon monoxide, and the effects of acid and carbon dioxide. Biochemical Journal, 15(5):577, 1921. [27] Edda Klipp, Wolfram Liebermeister, Christoph Wierling, Axel Kowald, Hans Lehrach and Ralf Herwig. Systems Biology. Wiley-Blackwell, 1 edition, 2009. [28] Jacque Monod, Jeffries Wyman and Jean-Pierre Changeux. On the nature of allosteric transitions: A plausible model. J. Mol. Biol, 12:88–118, 1965. [29] Albert Goldbeter. Biochemical Oscillations and Cellular Rhythms: The Molecular Bases of Periodic and Chaotic Behaviour. Cambridge University Press, 1997. [30] Mario Markus and Benno Hess. Control of metabolic oscillations: unpredictability, critical slowing down, optimal stability and hysteresis. In Control of Metabolic Processes, pages 303–313. Springer, 1990. [31] Mario Markus, Dietrich Kuschmitz and Benno Hess. Chaotic dynamics in yeast glycolysis under periodic substrate input flux. FEBS letters, 172(2):235–238, 1984. [32] Yves Termonia and John Ross. Oscillations and control features in glycolysis: numerical analysis of a comprehensive model. proceedings of the national academy of sciences, 78(5):2952–2956, 1981. [33] Finn Hynne, Sune Danø and Preben G. Sørensen. Full-scale model of glycolysis in saccharomyces cerevisiae. Biophysical chemistry, 94(1-2):121–163, 2001. [34] Francisco Montero y Federico Mor´an. Biof´ısica: Procesos de autoorganizaci´on en Biolog´ıa. Eudema, 1992. [35] Javier Buceta Fern´andez, Elka Radoslavova Koroutcheva y Juan Manuel Pastor Ruiz. Temas de biof´ısica. Editorial UNED, 2012. 50 Anexos: A. Reacciones con inhibici´on En esta secci´on se ampl´ıa el estudio para las enzimas analizando tres tipos diferentes de procesos de inhibici´on enzim´atica: la inhibici´on competitiva, la acompetitiva y la no competitiva. A.1. Inhibici´on competitiva Un inhibidor competitivo es aquel que compite con el sustrato por ocupar el centro catal´ıtico de la enzima. Si la enzima se liga al inhibidor, esta ya no puede ligarse al sustrato para dar lugar al producto hasta que quede liberada. Un ejemplo de este tipo de inhibici´on se encuentra en la enzima asparaginasa, la cual reacciona con la glutamina para producir glutamato y con el ´acido asparag´anico produciendo aspartato [35]. E+Sk1 k−1 ES k2 k−2 E+P + I k3 k−3 EI (170) En este caso las velocidades de cada etapa son: v1=k1[E][S]−k−1[ES] (171) v=v2=k2[ES]−k−2[E][P] (172) v3=k3[E][I]−k−3[EI] (173) De acuerdo con la f´ormula (27) la concentraci´on de complejo enzima-sustrato crece con velocidad v1y disminuye con velocidad v2. 51 d[ES] dt =v1−v2=k1[E][S]−k−1[ES]−k2[ES] + k−2[E][P] (174) De nuevo, se puede aplicar la condici´on de estado cuasiestacionario (18). (k1[S] + k−2[P])[E] = (k−1+k2)[ES] (175) Del mismo modo, puede asumirse que la reacci´on de formaci´on del complejo enzimainhibidor alcanza tambi´en un estado cuasiestacionario. d[EI] dt = 0 (176) Lo que equivale a decir que la velocidad (173) es nula y adem´as: k3[E][I] = k−3[EI] (177) Para reducir el n´umero de constantes, se define la constante de disociaci´on del inhibidor KI3como: KI3=k−3 k3 (178) Haciendo uso de la definici´on anterior y de la ecuaci´on (177), la expresi´on para la concentraci´on de enzima inhibida queda de la siguiente forma: [EI] = [E][I] KI3 (179) En este caso, la cantidad total de enzima ser´a la suma de la enzima libre m´as la que est´a formando parte de un complejo enzim´atico, ya sea unida a un sustrato o a un inhibidor. Despejando la cantidad de enzima libre se tiene: [E]=[Etotal]−[ES]−[EI] (180) Donde la cantidad [EI] viene dada por (179) [E]=[Etotal]−[ES]−[E][I] KI3 (181) As´ı, despejando [E]: 52 [E] = [Etotal]−[ES] 1 + [I] KI3 (182) Ahora es posible sustituir los t´erminos de concentraci´on de enzima libre tanto en la ecuaci´on de velocidad (172) como en la ecuaci´on (175) sustituyendo [E] por su expresi´on obtenida en (182). v=k2[ES]−k−2 1 + [I] KI3 [Etotal][P] + k−2 1 + [I] KI3 [ES][P] (183) k1[S][Etotal] 1 + [I] KI3 +k−2[P][Etotal] 1 + [I] KI3 −k−2[P][ES] 1 + [I] KI3 −k1[S][ES] 1 + [I] KI3 = (k−1+k2)[ES] (184) De la ´ultima f´ormula puede despejarse [ES]: [ES] = k1[S][E]total +k−2[P][E]total (k−1+k2)1 + [I] KI3+k1[S] + k−2[P] (185) Sustituyendo este valor en la ecuaci´on (183) se obtiene finalmente la ecuaci´on cin´etica. v= k2k1[S][Etotal] + k2k−2[P][Etotal] + k−2k1 1+ [I] KI3 [P][S][Etotal] + k2 −2 1+ [I] KI3 [P]2[Etotal] (k−1+k2)1 + [I] KI3+k1[S] + k−2[P] − − k−2k−1[P][Etotal] + k2k−2[P][Etotal] + k−2k1 1+ [I] KI3 [P][S][Etotal] + k2 −2 1+ [I] KI3 [P]2[Etotal] (k−1+k2)1 + [I] KI3+k1[S] + k−2[P] (186) Esta ecuaci´on se puede simplificar hasta obtener lo siguiente: v= [Etotal]k2k1[S]−k−2k−1[P] (k−1+k2)1 + [I] KI3+k1[S] + k−2[P] (187) Pero para hallar la forma m´as simple de esta ecuaci´on ser´a necesario, en primer lugar, dividir el numerador y el denominador por la cantidad k−1+k2y luego, aplicar las definiciones KM(26), K−M(38), Vmax yVmin introducidas en el cap´ıtulo 3. v= Vmax KM[S]−Vmin K−M[P] 1 + [I] KI3+[S] KM+[P] K−M (188) Como caso particular de esta ecuaci´on, cuando la reacci´on 2 del esquema (170) es irreversible, el valor de k−2es nulo y por tanto K−M, que es la ´unica constante que 53 depende de la anterior, tender´a al infinito. De esta forma se obtiene la ecuaci´on cin´etica para la inhibici´on competitiva irreversible. v=Vmax[S] KM1 + [I] KI3+ [S] (189) Figura 13: La gr´afica de la izquierda muestra la velocidad de reacci´on vfrente a la concentraci´on de sustrato [S] para distintas concentraciones de sustancia inhibidora. Las l´ıneas m´as azules (m´as hacia abajo) indican una mayor concentraci´on de inhibidor respecto a la constante de inhibici´on KI3. En concreto, los valores de [I]/KI3usados son 0, 1, 1,5, 2 y 2,5. La gr´afica de la derecha muestra el diagrama de Lineweaver-Burk correspondiente. Puede observarse que cuando la cantidad de inhibidor es mayor, el sistema requiere m´as sustrato para alcanzar la misma velocidad de reacci´on. Adem´as, puede verse c´omo peque˜nas variaciones en la cantidad de inhibidor pueden suponer cambios relevantes en el comportamiento de la parte lineal del sistema (aquella que se da para baja concentraci´on de sustrato). En esta representaci´on, la constante KMtoma el valor 0,2M. A.2. Inhibici´on acompetitiva En ocasiones, cuando una enzima se une a un sustrato, se producen en esta unos cambios conformacionales tales que hacen al complejo enzima-sustrato susceptible a la acci´on de ciertos inhibidores. Estos casos son conocidos como inhibiciones acompetitivas. E+Sk1 k−1 ES k2 k−2 E+P + I k4 k−4 ESI (190) Las velocidades para cada una de estas etapas son: 54 Figura 15: La gr´afica de la izquierda muestra la velocidad de reacci´on vfrente a la concentraci´on de sustrato [S] para distintas concentraciones de sustancia inhibidora. Las l´ıneas m´as azules (m´as hacia abajo) indican mayor concentraci´on de inhibidor respecto a la constante de inhibici´on. En concreto se ha recurrido a los valores 0, 1, 1,5,2y2,5 para [I]/KI. La gr´afica de la derecha muestra el diagrama de Lineweaver-Burk correspondiente. En este caso se observa tanto la influencia de la cantidad de inhibidor en la velocidad a niveles bajos de concentraci´on de sustrato presente en la inhibici´on competitiva, como la disminuci´on de la velocidad m´axima en condiciones de saturaci´on de sustrato que se observaba en la inhibici´on acompetitiva. KM= 0,2M 61 B. C´odigos de MatLab utilizados para realizar algunas simulaciones En esta secci´on se muestran los c´odigos para MatLab con los que se han obtenido los atractores extra˜nos de las figuras 8 y 12. Simulaci´on para la parte alta de la ruta glicol´ıtica: 1%% Parte alta de la gluc´olisis . 2% Modelo cin´etico para la parte alta de la gluc´olisis . 3clear all; 4close all; 5clc ; 6 7%% Par´ametros : 8%Constantes : 9 10 %Concentraci´on de glucosa en el medio : 11 Glu0 =6.4087; %mM . Amplitud de se~nal y m´ınima concentraci´on de glucosa . 12 w=0.05; %Frecuencia de oscilaci´on de la entrada de glucosa en el sistema. 13 14 %Hexoquinasa: 15 V1max =0.8621; %mM /s. Velocidad m´axima de la reacci´on . 16 KATP =0.37; %mM 17 KGlu =0.10; %mM 18 KII =0.037; %mM 19 20 %Fosfoglucoisomerasa: 21 V2max =2.338; %mM/s Velocidad m´axima a la que se puede generar producto . 22 V2min =2.338; %mM/s Velocidad m´axima a la que se puede generar sustrato . 23 KG6P =0.8; %mM 24 KF6P =0.15; %mM 25 26 %Fosfofructoquinasa : 27 V3max =0.74548; %mM /s. Velocidad m´axima de la reacci´on . 28 K_R =0.145; %mM. 29 k=0.15; 30 62 31 %Otras reacciones que tienen lugar en el entorno : 32 k4 =0.1008; %1/s 33 k5 =0.0377; %1/s*mM 34 %k6 =1.1413; %1/s 35 k6 =0.8828; %1/s 36 k7 =0.0535; %1/s 37 k8=7.215; %1/s*mM 38 k_8 =2.222; %1/s*mM 39 40 %Variables iniciales : 41 G6P0 =1; %mM . %Cantidad inicial de glucosa -6 - fosfato . 42 F6P0 =0; %mM. %Cantidad inicial de fructosa -6- fosfato . 43 FBP0 =0; %mM . %Cantidad inicial de fructosa -1 ,6- bisfosfato . 44 ATP0 =2.1; %mM. %Cantidad inicial de ATP. 45 ADP0 =1.4; %mM. %Cantidad inicial de ADP. 46 % AMP0 =0.1; %mM . %Cantidad inicial de AMP. 47 AMP0=0.09994199; %mM. %Cantidad inicial de AMP . 48 49 R0=[G6P0,F6P0,FBP0,ATP0,ADP0,AMP0]; %Vector de reactivos iniciales . 50 51 %Intervalo temporal : 52 tspan =[0 ,1500]; %s 53 54 %% Soluci´on del sistema : 55 options = odeset ( ’AbsTol’,1e -12 , ’RelTol’,1e -12) ; 56 Der =@(t,R) odeGluco (t,R,Glu0 ,w,V1max ,KATP ,KGlu ,KII ,V2max ,V2min ,KG6P , KF6P ,V3max ,K_R ,k,k4 ,k5 ,k6 ,k7 ,k8 ,k_8 ); 57 [t,R]= ode45 ( Der ,tspan ,R0 , options ); 58 G6P =R (: ,1); %Glucosa -6 - fosfato . 59 F6P =R (: ,2); %Fructosa -6- fosfato . 60 FBP =R (: ,3); %Fructosa -1 ,6 - bisfosfato . 61 ATP =R (: ,4); 62 ADP =R (: ,5); 63 AMP =R (: ,6); 64 65 %% Representaci´on temporal de la concentraci´on de glucosa . 66 figure(’Name ’,’Glucosa’); 67 Glucosa =Glu0 - Glu0* sin (w*t); 68 plot(t, Glucosa ,’ LineWidth ’ ,2,’Color ’ ,[1 ,0.5 ,0] ’) ; 63 69 70 %% Representaci´on temporal de las concentraciones : 71 figure(’Name ’,’Representaci´on temporal’); 72 plot(t,ATP ,’r’,’LineWidth ’ ,2); 73 axis([0 1500 0 3.5]) ; 74 title (’Representaci´on temporal .’); 75 xlabel(’ Tiempo (s)’); 76 ylabel(’ Concentraci´on (mM)’); 77 78 %% Representaci´on en el espacio de fases : 79 figure(’Name ’,’Espacio de fases ’); 80 plot( Glucosa , ATP ); %Concentraci´on de sustrato frente a concentraci´on de producto. 81 axis([0 14 0 3.5]) ; 82 title (’Espacio de fases .’); 83 xlabel(’[Glu ] (mM)’); 84 ylabel(’[ATP ] (mM)’); 85 86 %% Representaci´on en el espacio de fases completo : 87 figure(’Name ’,’Espacio de fases completo ’); 88 plot3 (Glucosa ,ATP , ADP); %Concentraci´on de sustrato frente a concentraci´on de producto . 89 axis([0 15 0 4 0 2]) ; 90 title (’Espacio de fases completo .’); 91 xlabel(’[Glu ] (mM)’); 92 ylabel(’[ATP ] (mM)’); 93 zlabel(’[ADP ] (mM)’); Donde la funci´on odeGluco es: 1function [ dRdt ]= odeGluco (t,R,Glu0 ,w,V1max ,KATP ,KGlu ,KII ,V2max , V2min , KG6P ,KF6P ,V3max ,K_R ,k,k4 ,k5 ,k6 ,k7 ,k8 ,k_8 ) 2% Funci´on para el script que modela la cin´etica de la parte alta de la gluc´olisis . 3 4%Variables del sistema : 5G6P =R (1); 6F6P =R (2); 7FBP =R (3); 8ATP =R (4); 64 9ADP =R (5); 10 AMP =R (6); 11 12 %Entrada de Glucosa: 13 Glu =Glu0 - Glu0 *sin (w*t); 14 15 %Velocidades de reacci´on : 16 v1 = V1max * ATP *Glu /( KII + KGlu * Glu + KATP * ATP + Glu * ATP ); %Hexoquinasa. Reacci´on bisustrato con inhibici´on no competitiva . 17 v2 =( V2max * G6P / KG6P - V2min * F6P / KF6P ) /(1+ G6P / KG6P + F6P / KF6P ); % Fosfoglucoisomerasa . Reacci´on tipo Michaelis - Menten reversible . 18 v3 = V3max * F6P ^2/( K_R ^2*(1+ k *( ATP / AMP ) ^2) + F6P ^2) ; %Fosfofructoquinasa . 19 v4=k4* FBP; %P´erdida de fructosa -1 ,6 - bisfosfato . 20 v5=k5* G6P* ATP; %Uso de glucosa -6 - phosfato en otras rutas . 21 v6=k6* ADP; %Generaci´on de ATP a partir de ADP . 22 v7=k7* ATP; %Generaci´on de ADP a partir de ATP . 23 v8=k8* ATP *AMP - k_8 *ADP ^2; %Adenilato quinasa . 24 25 %Sistema de ecuaciones diferenciales : 26 dG6Pdt =v1 -v2 -v5; 27 dF6Pdt =v2 - v3; 28 dFBPdt =v3 - v4; 29 dATPdt =-v1 -v3 -v5+v6 -v7 -v8; 30 dADPdt =v1+ v3+v5 -v6+ v7 +2* v8; 31 dAMPdt =- v8; 32 33 %Incremento de las concentraciones de reactivo : 34 dRdt=[dG6Pdt;dF6Pdt;dFBPdt;dATPdt;dADPdt;dAMPdt]; 35 end Simulaci´on para el modelo de Decroly y Goldbeter: 1%% Oscilador de Decroly y Goldbeter . 2% Modelo cin´etico de Decroly y Goldbeter de acuerdo a la cin´etica de MWC . 3clear all; 4close all; 5clc ; 6 7%% Par´ametros iniciales : 65 8v0 =1; %Velocidad de entrada de sustrato (->S). 9KR1 =2; %Constante de acci´on de Michaelis -Menten de la ecuaci´on 1 (S->P) . 10 sigma1 =10; 11 sigma2 =10; 12 L1 =5*10^8; %Constante alost´etica . 13 L2 =100; 14 q1 =50; 15 q2=0.02; 16 d=0.0001; 17 %k3 =2; %1/s. 18 k3 =2.5; %1/s. 19 %k3 =3; %1/s. 20 tspan =[0 ,2000]; %Intervalo de definici´on . 21 S0 =0; %Concentraci´on inicial de sustrato . 22 P10 =0; %Concentraci´on inicial de producto 1. 23 P20 =0; %Concentraci´on inicial de producto 2. 24 R0 =[ S0 ,P10 , P20 ]; %Vector de reactivos iniciales . 25 26 %% Soluci´on del sistema : 27 options = odeset ( ’AbsTol’,1e-5,’RelTol’,1e -5) ; 28 Der =@(t,R) odeOsciladorDG (t,R,v0 ,KR1 , sigma1 ,sigma2 ,L1 ,L2 ,q1 ,q2 ,d,k3); 29 [t,R]= ode45 ( Der ,tspan ,R0 , options ); 30 alpha =R (: ,1) ; 31 gamma1 =R (: ,2) ; 32 gamma2 =R (: ,3) ; 33 34 %% Representaci´on temporal de las concentraciones : 35 figure(’Name ’,’Representaci´on temporal’); 36 plot(t,alpha ,’r’); %Concentraci´on de sustrato frente al tiempo . 37 title (’Representaci´on temporal .’); 38 xlabel(’Tiempo.’); 39 ylabel(’ Concentraci´on .’); 40 hold on; 41 plot(t,gamma1 ,’g’); %Concentraci´on de producto 1 frente al tiempo . 42 plot(t,gamma2 ,’b’); %Concentraci´on de producto 2 frente al tiempo . 43 hold off ; 44 45 %% Representaci´on en el espacio de fases: 66 46 figure(’Name ’,’Espacio de fases ’); 47 plot( gamma1 , alpha ); %Concentraci´on de sustrato frente a concentraci´on de producto . 48 title (’Espacio de fases . ’); 49 xlabel(’\ gamma_1 ’); 50 ylabel(’\ alpha ’); 51 52 %% Representaci´on en el espacio de fases completo : 53 figure(’Name ’,’Espacio de fases completo ’); 54 plot3 (gamma1 , alpha , gamma2 ); %Concentraci´on de sustrato frente a concentraci´on de producto . 55 title (’Espacio de fases completo .’); 56 xlabel(’\ gamma_1 ’); 57 ylabel(’\ alpha ’); 58 zlabel(’\ gamma_2 ’); Donde la funci´on odeOsciladorDG es: 1function [ dRdt ]= odeOsciladorDG (t,R,v0 ,KR1 , sigma1 , sigma2 ,L1 ,L2 ,q1 ,q2 ,d, k3) 2% Funci´on para el script que modela la cin´etica del modelo de Decroly y Goldbeter . 3 4%Variables del sistema : 5alpha =R (1) ; %Concentraci´on de sustrato . 6gamma1 =R (2) ; %Concentraci´on de producto 1. 7gamma2 =R (3) ; %Concentraci´on de producto 1. 8 9%Sistema de ecuaciones diferenciales : 10 dSdt =v0 /KR1 - sigma1 * alpha *(1+ alpha ) *(1+ gamma1 ) ^2/( L1 +(1+ alpha ) ^2*(1+ gamma1 ) ^2) ; %Producci´on de sustrato normalizada . 11 dP1dt =q1 * sigma1 * alpha *(1+ alpha ) *(1+ gamma1 ) ^2/( L1 +(1+ alpha ) ^2*(1+ gamma1 ) ^2) - sigma2 * gamma1 *(1+ d* gamma1 ) *(1+ gamma2 ) ^2/( L2 +(1+ d* gamma1 ) ^2*(1+ gamma2 ) ^2) ; %Producci´on de producto 1 normalizada . 12 dP2dt =q2 * sigma2 * gamma1 *(1+ d* gamma1 ) *(1+ gamma2 ) ^2/( L2 +(1+ d* gamma1 ) ^2*(1+ gamma2 ) ^2) -k3 * gamma2 ; %Producci´on de producto 2 normalizada . 13 14 %Incremento de las concentraciones de reactivo : 15 dRdt =[ dSdt ; dP1dt ; dP2dt ]; 16 end 67