Modelado de filtros IIR usando un algoritmo inspirado en el electromagnetismo

Resumen

El filtraje de respuesta infinita al impulso (IIR) provee un enfoque poderoso para la solución de una variedad de problemas. Sin embargo, su diseño representa una tarea compleja, ya que su superficie de error es multi-modal con respecto a los coeficientes del filtro. Este hecho imposibilita la utilización de métodos iterativos de optimización. Para resolver este problema, se requiere la utilización de técnicas de optimización global que eviten quedar atrapados en mínimos locales. En este artículo se propone un método de identificación y modelado de filtros IIR basado en el enfoque con capacidades de optimización global llamado: Electromagnetism-Like Optimization Algorithm (EMO). EMO es un método heurístico que emplea los principios de la teoría del electromagnetismo para la definición de operadores evolutivos. En el algoritmo, las soluciones son construidas considerando la atracción y repulsión electromagnética entre las partículas cargadas; en donde la carga representa la afinidad que tiene cada partícula con la solución. Resultados de simulación demuestran que el enfoque propuesto tiene una mejor precisión y convergencia que los obtenidos por métodos reportados en la literatura para desempeñar la misma tarea.

Abstract

Infinite-impulse-response (IIR) filtering provides a powerful approach for solving a variety of problems. However, its design represents a very complicated task, since the error surface of IIR filters is generally multimodal, global optimization techniques are required in order to avoid local minima. In this paper, a new method based on the Electromagnetism-Like Optimization Algorithm (EMO) is proposed for IIR filter modeling. EMO originates from the electro-magnetism theory of physics by assuming potential solutions as electrically charged particles which spread around the solution space. The charge of each particle depends on its objective function value. This algorithm employs a collective attraction-repulsion mechanism to move the particles towards optimality. The experimental results confirm the high performance of the proposed method in solving various benchmark identification problems.

Introducción

Los filtros digitales son clasificados en dos principales categorías (Oppenheim et al., 1999): los filtros de respuesta finita al impulso (FIR) y los filtros de respuesta infinita al impulso (IIR). La salida de un filtro FIR es calculada con base al valor actual y previos de la entrada. Esto implica que cuando no existe una señal de entrada, la salida tenderá a desaparecer. Por otro lado, la salida de un filtro IIR es dominada no sólo por el valor presente y previo de la entrada, sino también por los valores pasados de su propia salida. Esta realimentación produce el hecho de que cuando la entrada deja de existir, las salidas pasadas del filtro seguirán excitándolo.

La identificación y modelado de filtros digitales IIR es de gran importancia en diversas áreas de la ingeniería, tales como procesamiento de señales, control, sistemas de comunicaciones, procesamiento de imágenes y eliminación de ruido (Howell y Gordon, 2001). El proceso de identificación consiste en alterar los coeficientes de un determinado filtro modelo para estimar los parámetros de un sistema desconocido, a partir del conocimiento de su entrada y salida. Para desempeñar esta tarea es necesario la minimización de una función objetivo, típicamente el error cuadrático medio (MSE, mean square error) existente entre la respuesta del sistema desconocido y la salida del filtro modelado.

Aún y tomando en cuenta la amplia gama de aplicaciones de los filtros IIR, su identificación y modelado se considera un proceso complejo que requiere de algoritmos de optimización globales capaces de operar sobre superficies de error multi-modales típicas de estos procesos. Para resolver este problema, se han desarrollado numerosos algoritmos que hacen uso de diversas técnicas heurísticas. Algunos ejemplos son: el trabajo realizado por Howell y Gordon en 2001, en el cual se emplea un algoritmo de identificación basado en autómatas de aprendizaje. La propuesta de Kalinli y Karaboga (2005), que introdujeron el uso de algoritmos inmunológicos (IA) como algoritmo de optimización para el diseño de filtros IIR. Recientemente, el enfoque presentado por Rashedi (2011), en el cual se propone el modelado de filtros digitales usando el método de búsqueda gravitational (GSA). Otra propuesta, que además en la literatura se usa como forma de comparación estándar, es la presentada en Lin et al. (2008), en la cual se propone un algoritmo de identificación basado en el popular método Particle Swarm Optimization (PSO). A pesar de la gran aceptación que tienen estas propuestas, existen algunas deficiencias que afectan significativamente su desempeño (Shilane et al., 2008). Entre ellas se destacan la necesidad de un gran número de iteraciones y de elementos de la población para encontrar una solución lo más parecida al mínimo global, repercutiendo directamente en el costo computacional del algoritmo (Elbeltagia et al., 2005). Además, para la mayoría de estos métodos, la convergencia no está comprobada, por lo tanto no existe garantía de que los valores encontrados como solución, puedan ser considerados como los mejores.

