Predicción de alfa-hélices en proteínas de

Transmembrana:
Análisis y propuesta de mejoramiento
del método HMMTOP con Algoritmos Genéticos

Ricardo Antonio Mendoza León

Asesor

Rafael A. García

Departamento de Ingeniería de Sistemas y Computación

Universidad de los Andes

Junio 2006
Introducción

Desde el descubrimiento histórico de la estructura de doble hélice del ADN

por los profesores James Watson y Francis Crick en 1953, la comunidad
científica inició una imparable carrera por el entendimiento de la vida y su
funcionamiento, la cual ha traído consigo grandes esfuerzos, retos, triunfos y
sobre todo, conocimiento. Pero también una creciente necesidad en el
desarrollo de herramientas que apoyen esta carrera, dando soporte a los
nuevos retos que se presenten en el futuro. Un caso en donde fue clara esta
necesidad, se presentó durante la realización del PGH o proyecto del
genoma humano, sin duda una de las más grandes empresas de
investigación en la historia, requiriendo la colaboración de un gigantesco
número de investigadores, recursos, tiempo y tecnología. Es evidente que la
realización de esta empresa no hubiese sido posible sin el desarrollo
tecnológico de la fecha, en la que las herramientas informáticas prestaron un
rol indispensable debido a la gigantesca cantidad de información que debía
ser analizada, categorízada y almacenada.

Lo expresado es un ejemplo de un nuevo reto. Cada día la cantidad de

secuencias proteicas identificadas crece a ritmos acelerados; sin embargo el
análisis detallado de su estructura y función no sigue el mismo ritmo, debido
a los altos costos de recursos requeridos para obtener un modelo
experimental de alta resolución para la estructura de estas proteínas,
información importante, para comprender a fondo su función. Debido a lo
anterior y a la importancia fundamental que representan las proteínas en la
mecánica de la vida, la utilización de herramientas computacionales
derivados de la disciplina de la bioinformática, como parte del proceso de
investigación, ha tomado cada día mas fuerza e importancia como parte de la
solución de estos problemas.

El tema de este trabajo de grado se concentra inicialmente en proveer de un

recuento de los métodos computacionales de la bioinformática, empleados
en la predicción de una clase muy importante de estructuras presentes en las
proteínas de trans-membrana, conocidas como Alfa-hélices y presentar una
propuesta de mejoramiento al método de predicción HMMTOP[1].

2
Definición del problema

Dado que los métodos a los cuales hace referencia este trabajo buscan
solucionar esencialmente el mismo problema y adicionalmente que también
se plantea una propuesta de mejoramiento del método de predicción
HMMTOP, el problema base que compete a este trabajo es esencialmente el
mismo, definido de la siguiente manera:

Dada una cadena de aminoácidos (residuos) correspondiente a una proteína

de trans-membrana, predecir computacionalmente con la mayor eficacia
(precisión y completitud) posible, las subsecuencias en la cadena tales que
presenten una topología estructural de tipo alfa-hélice en la proteína.

Objetivos

• Investigar los métodos disponibles (estado del arte) cuyo propósito es

solucionar el problema de la predicción de alfa-hélices.
• Realizar un análisis detallado al método HMMTOP.
• Proponer una modificación al método HMMTOP, de manera que el
nuevo método cuente con una eficacia comparable con respecto al
método analizado.
• Implementar el método modificado y evaluar la eficacia de este
mediante la utilización de un conjunto de prueba con topología
definida.
• Proveer una herramienta asequible, que permita su ajuste a casos
particulares del problema.
• Proveer de una referencia útil para futuras investigaciones no solo
relacionadas con el problema específico de este trabajo, sino a la
bioinformática y las ciencias de la computación.

3
1 Estado del Arte

Existe una gran diversidad de aproximaciones al problema de la

identificación de las Alfa-hélices. Cada método se fundamenta en
propiedades y características intrínsecas de las Alfa-hélices o de los
aminoácidos, como base para el diseño de la función de calificación. Las
características comúnmente usadas suelen ser químicas, físicas,
termodinámicas o estadísticas. Como segunda clasificación de los métodos
tenemos la técnica computacional utilizada para sacar provecho de estas
propiedades o características. A continuación se presenta una selección de
algunas técnicas desarrolladas en años anteriores y que para efectos de este
trabajo de grado tienen las siguientes características: a) Constituyen un
conjunto extenso de los métodos computacionales utilizados. b) Tienen una
importancia en la evolución histórica de los métodos desarrollados. c) En el
momento que fueron desarrollados hacían parte de los de mayor eficacia.

KD 5, 9 (Kyte y Doolittle, 1982), Eisenberg (Eisenberg, 1982), ALOM2 (Nakai

y Kanehisa, 1992), TopPred II (Claros y Heijne, 1994), DAS (Cserzo, 1997) y
Gonnet (Gonnet 2004) [2], son métodos que utilizan las propiedades locales
de la secuencia de aminoácidos, utilizando una aproximación tipo ventana
deslizante, en particular KD (5 y 9) y Eisenberg hacen uso de escalas
simples de hidrofobicidad para los aminoácidos. TopPred II aplica la ‘regla
positiva interior’ para evaluar la valides de modelos topológicos derivados de
análisis hidropáticos, para esto utiliza varias matrices de preferencia para
incrementar su efectividad. Por su parte, DAS está basado en características
obtenidas a partir de comparar la secuencia contra una colección de
proteínas de membrana no homólogas, utilizando una matriz de puntajes
especial derivada previamente. La aproximación de Gonnet se basa en
índices de hidrofobicidad de los aminoácidos y la interpretación geométrica
de estos con respecto a su localización en la cadena de aminoácidos. Este
método cuenta con gran eficacia, sin embargo aún no se tiene un buen
conocimiento en lo referente a cómo calcular estos índices de hidrofobicidad,
motivo por el cual, el algoritmo tiene poco uso práctico.

TMHMM 1.0, 2.0 (Sonnhammer, 1998), 2.0 re-entrenado (Krogh, 2001),

HMMTOP 2 (Tusnády y Simon, 1998), son aproximaciones globales, que
buscan determinar la topología estadísticamente más probable para toda la
proteína. Para esto utilizan modelos ocultos de Markov (HMM) que
representan varias restricciones en la proteína de membrana dadas por la bi-

4
capa lípida; históricamente estos han sido los métodos más efectivos a la
fecha. PHD (Rost, 1996), utiliza una red neuronal para realizar el análisis,
ésta es considerada una aproximación local.

MEMSAT 1.5 (Jones, 1994), TMAP (Persson y Argos, 1997), TMpred

(Hoffmann y Stoffel, 1993), SOSUI (Hirokawa, 1998), utilizan métodos
combinados, donde los resultados de nivel local son evaluados por
heurísticas globales como la “regla positiva interior” o otras diferencias en la
distribución de los aminoácidos. TMpred, usa preferencias estadísticas para
predecir Alpha-hélices, estas preferencias fueron obtenidas de una base de
datos de proteínas de membrana compilada por un experto.

BPROMPT (2003) [3], utiliza una red de creencia bayesiana para desarrollar
un modelo de consenso entre el siguiente conjunto de métodos de
predicción: TMpred, TopPred II, SOSUI, DAS, TMHMM y HMMTOP [1].
SOPMA (Georujon y Déleage, 1995) es otra aproximación basada en el
consenso de los métodos: GOR (Garnier, 1996) [4], Homologia de Levin
(Levin, 1986), El método de doble predicción (Déleage y Roux, 1987), PHD y
el método CNRS.

5
2. Marco Teórico

2.1 Biología Molecular

La biología molecular es la ciencia que estudia la estructura y función de las

proteínas, el acido nucleico y la interacción de estos, con el propósito de
entender los sistemas biológicos.

2.1.1 Información genética, Genes y el ADN

La información genética hace referencia a la manera como los organismos

vivos codifican la información que los diferencia de otros; esto es, el conjunto
de información que define la estructura, características y dinámica biológica
del organismo. Esta información se encuentra almacenada en estructuras
conocidas como Genes, los cuales son la menor estructura de información
independiente en un organismo, en otras palabras, codifican completamente
una función o característica simple de dicho organismo. Los genes son una
cadena de polínucleotidos en forma de doble hélice conocido como ácido
desoxirribonucleico o ADN (Figura 2.2.1). La composición de esta cadena
cuenta con un par de bases de acido nucleico (nucleótidos) unidos por
enlaces de hidrógeno; cada uno de estos pares se une al siguiente con una
columna de fosfato de carbono. Los nucleótidos cuentan con tres grupos, el
primero un fosfato que forma la columna que une a los pares de nucleótidos,
el segundo, un grupo desoxirribosa y el tercero, una base nitrogenada con la
cual se unen los pares de nucleótidos. En el ADN existen cuatro tipos de
nucleótidos, dos purínicas denominadas adenina y guanina y dos
pirimidínicas denominadas citosina y timina. En el ADN solo existen dos
combinaciones de pares de bases, estas son adenina-tiamina y citosina-
guanina. El aspecto importante sobre lo anterior, radica en que de esta
manera se encuentra codificada la información en el ADN; por este motivo
desde el punto de la información genética, lo relevante es la secuencia de los
nucleótidos en la cadena de ADN.

6
 Molécula de ADN y su estructura
(Figura 2.1.1)

2.1.2 Transcripción del ADN, Codones y ARN Polimerasa

El proceso de transcripción del ADN consiste en la separación, obtención y

recodificación de la información contenida en una sección del ADN. Este
proceso es realizado por el ARN polimerasa, el cual es una macromolécula
compuesta de proteínas enzimáticas, que inicialmente separa las uniones de
hidrógeno que unen a los pares de nucleótidos, de esta manera, separando
la doble hélice del ADN. Luego, genera una copia de la sección del ADN
correspondiente al gen ARN mensajero. Finalmente, se genera una cadena
de aminoácidos con esta información y se une de nuevo la doble hélice. La
información que codifica al aminoácido está organizada en grupos de tres
nucleótidos conocidos como Codones. Salvo raras excepciones, los codones
representan un total de 20 códigos, cada uno de los cuales representa un
aminoácido especifico. Es importante tener en cuenta que la información
presente en el ADN no solo consiste en codones, una parte importante de
ésta provee información al ARN Polimerasa de donde iniciar y terminar la
trascripción. Estas señales de inicio y final de trascripción son secuencias de
longitud variable denominadas promotores.

2.1.3 Proteínas

Son los componentes básicos (o como una analogía, ladrillos) de todo

sistema biológico vivo u organismo, estos le dan su estructura y forma física;
además, controlan los procesos físicos y químicos necesarios para captar,

7
almacenar y utilizar energía. Existen varios grupos de proteínas identificados
por su rol en la mecánica de la vida; entre los grupos más importantes, se
encuentran las proteínas estructurales y las enzimáticas.

Las proteínas estructurales actúan como los ladrillos con los cuales son
construidas las estructuras de los organismos como los son los tejidos,
membranas, etc. Las enzimas son proteínas que actúan como catalizadores
en las reacciones bioquímicas de los organismos vivos; los catalizadores son
substancias que aceleran las reacciones bioquímicas. Debido a esto, son
fundamentales, ya que algunas de las reacciones bioquímicas necesarias
para la vida no serian viables sin su acción. Es importante tener en cuenta
que las enzimas son muy especificas, por lo general cada una es útil
únicamente para una reacción particular, por lo cual hay una gran diversidad
de estas. También existen otros tipos de proteínas como las responsables de
transporte de oxígeno y anticuerpos de defensa entre otras.

2.1.4 Estructura de las proteínas

Las proteínas son cadenas de moléculas mas simples llamadas

aminoácidos. Cada aminoácido tiene un átomo central de carbono, conocido
como el carbón alfa, unido a un grupo amino (NH2), un grupo Carboxy
(COOH), y una cadena lateral. La cadena lateral es la que distingue un
aminoácido de otro. Se pueden encontrar cadenas laterales simples como un
átomo de hidrogeno (glycina) o complejas como dos anillos de carbono en el
caso de el tryptophan. En la naturaleza existen al menos 20 aminoácidos
comunes en la composición de las proteínas (Tabla 2.1), estos se dividen en
Hidrófobos (no solubles en agua) y Polares. Adicionalmente a estos, existen
otros aminoácidos no estándar, estos se encuentran principalmente en los
virus.

(Tabla 2.1)
Los 20 Aminoácidos principales

Letra Nombre Tipo

código

A Alanine Hidrófobo

C Cysteine Polar

D Aspartic Acid Acido

8
 E Glutamic Acid Acido

F Phenylalanime Hidrófobo

G Glycine Polar

H Histidina Básico

I Isoleucine Hidrófobo

K Lisien Básico

L Leucine Hidrófobo

M Methionine Hidrófobo

N Asparagine Polar

P Proline Hidrófobo

Q Glutamina Polar

R Arginine Básico

S Serine Polar

T Threonine Polar

V Valine Hidrófobo

W Tryptophan Hidrófobo

Y Tyrosine Polar

En una proteína los aminoácidos están unidos mediante uniones péptidas,

por lo cual también se conocen a las proteínas como cadenas polipéptidas.
En las uniones péptidas, el átomo de carbono perteneciente al grupo carboxy
(COOH) del aminoácido Ai, éste se une con el átomo de nitrógeno del grupo
amino del aminoácido Ai+1, liberando una molécula de agua; por esta razón,
se conoce como residuo a un aminoácido cuando hace parte de una
proteína. En general las proteínas tienen un tamaño entre 100 a 5,000
residuos.

9
La estructura de una proteína es una compleja configuración tridimensional.
En general, existen 4 etapas de la estructura de las proteínas: la primera o
estructura primaria, es la cadena de aminoácidos en sí, pero ésta estructura
primaria no es real, ya que la cadena de aminoácidos se dobla y solda
consigo misma dependiendo de los aminoácidos que la conforman. A esta
forma generada, se le conoce como la estructura secundaria y corresponde
al doblamiento inicial de esta, la que se encuentra fuertemente ligada a su
forma en tres dimensiones. La estructura terciaria, es la estructura que toma
la proteína luego de alcanzar el equilibrio termodinámico. Por ultimo, está la
estructura terciaria, la cual está relacionada con la topología que es
generada por la unión de conjuntos de proteínas creando una macromolécula
(El ARN polimerasa, es un ejemplo de esto). Una gran cantidad de
problemas de investigación en bioinformática están relacionados con la
identificación de la estructura secundaria y terciaria de una proteína y las
propiedades de está, dada su cadena de aminoácidos o estructura primaria.

2.1.5 Estructura secundaria de las proteínas

La estructura secundaria de una proteína hace referencia a la conformación

en columna de la cadena polipéptida. Existen distintas clases de tipos de
elementos en la estructura secundaria de las proteínas, como lo son las alfa-
hélices, los beta-sheets, los bucles y los dobleces; estos últimos reversan la
dirección de la cadena polipéptida.

2.1.6 Alfa-hélices y proteínas de transmembrana

Las Alfa-hélices juegan un papel crítico en la comprensión de la construcción

y funcionalidad de las proteínas, siendo el tipo de estructura más común
encontrado en éstas, razón por la cual se han propuesto a lo largo de las
últimas décadas, varios métodos que ayudan a su estudio. Existen varias
motivaciones relacionadas con el estudio de las Alfa-hélices, como el análisis
de la estructura de proteínas, el cual ha dado parámetros importantes para
realizar comparaciones entre proteínas de tipo estructural y geométrico entre
otros, lo que a su vez, ha permitido relacionar muchas proteínas en
determinadas familias, consiguiendo de esta manera, categorizar su
funcionalidad y propiedades.

Las proteínas de membrana (Figura 2.1.2), están clasificadas en dos