Ante estos problemas, en este artículo se presenta un nuevo algoritmo para la identificación y modelado de filtros IIR, basado en la técnica heurística Electromagnetism-Like Optimization EMO (İlker et al., 2003). EMO es un método estocástico poblacional basado en la teoría del electromagnetismo, cuyas propiedades de convergencia ya fueron probadas en (İlker et al., 2004; Rocha y Fernández, 2009a). Los operadores usados por EMO para el cálculo de nuevos puntos de búsqueda se asemejan a aquellos empleados por los algoritmos de enjambre de partículas (PSO) (Ying y Pei, 2010) y al grupo de algoritmos de optimización basados en colonias de hormigas (ACO) (Blum, 2005). El primer paso del algoritmo EMO es generar un grupo de soluciones aleatorias a partir de un dominio en el cual es factible encontrar las soluciones, suponiendo que cada solución emula a una partícula cargada. La carga de cada partícula la determina una función objetivo (función a optimizar); utilizando las evaluaciones de tal función el algoritmo EMO modifica la posición de cada partícula usando operadores que emulan la atracción o repulsión existente dentro de la población de partículas.

Dicho mecanismo de atracción-repulsión del algoritmo EMO corresponde a los operadores de mutación y crossover en los algoritmos evolutivos (Rocha y Fernández, 2009b).

De forma general, el algoritmo EMO se puede considerar como un algoritmo rápido y robusto que representa una alternativa real para resolver problemas de optimización complejos, no-lineales, no-diferenciables y no-convexos (Naderi et al., 2010; Yurtkuran y Emel, 2010; Jhen y Kun, 2009).

Las principales ventajas del algoritmo EMO radican en las siguientes características: no tiene operaciones de gradiente, se puede emplear directamente en el sistema decimal a diferencia de los algoritmos genéticos (AG), necesita pocas partículas para converger y se garantiza su convergencia (İlker et al., 2004; Rocha y Fernández, 2009a).

Este artículo presenta un nuevo enfoque para la identificación y modelado de filtros IIR, el cual está basado en EMO. Bajo este enfoque, la adaptación de los coeficientes es considerada como un problema de optimización, para la cual, se utiliza como función objetivo el error cuadrático medio (MSE) existente entre la respuesta del sistema desconocido y la estimada por el filtro, durante cierto periodo de tiempo. Guiados por los valores de la función objetivo el conjunto de coeficientes que modelan al filtro son operados usando el algoritmo EMO hasta encontrar la mejor configuración posible (mínimo global).

Los resultados de simulación obtenidos, demuestran que el enfoque propuesto tiene una mejor precisión y convergencia que los obtenidos por métodos reportados en la literatura para desempeñar la misma tarea.

La estructura del resto del artículo es la siguiente: la segunda sección ofrece un panorama general del algoritmo EMO. En la tercera sección se definen los filtros digitales y su diseño.

En la cuarta se formula el enfoque propuesto considerando la aplicación del algoritmo EMO para encontrar los coeficientes de los filtros digitales. La quinta sección presenta los resultados experimentales. Finalmente, en la última sección se discuten las conclusiones y se establece el trabajo futuro.

Algoritmo: Electromagnetism-Like Optimization (EMO)

El algoritmo EMO es un algoritmo de optimización, el cual está inspirado en principios físicos. El algoritmo permite la optimización global de funciones multi-modales. En comparación con los AG, no emplea operadores de mutación o crossover para explorar el espacio de búsqueda, ya que está basado en el fenómeno físico del electromagnetismo.

El algoritmo EMO permite resolver problemas de optimización, definidos de la siguiente forma:

donde [l, u] = {x ∈ ℜⁿ|l_d ≤ x_d ≤ u_d, d = 1, 2,…,n&}, siendo n la dimensión de la variable x, [l, u] ⊂ ℜⁿ, es un subconjunto no vacío, y f: [i, u] → ℜ es una función de valores reales. Por lo tanto, son conocidas las siguientes características del problema:

• n: es la dimensión del problema.

• u_d: es el límite superior de la dimensión.

• l_d: es el límite inferior de la dimensión.

• f (x): es la función que será minimizada.

Con estas condiciones, el algoritmo EMO emplea dos procesos básicos para la optimización: primero explora el espacio de búsqueda en forma aleatoria, el segundo proceso consiste en la explotación local de los puntos elegidos, para este propósito el algoritmo EMO usa los principios de la teoría del electromagnetismo. Utilizando ambos procesos se garantiza que el algoritmo converge en los mínimos de la función que son altamente atrayentes, y se aleja de los valores máximos. Partiendo de estos dos procesos fundamentales, el algoritmo EMO tiene cuatro fases para lograr la optimización global (ilker y Shu, 2003). Cada etapa es descrita a continuación:

• Inicialización:m partículas son tomadas aleatoriamente considerando el límite superior (u) y el límite inferior (l).

• Búsqueda local: se busca un mínimo en la vecindad de un punto x^p, donde p ∈ (1,…, m).

• Cálculo del vector de fuerza total: las cargas y fuerzas son calculadas para cada partícula.

• Movimiento: cada partícula es desplazada de acuerdo al vector de fuerza total calculada.

Inicialización

Se produce de manera aleatoria un grupo G de m soluciones n-dimensionales, el cual es considerado como la población de soluciones iniciales. Cada solución es considerada como una partícula cargada y se asume que todas las partículas están distribuidas uniformemente entre el límite superior (u) y el límite inferior (l). La mejor partícula (mejor solución) es encontrada por la evaluación de la función objetivo, la cual depende de cada problema de optimización. Este procedimiento termina cuando todas las m partículas son evaluadas, eligiendo la partícula que en relación a la función objetivo tiene un mejor resultado. Esta fase corresponde al proceso de explorar las regiones atrayentes en un espacio factible.

Búsqueda local

En esta etapa se intenta mejorar la solución ya encontrada, sin embargo, para algunos problemas puede resultar innecesaria. Considerando esto, es posible formular una clasificación de los algoritmos EMO: EMO sin búsqueda local, EMO con búsqueda local aplicada sólo a la mejor partícula actual y EMO con búsqueda local aplicada a todas las partículas, éste último es el caso de estudio de este artículo.

Considerando un determinado número de iteraciones, llamado LSITER y un parámetro de vecindad de búsqueda δ, el procedimiento para encontrar el valor óptimo local se lleva a cabo de la siguiente manera: el punto x^p es asignado a una variable temporal y para almacenar la información inicial. A continuación, para una coordenada dada d de una partícula, un número aleatorio (λ₁) es seleccionado y combinado con δ obteniendo la longitud de paso de búsqueda. El punto y es entonces desplazado en la dirección que la longitud de paso indique, el signo que esta dirección tendrá, también se calcula de manera aleatoria (λ₂). Si el valor obtenido tras evaluar y en la función a optimizar es mejor después de haber realizado LSITER iteraciones, el punto x^p es reemplazado por y, terminando así la búsqueda en la vecindad p, de otra forma x^p conserva su valor. Finalmente, el mejor punto actual es actualizado en la partícula. El pseudo-código de este algoritmo es descrito en la figura 1.

Figura 1.

Pseudo-código del método de búsqueda local

De forma general, la búsqueda local aplicada a todas las partículas puede reducir el riesgo de caer en un mínimo local, pero la desventaja es que tiende a aumentar el tiempo de cómputo. Mantener la búsqueda local centrada en la mejor partícula actual resulta más conveniente, ya que así es posible mantener la eficiencia y precisión computacional. En la búsqueda local la longitud de paso representa un importante factor que depende de los límites de cada dimensión y determina el desempeño del método de búsqueda local.

Cálculo del vector de fuerza total

El cálculo del vector de fuerza total se basa en el principio de superposición (figura 2) de la teoría del electromagnetismo, el cual establece que: “la fuerza ejercida en una partícula por medio de otra partícula es inversamente proporcional a la distancia entre los puntos y directamente proporcional al producto de sus cargas” (Cowan, 1968). Cada partícula se desplaza de acuerdo a la ley de Coulomb (figura 3), la cual emplea la fuerza producida entre las partículas, que depende del valor de carga poseída por cada una de ellas. Dicha carga es determinada por su valor de desempeño de la función objetivo y se calcula de la siguiente manera:

donde n denota la dimensión y m representa el tamaño de la población de partículas del algoritmo EMO, el vector x^h corresponde a la h-ésima partícula de la población. Un problema de optimización que se define con un número de dimensiones elevado generalmente requiere una mayor población. En la ecuación (2), la partícula con el mejor valor en la función objetivo x^mejor es llamada “la mejor partícula”, y tiene una mayor carga. La fuerza de atracción que ejerce la mejor partícula sobre otra partícula dada es inversamente proporcional a la distancia existente entre ellas. Por lo tanto, esta partícula atrae otras partículas que tengan peores resultados de afinidad y repele a las que tienen mejores valores de afinidad.