familias, las “Alpha-Helical” y “Beta-Barrels”, La familia “Alpha-Helical” cuenta
mayoritariamente con Alfa-hélices y representa generalmente la mayoría de
las proteínas de membrana encontradas en el genoma. Estas se caracterizan
por interactuar con ambas regiones de la membrana, ya que suelen
atravesarla en varias ocasiones, siendo las estructuras alfa-hélice en estas

10
proteínas, las que se encuentran en el interior de la membrana. Las
anteriores, participan en la traducción de señales a través de las membranas
celulares y controlan el flujo de moléculas que ingresan y salen de la célula.
Debido a esto, son de gran importancia en biomedicina, ya que sirven de
blanco para la mayoría de los agentes farmacéuticos de uso actual y futuro,
como es el caso del desarrollo de antibióticos.

Estructura terciaria de una proteína de transmembrana.

(Las flechas son beta-sheets, los bucles representan alfa-hélices.)
(Figura 2.1.2)

11
2.2 Métodos Computacionales

2.2.1 Modelos ocultos de Markov

Los modelos ocultos de Markov son utilizados para modelar y resolver

interrogantes concernientes a un proceso que no es observable
directamente, sino mediante la observación de otros eventos que están
ligados al proceso oculto.

Los modelos ocultos de Markov están caracterizados por:

Un numero de estados del modelo n, que conforman el proceso oculto.

Un conjunto de símbolos y el estado en el tiempo t es representado por:

S = {s1 , s2 ,..., sn }
qt ∈ S , 1 ≤ t < T

El numero de estados o símbolos observables m, corresponden al conjunto

de posibles tipos de observaciones que son provocados o que cuentan con
alguna relación en el comportamiento del proceso oculto. Este conjunto de
símbolos es representado por:

V = {v1 , v2 ,..., vn }
qt ∈ S , 1 ≤ t < T

Una matriz de distribución de estados de transición A = {ai , j } donde:

ai , j = P[qt +1 = S j | qt = Si ], 1 ≤ i, j < n

Una matriz de distribución de probabilidad de símbolos de observación o

matriz de difusión B = {bi ,k } donde:

bi , k = P[vk ∧ t = tc | qt = Si ], 1 ≤ i, j < n

El vector de distribución inicial de estados, indica la probabilidad de iniciar en

el estado si . El vector es representado por:

12
π i = P[q1 = Si ]

Por ultimo, definimos la secuencia de símbolos observados en el tiempo:

O = {O1,O2 ,..., OT }
Ot ∈ V , 1 ≤ t < T

Con lo anterior, un modelo oculto de Markov es representado por la tripla:

λ = (π , A, B) .

Los modelos de Markov son comúnmente utilizados para describir algún tipo
de proceso particular y con base en éste generar secuencias que simulan un
determinado proceso. Adicionalmente, estos cuentan con tres tipos de
problemas básicos de interés, los cuales se conocen como: el problema de la
evaluación, decodificación y aprendizaje.

El problema de la evaluación, en un modelo oculto de markov, consiste en

calcular la probabilidad de ocurrencia de determinada secuencia de estados
ocultos, dada una secuencia de observación. La solución a este problema la
provee el algoritmo Forward.

El problema de la descodificación plantea la búsqueda de la secuencia de

estados ocultos más probable para una secuencia de observaciones dada. El
algoritmo Viterbi provee la solución a este problema.

El problema mas complejo, es el problema del aprendizaje, que consiste en

que dada una secuencia de observaciones, determinar el modelo oculto, que
genere las secuencias de estados ocultos con mayor probabilidad, respecto
a las observaciones. El algoritmo de optimización Forward-Backward es la
propuesta de solución mas comúnmente usado en este problema.

[Link] Algoritmo Forward

Dada una secuencia de observaciones O y un HMM con parámetros

(π , A, B ) definidos, queremos determinan la probabilidad de la secuencia O .
Debido a que pueden existir varias secuencias de estados ocultos que den
como resultado una secuencia de observaciones O , la aproximación más
básica es calcular, para cada secuencia posible de estados ocultos, la

13
probabilidad de que genere la secuencia de observaciones O y luego sumar
los resultados.

Si m = numero de estados ocultos y t la longitud de la secuencia, entonces

tendríamos.

Pr(O | (π , A, B )) = ∑s∈C Pr(O | S )

t

Donde C t es el conjunto de las secuencias de estados ocultos de longitud t .

Nótese que si H es el conjunto de estados ocultos tenemos:

C1 = H , C 2 = H × H ,  , C t = C t −1 × H

De manera que:

C1 = H = n, C 2 = H × H = n 2 ,  , C t = C t −1 H = n t −1 × n = n t .

Con lo anterior es clara la imposibilidad práctica de utilizar esta

aproximación, con tamaños pequeños de secuencias, tamaños suficientes,
para sobrepasar el poder de cómputo de cualquier máquina, debido al
crecimiento geométrico de las posibles rutas a evaluar. Afortunadamente,
existe una aproximación más eficiente para este problema, conocida como el
algoritmo Forward.

El algoritmo Forward es un método de programación dinámica utilizado para

calcular la probabilidad de una secuencia de observaciones O , dado un
modelo oculto de Markov. La base de su efectividad radica en el supuesto de
que las matrices de probabilidad de transición A y de difusión B son
invariantes en el tiempo, de esta manera se puede hacer uso cálculos
realizados anteriormente.

[Link] Definición del algoritmo Forward

Sea k (t ) El estado observado en el tiempo t , bij la probabilidad de observar

el estado j dado que el estado oculto es i y aij la probabilidad de transición
del estado i al j . Definimos la probabilidad parcial del estado j , en el

14
tiempo t , como α t ( j ) , el cual es la suma de probabilidades para las sub-
secuencias de estados ocultos de longitud t , que terminan en el estado j :

n
α t +1 ( j ) = b jk ( t +1) ∑ α t ( i ) aij
i =1

Con lo anterior se define el algoritmo Forward de la siguiente manera:

Caso base: t = 1
Sea π j = probabilidad inicial del estado j.

α 1 ( j ) = π j b jk ( 1)

Caso recursivo: 0 < t ≤ T

n
α t +1 ( j ) = b jk ( t +1) ∑ α t (i )aij
i =1

Finalmente calculamos la suma de las probabilidades parciales al tiempo T ,

con lo cual se obtiene la probabilidad total.

n
Pr(O | (π , A, B )) = ∑ α T ( j )
j =1

El tiempo requerido por el algoritmo Forward es descrito mediante la

siguiente ecuación de diferencia:

a (0) =1
a (t ) = a (t − 1) + i + c

Solucionando la ecuación de diferencia tenemos:

a (t ) = t (i + c)

15
Utilizando notación asintótica [5] tenemos:

O(a (t )) = O(t (i + c ))
O(t (i + c)) = O(O(t )O(i + c))
O(O(t )O(i + c)) = O(O(t )O(1)) = O(t )

Como resultado, el algoritmo requiere de un tiempo del orden de O(t ) . El

cual es inferior en comparación a la aproximación de búsqueda exhaustiva
que se encuentra en el orden de O(nt ) .

[Link] Algoritmo Viterbi

Dada una secuencia de observaciones O y un HMM con parámetros ( π , A, B )

definidos, queremos determinar la secuencia de estados ocultos
E = { e1 , e2 ,..., eT } , que cuente con la mayor probabilidad de ocurrencia, esto
es:

max ( Pr ( E | O ∧ ( π , A, B ) ) )
i =1

De manera similar al caso del cálculo de la probabilidad total, para una

secuencia O , pueden existir varias secuencias de estados ocultos que den
como resultado una secuencia de observaciones O , por tanto, una primera
aproximación, es calcular para cada secuencia posible de estados ocultos la
probabilidad de que ésta genere a O y luego sumar los resultados. Este
cálculo es muy similar al visto anteriormente y como se mostró el número de
secuencias a analizar, crece geométricamente, por lo cual no es un método
de solución útil para tamaños de ejemplar grandes.

El algoritmo de Viterbi es un algoritmo de búsqueda que utiliza el método de

programación dinámica para encontrar la secuencia de estados ocultos con
mayor probabilidad de ocurrencia en un HMM, dada una secuencia de
observaciones. Como en el caso del algoritmo Forward se hace uso del
supuesto de invarianza en el tiempo de las matrices de transición de
probabilidad y de difusión.

16
Sea k (t ) El estado observado en el tiempo t , bij la probabilidad de observar
el estado j , dado que el estado oculto es i y aij la probabilidad de transición
del estado i al j . Definimos la probabilidad parcial del estado oculto j en el
tiempo t como δ t ( j ) , la cual es igual a probabilidad de la secuencia de
estados ocultos con longitud t , que termina en el estado j , dada una
secuencia de observaciones O = { o1 , o2 ,..., ot } , que cuenta con la máxima
probabilidad:

δ t ( j ) = max (δ t −1 ( i ) aij b jk ( t ) )
n

i =1

Sea arg max la función que retorna el subíndice j que maximiza δ t −1 ( j ) a ji .

Definimos φt ( i ) , como el estado oculto por el cual pasa la secuencia en el
tiempo t − 1 tal que esta cuenta con la máxima probabilidad al alcanzar el
estado i :
n
φt ( i ) = arg max (δ t −1 ( j ) a ji )
j =1

Con esto definimos el algoritmo de Viterbi de la siguiente manera:

Caso base: t = 1

δ1 ( j ) = π jb jk ( 1)

Caso recursivo:
δ t ( i ) = max (δ t −1 ( i ) aij b jk ( t ) )
n

j =1
n
φt ( i ) = arg max (δ t −1 ( j ) a ji )
j =1

17
Al finalizar se determina el estado final más probable s y se recupera la ruta
de máxima probabilidad Y T = { y1 , y2 ,..., yT } de la siguiente manera:
n
s = arg max δ T ( j )
j =1

yT = s

yt = φt +1 ( yt +1 )

El tiempo requerido por el algoritmo de Viterbi es descrito mediante la

siguiente ecuación de diferencia:

a (0) =1
a (t ) = a (t − 1) + i + c

Solucionando la ecuación de diferencia tenemos:

a (t ) = t (i + c)

Utilizando notación asintótica sobre a(t ) :

O(a (t )) = O(t (i + c ))
O(t (i + c)) = O(O(t )O(i + c))
O(O(t )O(i + c)) = O(O(t )O(1)) = O(t )

Este algoritmo, al igual que el Forward se encuentra en el orden de O(t ) ,

esto se debe a que el algoritmo Viterbi realiza de manera similar la
recurrencia.

[Link] Algoritmo Forward-Backward

Como se introdujo anteriormente, el algoritmo Forward-Backward es un

método de optimización de los parámetros de un modelo oculto dado, con
base en una secuencia de observaciones. De manera simplificada el
algoritmo busca determinar las probabilidades tanto en el vector π , como en
las matrices de distribución de transición y de difusión, de tal manera que

18
para una secuencia de observaciones, el modelo produzca una secuencia de
estados ocultos que cuente con la mayor probabilidad posible. La idea básica
de cómo determina los parámetros en cada iteración, es mediante el cálculo
de las probabilidades de las posibles secuencias de estados para la
secuencia dada, para finalmente realizar su ajuste.

[Link] Definición del algoritmo Forward-Backward

Primero iniciaremos repasando la definición de la probabilidad parcial

Forward:

n
α t +1 ( j ) = b jk ( t +1) ∑ α t ( i ) aij
i =1

Ahora definimos la función de probabilidad parcial Backward, como la

probabilidad para la secuencia de observación parcial O = {Ot +1,Ot + 2 ,..., OT } ,
partiendo en el tiempo t, hasta un tiempo final T, iniciando en el estado Si con
tiempo t como:

N
βT ( i ) = 1, β t ( i ) = ∑ aij b j (Ot +1 ) β t +1 ( j )
j =1

Con lo anterior, la probabilidad de encontrarse en el estado Si en el tiempo t,

dada la secuencia de observaciones O , corresponde a la función de medida
probabilística γ t ( i ) :

α t (i ) β t (i )
γ t (i ) = N

∑ α ( j )β ( j )
j =1
t t

La probabilidad de encontrarse en el estado Si en tiempo t, y de estar en el

estado Sj en t+1 esta dada por:

α t (i )aij b j (Ot +1 ) β t +1 ( j )
ξt ( i, j ) = N N

∑∑α (i)a b (O
i =1 j =1
t ij j t +1 ) β t +1 ( j )

19
Esta última corresponde a la probabilidad conjunta de 3 eventos. El primero,
que esté en el estado Si en el tiempo t dada por la probabilidad forward,
encontrarse en el estado Sj en el tiempo t+1 correspondiente a la
probabilidad Backward y finalmente que se produzca una transición de Si a Sj
aij b j (Ot +1 ) .

Con el resultado anterior tenemos:

T −1

∑ γ (i) : El numero esperado de transiciones desde Si.

t =1
t

T −1

∑ ξ (i, j ) : El numero esperado de transiciones de Si a Sj.

t =1
t

Finalmente utilizamos las expresiones anteriores para definir los estimadores

para los parámetros del modelo:

π i = γ 1 (i) : La frecuencia esperada para en estado Si en el tiempo t=1.

T −1

∑ ξ (i, j )
t
aij = t =1
T −1 : Numero de transiciones de Si a Sj, sobre las transiciones a Si.
∑ γ (i )
t =1
t

∑ γ ( j)
t
t =1, Ot =Vk
b j (k ) = T : Las veces en Sj, observando Vk, sobre veces en estado Sj.
∑ γ t ( j)
t =1

Con esto se obtendrá un nuevo modelo λ que cumple la condición:

P (O, λ ) > P(O, λ )

20
Esto asegura la convergencia del modelo mediante la aplicación de varias
optimizaciones sucesivas, sin embargo la debilidad de este método radica en
que no garantiza la convergencia del modelo al máximo global.

2.2.2 Algoritmos Genéticos

Los algoritmos genéticos son una clase de métodos de programación

evolutiva, comúnmente utilizados en problemas de optimización complejos.
Estos fueron inventados por John Holland (Holland, 1975) en los 60. Los
algoritmos genéticos son un método que mediante la “evolución” de una
población de “cromosomas” muestrea el espacio de búsqueda para un
problema, con el propósito de encontrar los parámetros óptimos para éste;
utilizando un conjunto de operadores inspirados en la selección natural,
como el cruce, mutación e inversión de la información representada en los
cromosomas. Comúnmente los cromosomas están representados en bits, los
cuales son denominados comúnmente como “genes”. Parte de las ventajas
con las que cuentan los algoritmos genéticos, radican en su capacidad para
explorar extensivamente el espacio de búsqueda del problema. Gracias a
esto logran obtener mejores resultados sobre otros métodos, en especial
cuando el problema cuenta con varios óptimos locales.

Como primera característica en la utilización de los algoritmos genéticos,

contamos con, un supuesto inicial sobre la representación de los parámetros,
que consiste en que son discretos; usualmente, los parámetros son
discretizados en potencias de 2, siguiendo el modelo canónico. Por lo
general la realización de esta discretización no es generalmente un problema
complejo cuando los parámetros son continuos. Sin embargo, ésto no es
trivial en los casos donde los parámetros son discretos, debido
principalmente a que la discretización puede generar cromosomas no válidos
y de esta manera dando lugar a un espacio de búsqueda no óptimo.

Los algoritmos genéticos comúnmente son descritos como un método de

optimización que no utiliza información de gradiente. Por tanto son aplicables
a problemas con funciones no diferenciables o con múltiples óptimos locales.
Por lo general si existe un método de optimización especializado para un
problema, este será mas eficiente que una aproximación basada en
algoritmos genéticos; por esta razón es común, el uso de una estrategia
hibrida que utiliza la sinergia entre los algoritmos genéticos con otro método.

2.3.1 Forma canónica de los algoritmos genéticos

La forma canónica representa el modelo general utilizado por Holland, para