Figura 2.

Principio de superposición

Figura 3.

La ley de Coulomb

La fuerza resultante que existe entre las partículas determina la modificación del valor de las partículas en el proceso de optimización. La fuerza de cada partícula se calcula por la ley de Coulomb y el principio de superposición, por medio de (3).

donde f(x^h) < f(x^p) representa la atracción y f(x^h) ≥ f(x^p) representa la repulsión (figura 3). La fuerza resultante de cada partícula es proporcional al producto de las cargas e inversamente proporcional a la distancia entre las partículas. Para que el proceso sea numéricamente consistente (3) debe ser normalizada como a continuación se muestra:

Movimiento

De acuerdo a la fuerza resultante cada partícula se desplaza como indica la siguiente ecuación:

En la ecuación (5), λ es un paso de búsqueda aleatorio, el cual es uniformemente distribuido entre cero y uno, por otro lado, u_d y l_d representa los límites superior e inferior de la dimensión d, respectivamente. Si la fuerza es positiva, la partícula se desplaza hacia el límite superior, en caso contrario se desplaza hacia el límite inferior, en ambos casos utilizando una longitud de paso aleatoria. La mejor de las partículas no se mueve, debido a que es la que tiene una atracción absoluta y repele o atrae a los demás elementos de la población.

El proceso termina cuando se alcanza un número máximo de iteraciones o cuando un valor f(x^mejor) es óptimo en algún sentido. Estas tres fases de EMO (búsqueda local, cálculo del vector de fuerza total y movimiento), representan el proceso de explotación para encontrar el valor óptimo.

Modelado de filtros digitales

El proceso de identificación de un filtro es un proceso que consiste en modificar los coeficientes del filtro, hasta que se asemeje a la función de transferencia de un sistema desconocido (figura 4).

Figura 4.

Proceso de identificación usando algoritmos meta-heurísticos

Dadas estas características, es posible definir el modelado de filtros como un problema de optimización, en el cual se busca un conjunto de soluciones, que den como resultado un mínimo error entre la salida del filtro y la salida del sistema desconocido. Para dar solución a esto, es posible hacer uso de las técnicas de optimización meta-heurísticas.

Filtros FIR e IIR

Dentro del procesamiento digital de señales los filtros son representados por medio de una función de transferencia. De esta forma, los filtros FIR son definidos de forma general por (6).

donde, B = {b₀, b₁, b₂, …, b_r} es un vector R-dimensional, que contiene los parámetros del numerador. Del mismo modo, Y(z) corresponde a la respuesta o salida del filtro, mientras que X(z) es la señal de entrada del filtro FIR.

Por otra parte, un filtro IIR se define por la siguiente función de transferencia:

donde Y (z) y X (z) son la salida y la entrada del filtro, respectivamente, el vector de parámetros está representado por A ∪ B = {b₀, b₁, b₂, …, b_r, a₀, a₁, a₂, …, a_c} de forma que sus elementos se encuentran en las dimensiones {R + C}. Estos dos tipos de filtros son comúnmente usados dentro del procesamiento digital de señales.

El problema de diseño de estos filtros se puede resumir en la búsqueda de los valores de los vectores B y A ∪ B, que mejor aproximen la salida del filtro con la salida de un modelo dado.

Identificación de filtros digitales usando el algoritmo EMO

Representación de partículas

Bajo el enfoque EMO, las partículas son generadas de forma aleatoria dentro de un despacio acotado d-dimensional. Para la identificación de filtros IIR, la dimensión es de d = R + C, que corresponde a la longitud de los vectores del numerador y denominador de dicho filtro. De esta manera, las partículas son almacenadas en un vector de soluciones S, el cual contendrá m partículas de d dimensiones cada una. De manera tal que S = {s₁, s₂, s₃, …, s_m} } contiene partículas de la forma sm=sm1,sm2,sm3,…,smdT donde T se refiere al vector s_m transpuesto.

Función objetivo

La calidad de dichas partículas es evaluada con respecto a la función objetivo, que mide la afinidad de las respuestas entre el sistema desconocido y el filtro estimado. Al evaluar cada partícula en el filtro digital IIR, éste genera una salida Y_fil, la cual debe ser comparada con la salida del sistema desconocido Y_p. Esta diferencia es medida empleando el error cuadrático medio (MSE) (Howell et al., 2001). Este valor es considerado en este trabajo como la función objetivo J, la cual se define formalmente como.

Para la medición del MSE tanto el sistema como el filtro son simulados durante un número de iteraciones. Por objetivos de comparación, este parámetro se establece como (Howell y Gordon, 2001), que es el valor usado por las demás propuestas presentadas en la literatura. El valor de J debe ser minimizado por el algoritmo EMO, tras cada iteración los valores serán simulados y un nuevo valor de afinidad será asignado a cada partícula. Finalmente, tras un número definido de iteraciones o un valor J permisible se encuentran los parámetros óptimos del filtro.

Implementación

La implementación del algoritmo propuesto se puede resumir en los siguientes pasos:

• Paso 1:

  Se generan m partículas iniciales. Cada partícula s_p es evaluada en la función objetivo J. La mejor partícula s^mejr se determina a partir de los valores de J, donde s^mejor = arg min {J (s_p), ∀ p}.

• Paso 2:

  La partícula s_p es almacenada temporalmente en y. Se selecciona un número aleatorio que es combinado con δ para obtener la longitud de paso para una dimensión d dada. Por lo tanto, la partícula s_p es desplazada a lo largo de esa dirección. Si Js_p es minimizada, la partícula s_p es remplazada por su nuevo valor, de lo contrario, el valor almacenado temporalmente se mantendrá.

• Paso 3:

  La carga entre las partículas es calculada usando (2), y su vector de fuerza se calcula con (3). La partícula s^mejpr, con un mejor valor en la función objetivo, mantiene una carga superior, y por lo tanto, una mayor fuerza de atracción o repulsión.

• Paso 4:

  Las partículas son desplazadas de acuerdo a la magnitud de su fuerza. La nueva posición de la partícula es calculada por (4). s^mejor no es desplazada porque tiene la fuerza más grande y ésta atrae las otras partículas hacia ella.

• Paso 5:

  El índice n es incrementado. Si n = MAXITER o si el valor de J(S) es menor que un valor de umbral predefinido entonces el algoritmo se detiene y el flujo continua en el paso 6. De lo contrario, se regresa al paso 3.

• Paso 6:

  La mejor particular s^mejr (solución) es seleccionada de la última iteración.

Resultados experimentales

Para evaluar el desempeño del algoritmo propuesto en la identificación de filtros IIR, se empleó un conjunto de problemas propuestos originalmente por Johnson y Larimor (1977), y que son utilizados como pruebas estándar para calificar la calidad de identificación de una determinada propuesta. Estos problemas en general incluyen la identificación de sistemas, considerando la utilización de filtros IIR de diferente orden. Con el objetivo de verificar la eficiencia del algoritmo propuesto, los resultados se compararon con el enfoque de identificación propuesto en (Lin et al., 2008), basado en el algoritmo de optimización heurístico PSO.

La descripción del trabajo experimental se presenta en el siguiente orden:

• •

  Planta de segundo orden con un modelo de primer orden.

• •

  Planta de segundo orden con un modelo de segundo orden.

• •

  Planta de tercer orden con un modelo de segundo orden.

• •

  Planta y filtro de orden superior.

• •

  Comparaciones con el método basado en PSO.

Cada experimento conducido fue llevado a cabo con un conjunto de 10 partículas, un máximo de iteración para la búsqueda local de LSITER = 4, la longitud de paso para la búsqueda local es δ = 0.025, y un valor máximo de iteraciones n = 3000. El numero de iteraciones n fue elegido con fines comparativos con otros trabajos similares (Howell y Gordon, 2001; Lin et al., 2008; Rashedi et al., 2011).

Planta de segundo orden con un modelo de primer orden

Para este problema, la función de transferencia que define a la planta a modelar, está dada por la siguiente ecuación:

Dicha planta se pretende modelar por un filtro IIR con la siguiente función de transferencia:

donde H_p (z⁻¹) y H_m (z⁻¹) corresponden a la relación entrada-salida de la planta y del modelo propuesto, respectivamente. u (z⁻¹) es la entrada, que es una señal de ruido blanco, la cual es aplicada tanto a la planta como al modelo. Los valores a y b son modificados por el algoritmo EMO, este proceso se realiza haciendo modificaciones en base a los valores que toma la función objetivo J. Para este problema, el espacio de búsqueda se encuentra definido por los límites superior u = [1, 1] e inferior l = [−1, −1].

La figura 5 muestra la evolución del MSE a lo largo de las 3000 iteraciones, se puede apreciar que dentro de las primeras 500 iteraciones el algoritmo logra converger a un valor cercano al optimo. Por otra parte, la manera en que los valores a y b van evolucionando a lo largo de las iteraciones se muestra en la figura 6, donde se percibe que de forma similar al MSE, ambos parámetros logran establecer sus valores dentro de las primeras 500 iteraciones y después sólo se realizan ajustes mínimos para llegar a la convergencia.