representar un problema con algoritmos genéticos. El modelo de un

21
algoritmo genético canónico consta de 4 objetos: Una población, individuos,
una función de aptitud o calificación y un conjunto de operaciones básicas
sobre los cromosomas denominadas selección, cruce y mutación.

La población es una colección de individuos que realizan un muestreo al

espacio de búsqueda; su tamaño está relacionado con las características
especificas del problema tales como el tamaño del espacio de búsqueda y la
función de calificación. En el modelo canónico este tamaño es constante en
el tiempo.

Los individuos o cromosomas, son cadenas binarias que representan

instancias de valores pertenecientes al conjunto de parámetros a optimizar. A
cada valor binario en la cadena se le conoce como un gen; de esta manera
un cromosoma representa una posición en el hipercubo que representa al
espacio de búsqueda. Para el caso canónico la longitud de la cadena es
constante.

La función de aptitud, objetivo o de calificación, es comúnmente una función

de B k → R . El conjunto de k parámetros es representado por los
cromosomas cuyo valor se desea optimizar. Esta función generalmente
relaciona y/o pondera un conjunto de características abstractas del problema
particular que se desea solucionar.

Finalmente se tiene un conjunto de operaciones sobre la población y los

cromosomas conocidos como selección, cruce y mutación que realizan
transformaciones sobre los genes:

La operación de selección escoge entre la población, los cromosomas que

cuentan con el menor o mayor valor dado por la función de aptitud, esto
dependiendo de si se busca minimizar o maximizar la función, De esta
manera, se crea una nueva colección de individuos conocida como la
población intermedia.

La operación de cruce intercambia los genes de dos individuos “padres”,

creando dos nuevos individuos o cromosomas “hijos”, esto con el propósito
de que alguno de los hijos cuente con las características fuertes de ambos
padres.

Finalmente la operación de mutación, modifica los genes de los individuos

con el propósito de agregar diversidad a la población y de esta manera
incrementar la efectividad del muestreo del espacio de búsqueda.

22
[Link] Etapas en un algoritmo genético

La mayoría de los algoritmos genéticos cuentan 3 etapas: La primera, es la

generación de la población inicial; aquí se determinan de manera aleatoria
los genes de los individuos en la población y con esto se busca obtener una
dispersión uniforme de los individuos en el espacio de búsqueda. En la
segunda etapa, los individuos en la población son evaluados por medio de la
función de aptitud y se aplica la operación de selección con base en los
resultados generando la población intermedia. En la última etapa, se aplican
los operadores de cruce y mutación sobre la población intermedia, utilizando
generalmente como parámetros la calificación de los individuos, de esta
manera se crea una nueva población. Al ciclo de las etapas de evaluación,
selección y aplicación de los operadores de cruce y mutación se le conoce
como una “generación”. Este ciclo puede ser realizado un número
establecido de veces o hasta que se cumpla con una condición establecida
sobre las calificaciones, bien sea de la población o de un individuo particular.

[Link] Teoría del Esquema

La teoría del esquema es la aproximación mas comúnmente utilizada en la

teoría de los algoritmos genéticos canónicos. Esta se basa en el análisis de
como los algoritmos genéticos muestrean el espacio de búsqueda y la taza a
la cual lo realizan. Como preliminar se presentará el espacio de búsqueda
como un hipercubo y una de las propiedades más interesantes de los
algoritmos genéticos conocida como el “paralelismo implícito”.

[Link] El espacio de búsqueda como un hipercubo y esquemas

Si se considera la cadena binaria que modela a los individuos en un

algoritmo genético, como un vector que indica una posición en un espacio
B L , de esto podemos interpretar geométricamente el espacio B L como un
hipercubo H de L dimensiones. Ahora si utilizamos los símbolos {0,1,*}
para representar los valores de la cadena en donde ‘*’ significa indiferente (0
o bien 1), entonces los posibles cromosomas que una de estas cadenas
puede representar, genera una partición o hiperplano h ; al conjunto de las
particiones de H que se pueden generar utilizando el símbolo ‘*’ se conoce
como el conjunto de esquemas de H . Un detalle importante de lo anterior
radica en que una cadena binaria arbitraria puede ser representada por
varios esquemas y por tanto pertenecer a varios hiperplanos al mismo
tiempo.

23
[Link] Paralelismo implícito

El paralelismo implícito hace referencia a la propiedad de los algoritmos

genéticos de muestrear varios hiperplanos de manera no explícita, esto es,
cada ocasión en que es evaluada una cadena, implícitamente se está
realizando una evaluación de los hiperplanos a los cuales pertenece. La
razón de contar con poblaciones obedece al hecho de que la evaluación de
una sola cadena, provee a manera de analogía, de “un grano de arena” en
información sobre las propiedades (aptitud) del conjunto de puntos
pertenecientes a cada uno de los hiperplanos, que corresponden a los
esquemas representables en la cadena. Cuando esta información sobre los
hiperplanos es acumulada mediante la evaluación de toda una población, se
obtiene un muestreo aproximado de las propiedades de estos hiperplanos (al
contar con un número superior de individuos distintos realizando el muestreo
de un hiperplano dado, la aproximación es mejor). De esta manera, el
algoritmo genético potencializa la exploración del espacio de búsqueda.

[Link] Esquemas y disrupciones

Como ya se introdujo brevemente, los esquemas representan particiones o

hiperplanos de H , ahora realizaremos algunas definiciones notacionales
concernientes a los esquemas y al concepto de disrupción de esquemas.
Sea S = [ s1 , s2 ,..., sL ] ∧ s1 , s2 ,..., sL ∈{0,1,*} , un esquema para una cadena binaria
de longitud L . Definimos o( S ) = ( + i | 1 < i ≤ L ∧ ( si ='0'∨ si ='1') : 1) , como el orden
del esquema S . Con esto definimos a una disrupción como una modificación
del esquema S . Es importante tener en cuenta que esto no necesariamente
implica un cambio en el orden del esquema, sino un cambio en el hiperplano
representado por este, como se ilustra en (figura 2.3.1), al realizar una
modificación a la cadena C , el esquema S 2 a diferencia de S1 sufre una
disrupción cambiando a S '2 .

[Link] Efecto en el muestreo del operador de selección

La aplicación del operador de selección, provoca una disminución del

número de cadenas que muestrean los hiperplanos del espacio de
búsqueda, como consecuencia, de la pérdida de cromosomas en la
población que realizaban un muestreo a los hiperplanos. El cambio en el
número de cadenas que muestrean un hiperplano dado h , luego de la
selección se define de la siguiente manera: Sea M (h, t ) , el número de
cadenas que muestrean el hiperplano en la generación t . Sea (t + ∆) el

24
índice para la generación t luego de la operación de selección (población
intermedia de t ). Sea f (h, t ) el promedio de aptitud para las cadenas que
muestrean a h en la generación t y f el promedio de aptitud de la población
en la generación t . La relación del cambio en el muestreo del hiperplano h ,
entre la población en la generación t y la población intermedia es:

f (h, t )
M (h, t + ∆) = M ( h, t )
f

C = 1001001110111
S1 = 10 * * * * * * * * *11
S 2 = * * 01 * * * * * * * *1
o( S1 ) = 4
o( S 2 ) = 3