Figura 5.

Evolución del error cuadrático medio para la planta de segundo orden con un modelo de primer orden

Figura 6.

Parámetros a y b a lo largo de las iteraciones

Los valores finales que calcula el algoritmo EMO comparados con los que dan solución al problema son presentados en la tabla 1.

Tabla 1.

Valores encontrados por el algoritmo EMO y valores reales del modelo

Valores	a	b	J
EMO	0.9034	−0.3030	0.0107
Reales	0.906	−0.311	-

Planta de segundo orden con un modelo de segundo orden

Para este problema, la planta se define por la siguiente función de transferencia:

Mientras que (12) muestra la función de transferencia del filtro IIR propuesto.

En el modelo, los parámetros a calibrar son a₁a₂y b, dando como resultado un problema de optimización de 3 dimensiones. La entrada u (z⁻¹) en este caso, corresponde a una señal de ruido “coloreado”, el cual se genera a partir del filtro FIR que se muestra en (13). Los límites de cada una de las dimensiones son: u = [2, 1, 1] e l = [−2, −1, −1].

La figura 7 muestra como el MSE se minimiza drásticamente dentro de las primeras 100 iteraciones, finalmente se estaciona en un valor cercano a cero.

Figura 7.

Evolución del error cuadrático medio para la planta de segundo orden con un modelo de primer orden

La evolución de los tres parámetros durante el proceso de identificación es mostrada en la figura 8. Estos parámetros toman diversos valores iniciales que son modificados en base a la función objetivo J, hasta que finalmente logran converger a un valor que ha minimizado el MSE. La tabla 2 presenta los valores finales de a₁a₂ y b encontrados por el algoritmo EMO.

Figura 8.

Parámetros a₁, a₂ y b a lo largo de las iteraciones

Tabla 2.

Valores encontrados por el algoritmo EMO y valores reales del modelo

Valores	a1	a2	b	J
EMO	−1.3936	0.4803	1.0005	2.463e-5
Reales	−1.4000	0.4900	1.0000	-

Planta de tercer orden con un modelo de segundo orden

El tercer problema propuesto por el conjunto de pruebas elegido (Jonson y Larimore, 1977), corresponde a una planta de tercer orden que se pretende modelar empleando un filtro IIR con una función de transferencia de segundo orden. La ecuación (14) define la estructura planta propuesta.

Por otra parte, en (15) se propone el modelo IIR de segundo orden que aproxima a la planta.

Al igual que en las sub-secciones anteriores, el problema de optimización se define como encontrar los valores de los parámetros a₁a₂ y b, que minimicen el MSE existente entre las respuestas de y_p y y_m. Las cuales dependen de una señal de entrada u(z⁻¹), que contiene el ruido añadido por medio del filtro FIR, que se muestra a continuación:

El espacio de búsqueda de este problema consta de tres dimensiones, las cuales están acotadas por los siguientes límites: superior u = [1, 1, 1] e inferior l = [−1, −1, −1]. En la figura 9a se muestra la gráfica del MSE a lo largo de 3000 iteraciones. Es importante mencionar que en este tipo de problemas el MSE no logra llegar a valores cercanos a cero, ya que el filtro con el que se está modelando no tiene el orden apropiado, por lo tanto, el algoritmo EMO trata de buscar los mejores valores que minimicen ese valor de error; sin embargo, esto no significa que los parámetros encontrados modelen con exactitud a la planta. En la figura 9b se puede ver como los valores de los coeficientes evolucionan, se aprecia en la gráfica que aunque llegan a converger en un valor fijo, estos valores son los que corresponden al mínimo global de la superficie de error, sin que ello signifique que el modelo, debido a su bajo orden, identifique plenamente a la planta.

Figura 9.

a) evolución del error cuadrático medio para la planta de tercer orden con un modelo de segundo orden, b) parámetros a₁a₂ y b a lo largo de las iteraciones

Planta y filtro de orden superior

Hasta ahora, los filtros y las plantas propuestas han sido problemas en los cuales se han calibrado pocos parámetros. El algoritmo EMO ha arrojado resultados competitivos con esos problemas que pueden ser considerados como sencillos. Sin embargo, la mayor complicación en la identificación es encontrar los valores de los parámetros en filtros de orden mayor a tres. Para ejemplificar estos filtros se propone una planta de orden 6 (17), que debe ser modelada por un filtro IIR de orden 4 (18).

Los parámetros que deben calibrarse son a₁, a₂, a₃, a₄ así mismo b₀, b₁, b₂, b₃, b₄. Por lo tanto, como se explicó en la sección de identificación de filtros digitales usando el algoritmo EMO, el espacio de búsqueda es de nueve dimensiones y cada una tiene como limites superior u = 1.5 e inferior l = −1.5. El algoritmo EMO presenta un buen rendimiento al momento de trabajar con grandes dimensiones (İlker et al., 2003), por tal motivo resulta idónea su aplicación para encontrar los múltiples valores que den solución a la identificación con filtros de orden superior.

Los resultados experimentales muestran la eficiencia del algoritmo para minimizar el MSE. En la figura 10a se puede ver la gráfica de los valores que va tomando el error cuadrático medio durante las iteraciones. Este valor va tendiendo a cero, hasta llegar a estabilizarse; es de esperarse que al tener más dimensiones, tiene un grado de dificultad mayor. Por otra parte, la figura 10b presenta los valores que cada coeficiente IIR del denominador (a₁, a₂, a₃, a₄) toma a lo largo de las iteraciones.

Figura 10.

a) grafica de la evolución del error cuadrático medio para un filtro de orden 4, b) evolución de los parámetros del denominador de un filtro de orden 4, c) evolución de los parámetros del numerador de un filtro de orden 4

Dentro de las primeras 500 iteraciones el algoritmo basado en EMO realiza una exploración exhaustiva a lo largo del espacio de búsqueda hasta que se estabiliza o converge a los valores finales. Después, se realizan algunos ajustes finos sobre los valores de los parámetros. Ambos efectos corresponden respectivamente a la exploración y explotación del algoritmo EMO. Los parámetros del numerador (b₀, b₁, b₂, b₃, b₄), de sus valores iniciales aleatorios cambian abruptamente, pero después de transcurrir cierto numero de iteraciones logran converger para finalmente realizar los ajustes finales que minimicen el error de la manera más óptima posible. La gráfica de dicho procedimiento se muestra en la figura 10b. Los resultados que se presentan en las figuras 10a, 10b y 10c exponen las capacidades del algoritmo EMO para la optimización multi-dimensional. Los valores finales encontrados por el algoritmo después de las 3000 iteraciones se expresan en la tabla 3, mientras el valor final del MSE es de 0.0069.

Tabla 3.

Valores encontrados por el algoritmo EMO para los parámetros del numerador y los del denominado

	Valores a		Valores b
a1	−0.4950	b0	1.0335
a2	−0.7049	b1	−0.6670
a3	0.5656	b2	−0.4682
a4	−0.2691	b3	0.6961
		b4	−0.2100

Comparaciones con el método de identificación basado en PSO

Se realizó la comparación del método presentado en este artículo con el algoritmo de identificación propuesto en Lin et al. (2008), basado en el popular método PSO (Kennedy y Eberhart, 1995). La comparación se llevó a cabo con cada uno de los problemas propuestos en las sub-secciones: planta de segundo orden con un modelo de primer orden, planta de segundo orden con un modelo de segundo orden y planta y filtro de orden superior, del mismo modo cada prueba se efectuó con 3000 iteraciones. Los resultados se describen a continuación.

La figura 11a muestra las gráficas del MSE, tanto del algoritmo basado en EMO, como del algoritmo PSO, para el caso del experimento descrito en la planta de segundo orden con un modelo de primer orden. Se puede ver claramente que ambos algoritmos minimizan el error en un número de iteraciones reducido, sin embargo, el algoritmo PSO mantiene una serie de iteraciones a lo largo de un valor cercano a cero. Mientras que el algoritmo EMO tiene un comportamiento estable, realizando sólo un pequeño ajuste para minimizar el valor. La tabla 4 presenta los resultados de ambos algoritmos, asimismo, los valores óptimos reales que dan solución al problema. La evolución de los valores que toma cada uno de los parámetros durante las iteraciones realizadas es presentada en la figura 11b.

Figura 11.

a) comparación de los valores MSE del algoritmo EMO y del PSO para una planta de segundo orden modelada con un filtro de primer orden, b) comparación de los parámetros a y b, encontrados por el algoritmo EMO y PSO durante 3000 iteraciones y c) comparación de los valores MSE del algoritmo EMO y del PSO para una planta de segundo orden modelada con un filtro de segundo

Tabla 4.

Valores comparativos de los parámetros a, b y J encontrados tanto por el algoritmo EMO como el PSO

Valores	a	b	J
EMO	0.9034	−0.3030	0.0107
PSO	0.8950	−0.3142	0.0121
Reales	0.906	−0.311	–