C ' = 1011001110111
S1 = 10 * * * * * * * * *11
S '2 = * *11 * * * * * * * *1
o( S1 ) = 4
o( S '2 ) = 3

Esquemas y el efecto de la disrupción.

Figura 2.2.1

[Link] Efecto en el muestreo del operador de cruce

El operador de cruce es aplicado con probabilidad p sobre la población, de

manera que solo algunos de los cromosomas se vean afectados. Desde la
perspectiva del muestreo del hiperplano h , los nuevos cromosomas
generados mediante cruce, generarán pérdidas o ganancias a la cantidad de
cadenas que muestrean a h . Lo anterior es consecuencia de la disrupción
(separación de las posiciones que definen a los esquemas) causado por el
operador de selección. Con esto la expresión correspondiente al cambio en
el muestreo para el hiperplano h , luego de aplicar el operador de selección y
cruce a la población es:

25
 f (h, t )  f (h, t ) 
M (h, t + 1) = (1 − p ) M (h, t ) + p M ( h, t ) (1 − perdidas ) + ganancias 
f  f 

En este punto se aplica un supuesto conservativo, el cual consiste en asumir

que el operador de cruce siempre genera disrupciones en los esquemas
representando a h , de esta manera el término de ganancias es eliminado de
la expresión, obteniendo de esta manera la siguiente desigualdad:

f (h, t )  f (h, t ) 
M (h, t + 1) ≥ (1 − p ) M (h, t ) + p M (h, t ) (1 − disrupciones) 
f  f 

Ahora para el término de disrupciones tenemos: Sea P (h, t ) la representación

proporcional de h , la cual se obtiene dividiendo M (h, t ) por el tamaño de la
población, ∆(h) L − 1 es la probabilidad de crear una disrupción del esquema
que define a h al aplicar la operación de cruce de un punto (ver [6]). La
probabilidad de escoger una pareja en una población tal, que éstas realicen
un muestreo a h es P (h, t ) . Finalmente, la probabilidad de que la operación
de cruce genere una disrupción en el esquema que define a h es definida
como:

∆ ( h)
(1 − P(h, t ) )
L −1

El cambio en el muestreo del hiperplano h , de la generación t , a la t + 1

luego de aplicar las operaciones de selección y cruce se define como:

f (h, t )  ∆ ( h)  f (h, t )  
M (h, t + 1) ≥ P (h, t ) 1− p 1 − P (h, t ) 
f  L − 1  f  

26
[Link] Efecto en el muestreo del operador de mutación

Como se había comentado anteriormente la mutación genera cambios

aleatorios en los genes del cromosoma, generando con esto disrupciones en
los esquemas y modificando el muestreo de los hiperplanos. Al incluir el
efecto de la mutación en la ecuación que describe el muestreo de un
hiperplano tenemos: Sea o(h) = o( S ) , donde S es el esquema que describe
al hiperplano h . Sea q la probabilidad de mutación para cada gen del
cromosoma. Sea (1 − q) o ( h ) la probabilidad de que la operación de mutación
no afecte al esquema S . El teorema del esquema, que es la desigualdad que
describe el cambio en el muestreo del hiperplano h luego de realizar las
operaciones de selección, cruce y mutación se define como:

f (h, t )  ∆ ( h)  f (h, t )  
M (h, t + 1) ≥ P (h, t ) 1− p 1 − P (h, t )  (1 − q ) o ( h )
f   L −1   f  

[Link] Comentarios sobre la teoría del esquema

El teorema del esquema se concentra en el comportamiento del muestreo de

hiperplanos para un algoritmo genético y por tanto de su preservación. La
debilidad mas grande del teorema radica en que explica muy poco sobre la
dinámica de convergencia del algoritmo genético a un óptimo, esto es, no
provee de información con respecto a los parámetros a utilizar en una
determinada implementación o problema, así como tampoco con respecto a
la complejidad computacional del mismo. Como ejemplo de esta debilidad se
tiene, bajo el teorema del esquema que la operación de mutación, debido a
sus propiedades disruptivas es un operador indeseable. Esto como
consecuencia de la incapacidad del teorema de proveer de un modelo que
concentre el comportamiento del algoritmo genético en conjunto para todos
los hiperplanos del espacio de búsqueda (si bien la mutación puede generar
una disrupción en un determinado esquema, este también creará una
ganancia en el muestreo de otro hiperplano).

27
3. Análisis método HMMTOP

La metodología para el análisis del método HMMTOP tuvo como primera

parte la realización de un análisis tanto de la técnicas utilizadas en éste como
la manera en que fueron empleadas; esto incluye el modelo de estados,
parámetros iniciales y resultados. Con base en lo anterior, se procedió a
realizar una implementación de éste, con la información proveída sobre el
método con el ánimo de reproducir los resultados y observar otras
debilidades y fortalezas no documentadas.

3.1 Introducción

HMMTOP [1] es un método para la predicción de la estructura secundaria en

proteínas de transmembrana el cual utiliza un modelo oculto de Markov
(HMM) con este propósito. El modelo de Markov cuenta con 5 estados
ocultos y 20 estados observables que representan a los 20 distintos
aminoácidos que conforman las secuencias proteicas. La base teórica del
método, radica en la hipótesis de los autores, la cual propone que la
localización de los segmentos y la topología de la proteína, son determinadas
por las diferencias en las distribuciones de aminoácidos entre los segmentos
y la distribución general de la secuencia. Esta base teórica es contraria a una
composición específica en los aminoácidos de la secuencia.

3.2 El Método HMMTOP

3.2.1 Base Teórica de HMMTOP

La hipótesis propuesta por los autores postula que es la divergencia entre la

distribución de aminoácidos de un segmento estructural y la distribución de
toda la secuencia y no una composición especifica de aminoácidos en el
segmento, lo que determina el tipo de estructura (Bucle, Hélice o cola de
hélice) y la localización en la membrana celular que éste tendrá (Interior o
exterior). Este tipo de hipótesis es muy similar a la usada en otro tipo de
problemas de identificación de patrones como en visión artificial, en el cual
se identifican diferencias con base en el cambio de la distribución de los
píxeles en un área determinada con respecto a la imagen completa.

28
Para realizar la identificación del tipo de estructura y su localización, para
cada segmento, los autores postulan que para cada tipo de estructura y
localización en la membrana, se presenta una distribución específica, es
decir la distribución presentada por un bucle interno difiere de la observada
para un bucle externo o una Alfa-hélice.

Con base en esta hipótesis, los autores proponen que la estructura

presentada por la secuencia proteica será la que presente la máxima suma
de las divergencias entre las distribuciones de los segmentos y la distribución
de la secuencia completa, motivo por el cual se puede decir que HMMTOP
resuelve un problema de optimización.

3.2.2 La función de divergencia

Los autores utilizan la función de divergencia de Lullback-Leibler, como

norma para medir la diferencia entre las distribuciones del segmento y la
secuencia. La función de divergencia también es llamada algunas veces la
función de ganancia de información: Si X es un conjunto de posibilidades, P y
Q son distribuciones para X, entonces la función de divergencia retorna la
ganancia de información sobre X obtenida si P es usada en vez de Q.

 P (i ) 
Dkl = ∑ P (i ) log 
i  Q (i ) 
Función de divergencia

En el caso particular de HMMTOP, Pj es la distribución del j-esimo

segmento de aminoácidos y Q corresponde a la secuencia completa,
entonces HMMTOP busca encontrar una combinación de segmentos Dk
para la secuencia tal que:

n
 P (i ) 
Dk = ∑∑ Pj (i ) log j  y Dk es máximo
j i  Q (i ) 

29
3.2.3 Modelo HMMTOP de la estructura para proteínas de
transmembrana

El modelo que crean los autores de la secuencia proteica está basado en la

definición de 3 tipos de estructura, el bucle (loop), la cola de hélice (tail) y las
hélices.

Las hélices son estructuras que se encuentran embebidas en el interior de la

membrana. Los bucles son segmentos de la secuencia que se extienden
hacia el exterior de la membrana celular, estos cuentan con una longitud
arbitraria y pueden formar un dominio o una estructura simple. Las colas de
hélice son elongaciones de las hélices de membrana localizadas al exterior
de la misma, las cuales pueden estar seguidas por un bucle u otra cola
formando un bucle corto que interactúa con las áreas interior o exterior de la
membrana.

Adicionalmente, los autores también se refieren a 2 áreas separadas por la

membrana, llamadas interior (citoplasma) y exterior (el medio de la célula).
Para el caso de los bucles y las colas, se diferencian si se encuentran en el
área interior o exterior. De esta manera los autores definieron 5 estados base
para el modelo oculto de Markov (figura 3.1), denominados cola interna
(inside tail), cola externa (outside tail), bucle interno (inside loop), bucle
externo (outside loop) y hélice (helix).

La importancia de distinguir las estructuras internas y externas se encuentra

en el hecho de que la composición del medio citoplasmático y extra-
citoplasmático son distintos, razón por la cual se considera que la
composición de las estructuras en cada medio debe diferir.

3.2.4 El modelo oculto de Markov de HMMTOP

El HMM utilizado por HMMTOP, utiliza 5 estados base, que representan a

cada una de las estructuras definidas anteriormente. Adicionalmente, define
dos tipos de estados, los FL o de longitud fija, que cuentán con una cota
superior e inferior para la longitud de secuencia en la estructura y los NFL, de
longitud no fija para las estructuras con longitud arbitraria. Los estados que
caracterizan a las estructuras hélice, y cola interior y superior son de tipo FL,
los estados Bucle son de tipo NFL (figura 3.2).

30
 Estructuras de la proteína y su localización en la membrana celular.
Figura 3.1
(Tomada de [1]).

Los estados FL cuentan con dos cotas una inferior (MINL) y una superior
(MAXL), que limitan la longitud de la estructura que representan; por tanto,
cada estado FL cuenta con MAXL sub-estados. Para garantizar esto, los
primeros MINL-1 sub-estados cuentan con una única transición al siguiente
sub-estado de manera incremental, es decir, el i-esimo sub-estado cuenta
con una única transición al estado i+1. Los sub-estados MINL a MAXL-1
cuentan con 3 transiciones posibles para las colas iniciales (saliendo de la
hélice), 2 para las colas terminales y hélices y por ultimo el sub-estado
MAXL solo cuenta con 1 o 2 transiciones posibles dependiendo si es cola o
hélice.

Los estados NFL cuentan con 2 transiciones, la primera para el siguiente

estado y otra hacia si mismo, de esta manera define estructuras (Bucles) de
longitud arbitraria.

31
 Modelo de estados ocultos de HMMTOP.
Figura 3.2.
(Tomada de [1])

Los valores de los parámetros MAXL y MINL para los estados tipo FL fueron
determinados en base a observaciones de la distribución de longitud para
cada estructura. En el caso de las colas internas y externas MINL y MAXL los
valores son 1 y 15, y para las hélices los valores son 17 y 25
respectivamente. Con estos parámetros se determinó el número total de
estados ocultos del modelo y utilizando las transiciones entre los estados, se
definieron las restricciones que la matriz de transición y el vector π deben
cumplir.

32
3.2.5 Método de predicción de HMMTOP

HMMTOP consta de tres pasos para realizar este método: definición de los
parámetros iniciales, optimización del modelo y predicción.

En la definición de los parámetros iniciales del modelo oculto, se define el

vector de probabilidad inicial π , determinando los estados en donde es
valido iniciar la secuencia, debido a que las secuencias no inician con una
hélice, solo los estados iniciales para colas y bucles tienen una probabilidad
de inicio mayor a cero. En el caso de la matriz de transición, las no definidas
en el modelo, tienen una probabilidad de cero y los estados con una única
transición tienen una probabilidad de 1, este es el caso para los sub-estados
1 a MINL-1 en el tipo FL. Los autores aclaran que las probabilidades de los
parámetros no restringidos por el modelo, pueden ser definidos con valores
fijos o aleatorios. En este caso particular, los autores los definen con valores
fijos, obtenidos a partir de las distribuciones en proteínas del conjunto de
entrenamiento con estructuras bien identificadas experimentalmente.

La optimización del modelo se realiza mediante utilización del algoritmo

Forward-Backward y un conjunto de secuencias de entrenamiento con el
propósito que el modelo “aprenda” el patrón estructural particular presente en
las secuencias. Debido al problema de la convergencia en máximos locales,
los autores utilizaron adicionalmente el método de pseudo-conteo sobre la
información estadística de las secuencias de aprendizaje, con el propósito de
ubicar en un “buen” lugar del espacio de búsqueda, al modelo inicial. El
método de pseudo-conteo es un caso especial de mixtura de Dirichlet, el cual
ayuda a incrementar la efectividad de la información estadística a la hora de
definir las probabilidades iniciales. De esta manera incrementando la
probabilidad de convergencia en el máximo global del algoritmo Forward-
Backward, el arreglo utilizado por el método de pseudo-conteo fue derivado,
al igual que los parámetros iniciales de un conjunto de proteínas con
topología bien definida.

El último paso de predicción se realiza mediante el algoritmo de Viterbi sobre

el modelo generado, luego del proceso de optimización y la secuencia a
predecir. Viterbi retorna a la secuencia de estados mas probable dada la
secuencia y el modelo para luego codificar los estados ocultos respecto a su
estado estructural correspondiente (figura 3.3).

33
 Predicción de estructura secundaria de HMMTOP.
Figura 3.3.
(Tomado de [1])

3.2.6 Conjunto de pruebas y resultados de HMMTOP

Los autores utilizaron tres conjuntos de secuencias (83TMP [7], 48TMP [8] y
prokTMP [9]), previamente utilizadas en otros métodos de predicción y
recolectados de la base de datos SwissProt [10], para la obtención de los
parámetros del modelo y para la realización de las pruebas. Las información
de secuencia y estructura del conjunto 83TMP, fue utilizado en la obtención
del vector de pseudo-conteo y en la optimización mediante el algoritmo
Forward-Backward.

Los conjuntos 48TMP y prokTMP y 83TMP fueron utilizados en la realización

de las pruebas. Previamente fueron descartados del conjunto 83TMP las
secuencias que no contaban con una identificación experimental confiable y
que contaban con una longitud superior a 500 residuos.

El conjunto de pruebas lo conformaron 83 proteínas para 83TMP, 47 para

48TMP y 44 en el caso de prokTMP, para un total de 174 proteínas y 698
segmentos de transmembrana. Con esto HMMTOP predijo un total de 709
segmentos de transmembrana de los cuales 694 fueron correctos. Todos los
segmentos de transmembrana fueron predichos correctamente en 74 (89%)
proteínas del conjunto 83TMP, 45 (96%) en 48TMP y 38 (87%) para
prokTMP para un total de 157 (90%).

3.3. Reproducción del experimento HMMTOP

Como se indicó en la introducción de esta sección, luego de realizar el

estudio sobre el método HMMTOP se procedió a realizar una reproducción
del método, con base en la información estudiada. A continuación se
presentarán los pasos que se siguieron para realizar esta tarea, los
resultados obtenidos y las conclusiones derivadas de esta implementación.

34
3.3.1 Definición del modelo oculto de Markov

Como primer paso, se procedió a definir el modelo oculto de Markov a ser

utilizado en la implementación. En la definición se siguió el mismo esquema
de estados básicos y condiciones para los sub-estados (parámetros MINL y
MAXL) definidos por los autores de HMMTOP. Debido a que estos sub-
estados cuentan con parámetros a ser optimizados, se tomó la decisión de
incluirlos en el modelo. De esta manera se definieron 60 estados ocultos; la
razón de que se definieran tantos estados, radica en las condiciones de
longitud para los sub-estados de tipo FL. El diagrama de estados (figura 3.4)
muestra en detalle los estados definidos y sus transiciones. Primero, notar
que tenemos unos estados iniciales X, x, Z y z; éstos tienen como función
primaria indicar si la secuencia inicia con una cola (x y z) o un bucle (X, Z) y
segundo, indican si la secuencia inicia en la zona interior (X y x) o exterior (Z
y z) de la membrana.

Los estados “i” y “o” tienen una semántica distinta aquí, el estado “i” indica
una cola que entra a la hélice, “o” indica una cola que sale de la hélice y
finalmente “I” define a los bucles (interiores y exteriores). Lo anterior con el
propósito de garantizar secuencias correctas (cola entrando a hélice  hélice
 cola saliendo de la hélice) al utilizar estos estados en conjunto con la
información de los estados iniciales.

Con respecto a los estados observables, no hay modificación en relación a

los definidos en HMMTOP, contando con 20 estados, correspondientes a los
20 distintos aminoácidos. Con esto se obtuvo una matriz de transición de
60X60, una matriz de difusión de 60x20 y un vector π de 60 posiciones.

2.5 Método de predicción

El método de predicción en la implementación siguió un esquema similar al

utilizado por HMMTOP. Como diferencias tenemos que se omitió el uso del
método de pseudo conteo y segundo, la forma como se realiza la
optimización con el algoritmo Forward-Backward, utilizando múltiples
secuencias. Debido a que no se utilizó el método de pseudo conteo y
adicionalmente con la intención de probar el comportamiento del algoritmo
Forward-Backward, se realizaron pruebas utilizando 8 modelos iniciales que
contaban con las restricciones del modelo de estados definido anteriormente.

35
 Z x
z

X
I

H(25) i(1)
MAXL MINL

H(18) i(2)

o(1) o(15)
MINL MAXL

H(17)
i(j)
MINL
o(2)

i(15)
H(2)
MAXL

H(1)

Diagrama del modelo oculto de Markov de 60 estados.

(Figura 3.4)

El uso de múltiples secuencias de entrenamiento involucro la toma de una

decisión con respecto a qué método utilizar para la optimización. Esto,
debido al hecho que el algoritmo Forward-Backward está diseñado para
realizar la optimización con respecto a una única secuencia de
observaciones. Adicionalmente, no se encontró una referencia con respecto
a cómo atacar este problema y los autores de HMMTOP tampoco
documentaron el método utilizado. Debido a esto, se propusieron dos
posibles aproximaciones definidas de la siguiente forma:

Aproximación 1: Sea Si la i-esima secuencia de entrenamiento, λ0 el

modelo inicial, λi el i-esimo modelo generado y FB (λ , S ) → λ la función

36
iterativa de optimización Forward-Backward. La aproximación consiste en la
ejecución de los siguientes pasos de manera iterativa:

Inicialización:

λ1 = FB(λ0 , S1 )

Iteración:

λi = FB(λi −1 , Si ),1 < i < T

Aproximación 2. Sea Si la i-esima secuencia de entrenamiento, λ0 el

modelo inicial, λi el i-esimo modelo generado, λT el modelo final y
FB (λ , S ) → λ la función iterativa de optimización Forward-Backward. La
aproximación consiste en la ejecución de los siguientes pasos de manera
iterativa:

Iteración:

λi = FB(λ0 , Si )

Finalización:

λT = ∑
λi
T

Se realizó un análisis de las propiedades de cada aproximación, para

seleccionar la aproximación mas adecuada. De este análisis se determinó
que la segunda aproximación era la más adecuada. El factor determinante
fue la convergencia, ya que en el caso de la primera aproximación, no se
encontraron argumentos que la garantizaran al contrario de la aproximación
2. La (Figura 3.5) ilustra simplificadamente el comportamiento del modelo en
el espacio de búsqueda, para las dos aproximaciones a medida que se
realizan iteraciones. Los rombos indican el máximo local para cada
secuencia, las flechas rotuladas son las iteraciones y los círculos indican la
posición de un modelo en el espacio. Se observa que para el primer caso
existe la posibilidad que el modelo ingrese en un ciclo, ya que en cada
iteración busca aproximarse al máximo de la secuencia correspondiente. La
aproximación 2, por el contrario, busca el modelo que mejor se ajuste a las
distintas secuencias; sin embargo, este modelo no necesariamente es el
máximo global.

37
 Aproximación 1

Aproximación 2, dos iteraciones.

(Figura 3.5)

Con esto, finalmente el método realizará iteraciones hasta alcanzar un

umbral de convergencia dado. Este umbral representa la norma máxima
entre el modelo λn −1 y λn , siendo esta norma definida de la siguiente manera:

38
 ε= ∑ (π ( i ) −π ( i ) ) + ∑∑ ( A ( i, j ) −A ( i, j ) ) + ∑∑ ( B ( i, j ) −B ( i, j ) )
n n −1
2
n n −1
2
n n −1
2

3.2.6 Conjunto de pruebas y resultados de la implementación

En las pruebas se emplearon las secuencias de aprendizaje definidas en el

capitulo 5. Se realizaron 8 pruebas para evaluar el poder de predicción de la
implementación (Tabla 3.1) y otras características de ésta. El umbral de
convergencia utilizado en los modelos correspondientes a las 8 pruebas fue
de 0.001.

(Tabla 3.1)
Resultados de las pruebas para la implementación HMMTOP.

M1 M2 M3 M4 M5 M6 M7 M8
Secuencias E. 100 100 100 100 100 100 100 100
Secuencias P. 75 81 69 65 46 40 61 55
Verdaderos Positivos 797 850 756 674 538 314 626 614
Verdaderos Negativos 914 868 973 967 1197 1166 1047 1097
Falsos Positivos 619 665 559 565 335 367 486 436
Falsos Negativos 1063 1010 1104 1186 1322 1546 1234 1246
Completitud 50.42% 50.63% 50.95% 48.36% 51.13% 43.61% 49.3% 50.42%
Precisión 25.44% 25.63% 25.98% 23.4% 26.16% 19.02% 24.31% 25.42%

Los resultados indicaron una alta variabilidad de predicción, en la cantidad de

secuencias alfa-hélice, así como también en lo que respecta a los
verdaderos positivos. Esto indicó que la convergencia de los modelos se
localizó en distintos óptimos locales. Algunos de éstos, más eficientes para
predecir residuos de hélice o verdaderos positivos y otros mas eficientes
para predecir residuos no-hélice o verdaderos negativos.

Se detectó adicionalmente un patrón en los resultados para algunas

predicciones, en las cuales ninguno de los modelos realizó predicción alguna
de una secuencia de alpha-hélice; éste fue el caso para las secuencias de
prueba 1JFP, 1E12 y 1UAZ. En particular, las secuencias de predicción se
caracterizaron por permanecer en un estado de bucle interior “I”, esto es,
atrapado en un estado; este comportamiento indicó que el método cuenta
con un “sesgo” o tendencia natural de converger en cierta zona del espacio
de búsqueda. Por último, es notoria la efectividad del método para algunas
secuencias, lo que evidencia la potencia que aporta el modelo de estados,
haciendo a éste responsable de buena parte de la efectividad del método de
predicción.

39
Como conclusión, los resultados de las pruebas fueron buenos, esto con
base a la expectativa que se tenia inicialmente para el umbral utilizado. Sin
embargo, estos no lo fueron en promedio, ya que se contó con algunos
casos en donde la predicción sencillamente fracasó totalmente como se
explicó anteriormente; por tanto, es evidente la importancia de contar con un
buen punto de partida inicial de optimización con el fin de evadir las
debilidades planteadas; de esta manera, generar una sinergia mayor entre el
modelo de estados y los parámetros.

En la realización de las pruebas fue necesaria la utilización de una librería

para extender la precisión de punto flotante, Dfp [11]. Esto debido
principalmente a que los cálculos realizados por los algoritmos Forward-
Backward y Viterbi utilizan valores muy pequeños, provocando asi que los
resultados decaigan a 0, con una tasa exponencial relacionada directamente
al tamaño de la secuencia, superando el límite para el tipo de datos (double)
de 64bits. Esto generaba operaciones indeterminadas como el caso 0/0.

3.2 Conclusiones sobre el método HMMTOP

3.2.1 Fortalezas

Como fortalezas se destaca el modelo de estados. Con base en los

resultados se observó la gran influencia que éste tiene en la efectividad del
método, ya que estas restricciones limitan de manera muy efectiva el número
posible de secuencias de predicción que el algoritmo de Viterbi puede
realizar. Adicionalmente de este conjunto posible de predicciones, la gran
mayoría son consistentes con las estructuras determinadas
experimentalmente, esto reduce notoriamente el espacio de búsqueda y
adicionalmente incrementa la probabilidad de convergencia a parámetros
cercanos al óptimo global.

Los resultados de los autores de HMMTOP, demuestran la efectividad de

usar el método de pseudo conteo para incrementar la probabilidad de
convergencia en óptimos locales, posicionando el modelo inicial en la
periferia de estos óptimos y en consecuencia evitando al algoritmo Forward-
Backward quedar atrapado en un valor local no óptimo. Los resultados de la
reproducción del experimento confirmaron esto, ya que se evidenció una
tendencia a converger en distintos máximos locales no óptimos, lo que indicó
la existencia de un número elevado de éstos.

Finalmente los resultados de los autores demuestran la gran efectividad del

Método de predicción, el cual es el resultado de una buena sinergia entre el
modelo de estados y el método de pseudo conteo.

40
3.2.2 Limitaciones

Como limitaciones del método, tenemos el problema de ajustarlo a sub-

problemas de predicción; esto se debe a que la responsabilidad en la eficacia
del método, no recae completamente en las manos del modelo de Markov;
por el contrario, está asociada a un gran trabajo y cuidado por parte de los
autores en la selección y afinación de los parámetros iniciales y secuencias
de entrenamiento. A causa de ello, la tarea de entrenar un nuevo modelo con
el propósito de realizar predicciones en un subconjunto de proteínas,
requiere de un gran esfuerzo y tiempo, lo que limita su aplicabilidad práctica.

Otra limitante bastante fuerte recae en el limitado tamaño de ventana del

Método, reduciendo la cantidad de secuencias en las cuales puede realizar
predicciones con un buen grado de precisión. En buena parte este
inconveniente es derivado de la demandante precisión de punto flotante que
los algoritmos Forward-Backward y Viterbi requieren para realizar una
optimización y predicción confiable, lo cual se evidenció durante la
implementación y pruebas realizadas al reproducir el experimento, que forzó
la utilización de una librería especial Dfp [11], para incrementar esta
precisión. Como consecuencia, los requerimientos temporales y espaciales
del método se incrementaron de manera notoria.

El tamaño del modelo de estados también representa una debilidad a la hora

de realizar la optimización, debido al orden temporal de cada iteración en el
algoritmo Forward-Backward, el cual es O(n 2 ) con respecto al número de
estados y O(n) con respecto al tamaño de la secuencia. Para el caso de la
reproducción del experimento que contaba con un modelo de 60 estados, los
requerimientos temporales fueron importantes.

41
4. Propuesta de mejoramiento a HMMTOP

Con base en el análisis del método utilizado por HMMTOP se determinó que
el factor que mas afectaba su desempeño está relacionado con la limitante
del algoritmo Forward-Backward para encontrar el máximo global. Aunque el
método de pseudo-conteo logra conseguir un buen punto de partida para el
algoritmo, esto tampoco garantiza la convergencia al máximo global; en sí,
solo incrementa la probabilidad que el algoritmo logre convergencia en un
máximo global relativamente óptimo. Esto es, como lo describen los autores,
una reducción de los óptimos locales deficientes. Por este motivo, se decidió
que la propuesta de mejoramiento se enfocaría en proveer una mejor
alternativa para esta limitación y para ello se revisaron otras propuestas
alternativas como las redes neuronales, redes Bayesianas y los algoritmos
genéticos. Como resultado se determinó que los algoritmos genéticos
proveían de una poderosa alternativa para este caso particular, dada su
capacidad para explorar distintas regiones del espacio de búsqueda de
manera simultánea, característica ésta de la cual carecen el método
Forward-Backward y en general, las otras técnicas empleadas.

4.1 Modelo del algoritmo genético para el problema

El diseño del algoritmo genético está ligado a la toma de varias decisiones

con respecto a sus parámetros, tales como el modelo de la población, el
método de selección, cruzamiento, mutación y el diseño de la función de
aptitud o calificación de la población. Las decisiones correspondientes a
estos parámetros se presentan a continuación.

Como punto importante se tomó la dedición de reducir el modelo de estados

para el modelo oculto de Markov, pasando de los 60 estados del modelo
inicial a 24, mediante la eliminación de las restricciones MAXL en los estados
FL (figura 4.1), de manera que el modelo aun fuese válido para representar
la estructura secundaria de la proteína.

El motivo de esta reducción se fundamentó en los recursos que requeriría el

algoritmo genético para efectuar la optimización del problema, ya que implica
mas recursos espaciales para la representación de los individuos y la
población; un tamaño de la población mayor y finalmente un numero de
generaciones muy elevado y por ende, un incremento del tiempo de

42
cómputo, para obtener resultados comparables con la propuesta inicial. Todo
esto, finalmente hace inviable la utilización del método en un equipo Desktop
convencional.

X
I

o(1)
x
MINL

H(17)
MINL

i(1)
MINL

H(2) H(1)

Modelo de estados y transiciones para el Modelo de Markov de HMMTOP-

GA.
(Figura 4.1)

43
 Esquema de los 5 pasos del algoritmo genético.
(Figura 4.2)

En general, los 5 pasos (Figura 4.2) que sigue la implementación del método
son los siguientes:

• Primero, generar una población inicial de manera aleatoria.

• Segundo, utilizar cada uno de los individuos en la población, como
modelo oculto para el algoritmo Viterbi. Utilizando este modelo, se
realiza la predicción a todo el conjunto de prueba.
• Tercero, calificar todas las secuencias predichas contra las
determinadas experimentalmente y consolidar el resultado para el
individuo particular.
• Cuarto, seleccionar los individuos que sobrevivirán a la siguiente
generación; la selección se realiza acorde con la calificación del
grupo, en este caso los individuos con la menor puntuación.

44
 • El quinto paso, consiste en cruzar y aplicar mutaciones a los
cromosomas de los distintos individuos seleccionados de esta
manera, generando en algunos casos, nuevos individuos.
• Por ultimo se completa la población mediante el cruce adicional entre
los individuos seleccionados y esto se realiza de manera aleatoria.

Los pasos del 2 al 5 se repiten de manera indefinida hasta que un individuo

en la población, alcance la puntuación umbral que es definida al inicio del
método.

4.1.1 Modelaje de la población

Los individuos de la población fueron modelados con base a los parámetros

del modelo oculto de Markov de 24 estados, Como primera aclaración, en
este contexto el término “cromosoma” no es sinónimo de “individuo”, como si
lo es bajo el modelo canónico de algoritmos genéticos. Esto debido al orden
(tamaño) de la representación de los individuos para este problema
particular. Por el contrario y en analogía con la organización de la
información genética en organismos más complejos, el genoma se compone
de un conjunto de cromosomas. Con ello tenemos cromosomas que
representan al vector π y las filas de la matriz de transición y difusión; la
razón de esta elección se debe a la propiedad estocástica de estos
cromosomas, ya que la suma de sus componentes siempre es igual a 1 y
cada individuo en la población cuenta con 49 cromosomas.

Un aspecto importante a decidir sobre el modelo de la población, es el

método de generación aleatoria de cromosomas, ya que debido a la
propiedad invariante que estos cumplen y al hecho que el espacio de
búsqueda es continuo y por lo tanto existe un número infinito de posibles
cromosomas, fue necesario utilizar un método para subdividir de manera
discreta el espacio de búsqueda y de asi obtener un conjunto finito de
posibles cromosomas. La definición del método utilizado para este propósito
se presenta a continuación.

Definición 4.1: Sea C el conjunto de los 49 cromosomas (vectores)

genéricos del genoma para los individuos en una población. Sea k ∈ C un
cromosoma genérico particular y i = 1,2,.., k , el i-esimo componente del
cromosoma. Definimos a R(k , i ) : C × Ν → {0,1} , como la función que retorna si
el i-esimo termino de k cuenta con una restricción o no.

Las restricciones hacen referencia a las restricciones del modelo de estados,

el caso R (k , i ) = 0 , indica que ki tiene un valor de 0 estático.

45
Definición 4.2: Sea χ , α ∈ Ν + , el coeficiente de discretización y coeficiente
de dispersión respectivamente. Sea rand ( x) ∈{1,2,..., x} , una función aleatoria
con distribución constante. Definimos como la función de generación a
y = rand ( χ )α R(π , i ) .

Utilizando como base la función de generación, generamos los valores

iniciales de los parámetros para los 49 cromosomas del individuo y luego
realizamos una normalización para cumplir con la condición estocástica de
estos de la siguiente manera:

Generación:
π i = rand ( χ )α R(π , i )
Ai , j = rand ( χ )α R( Ai , j )
Bi , j = rand ( χ )α R( Bi , j )

Normalización:
1
πi = πi
∑π i
1
Ai , j = Ai , j
∑ Ai , j
j

1
Bi , j = Bi , j
∑ Bi, j
j

Definición 4.3: Sea Υ j , la representación de las restricciones del j-esimo

cromosoma, c un cromosoma y U (c, Υ j ) una función boleana valida cuando
c cumple las restricciones Υ j . Definimos el conjunto de j-cromosomas K j ,
como el conjunto de los cromosomas tales que (∀c |: U (c, Υ j )) .

El coeficiente de discretización establece el tamaño del conjunto de

cromosomas posibles; de esta manera, en una dimensión del espacio de
búsqueda, el coeficiente de dispersión asegura que los individuos generados
se dispersen de la manera mas uniforme por dicho espacio. A un menor valor
de coeficiente, los individuos se generarán con mayor tendencia a una

46
media. Finalmente, este método busca obtener una discretización uniforme
del espacio sin generar espacios no válidos (individuos que no cumplan con
la propiedad estocástica del modelo de Markov).

4.1.2 Función de calificación

La función se definió como una medida de la diferencia entre dos secuencias

estructurales; esta diferencia solo toma en cuenta las diferencias de
estructura para cada posición de la secuencia. El resultado final es una
densidad normalizada de diferencia, por residuo en la secuencia. Debido a
que en este caso, la prioridad es la adecuada identificación de las alfa-
hélices, las diferencias entre la estructura secundaria experimental y la
predicha, cuando alguno de estos representa una sección de alfa-hélice,
tienen un mayor peso. La (tabla 4.1) lista los diferentes valores de diferencia.

(Tabla 4.1)
Valores para cuantificación de diferencia en la función de calificación.

Estructura a Estructura b Valor de diferencia

S S 0.0
( I | O) ( O | I ) 0.01
( i | o ) ( o | i ) 0.01
( I | i ) ( i | I ) 0.0025
( O | o ) ( o | O ) 0.0025
( I | o ) ( o | I ) 0.04
( O | i ) ( i | O ) 0.04
H ( o | i ) 0.6
( o | i ) H 0.6
H ( O | I ) 1.0
( O | I ) H 1.0

Es importante hacer la aclaración que existe gran diversidad de propuestas

posibles y metodologías válidas para diseñar la función de calificación. Sin
embargo en común, éstas requieren de una afinación de sus parámetros
óptimos lo cual es muy costoso. Existe un método propuesto que promete
eliminar estos problemas, el cual se fundamenta en la utilización de una
norma universal de información basada en el concepto de complejidad de
Kolmogorov [12].

47
4.1.3 Método de selección, cruce y mutación

La principal dificultad con los métodos de selección y cruce es la pérdida de

diversidad en la población de individuos al paso de las generaciones lo que
normalmente genera en el algoritmo genético una reducción sustancial de su
poder para explorar el espacio de búsqueda y por tanto, se incrementa la
probabilidad de convergencia en un óptimo local (convergencia prematura),
esto explicado mediante el teorema del esquema. El muestreo de los
hiperplanos decrece hasta el punto en que todos los individuos muestrean el
mismo hiperplano de la misma manera, eliminando la propiedad emergente
del paralelismo implícito. Cuando adicionalmente se toma en cuenta la
existencia de múltiples óptimos locales, el problema se torna mas crítico.
Debido a esto, se propone en el presente trabajo un método para solventar
ambas dificultades que denomino como “Método de diversidad controlada”
con el propósito de evitar esta dificultad.

4.1.4 Introducción al Método de diversidad controlada

La idea básica del método radica en controlar la pérdida de la diversidad en

la población, mediante la realización de una selección mínima, esto es,
eliminando únicamente al peor individuo de la población. Se realizaran
cruces entre dos individuos de tal forma que se conserve una buena
proporción de sus cromosomas originales y posteriormente inyectando en la
población nuevos cromosomas, mediante mutación a una tasa similar en la
que pierde diversidad, debido al cruzamiento y selección. Lo anterior, busca
en las primeras generaciones maximizar la exploración del espacio de
búsqueda, utilizando poblaciones relativamente pequeñas en comparación
con éste; en generaciones subsiguientes, permite que los mejores individuos
tengan mayores probabilidades de sobrevivir en la población y transferir sus
características útiles a otros, para que la población converja de una manera
bastante estable.

4.1.5 Operador de Selección

Como se describió anteriormente, el operador de selección busca realizar

una selección mínima en la población, y solo se pierde una cantidad pequeña
de cromosomas únicos. Los individuos que cuentan con buenas
características, permanecerán una mayor cantidad de generaciones en la
población y contarán con una mayor probabilidad de ser mejorados. De lo
anterior surge, la necesidad de modelar los efectos del operador de selección
en la diversidad de la población al paso de las generaciones.

48
Definición 4.4: Sea n el tamaño de la población de individuos G . Sea
ϕT = 49n , el total de cromosomas presentes en G y ϕ D el número de
cromosomas distintos presentes en G . Definimos el factor de diversidad
ρ = ϕ D ϕT , como la probabilidad de seleccionar un cromosoma θ en G , tal
que θ es único.

El factor de diversidad ρ , también se puede interpretar como la densidad de

cromosomas distintos en la población. En este punto es importante
plantearse la pregunta: ¿Cuál es el efecto del operador de selección en la
diversidad de la población?. El teorema del esquema nos indica, que la
selección tiene un efecto importante en el muestreo de un hiperplano dado,
lo se puede interpretar como una pérdida en la diversidad de la población. La
aplicación del operador de selección y el factor de diversidad están
relacionados de la siguiente manera:

Definición 4.5: Sea n y ρi el tamaño de la población y el factor de

diversidad en la generación i . El valor esperado del factor de diversidad para
la población intermedia ρ I , luego de aplicar el operador de selección es
representado por.

n −1
ρ I = ρi
n

Nótese que el efecto del operador de selección en la diversidad disminuye a

medida que la diversidad es menor, esto debido a que un número mayor de
los cromosomas presentes en el individuo no-seleccionado, se encuentran
replicados en la población a medida que la diversidad disminuye. Ahora bien:
si con el propósito de compensar la pérdida de un individuo y por tanto la
disminución del tamaño de la población, seleccionamos aleatoriamente a un
individuo de la población intermedia y lo clonamos; entonces el efecto en la
diversidad, dada la aplicación sistemática del operador de selección en cada
generación, se define como:

Definición 4.6: Sea ρ0 el factor de diversidad en la generación inicial e i la

generación de la población. El valor esperado del factor de diversidad para la
generación i , aplicando el operador de selección esta dado por la siguiente
función recursiva ρ (i ) :

ρ (0) = ρ 0

49
 n −1
ρ (i) = ρ (i − 1)
n

4.1.6 Operador de cruce

El método de cruce busca modificar la información de un individuo en la

población, agregando características de otro. El propósito de ésto radica en
obtener un mecanismo para distribuir las características útiles presentes en
la población. La siguiente definición describe al operador de cruce:

Definición 4.7: Sean λP ,i , λD , k ∈ G , el individuo principal, el individuo donante

e i como la generación actual de la población. Sea j el índice del j-esimo
cromosoma de un individuo, λ ( j ),1 ≤ j ≤ 49 el j-esimo cromosoma del
individuo y X una variable aleatoria de distribución Bernoulli con
probabilidad p . Definimos Cross ( X , j , λP ,i , λD , k ) , como la función que
determina el j-esimo cromosoma del individuo objetivo λO ,i +1 :

 X = 0, λP ,i ( j ) 
Cross ( X , j , λP ,i , λD , k ) =  
 X = 1, λD , k ( j ) 

λO ,i ( j ) = Cross ( X , j, λP ,i , λD , k )

Un punto importante sobre este operador radica en la pérdida de una parte

de la diversidad cuando X = 0 ; en este evento, el j-esimo cromosoma del
donante se replica en el individuo objetivo, y por tanto, tenemos la pérdida de
una copia del el j-esimo cromosoma (esta copia puede ser la única) del
principal, para la siguiente generación. Como los eventos X j = 0 son
independientes e idénticamente distribuidos, tenemos que la probabilidad de
ocurrencia de k cruces en la población, en una generación, sigue una
distribución binomial, cuyo valor esperado corresponde a 49np . Teniendo en
cuenta una densidad de diversidad ρ x , el valor esperado de pérdida de
cromosomas distintos ésta dado por 49npρ x . Replanteando la expresión para
el factor de diversidad (Definición 4.6) de manera que incluya el efecto del
operador de cruce, tenemos:

50
Definición 4.8: Sea ρ0 el factor de diversidad en la generación inicial, p la
probabilidad de cruce, 49np el valor esperado de cruces en la población. El
valor esperado del factor de diversidad para la generación i , luego de aplicar
el operador de selección y cruce esta dado por la siguiente función recursiva
ρ (i) :

ρ (0) = ρ 0
(n − 1) 49(n − 1) p  n − 1 − (n − 1) p 
ρ (i) = ρ (i − 1) − ρ (i − 1) = ρ (i − 1) 
n 49n  n 

Si notamos el efecto del cruce en la población de manera análoga a la

disrupciones generadas por el operador de cruce en el teorema del
esquema, tenemos una pérdida de diversidad, en consecuencia, a medida
que un individuo sobrevive un mayor número de generaciones, los demás
individuos en la población adquieren eventualmente más características de
éste. Desde la perspectiva del teorema del esquema, la población muestrea
un conjunto de hiperplanos mas reducido y de mayor orden.

4.1.7 Método de mutación

El método de mutación sigue un esquema muy similar al cruce de

cromosomas. La diferencia radica en el cromosoma donante, el cual en esta
oportunidad es tomado del conjunto de cromosomas posibles de manera
aleatoria. Adicionalmente, la probabilidad de mutación es dependiente del
evento de cruce; si no se realiza un cruce, entonces se tiene una
probabilidad de realizar una mutación.

Definición 4.9: Sea λP ,i ∈ G , el individuo principal e i la generación actual de

la población. Sea j el índice del j-esimo cromosoma de un individuo,
λ ( j ),1 ≤ j ≤ 49 el j-esimo cromosoma de un individuo, Rcrom( j ) = c ∧ c ∈ K j
una función que selecciona un cromosoma de manera aleatoria en K j y X
una variable aleatoria de distribución Bernoulli con probabilidad q . Definimos
Mutate( X , j , λP ,i ) , como la función que determina el j-esimo cromosoma del
individuo objetivo λO ,i +1 :

51
  X = 0, λP ,i ( j ) 
Mutate( X , j ) =  
 X = 1, Rcrom( j )

λO ,i ( j ) = Mutate( X , j )

La operación de mutación contrario a los operadores de selección y cruce,

agrega diversidad a la población (la probabilidad de agregar un cromosoma
de manera aleatoria tal que éste se encuentre en la población, es muy
cercana a 0). Adicionalmente, debido a que los eventos de mutación son
independientes entre si e idénticamente distribuidos y la probabilidad de
mutación es dependiente del evento de cruce (solo se muta si no se realiza
el cruce), tenemos con base en la vease (Definición 4.8).

Definición 4.10: Sea ρ0 el factor de diversidad en la generación inicial, q la

probabilidad de cruce, q la probabilidad de mutación y 49n(1 − p )q el valor
esperado de mutaciones en la población. El valor esperado del factor de
diversidad para la generación i , luego de aplicar el operador de selección,
cruce y mutación esta dado por la siguiente función recursiva ρ (i ) :

ρ (0) = ρ 0
49(n − 1) 49( n − 1) p  49(n − 1) 
ρ (i ) = ρ (i − 1) − ρ (i − 1) + 1 − ρ (i − 1) (1 − p) q
49n 49n  49n 
 (n − 1) − (n − 1) p − nq(1 − p ) 
= ρ (i − 1)  + q (1 − p )
 n 

La anterior, no es más que una ecuación de diferencia no homogénea y su

solución nos provee del comportamiento dinámico, para la diversidad de la
población con respecto a las generaciones, la cual esta dada por:

Definición 4.11: Dado u = ((n − 1) − (n − 1) p − nq (1 − p )) n y z = q (1 − p ) ,

Definimos la función de densidad de diversidad para la generación i como:

 ρ u + z − ρ 0   ρ 0u + z − ρ 0  i
ρ (i ) =  ρ 0 − 0 + (u )
 (1 + u )   (1 + u ) 

52
Notar que u < 1 y por tanto (u )i → 0 cuando i → ∞ . De este resultado
tenemos:

Definición 4.12: Sea p(i) la función de densidad de diversidad, definimos la

densidad de diversidad estable ρ ∞ como:

 ρ 0u + z − ρ0 
lim ρ (i) =  ρ
i →∞
0 −
(1 + u ) 
 = ρ∞

La densidad de diversidad estable, finalmente provee el fundamento principal

del método de diversidad controlada, ya que solo requerimos definir los
valores de n, p y q para obtener una densidad de diversidad deseada durante
la ejecución del algoritmo genético (dado un numero suficientemente grande
de generaciones); de esta manera, logramos una diversidad consistente con
el tamaño del espacio de búsqueda y el umbral de la función de calificación
que se desee obtener. Esto también permite controlar la supervivencia media
de un cromosoma, lo cual es importante cuando en la población se cuenta
con más de un individuo con buena puntuación, pero pueden encontrarse
también, en regiones del espacio de búsqueda pertenecientes a diferentes
óptimos locales; el permitir a estos sobrevivir un numero adecuado de
generaciones, incrementa la probabilidad de finalmente escoger al candidato
mas adecuado (figura 4.3).

4.2. Elección de los parámetros del algoritmo genético

Uno de los inconvenientes más grandes al utilizar el método de algoritmos

genéticos recae en la selección de los parámetros óptimos, tales, como el
tamaño de la población y las probabilidades de cruce y mutación. La
selección de estos parámetros depende principalmente de factores como el
tamaño del espacio de búsqueda, la función de calificación, la complejidad
de la información usada en la función de calificación y el umbral de
optimización deseado o las generaciones a realizar.

De los anteriores factores, la complejidad de la información es la más

importante, debido a que todos los otros factores finalmente dependen de
este. El problema principal de lo anterior es el como cuantificar esta
complejidad y cual es su relación con los demás parámetros. Para ello
tomaremos una aproximación conceptual “a la Shannon” [13], esto es,
basada en la entropía de la información.

53
 (Figura 4.3)

La entropía en este caso particular, corresponde a la información del

conjunto de aprendizaje y la transformación que realiza la función de
calificación a la información. El efecto general del incremento en la entropía,
consiste en una dificultad cada vez mayor en la convergencia de la población
a un óptimo, en otras palabras, gradientes más fuertes. Para ilustrar esto
definiremos inicialmente el concepto de región de calificación:

Definición 4.13: Sea Ξ , el espacio de búsqueda del algoritmo genético,

ε ∈ R[ 0,1] un valor umbral de calificación y fΞ (λ ) una función de aptitud o
calificación tal que f Ξ (λ ) ∈ R[ 0,1] ∧ λ ∈ Ξ . Definimos Wε = (λ | λ ∈ Ξ ∧ f Ξ (λ ) ≤ ε )
como la región de calificación para el umbral ε .

Como se anotó antes, al contar con mayor entropía en la información, la

convergencia de la población es cada vez mas lenta, de esto tenemos que el
tamaño de una región de calificación arbitraria Wε está ligada a la entropía,
esto se ilustra en la figura (Figura 4.4) donde se presenta el tamaño de
distintas regiones Wε con respecto a la entropía del conjunto de aprendizaje.

54
 (Figura 4.4)

Desde el punto de vista de los algoritmos genéticos, se tiene que la

probabilidad de contar con un individuo en la población que se encuentre en
un área del espacio de búsqueda con una calificación inferior a ε , decrece al
aumentar la entropía. Al no contar con una población de un tamaño
adecuado, se corre el riesgo de no lograr un umbral deseado de calificación
(Figura 4.5). Finalmente, definimos la relación entre la región Wε y la
entropía de la información como:

Definición 4.14: Sea Wε la región de calificación para el umbral ε en un

espacio de búsqueda Ξ continuo y ξ la entropía en bits del conjunto de
aprendizaje. Definimos la proporción del espacio de búsqueda
correspondiente a la región de calificación Wε , como:

2ξε 1 W
S (ε , ξ ) = ξ
= ξ (1−ε ) = ε
2 2 Ξ

55
 (Figura 4.5)

4.2.1 Factor de Discretización del espacio de búsqueda

Con respecto al factor de discretización χ , se tiene un argumento similar al

anterior, si la subdivisión del espacio de búsqueda continuo es reducido, la
probabilidad de contar con un individuo en la población que se encuentre en
un área del espacio de búsqueda con una calificación inferior a un valor ε ε

decrece e incluso, puede no existir. Por el contrario, si éste se incrementa, se

cuenta con un mayor número de individuos posibles en el área. Sin embargo,
esto a costa de un incremento del espacio de búsqueda discretizado.

56
Definición 4.15: Sea χ el factor de discretización del espacio de búsqueda
Ξ . Definimos el tamaño del espacio de búsqueda discretizado como
κ = ( χ + 1) 49 .

El valor esperado de individuos contenidos en la región del espacio de

búsqueda, cuyo umbral es menor a χ , se define como:

Definición 4.16: Sea ξ la entropía de la información, ε un umbral de

calificación, Wε la región de calificación y κ el tamaño del espacio de
búsqueda discretizado. Definimos el valor esperado de individuos que
pertenecen a Wε para el espacio discretizado Ξ χ como:

Κ (κ , ε , ξ ) = S (ε , ξ )κ

El valor de Κ (κ , ε , ξ ) a utilizar está ligado a la distribución en el espacio de

búsqueda de la región de calificación Wε y por tanto, dependiente de la
función de calificación utilizada. La recomendación es utilizar un valor alto (
Κ (κ , ε , ξ ) > 10100 ), con el propósito de cubrir de manera adecuada los óptimos.
Con el valor de Κ (κ , ε , ξ ) , determinamos el valor para χ , utilizando:

 1

  Κ (κ , ε , ξ )  49
χ =   − 1
 S (ε , ξ )  
 

4.2.2 Aproximación para S (ε , ξ )

Debido a que calcular S (ε , ξ ) , no es trivial, se puede aproximar tomando una

muestra grande de individuos m en el espacio de búsqueda y realizar una
regresión lineal, con base en el conteo qk de individuos en la muestra con
calificaciones superiores a varios valores ek . Con el cálculo de un conjunto
de valores Qk = qk m , obtenemos la entropía experimental ξ R linealizando
con la ecuación:

− log 2 Qe = ξ R (1 − ε e )

57
Finalmente, con el valor ξ R obtenido mediante regresión se tiene:

2ξ R ε
S (ε , ξ ) ≈ ξ R
2

4.2.3 Probabilidades para la mutación y cruce

Existe un gran número de opciones a escoger con respecto a las

probabilidades de los operadores de mutación y cruce. La definición 4.12
nos proporciona un primer parámetro en esta selección; sin embargo, no
contamos con argumentos respecto al nivel de diversidad mas adecuado a
seleccionar. Por este motivo, a continuación indagaremos con mayor
profundidad las propiedades de las operaciones de mutación y cruce.

4.2.4 Efecto de la Mutación

Ahora nos interesa conocer el efecto de la mutación sobre un individuo, esto

es, el comportamiento del individuo con respecto al espacio de búsqueda. En
términos generales una mutación genera un desplazamiento del individuo en
el espacio de búsqueda; qué tan abrupto sea este movimiento, estará dado
por el factor de discretización. Lo primero a tener en cuenta es la naturaleza
de este desplazamiento; cuando ocurre una mutación, uno de sus
cromosomas es reemplazado por uno aleatorio. Si pensamos en el espacio
de búsqueda de manera análoga a la teoría del esquema, como un
hipercubo de m dimensiones y el valor para cada dimensión es representado
por un cromosoma, entonces el operador mutación genera un
desplazamiento aleatorio en la dimensión correspondiente al cromosoma
modificado. Ahora, si el individuo se encuentra en una región, de calificación
Wε k antes de la mutación, nos interesa conocer cual es la probabilidad de que
el individuo permanezca en la región luego de sufrir una mutación. Utilizando
el concepto de conteo de eventos, tenemos que esta probabilidad
corresponderá a los posibles sitios para moverse, dentro de la región
definida, frente a los posibles sitios fuera de ésta, pero teniendo en cuenta
que el desplazamiento, es en una sola dimensión:

Definición 4.17: Sea Wε k una región de calificación, m el numero de

cromosomas de los individuos o las dimensiones de Ξ y S (ε , ξ ) = Wε Ξ la
proporción de la región de calificación. Definimos la proporción
unidimensional media para la región Wε k como:

58
 1 1
Su (Wε ) = S (ε , ξ ) m −1
= S (ε , ξ ) 48

Un aspecto interesante de Su (Wε ) radica en que es el valor esperado de la

proporción que cubre Wε para cualquier dimensión (cromosoma), esto es
importante ya que no tenemos idea de cual es la forma de la región Wε . Con
base en la definición anterior tenemos:

Definición 4.18: Sea Wε una región de calificación y λ ∈Wε un individuo en

la región de calificación. La probabilidad de que el individuo λ permanezca
en Wε luego de una mutación o probabilidad de permanencia esta dada por:

1
PM ,ε (λ '∈Wε | λ ∈Wε ) = Su (Wε ) = S (ε , ξ ) 48

De manera análoga a la definición anterior tenemos:

Definición 4.19: Sea in : λ ∈Wε y out : λ ∉Wε estados que definen el atributo
de pertenencia a Wε de λ . Definimos la probabilidad de no permanencia
PM ,ε (in | out ) , la probabilidad de retorno PM ,ε (out | in) y la probabilidad de
permanecía fuera de Wε o PM ,ε (out | out ) como:

1
PM ,ε (in | out ) = 1 − S (ε , ξ ) 48

1
PM ,ε (out | in) = 1 − (1 − S (ε , ξ )) 48
1
PM ,ε (out | out ) = (1 − S (ε , ξ )) 48

Si contamos con un individuo en Wε y le aplicamos un numero k de

mutaciones, el individuo realizará una serie de transiciones entre los estados
in y out , con probabilidades de transición dadas en definición 4.18 y 4.19.
El comportamiento de este sistema estocástico se describe mediante una
cadena de Markov [14] y su análisis nos provee del comportamiento del
individuo bajo un número arbitrario de mutaciones. Nuestro interés sin

59
embargo, radica en conocer la probabilidad de permanencia luego de aplicar
k mutaciones. Con esto en mente, si suponemos valores ε pequeños
tenemos que P(out | in) ≈ 0 , lo que indica que luego que una mutación
desplace a un individuo fuera de la región Wε , es muy improbable que otra
mutación lo retorne a esta. En base a este supuesto definimos la
probabilidad de permanencia luego de k mutaciones:

Definición 4.20: Sea Wε una región de calificación y λ ∈Wε un individuo en

dicha región. Definimos probabilidad de permanencia para k mutaciones,
como la probabilidad de que el individuo λ permanezca en Wε luego de la
aplicación de k mutaciones, la cual es aproximada por:

k
PM ,ε (in | in, k ) ≈ Pε (in | in) = S (ε , ξ )
k 48

Por ultimo bajo el modelo de convergencia controlada, para una probabilidad

de mutación q, y de cruce p, tenemos un valor esperado de mutaciones por
generación en un individuo de d mut = 49(1 − p )q . Con base en lo anterior y la
definición 4.20 tenemos:

d mut 49 (1− p ) q
PM ,ε (in | in, k ) ≈ S (ε , ξ ) 48
= S (ε , ξ ) 48

4.2.5 Mutación y convergencia

Ahora buscamos una relación para el desplazamiento de un individuo bajo el

efecto de d mut mutaciones y cómo esto se relaciona con la convergencia;
para esto, definimos la probabilidad de que un individuo en Wε , ingrese a una
región objetivo Wo :

Definición 4.21: Sea PM ,ε a (in | in, d mut ) la probabilidad de permanencia para λ

y ε a > ε o el umbral de calificación para Wε a y Wε o respectivamente. La
probabilidad de ingreso a Wε o para λ por efectos de mutación en una
generación, dado que λ ∈Wε a y Wε o ⊆ Wε a es:

60
 S (ε o , ξ )
PM ,i (λi ∈Wε a | λi +1 ∈Wε o , d mut ) = PM ,ε a (in | in, d mut )
S (ε a , ξ )
1
PM ,i (λi ∈Wε a | λi +1 ∈Wε o , d mut ) = ξ ( ∆ε ) PM ,ε a (in | in, d mut )
2

La última expresión se puede interpretar como la probabilidad de que un

individuo avance ∆ε = ε a − ε o en su calificación, estando en la región Wε a del
espacio de búsqueda, es decir, que converja.

4.2.6 Efecto del operador cruce

De manera análoga al operador de mutación, nos interesa conocer el efecto

del operador de cruce o intercambio de cromosomas sobre un individuo dado
λ . Para ello inicialmente definimos la región que contiene a la población G .

Definición 4.22: Sea G la población de individuos en el algoritmo genético y

n su tamaño.
Se define la región de población Wε G , como la mínima región del espacio de
búsqueda que contiene a todos los individuos en G y cumple la condición:

(∀λ ∈ G |: λ ∈Wε G ) ∧ (∀ε k | 0 < ε k ≤ 1 : (∀λ ∈ G |: λ ∈Wε k ) → ε G ≤ ε k )

Al igual que en el caso de una mutación, el cruce de cromosomas produce

un desplazamiento del individuo en el espacio de búsqueda y la diferencia
radica en la probabilidad de salir de la región Wε G , la cual es notoriamente
inferior en este caso, comparado con la mutación. Debido a que existe una
cantidad inferior de cromosomas posibles en los individuos contenidos en
Wε G y por tanto una menor cantidad de cromosomas que desplacen al
individuo fuera de Wε G , la probabilidad esta relacionada con la densidad de
diversidad ρ (i ) en la población G :

Definición 4.23: Sea Wε G la región de calificación mínima para la población

G , ρ (i ) la densidad de diversidad de la población en la generación actual y

61
 λ ∈Wε G un individuo en la región de calificación. La probabilidad de que el
individuo λ permanezca en Wε G luego de un cruce esta dada por:

1
 S (ε G , ξ )(1 − ρ (i ))  48 1
PC (in | in) ≈   = (1 − ρ (i )) 48
 S (ε G , ξ ) 

Tomando como supuesto que la población ha sido sometida a un número

importante de generaciones, aproximamos ρ (i) ≈ ρ∞ (Definición 4.12)
obteniendo:

1
 S (ε G , ξ )(1 − ρ∞ )  48 1

PC (in | in) ≈  
 = (1 − ρ∞ )
48

 S (ε G , ξ ) 

De manera análoga en el análisis del operador mutación, tenemos, la

probabilidad de no permanencia PC (in | out ) , la probabilidad de retorno
PC (out | in) y la probabilidad de permanecía fuera de Wε G o PC (out | out ) :

1
PC (in | out ) = 1 − (1 − ρ ∞ ) 48

1
PC (out | out ) = ( ρ ∞ ) 48

1
PC (out | in) = 1 − ( ρ ∞ ) 48

La probabilidad de que un individuo contenido en Wε G permanezca en éste,

luego de d cross =  49 p  , cruces, también puede ser aproximada. Bajo el
supuesto PC (out | out ) ≈ 1 , tenemos PC (out | in) ≈ 0 , cuando la población ha
pasado por un numero importante de generaciones y obtenemos la siguiente
definición:

Definición 4.24: Sea Wε G una región mínima de calificación para la

población G , λ ∈Wε un individuo en la región de calificación. Definimos

62
probabilidad de permanencia para k = d cross =  49 p  cruces, como la
probabilidad de que el individuo λ permanezca en Wε G luego de la aplicación
de k cruces y es aproximada por:

dcross
k
PC (in | in, k ) ≈ PC (in | in) 48
= PC (in | in) 48

4.2.7 Cruce y convergencia

La relación entre el cruce y la convergencia es análoga al caso de la

mutación. Sin embargo se tiene un inconveniente: a menor diversidad, el
cubrimiento del espacio Wε G pierde uniformidad, lo que reduce la
probabilidad de que existan individuos en un área de calificación inferior
Wε o ⊆ Wε G . A diferencia del caso de la mutación. Con base en la definición
4.12, para valores de i suficientemente grandes, el cambio en la uniformidad
esta dado por:

ρ (i + 1)
∆u = ≈1
ρ (i )

Con esto último la probabilidad de que un individuo en Wε G , ingrese en el

área Wε o ⊆ Wε G , luego de la aplicación de dcross =  49 p  cruces se define
como:

Definición 4.25: Sea PC (in | in, d cross ) la probabilidad de permanencia para λ

y ε G > ε o el umbral de calificación para Wε G y Wε o respectivamente. La
probabilidad de ingreso a Wε o para λ , causada por efectos de dcross =  49 p 
cruces en una generación, dado que λ ∈Wε G y Wε o ⊆ Wε G es:

∆uS (ε o , ξ )
PC ,i (λi +1 ∈Wε o | λi ∈Wε G , dcross ) = PC (in | in, dcross )
S (ε G , ξ )
1
PC ,i (λi +1 ∈Wε o | λi ∈Wε G , d cross ) = ξ ( ∆ε ) PC (in | in, dcross )
2

63
 d cross
1
PC ,i (λi +1 ∈Wε o | λi ∈Wε G , d cross ) = (1 − ρ∞ ) 48
2ξ ( ∆ε )

4.2.8 Convergencia con Mutación y cruce

Como parte final, queremos modelar el efecto combinado de los operadores

de selección, mutación, cruce y la manera como estos se relacionan con la
convergencia de la población. Para esto tomaremos como base el avance en
la convergencia ∆ε para la población.

Definición 4.26: Sea G la población del algoritmo genético, m el numero de

individuos en G, Wε G la región mínima de calificación para la población en la
generación i, Wε 'G la región mínima de calificación para la generación i+1 y
∆ε = ε G − ε 'G el avance en la convergencia para la población. La probabilidad
de que la aplicación de los operadores cruce y mutación desplacen a un
individuo en la región Wε 'G como:

Ψi (∆ε ) = PC ,i (λi +1 ∈Wε o | λi ∈Wε G , d cross ) + PM ,i (λi +1 ∈Wε o | λi ∈Wε G , d Mut )

Desarrollando la expresión anterior tenemos:

d
1 cross
1
Ψi (∆ε ) = ξ ( ∆ε )
(1 − ρ ∞ ) 48
+ ξ ( ∆ε ) PM ,ε G (in | in, d mut )
2 2
1  
d cross d mut
Ψi (∆ε ) = ξ ( ∆ε )  (1 − ρ ∞ ) 48 + S (ε G , ξ ) 48 
2  
1  
49 (1− p ) 49 (1− p ) q
Ψi (∆ε ) = ξ ( ∆ε )  (1 − ρ ∞ ) 48 + S (ε G , ξ ) 48 
2  
 49 (1− p ) 49 (1− p ) q

α =  (1 − ρ ∞ ) 48 + S (ε G , ξ ) 48 
 
α
Ψi (∆ε ) = ξ ( ∆ε )
2

Tomando en cuenta el efecto del operador de selección para lograr un

avance ∆ε para toda la población, requerimos que exactamente m-1
individuos de la población ingresen a la región Wε 'G . Si interpretamos la

64
probabilidad Ψi (∆ε ) como la probabilidad de un evento positivo y 1 − Ψi (∆ε )
para uno negativo, tenemos que la probabilidad de que ingresen m-1
individuos a Wε 'G esta dada por la distribución binomial así:

 m  m −1
Φ i (∆ε ) =  Ψi (∆ε )(1 − Ψi (∆ε ))
 m − 1

Con base en lo anterior, el valor esperado de la distribución binomial y

nuestro valor esperado objetivo m-1, obtenemos la definición del avance de la
convergencia para la población:

Definición 4.27: Sea G la población del algoritmo genético, m el numero de

individuos en G, Wε G la región mínima de calificación para la población en la
generación i, Wε 'G la región mínima de calificación para la generación i+1 y
Φ i (∆ε ) la probabilidad de que ingresen m − 1 individuos a la región Wε 'G . El
avance esperado de convergencia para la generación i, es el valor ∆ε tal
que:

mΦ i (∆ε ) = m − 1

Utilizando la Definición 4.27 y despejando ∆ε obtenemos:

m −1  m  m −1
= Φ i (∆ε ) =  Ψi (∆ε )(1 − Ψi (∆ε ))
m  m − 1

m −1 m −1
= Φ i ( ∆ε ) = m(Ψi (∆ε )(1 − Ψi (∆ε )))
m
m −1 αm α m −1
m − 1 = (Ψi (∆ε )(1 − Ψi (∆ε ))) = ξ ( ∆ε ) − ξ ( ∆ε )
2 m 2 m
ξ ( ∆ε )
(α m − α m −1 )
=2 m
m −1

∆ε =
( )
log 2 α m − α m −1 − log 2 ( m − 1)
ξm

65
A parte de este resultado, contamos con una herramienta para evaluar el
valor esperado de generaciones que requiere el algoritmo genético para
alcanzar un umbral de calificación dado.

Definición 4.28: Sea G la población del algoritmo genético, m el numero de

individuos en G, ξ la entropía del conjunto de aprendizaje, ε F el umbral de
calificación objetivo para el problema y ∆ε u , v = ε u − ε v . Definimos la función
acumulativa de convergencia Λ(i ) :

Λ(0) = 1
Λ(i ) = Λ(i − 1) − ∆ε Λ ( i −1),v

Finalmente:

Definición 4.29: Sea Λ(i ) la función acumulativa de convergencia. El valor

esperado de generaciones requerido para alcanzar el umbral objetivo ε F , es
el valor I tal que:

(∀i | 0 < i : ε F ≤ Λ (i ) ∧ i ≥ I )

4.3 Parámetros seleccionados

Debido a que encontrar los parámetros m, p, q tales que el valor esperado de

generaciones I sea mínimo, es un problema de optimización, para el cual no
contamos con una solución, se determinaron los parámetros con base en
pruebas realizadas con el algoritmo. Se escogieron los siguientes
parámetros:

El valor de ξ , fue obtenido mediante regresión sobre los mejores 100

individuos en una población de tamaño m=1000 (Figura 4.6). De esto se
obtuvo ξ = 118,95 .

Se seleccionó un umbral de calificación objetivo ε F = 0,2 y un valor

K (κ , ε , ξ ) = 10118 . Obteniendo un valor para el coeficiente de discretización
χ = 1029 . Se definió un tamaño de población m=50 y las respectivas
probabilidades de mutación y cruce como p = 1 150 y q = 1 147 , de manera

66
que los individuos con buena calificación tuvieran en promedio 3
generaciones para donar sus cromosomas a la población.

Regresion de parametro Xi.

Selección mejores 100 de poblacion de 1500 individuos
Con 6 Secuencias de entrenamiento.
12

10

8
(1-Eval)

6 Serie1
y = 118,85x - 57,844 Tendencia Lineal
2
R = 0,9841
4

0
0,51 0,52 0,53 0,54 0,55 0,56 0,57
-log2(Qval)

Regresión de la información experimental para el parámetro ξ

(Figura 4.6)

4.4 Resultados de las pruebas

Las pruebas se efectuaron sobre los mejores 8 individuos, obtenidos luego

del paso de 300 generaciones sobre la población de 50 individuos. Los datos
registrados del comportamiento de la población se encuentran en (Figura
4.7), (Figura 4.8) y (Figura 4.9). En ellas se puede observar cómo el
promedio de calificación para el grupo de individuos decae de manera
uniforme. A diferencia de este resultado, el comportamiento del máximo y el
mínimo para la población, presentó una alta variabilidad en la calificación
para el peor y mejor de los individuos en la población, dando así, un buen
ejemplo del efecto de las operaciones de cruce y mutación.

Los resultados obtenidos para los 8 individuos fueron muy similares,

contando con una fuerte similitud entre ellos, indicando que una buena
proporción de cromosomas eran iguales, lo que implica que se encontraban

67
 en una región similar del espacio de búsqueda, la cual era muy
posiblemente, un mínimo local.

En los resultados de las tablas 4.2 y 4.3, los individuos presentaron un alto
número de verdaderos positivos y falsos positivos, frente a un número muy
bajo de verdaderos negativos como de falsos negativos. Esto es
consecuencia del sesgo particular del óptimo local en donde se encontraban
los individuos utilizados en la prueba. Este sesgo provocó que el algoritmo
realizara predicciones con un alto número de hélices, las cuales
adicionalmente eran de gran longitud y cubriendo secuencias experimentales
de hélices como de colas y loops, esto ilustrado por el Ejemplo 4.1,
generando así una gran cantidad de verdaderos positivos y falsos positivos.

(Ejemplo 4.1)
Diferencias en predicción y secuencia experimental típicas en HMMTOP-GA

Estructura Experimental: iiiiIIIIiiihhhhhhhoooohhhhhhhiii

Predicción HMMTOP-GA típica: iihhhhhhhhhhhhhhhhhhhhhhhhhhhhoo

Estas diferencias en predicción, fueron la razón por la cual, el numero de

secuencias alfa-hélice predichas correctamente (Secuencias P.), fue bajo,
para las pruebas realizadas.

(Tabla 4.2)
Resultados para los 8 individuos de prueba seleccionados.

M1 M2 M3 M4 M5 M6 M7 M8
Secuencias E. 100 100 100 100 100 100 100 100
Secuencias P. 40 79 39 39 44 38 52 71
Verdaderos Positivos 1824 1632 1853 1853 1848 1852 1833 1781
Verdaderos Negativos 129 369 85 86 94 82 99 189
Falsos Positivos 1404 1162 1447 1446 1439 1451 1434 1340
Falsos Negativos 36 228 7 7 12 8 27 79
Completitud 57,55% 58,97% 57,11% 57,14% 57,23% 56,99% 56,94% 58,06%
Precisión 33,13% 34,81% 32,64% 32,67% 32,75% 32,48% 32,42% 33,77%

Datos del Modelo

Score 0,4315 0,426 0,4324 0,4331 0,4331 0,4323 0,4357 0,4287
Generación 158 183 232 233 242 243 246 298

68
 (Tabla 4.3)
Comparación de los resultados promedios de los métodos HMMTOP y
HMMTOP-GA.

Media HMMTOP-GA Media HMMTOP

Secuencias E. 100 100
Secuencias P. 50,25 61,5
Verdaderos Positivos 1809,5 646,125
Verdaderos Negativos 141,625 1028,625
Falsos Positivos 1390,375 504
Falsos Negativos 50,5 1213,875
Completitud 0,5749875 0,493525
Precisión 0,3308375 0,2442

En comparación con los resultados del experimento HMMTOP, las métricas

de precisión y completitud fueron mas altas en los algoritmos genéticos, pero
el número de segmentos alfa-hélice predichas correctamente, fue superior en
el caso de HMMTOP. Este análisis permite establecer que los resultados
obtenidos por el experimento HMMTOP fueron más favorables frente a los
obtenidos mediante HMMTOP-GA, es decir, que HMMTOP realiza
predicciones más confiables.

Lo anterior se atribuye a dos factores: el tamaño de la población utilizado y el

número de generaciones al cual fue sometida la población. El tamaño de la
población empleada no fue lo suficientemente grande para cubrir de manera
adecuada el espacio de búsqueda, provocando que la población se
estancara en una región local poco óptima. Con respecto al número de
generaciones, ésta fue reducida y por tanto no permitió que la población
desarrollara un mejor umbral dentro de la región local y así generar mejores
resultados.

Por tal motivo, es necesario incrementar la población y el número de

generaciones con el propósito de minimizar la probabilidad de que dicha
población termine estancada en una región no-óptima y adicionalmente,
cuente con el suficiente número de generaciones para que la población
desarrolle mejores resultados en la región optima. Los resultados Tabla 4.2,
(columna Generación) son evidencia de lo anterior. Teóricamente los
mejores individuos deberían tener una diferencia generacional mínima (La
diferencia de generación entre los individuos)

69
 Máxima puntuación en Población por generación.

0,9

0,8

0,7

0,6

0,5 Máximo Población

0,4

0,3

0,2

0,1

0
0 50 100 150 200 250 300 350

Puntuación del peor individuo en la población en cada generación.

(Figura 4.7)

70
5 Conjunto de pruebas y aprendizaje

El conjunto de pruebas y aprendizaje se recopiló con base en el conjunto de

secuencias 83TMP, la cual fue utilizada por los autores de HMMTOP como
fuente para generar los parámetros en el método de pseudo conteo, el
proceso de optimización en el algoritmo Forward-Backward y el conjunto de
pruebas.

De este se conjunto se seleccionaron 15 secuencias para la realización de

pruebas, tanto en la reproducción del experimento HMMTOP, como para la
propuesta de mejoramiento con algoritmos genéticos. De las 15 secuencias,
6 fueron utilizadas para el proceso de optimización en la reproducción del
experimento HMMTOP y el entrenamiento de la propuesta de mejoramiento
con algoritmos genéticos.

A continuación se procedió a recopilar las secuencias con sus respectivas

estructuras experimentales. La metodología utilizada para la recopilación y
las secuencias correspondientes se presentan a continuación.

5.1 Obtención del conjunto de pruebas

La metodología de obtención de la información tanto de las secuencias y su

estructura consistió en 5 pasos: selección de proteínas referencia, obtención
de la secuencia, búsqueda de referencia cruzada, selección de segmento y
obtención de estructura secundaria.

5.1.1 Selección de proteínas referencia

Las proteínas referencia se utilizaron como base para la búsqueda de

proteínas que contaban con su estructura secundaria determinada
experimentalmente. Estas proteínas referencia se escogieron con base en el
subconjunto utilizado por los autores de HMMTOP, recopilados del conjunto
83TMP, cuyas características incluían el contar con una longitud de
secuencia inferior a 500 residuos y una identidad inferior al 25% con respecto
a las proteínas del conjunto de prueba.

71
5.1.2 Obtención de la secuencia

Con el identificador de la proteína, la obtención de la secuencia se realizó

mediante búsqueda en la base de datos Swissprot, descargando la
secuencia de referencia correspondiente en formato Fasta.

5.1.3 Búsqueda de referencia cruzada

Las secuencias de aminoácidos utilizadas en el conjunto de pruebas, fueron

obtenidas mediante búsqueda de la proteína con estructura terciaria,
determinada experimentalmente con mayor identidad en la proteína
referencia, mediante la utilización del algoritmo Blast-P, en la base de datos
ModBase[15]. Así se obtuvo el ID del “protein database” o PDB y la cadena
correspondiente al resultado de la búsqueda, de la cual se extrajo la
secuencia.

5.1.4 Selección de segmento

La selección del segmento a ser finalmente utilizado, se realizó mediante la

selección del primer segmento que no contara con residuos sin estructura
secundaria determinada y utilizando la información contenida en el archivo
PDB y el visor de proteínas Cn3D.

5.1.5 Obtención de estructura secundaria

La estructura secundaria fue obtenida a partir del resultado experimental

registrado por los autores del modelo tridimensional de la secuencia, el cual
se obtuvo de la base de datos del RCSB [16].

5.2 Conjunto de aprendizaje y pruebas

En las Tablas 5.1 y 5.2 se compila el conjunto de pruebas y aprendizaje

utilizado tanto en la reproducción del experimento HMMTOP y la propuesta
de mejoramiento con algoritmos genéticos. La columna “ID Swissprot”
corresponde a la identificación de la proteína referencia en la base de datos
Swissprot. La columna “Referencia Swissprot” corresponde al identificador de
la proteína en el conjunto 83TMP. La columna “PDB” corresponde al
identificador “protein database” de estructura terciaria para la secuencia
encontrada mediante Blast-P. La columna “cadena” corresponde a las
cadenas se aminoácidos en el archivo PDB que corresponde al resultado. La
columna “longitud” corresponde al tamaño del segmento seleccionado en la
cadena. Finalmente, las columnas “E-Value” e “identidad” corresponden a las

72
métricas de error en identidad y completitud de secuencia, calculadas por
Blast-P para la cadena PDB y la secuencia de la proteína referencia.

(Tabla 5.1)
Información del conjunto de pruebas.

ID Referencia
Swissprot Swissprot PDB Cadena longitud E-value Identidad
P06654 spg_strsp 1IGD A 65 2e-18 98.0%
P0ABI8 cyoa_ecoli 1FFT A 501 0 100.0%
P0ABJ1 cyob_ecoli 1FFT B 257 9e-73 100.0%
P0ABJ3 cyoc_ecoli 1FFT C 185 2e-37 100.0%
P02725 glp_pig 1AFO A 40 2e-4 65.0%
P02921 melb_ecoli 1PW4 A 222 4e-4 13.0%
P03999 opsb_human 1U19 A 348 4e-79 44.0%
P04233 hg2a_human 1IIE A 75 3e-32 100.0%
P06009 rcel_rhovi 1DXR L 273 3e-91 100.0%
P06010 rcem_rhovi 1DXR M 323 1e-109 100.0%
P08099 ops2_drome 1L9H A 235 1e-44 27.0%
P08138 ngfr_human 1NGR * 85 3e-34 93.0%
P08908 5h1a_human 1JFP A 309 2e-48 22.0%
P33742 bach_halss 1E12 A 239 2e-58 69.00%
P69052 bac1_hals1 1UAZ Z 236 1e-74 100.0%

(Tabla 5.2)
Información del conjunto de entrenamiento.

ID Referencia
Swissprot Swissprot PDB Cadena longitud E-value Identidad
P0ABI8 cyoa_ecoli 1FFT A 501 0 100.0%
P04233 hg2a_human 1IIE A 75 3e-32 100.0%
P06009 rcel_rhovi 1DXR L 273 3e-91 100.0%
P06654 spg_strsp 1IGD A 65 2e-18 98.0%
P08099 ops2_drome 1L9H A 235 1e-44 27.0%
P08138 ngfr_human 1NGR * 85 3e-34 93.0%

73
6. Conclusiones

Los resultados obtenidos mediante la implementación realizada de

HMMTOP, demostraron que es un método muy poderoso para predicción de
la estructura secundaria en proteínas de transmembrana. Incluso, sin el uso
del método de Pseudo conteo, dada la flexibilidad que proveen los modelos
ocultos de Markov para modelar el comportamiento de la estructura y el
poderoso método de optimización Forward-Backward. Sin embargo, es
importante anotar que la debilidad de esta técnica al realizar la optimización
en problemas con múltiples óptimos locales radica en la realización de un
gran esfuerzo para su afinación. Esta conclusión prueba el cumplimiento de
uno de los objetivos de este trabajo, consistente en la realización de un
análisis detallado del método HMMTOP.

Por otra parte, los resultados obtenidos con la propuesta de mejoramiento,

aportaron evidencia para concluir que si bien el hecho que los resultados en
las pruebas fueron inferiores a los obtenidos por HMMTOP, es también un
método valido y aplicable al problema de predicción de la estructura
secundaria en proteínas de transmembrana, por contar con la capacidad de
solventar algunas de las debilidades de HMMTOP. Al método HMMTOP-GA
le queda todavía un gran potencial por ser investigado para que sea
aprovechada su gran utilidad práctica. El método propuesto HMMTOP-GA,
provee además, la solución de una primera aproximación sobre cómo
desarrollar estrategias que usen la sinergia de los algoritmos genéticos y los
modelos ocultos de Markov, aproximación útil no solo en el contexto
particular de la predicción de estructura secundaria, sino también en gran
cantidad de problemas de optimización multivariada no lineal en diferentes
áreas y aplicaciones. A su vez, también constituye el logro de importantes
objetivos de este trabajo de grado, como lo era el proveer de una referencia
útil en futuras investigaciones.

Como punto final, es importante resaltar algunas de las propuestas más

interesantes de este trabajo de grado. La proposición del método de
diversidad controlada y el modelo para la dinámica de población y
convergencia en algoritmos genéticos, que constituyen dos de los logros de
mayor valor; por su aporte al área de optimización evolutiva y las propuestas
subsiguientes que puedan gestarse a raíz de estas. Por tal motivo considero
importante que sean estudiadas formalmente, en un futuro próximo.

74
Bibliografía y Referencias

[1] Gábor E. Tusnády and István Simon. (1998), Principles Governing Amino
Acid Composition of Integral Membrane Proteins: Application to Topology
Prediction. J. Mol. Biol. 283, 489-506.

[2] Gonnet H. Gaston. (2004), Some string matching problems from

Bioinformatics which still need better solutions. Journal of Discrete
Algorithms. No 2, pp. 3-15.

[3] Paul D. Taylor, Teresa K. Attwood and Darren R. Flower. (2003),

BPROMPT: A consensus server for membrane protein prediction. Nucleic
Acids Research. Vol 31, No. 13.

[4] Sen Z. Taner, Jernigan L. Robert, Garnier Jean and Kloczkowski Andrzej.
(2005), GOR V server for protein secondary structure prediction.,
Bioinformatics. Vol 21, No 11, pp. 2787-2788.

[5] Cormen. H. (2001), Introduction to Algoritms, Second Edition. MIT Press.

[6] Whitley Darrell. A Genetic Algorithm Tutorial. Computer Science

Department. Colorado State University.

[7] Jones, D. T., Taylor, W. R. & Thorton, J. M. (1994). A model recognition

approach to the prediction of all helical membrane protein structure and
topology. Biochemistry, 33, 3038-3049.

[8] Rost, B., Fariselli, P. & Casadio, R. (1996). Topology prediction for helical
transmembrane proteins at 86% accuracy. Protein Sci. 5, 1704-1718.

75
[9] Cserzo , M., Wallin, E., Simon, I., von Heijne, G. & Elofsson, A. (1997).
Prediction of transmembrane alpha-helices in prokariotic membrane proteins:
the dense alignment surface method. Protein Eng. 10, 673-676.

[10] Swiss Institute of Bioinformatics, Swiss-Prot. Release 50.1, ExPASy

Proteomics Server. [Link]

[11] Rossi, William. (2001). Dfp floating point Library.

[Link]

[12] Ming li, Xin Chen, Xin Li, Bin Ma, Paul Vitányi. The Similarity Metric.

[13] Shannon, C.E. (1948), A mathematical theory of communication. Bell

system Technical Journal, Vol 27, pp. 379-423, 623-656.

[14] Freedman, David. H. (1983), Markov chains. Springer-Verlag.

[15] University of California, San Francisco. ModBase., Database of

Comparative Protein Structure Models., [Link]

[16] RCSB protein data bank. [Link]

[17] Baxevanis,A.D. and Ouellette,B.F.F.(2001), Bioinformatics: A Practical

guide to the Analysis of Genes and Proteins. John Wiley & Sons, Inc.

[18] Moller Steffen, Croning D. R. Michael and Apweiler Rolf. (2001),

Evaluation of methods for the prediction of membrane spanning regions.
Bioinformatics. Vol 17, No. 7, pp. 646-653.

[19] Cuthbertson M. Jonathan, Doyle A. (2005) Declan and Sansom S. P.

Mark. Transmembrane helix prediction: a comparative evaluation and
analysis. Protein Engineering, Design and Selection. Vol 18, No. 6,
pp.295-308.

76
[20] Lesk M. Arthur. (2002), Introduction to Bioinformatics. Second ed.
Oxford University Press.

[21] Baldi. Pierre, Brunak. Soren. (2001), Bioinformatics. The machine

learning Approach. Second ed. The MIT press.

[22] A. Raval, Ghahramani Z y D. L Wild. (2002) A bayesian network model

for protein fold and remote homologue recognition. Bioinformatics. Vol 18,
No 6, pp. 788-801.

[23] Ouzounis A. Christos and Valencia, Alfonzo. (2003), Early bioinformatics:

the birth of a discipline — a personal view. Bioinformatics. Vol 19, No 17,
pp. 2176-2190.

[24] Moti Zviling, Leonov. Hadas y Arkin T. Isaiah. (2005), Genetic algorithm-
based optimization of hydrophobicity tables. Bioinformatics. Vol 21, No 11,
pp. 2651-2656.

[25] Kuang Lin, Simossis A. Victor, Taylor R. Willam and Jaap Heringa.
(2005), A simple and fast secondary structure prediction method using hidden
neural networks. Bioinformatics. Vol 21, No 2, pp. 152-159.

[26] Ingrid Y. Y. Koh1, Volker A. Eyrich, Marc A. Marti-Renom, Dariusz

Przybylski, Mallur S. Madhusudhan, Narayanan Eswar, Osvaldo Graña,
Florencio Pazos, Alfonso Valencia, Andrej Sali and Burkhard Rost. (2003),
EVA: evaluation of protein structure prediction servers. Nucleic Acids
Research. Vol 31, No13, pp 3311-3315.

[27] Junwen Wang and Jin-An Feng. (2003), Exploring the secuence patterns
in the alpha-helices of proteins. Protein Engineering. Vol 16, No 11, pp 799-
807.

[28] H. Mamitsuka y K Yamanishi. (1995), Alpha-Helix region prediction with

stochastic rule learning. Comput. Appl. Biosci., Vol 11, pp 399 - 411.

77
[29] Okamoto, M Fukugita, T Nakazawa y H Hawai. Alpha-helix folding by
Monte Carlo simulated annealing in isolated C-peptide of ribonuclease.
Department of Physics. Nara Women's University

[30] Oppon. Ekow. (2000), Synergistic use of promoter prediction algorithms:

A choice form small training dataset? Department of computer Science,
University of the West Cape.

[31] Vivas. Sonia. (2005), Reconocimiento de promotores: Algoritmos y

sinergia.
Departamento de Ingeniería de Sistemas y Computación. Universidad
de los Andes.

[32] Niño. Marcos y García. Rafael. Predicción de Pseudonodos en la

estructura secundaria de ARN vía Gramáticas Estocásticas Independientes
del contexto. Departamento de Ingeniería de Sistemas y Computación.
Universidad de los Andes.

[33] Rabiner R. Lawrence. (1989), A tutorial on Hidden Markov Models and

Selected Applications in Speech Recognition. Proceedings of the IEEE, Vol
77 No 2

[34] Mitchell, Melanie.(1999), An Introduction To Genetic Algorithms. MIT

Press.

78