Continuando con la comparación, se toma el segundo ejemplo, que consiste en una planta de segundo orden modelada con un filtro IIR de segundo orden, para este caso la entrada corresponde a una señal a la que se le ha añadido ruido. La figura 11c muestra los valores del MSE que ambos algoritmos (EMO y PSO) generaron durante las iteraciones.

Para este problema, ambos algoritmos presentan un rendimiento similar, sin embargo, el algoritmo EMO logra converger en un menor número de iteraciones y mantener un valor MSE menor que el algoritmo PSO. La evolución de los valores de ambas propuestas se presenta en la figura 12, en la cual se aprecia que tanto el método basado en EMO como el PSO, requieren un número reducido de iteraciones para encontrar los valores de los parámetros a₁, a₂ y b, con los cuales el MSE tiende a un valor cercano a cero.

Figura 12.

Comparación de los parámetros a y b encontrados por el algoritmo EMO durante 3000 iteraciones

Cabe señalar, que a pesar de que las dos técnicas muestran comportamientos similares, el método propuesto es el que mejor aproxima los valores óptimos, tal como se puede apreciar en la tabla 5.

Tabla 5.

Valores comparativos de los parámetros a,b y J encontrados tanto por el algoritmo EMO como el PSO

Valores			b	J
EMO	−1.3901	0.4802	1.0091	2.6027e-5
PSO	−1.4024	0.4925	0.9706	1.4113e-4
Reales	−1.4000	0.4900	1.0000	–

Finalmente se realizó la comparación de ambas técnicas de optimización con el problema que propone modelar una planta de orden 6 con un filtro IIR de orden 4. Del mismo modo que en las otras comparaciones, ambos algoritmos se ejecutan durante un número definido de iteraciones, durante las cuales se buscan los mejores valores que minimizan la función de error MSE. En este caso, el problema consta de 9 dimensiones.

La figura 13 presenta las evoluciones de cada parámetro del denominador del filtro IIR identificado, considerando ambos algoritmos.

Figura 13.

Gráfica de los parámetros del denominador del filtro de los algoritmos: a) EMO, b) PSO

Por otro lado, los parámetros del numerador del filtro, que tanto el algoritmo EMO como el PSO identificaron durante las iteraciones, se muestran en la figura 14.

Figura 14.

Gráfica de los parámetros del numerador del filtro del algoritmo: a) EMO, b) PSO

El error cuadrático medio de cada algoritmo es mostrado en la figura 15. Los valores encontrados en la última iteración se presentan en la tabla 6.

Figura 15.

Comparación de los valores de MSE obtenidos por los algoritmos EMO y PSO

Tabla 6.

Valores encontrados por el algoritmo EMO y PSO para los parámetros del numerador y los del denominado

	Valores EMO	Valores PSO		Valores EMO	Valores PSO
a1	−0.4950	−0.2327	b0	1.0335	0.7254
a2	−0.7049	−0.1496	b1	−0.6670	−0.1618
a3	0.5656	−0.0469	b2	−0.4682	0.3113
a4	−0.2691	−0.3720	b3	0.6961	0.0601

Conclusiones

Este artículo presentó un nuevo algoritmo para la identificación y modelado de filtros IIR, el cual está basado en EMO. Bajo este enfoque, la adaptación de los coeficientes es considerada como un problema de optimización, para la cual se utiliza como función objetivo el error cuadrático medio (MSE) existente entre la respuesta del sistema desconocido y la estimada por el filtro durante cierto periodo de tiempo. Guiados por los valores de la función objetivo, el conjunto de coeficientes que modelan al filtro son operados usando el algoritmo EMO hasta encontrar la mejor configuración posible (mínimo global).

El desempeño del algoritmo propuesto, fue comparado experimentalmente contra el algoritmo de identificación IIR basado en la técnica PSO (Lin et al., 2008), usando un conjunto de pruebas estándar (Jonson y Larimore, 1977). De los resultados se desprende que el algoritmo basado en EMO, entrega resultados más precisos en prácticamente el mismo número de iteraciones.

Aunque el algoritmo propuesto ofrece evidencia de un mejor desempeño, el objetivo de este artículo no es presentar un algoritmo que demerite las técnicas que actualmente son utilizadas para este fin, sino que tiene como propósito mostrar que los sistemas basados en principios físicos como el “Electromagnetism-Like”, se pueden considerar como una alternativa sumamente atractiva para la identificación de los filtros IIR.

Agradecimientos

El segundo autor agradece al Consejo Nacional de Ciencia y Tecnología (CONACyT), por la beca doctoral número 215517 y al Instituto Jalisciense de la Juventud por el apoyo brindado para realizar esta investigación.