Clasificación de fallos en motores en estado transitorio mediante redes neuronales
Abstract
Grado en Estadística
Full text
Facultad de Ciencias Clasificaci´on de fallos en motores en estado transitorio mediante redes neuronales Autor: Miguel Toquero Bar´on Tutores: Miguel Alejandro Fern´andez Temprano Alejandro Bar´on Garc´ıa Trabajo de fin de grado del Grado en Estad´ıstica Universidad de Valladolid Junio 2021
No hay nada que sea m´as grato Que ver como lo has logrado Despu´es de haber trabajado A˜nos de sol a sol Asi que sin ofenderte Mejor no me desees suerte Desea que sea persistente En mi lucha por ser mejor. El Chojin. Agradecimientos Quiero agradecer a todos los que me hab´eis apoyado y motivado para llegar hasta aqu´ı, los que d´ıa a d´ıa me anim´ais y dais fuerza para sacar lo mejor de m´ı. 1
Resumen Los motores el´ectricos de inducci´on son una pieza clave en el desarrollo industrial en la actualidad. Cuando un motor falla requiere la detenci´on de toda la producci´on para repararlo o cambiarlo. De ello, nace la necesidad de detectar los fallos antes de que ocurran, a ser posible mediante t´ecnicas no invasivas que permitan la monitorizaci´on del motor mientras est´a funcionando. En este trabajo de fin de grado se considera un problema de clasificaci´on de fallos en motores mediante los datos del arm´onico de la onda superior e inferior de la corriente de alimentaci´on del motor en una situaci´on que, seg´un los expertos, es especialmente compleja ya que se considera la corriente durante el estado transitorio del motor, motores alimentados por distintos inversores, lo cual tambi´en complica la clasificaci´on, y se efect´ua una clasificaci´on multiestados, ya que se clasifica el estado del motor en cinco estados diferentes seg´un el nivel de degradaci´on del motor. Se utilizan t´ecnicas de aprendizaje profundo, como son las redes neuronales, comprobando su desempe˜no en el problema. Se obtiene que, en este problema, las redes neuronales necesitan, al igual que suced´ıa con las t´ecnicas boosting, un preprocesado de los datos para una mejor soluci´on del problema. La precisi´on de la clasificaci´on obtenida es considerablemente buena, con valores similares a las t´ecnicas boosting. Se obtienen adem´as conclusiones interesantes desde el punto de vista industrial, como la confirmaci´on de que el nivel de carga alto permite mayor precisi´on en la clasificaci´on, y que el arm´onico de la banda inferior resulta m´as informativo que el arm´onico de la banda superior. Tambi´en se hace uso de las t´ecnicas recientemente desarrolladas en la literatura que permiten la interpretabilidad de este tipo de modelos de caja negra, posibilitando extraer informaci´on sobre c´omo se realizaron las predicciones y cuantificar de esta manera la intuici´on previa que se ten´ıa sobre cu´ales son las variables m´as interesantes en el problema. 2
Abstract Nowadays, electric engines are a key factor in industrial development. Whenever there is an engine failure, the whole production needs to stop in order to repair or change the piece. Hence, there is an arising need to detect faults before they occur. If possible, non-invasive techniques are preferred since they allow monitoring the engine while it is running. In this Bachelor thesis, we evaluate an engine failure classification problem through both upper-side and lower-side harmonic band power supply. According to the experts, the considered situation is especially complex due to several factors: power supply during the transitory state of the engine is assessed; inverter-fed induction engines are studied, a fact that complicates the analysis; and a multi-state classification takes place, which means the engine is categorized in five different states according to its degradation level. Deep Learning techniques, such as neural networks, are applied, and their performance solving the problem is evaluated. Obtained results show neural networks need, same as boosting techniques, pre-processed data to improve the solution to the problem. Accuracy of the obtained classification is remarkable, with similar values to those obtained with boosting techniques. Furthermore, we also get interesting conclusions from an industrial point of view, such as the confirmation that a high charge level allows better precision in the classification task and that the lower-side harmonic band turns to be less informative than the upper-side harmonic band. In addition, recently developed techniques that have been previously described in the literature are adopted to interpret the black box model, enabling the extraction of the information on how predictions were made. In such manner, our thoughts on which were the most useful variables to solve the problem can be quantified. 3
´ Indice general 1 Introducci´on 6 1.1 Descripci´on del problema . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 1.2 Objetivosdeltrabajo............................... 7 1.3 Estructura del documento . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 1.4 Asignaturas relacionadas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 1.4.1 Ampliaci´on de materia . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2 Metodolog´ıa 10 2.1 El problema de la clasificaci´on . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.1.1 Discriminante Lineal de Fisher . . . . . . . . . . . . . . . . . . . . . . 10 2.1.2 Regresi´on log´ıstica . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.1.3 ´ Arbolesdedecisi´on............................ 12 2.1.4 Redesneuronales ............................. 13 2.2 Redesneuronales ................................. 13 2.2.1 Neurona.................................. 13 2.2.2 Funciones de activaci´on . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.2.3 Arquitectura o topolog´ıa . . . . . . . . . . . . . . . . . . . . . . . . . 23 2.2.4 Entrenamiento .............................. 25 2.2.5 Funci´on de coste o p´erdida . . . . . . . . . . . . . . . . . . . . . . . . 25 2.2.6 Descenso del gradiente . . . . . . . . . . . . . . . . . . . . . . . . . . 26 2.2.7 Backpropagation-Retropropagaci´on del error . . . . . . . . . . . . . . 28 2.2.8 Regularizaci´on............................... 31 2.2.9 Hiperpar´ametros ............................. 36 2.3 Estimaci´ondelerror ............................... 37 2.3.1 Tarea de inducci´on de un clasificador . . . . . . . . . . . . . . . . . . 38 2.3.2 Estratificaci´on............................... 40 2.3.3 Hold-out.................................. 40 2.3.4 K-fold ................................... 41 2.4 Interpretabilidad de modelos complejos . . . . . . . . . . . . . . . . . . . . . 42 2.4.1 LIME ................................... 43 2.4.2 Shapleyvalues............................... 45 2.4.3 SHAP ................................... 47 4
3 Datos 50 3.1 Descripci´on .................................... 50 3.2 Procesamiento................................... 52 3.2.1 Conjunto de datos alternativo . . . . . . . . . . . . . . . . . . . . . . 54 3.2.2 Subconjuntos de datos . . . . . . . . . . . . . . . . . . . . . . . . . . 54 3.3 Terminolog´ıa ................................... 55 4 An´alisis y resultados 56 4.1 Experimentorealizado .............................. 56 4.2 Tasasdeerror................................... 57 4.3 An´alisisderesultados............................... 61 4.4 Clasificaci´on seg´un modelo de inversor . . . . . . . . . . . . . . . . . . . . . 65 4.4.1 InversorAB................................ 65 4.4.2 InversorABB............................... 70 4.4.3 InversorTM................................ 76 4.5 Comparaci´on de m´etodos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 5 Conclusiones y trabajo futuro 82 5.1 Conclusiones.................................... 82 5.2 Trabajo futuro y posibles mejoras . . . . . . . . . . . . . . . . . . . . . . . . 83 Bibliograf´ıa 85 A Anexos 89 A.1 C´odigo....................................... 89 ´ Indice de figuras 108 ´ Indice de tablas 110 5
Cap´ıtulo 1 Introducci´on 1.1 Descripci´on del problema Los motores de inducci´on son m´aquinas el´ectricas que se encargan de obtener energ´ıa mec´anica a partir de energ´ıa el´ectrica por medio de la rotaci´on de los campos magn´eticos generados en sus bobinas. En 2015 el 80 % de los motores utilizados en la industria eran motores de inducci´on [1]. La rotura de uno de los componentes de un motor de inducci´on resulta crucial en el proceso de producci´on, ya que este suceso conlleva el cese de la actividad de producci´on para reparar la pieza, lo que puede suponer un gran coste en la evaluaci´on de p´erdidas. Por todo esto, est´a en alza el conocido como mantenimiento predictivo, es decir la investigaci´on sobre t´ecnicas y m´etodos de diagn´ostico previos al fallo. Una de estas t´ecnicas es el denominado Motor Current Signature Analisis, MCSA [2], que consiste en la evaluaci´on del estado del motor, a partir de sus corrientes de alimentaci´on. Su principal ventaja es ser una t´ecnica no invasiva, es decir que no requiere interrumpir el funcionamiento para llevar a cabo el an´alisis. Dentro de este an´alisis de corriente del motor, en este trabajo de fin de Grado, se utilizan las frecuencias de la tensi´on de alimentaci´on del motor upper-sideband harmonic, USH, y lower-sideband harmonic, LSH, se llevar´a a cabo el proceso de detecci´on y diagn´ostico. El problema que se considera en este trabajo es particularmente interesante puesto que tiene varias caracter´ısticas que lo complican. En primer lugar, se consideran motores alimentados por inversores, que es la situaci´on habitual en la pr´actica, y no motores alimentados directamente de la red. Esto hace que la se˜nal sea menos limpia y m´as dif´ıcil de clasificar. En segundo lugar, se considera la se˜nal emitida por el motor en estado transitorio, es decir, cuando est´a comenzando a funcionar y no cuando ya est´a estabilizada (estado estacionario) lo que requiere una transformaci´on funcional de los datos m´as compleja y tambi´en complica el problema de la clasificaci´on. Adicionalmente se considera la clasificaci´on en cinco estados diferentes de deterioro y no simplemente la clasificaci´on en motor sano o motor da˜nado, que 6
ser´ıa mucho m´as sencilla de conseguir. Tambi´en se consideran varias covariables en el problema, como tipo de inversor o el nivel de carga, cuya influencia en la calidad de la clasificaci´on de los fallos se quiere valorar para poder determinar si puede establecerse un ´unico modelo de clasificaci´on o deben considerarse varios en funci´on de los valores de esas covariables. Tras observar en trabajos previos el buen desempe˜no de las t´ecnicas boosting en la detecci´on y clasificaci´on de este tipo de fallos [3], en este trabajo se quiere comprobar si las redes neuronales son capaces de mejorar los resultados obtenidos con esa t´ecnica. Se tratar´a de entrenar y obtener modelos v´alidos con los que llevar a cabo el diagn´ostico con un alto grado de acierto, explorando, en primer lugar, la metodolog´ıa de redes neuronales y su comportamiento ante este problema, para, a continuaci´on, comprobar si en este problema las redes neuronales son capaces de seleccionar las variables m´as relevantes para clasificar los fallos o si les es conveniente la ayuda proporcionada por una reducci´on de dimensi´on en los datos del problema. A continuaci´on se establecer´an uno o varios modelos de clasificaci´on que permitan obtener un diagn´ostico preciso y comprobar hip´otesis anteriores de expertos y otros trabajos (como que niveles de carga superiores mejoran la precisi´on en el diagn´ostico, que cada modelo de inversor necesita un clasificador propio o si la informaci´on que porporcionan las bandas LSH y USH es o no redundante). Hay alguna referencia en la literatura cient´ıfica del uso de redes neuronales en la clasificaci´on de fallos en motores de inducci´on. Tal vez la referencia m´as interesante sea [4] donde utilizando un m´etodo basado en wavelets y redes neuronales se analiza la corriente en estado estacionario de motores conectados directamente a la red y se clasifican dichos motores en dos estados (motor da˜nado o motor sano) logrando un alto porcentaje de acierto. En este trabajo se va m´as all´a puesto que se consideran motores alimentados por distintos inversores y se analiza la se˜nal de su estado transitorio mediante una transformada funcional m´as compleja y se clasifican los motores en cinco estados diferentes dependiendo del grado de deterioro del motor. No conocemos referencias en la literatura en las que se trate este tipo de problema mediante el procedimiento que aqu´ı se desarrolla. Empero, la modelizaci´on estad´ıstica no se realiza ´unicamente con objeto de predicci´on, si no de adquirir conocimiento en base a los datos en base a una representaci´on (modelo) de la realidad. Es por esto que, mediante t´ecnicas de interpretaci´on, se pueden extraer los conceptos aprendidos por el clasificador de forma comprensible por el ser humano. Esto permite no solo verificar un correcto aprendizaje del clasificador si no que permite extender el conocimiento previo que se tiene sobre un dominio de aplicaci´on, en este caso, del MCSA. 1.2 Objetivos del trabajo Las redes neuronales actualmente est´an presentes en casi todos los ´ambitos hasta donde llega la tecnolog´ıa. Su auge en la actualidad nos motiva a tratar de resolver este problema con su utilizaci´on, esperando obtener unos buenos resultados. Los objetivos que se esperan alcanzar en este trabajo son los siguientes: 7
•Valorar la capacidad de las redes neuronales para abordar el problema de clasificaci´on de fallos de motores anteriormente descrito. •Comprobar si las redes neuronales son capaces de seleccionar las variables m´as interesantes de los datos para la soluci´on del problema o si necesitan de un preprocesado de los datos que reduzca su dimensionalidad. •Comprobar si las redes neuronales mejoran a los m´etodos boosting. •Conocer qu´e factores o covariables influyen en la clasificaci´on. •Obtenci´on de modelos interpretables sobre modelos opacos o de caja negra a trav´es de t´ecnicas punteras como SHAP. A trav´es de estos, verificar qu´e conceptos ha aprendido el clasificador. 1.3 Estructura del documento Esta memoria se compone de los siguientes cap´ıtulos: •Metodolog´ıa. Se exponen los conceptos te´oricos y el fundamento de los m´etodos utilizados en este TFG. •Datos. Descripci´on, obtenci´on y procesamiento del conjunto de datos. •An´alisis y resultados. Se desarrolla un an´alisis de los resultados obtenidos y se comparan los resultados obtenidos en este TFG con otros resultados previos. •Conclusiones. Se resumen los resultados obtenidos en este trabajo y se proponen nuevas lineas de trabajo o mejora. 1.4 Asignaturas relacionadas En esta secci´on veremos la relaci´on de t´ecnicas utilizadas en este trabajo y el aprendizaje de estas en diferentes asignaturas del grado. •An´alisis de Datos, An´alisis Multivariante y An´alisis de Datos Categ´oricos: en ellas se estudia el fundamento y las bases de los clasificadores y modelos predictivos, as´ı como las m´etricas para la evaluaci´on de los clasificadores. •Fundamentos de Programaci´on, Programaci´on Orientada a Objetos y Paradigmas de Programaci´on: se asientan las bases sobre la programaci´on en distintos lenguajes, entre ellos el utilizado en este trabajo (Python). •Regresi´on y ANOVA, Modelos Lineales y Modelos Estad´ısticos Avanzados: en los que se estudia las t´ecnicas de an´alisis de la varianza. 8
Figura 2.6: Puerta XOR. Podemos comprobar que una sola neurona artificial no es capaz de representar la funci´on boleana XOR ni ninguna otra que no sea linealmente separable. Debemos fijarnos en que los pesos de cada variable de entrada son iguales entre ellos e iguales a uno. Partiendo de esta idea, Frank Rosenblantt invent´o el perceptr´on simple [11] o la neurona artificial tal y como la conocemos hoy en d´ıa. El fundamento es similar al empleado por McCulloch y Pitts, sin embargo, ahora las restricciones son menos estrictas: las entradas tienen unos pesos, de forma que juegan el papel de una regresi´on en la que no todas las entradas tienen la misma contribuci´on, el valor de la inhibici´on deja de ser absoluto. w1·x1+w2·x2+... +wN·xN+b(2.6) Encontramos un s´ımbolo a˜nadido a la suma b, este es el bias o lo que hemos llamado anteriormente Θ, que como ya se ha mencionado es el l´ımite de la activaci´on. Estas modificaciones dieron paso al uso de la neurona artificial para aprendizaje supervisado, permitiendo a la neurona aprender los pesos de las variables de entrada por s´ı misma. Figura 2.7: Neurona artificial 15
Como cada variable xide entrada tiene un peso asociado, nos permite modificar su influencia en el resultado final libremente. Tras realizar la suma la se˜nal o resultado se propaga a la funci´on de activaci´on, la cual veremos a continuaci´on con m´as detalle. Por todo esto, podemos definir la neurona artificial como nodo computacional, basado en el comportamiento de las neuronas biol´ogicas humanas, que a partir de unas entradas num´ericas procesadas con una funci´on matem´atica devuelven un valor de salida tambi´en num´erico. La uni´on de varias neuronas artificiales da lugar a las redes neuronales, como veremos m´as adelante. 2.2.2 Funciones de activaci´on Hasta ahora nada diferencia el modelo de la neurona artificial de una regresi´on m´ultiple. Sin embargo, la funci´on de activaci´on lo cambia todo. Esta es la encargada de introducir no-linealidades en el modelo. De esta manera, se permite aprender funciones mucho m´as complejas y adaptarse a aquellas que no sean linealmente separables. Podemos decir que una funci´on de activaci´on es una funci´on matem´atica aplicada a la salida de la neurona para modificar y a˜nadir deformaciones no lineales [12]. Estas son algunas de las propiedades que buscamos en una funci´on para poder utilizarla como funci´on de activaci´on: •No lineal: introduce irregularidades en el modelo para poder aproximar una gran variedad de funciones. •Diferenciable: interesa que sea diferenciable en todos los puntos para poder aplicar m´etodos de optimizaci´on basados en las derivadas, descenso del gradiente. •Rango: acota los valores de salida o los transforma dentro de un rango de valores. •Mon´otona: cuando la funci´on de activaci´on es mon´otona se garantiza que la superficie de error asociada a un modelo monocapa es convexa y, por tanto, el ´optimo local coincide con el ´optimo global. •Simplicidad: deben ser simples para minimizar el coste computacional. •Se aproximan a la identidad en el origen: necesario para evitar la influencia de los pesos iniciales aleatorios. Funci´on identidad Es una funci´on lineal que obtiene el mismo valor que su entrada, retorna un valor id´entico. f(x) = x(2.7) 16
Figura 2.8: Funci´on identidad. Esta funci´on, aunque puede ser ´util en otros problemas, no se adapta a las redes neuronales, ya que no introduce no-linearidades. Por tanto, no nos permitir´ıa adaptarnos a problemas que no sean linealmente separables. Si tuvi´esemos un modelo de red neuronal muy complejo con muchas capas y neuronas pero todas ellas tuviesen esta funci´on de activaci´on, ser´ıa equivalente a tener un modelo con una sola neurona. Funci´on signo En las neuronas presentadas previamente se ha utilizado como funci´on de activaci´on la funci´on signo. sgn(x) = +1 si x ≥0 −1en otro caso. (2.8) Transforma la entrada en una salida binaria en funci´on del signo. Esto nos da libertad ya que el bias, que hemos presentado con anterioridad, ser´ıa el par´ametro que nos permitir´ıa ajustar el limite de los dos valores de salida. 17
Figura 2.9: Funci´on signo. Funci´on ReLU Se conoce con el nombre de rectificador lineal(Rectified Linear Unit). Es la funci´on de activaci´on m´as usada [13]. Tiene una salida igual a cero si su entrada es negativa e igual a la entrada si esta es positiva, es estrictamente creciente. f(x) = x si x ≥0 0en otro caso. (2.9) Tambi´en se puede expresar como: f(x) = max(0, x) (2.10) La velocidad de construcci´on, entrenamiento, de modelos con esta funci´on de activaci´on es mayor que otras como la sigmoide. Tiene activaci´on dispersa, es decir, si las neuronas toman valores iniciales aleatorios solo la mitad de ellas se activar´ıan. Tambi´en tiene desventajas, como que no est´a centrada en cero y no es diferenciable en 0, pero s´ı en todos los dem´as puntos. Tambi´en encontramos la dificultad de que todos los puntos a la izquierda del cero toman el mismo valor, lo cual nos impide encontrar diferencias. Por desgracia, las unidades ReLU pueden ser fr´agiles durante el entrenamiento y pueden ”morir”. Por ejemplo, un gran gradiente que fluya a trav´es de una neurona ReLU puede hacer que los pesos se actualicen de tal manera que la neurona no vuelva a activarse en ning´un punto de datos. Si esto ocurre, 18
el gradiente que fluye a trav´es de la unidad ser´a siempre cero a partir de ese punto. Dying ReLu se refiere al problema cuando las neuronas ReLU se vuelven inactivas y s´olo dan salida a 0 para cualquier entrada [14]. Podemos observar que no es una funci´on lineal, es lineal a trozos. Figura 2.10: Funci´on ReLU Una de sus mayores desventajas es que no est´a acotada. Suele ser ´util en el reconocimiento de im´agenes y redes convolucionales. Funci´on LeakyReLU Utilizando la idea anterior del rectificador lineal nace esta variante. Esta funci´on proporciona una correcci´on para los valores negativos de 2.9, de forma que su valor no sea cero sino un valor muy peque˜no que depende de la entrada. f(x) = x si x ≥0 0.01x en otro caso. (2.11) Esto corrige el problema de que todos los valores negativos tomasen el mismo valor igual a cero y el problema de Dying ReLu. 19
Figura 2.11: Funci´on Leaky ReLU Existe una versi´on param´etrica en la que el coeficiente para la x, cuando es negativa, es un valor muy peque˜no denominado α. Al igual que la funci´on ReLU, no est´a acotada y se comporta bien con im´agenes y redes convolucionales. Funci´on sigmoide Tambi´en conocida como funci´on log´ıstica. Es especialmente ´util en problemas de clasificaci´on binaria ya que su rango de valores es (0,1), siempre toma valores positivos. Por tanto, el resultado puede interpretarse como una probabilidad. Es una funci´on estrictamente creciente y derivable, tiene una derivada no nula en cada punto, es diferenciable. Por esto, el m´etodo del descenso del gradiente puede lograr un mejor resultado en cada paso de la optimizaci´on. f(x) = 1 1 + e−x(2.12) 20
Figura 2.12: Funci´on sigmoide Esta funci´on tiene lenta convergencia, satura o mata el gradiente, sin embargo, tiene un buen rendimiento en la ´ultima capa por su alta interpretabilidad. Funci´on softmax Cuando tratamos con problemas de clasificaci´on multiclase utilizamos la funci´on softmax. Con una base matem´atica muy similar a la de la funci´on sigmoide vista en 2.12, esta funci´on puede tener tantas salidas como se deseen. Por lo general, esto se iguala al n´umero de clases. De esta manera, cada salida con una rango de valores (0,1) se interpreta como probabilidad de pertenencia a esa clase. Esta funci´on est´a definida para que la suma de todas sus salidas se igual a uno. Para kclases tenemos la siguiente formulaci´on: fi(yi) = eyi Pk j=0 eyi (2.13) 21
Figura 2.13: Funci´on de activaci´on softmax. Imagen de [15] Funci´on tangente hiperb´olica Tambi´en se parece a la funci´on sigmoide. Esta toma valores en el rango (−1,1). Esta centrada en torno al cero y es estrictamente creciente. En la pr´actica la optimizaci´on es m´as f´acil que para la funci´on log´ıstica pero sufre del problema del desvanecimiento del gradiente, problema que se comentar´a m´as adelante. f(x) = tanh(x) = ex−e−x ex+e−x(2.14) En ocasiones se usa esta otra reformulaci´on, que es equivalente, pero m´as eficiente computacionalmente al reducir el calculo de exponenciales f(x) = tanh(x) = e2x−1 e2x+ 1 (2.15) 22
Figura 2.14: Funci´on tangente hiperb´olica Al igual que la funci´on sigmoide tiene lenta convergencia y satura o mata el gradiente. Se utiliza para decidir entre una opci´on y la contraria por su rango de valores. Tiene buen desempe˜no en redes recurrentes. Existen muchas otras funciones de activaci´on y variantes de las aqu´ı presentadas. Sin embargo, en este trabajo solo se usar´an las funciones ya mencionadas por su suficiencia, utilidad y popularidad. Las funciones de activaci´on consituyen uno de los hiperpar´ametros que se deben ajustar en el modelo. 2.2.3 Arquitectura o topolog´ıa Como ya hemos visto, las redes neuronales son combinaciones de neuronas. Nos referimos a topolog´ıa o arquitectura cuando hablamos de la forma en que se organizan las neuronas dentro de una red neuronal. En este trabajo se desarrollar´an las redes densamente conectadas, con propagaci´on hacia delante. En este tipo de redes, cada neurona toma como entrada las salidas de las neuronas de la capa inmediatamente anterior y su salida se propaga a las todas las neuronas de la capa inmediatamente posterior. Las neuronas se agrupan en capas seg´un su funci´on y su situaci´on dentro de la red. Podemos diferenciar tres capas: •Capa de entrada (input layer): Es la capa que conecta nuestros datos con la red neuronal. Cada neurona perteneciente a esta capa tiene tantas entradas como variables tiene el conjunto de datos utilizado. •Capa de salida (output layer): Recibe como entrada las salidas de las neuronas de la capa anterior y tiene como salida la salida final del sistema. Juega un papel crucial ya 23
que el n´umero de neuronas en esta capa y el tipo de salida que se desea obtener tiene que adaptarse al problema tratado. •Capas ocultas (hidden layers): Todas las capas que se sit´uan entre la capa de entrada y la capa de salida se denominan capas ocultas. Suelen llevar el peso computacional por su gran cantidad de neuronas y porque el n´umero de capas es ilimitado. Figura 2.15: Arquitectura de una red neuronal En la imagen se muestra la topolog´ıa de una red neuronal con una ´unica capa oculta pero podr´ıa tener tantas como se deseen, es decir, puede tener un n´umero arbitrario de capas y un n´umero arbitrario de neuronas dentro de cada capa. El n´umero de neuronas en cada capa no est´a limitado y se puede adaptar a las condiciones del problema. La posibilidad de encadenar capas con grandes n´umeros de neuronas proporciona la posibilidad de representar funciones muy complejas, de esto nace el t´ermino de aprendizaje profundo [16]. Encontramos aqu´ı dos hiperpar´ametros muy interesantes y ´utiles para modelar diferentes problemas. El teorema de aproximaci´on universal de Hornik nos dice que una red neuronal de una sola capa oculta puede aproximar arbitrariamente bien a cualquier funci´on continua si se le da un n´umero de neuronas suficientes [17]. Como ya hemos comentado, es crucial que el n´umero de neuronas y funciones de activaci´on en la capa de salida se elija cuidadosamente para que se adapte a nuestro problema. Por ejemplo, si trat´asemos con una predicci´on de un valor que juega el papel de tiempo podr´ıamos usar una sola neurona con una funci´on de activaci´on ReLU, que no tiene valores negativos. Sin embargo, en un problema de clasificaci´on binaria deber´ıamos usar una neurona pero con funci´on de activaci´on sigmoide para poder interpretar la salida como probabilidad de pertenencia a cada clase. Como ´ultimo ejemplo, si nos encontramos con un problema de clasificaci´on multi-clase debemos utilizar tantas neuronas como clases haya y la funci´on de activaci´on softmax, de forma que tambi´en se pueda interpretar la salida como probabilidades de pertenencia. 24
Se denomina gradiente explosivo cuando tenemos la situaci´on contraria, valores excesivamente grandes de la funci´on gradiente. Optimizadores Con el conocimiento obtenido a cerca de la tarea de optimizaci´on de los par´ametros en la etapa del entrenamiento, veremos a continuaci´on que t´ecnicas se van a utilizar en la pr´actica. No debemos olvidar que todas estas t´ecnicas hacen uso del descenso del gradiente y del algoritmo de retropropagaci´on. •Descenso del gradiente estoc´astico o SGD [22]: optimizador con descenso de gradiente y momento. Puede incluirse la aceleraci´on de Nesterov [23]. •RMSprop [24]: mantiene una media m´ovil del cuadrado de los gradientes y divide el gradiente por la ra´ız de esta media. Esta implementaci´on de RMSprop utiliza el impulso simple, no el impulso Nesterov. Utiliza esa media m´ovil para estimar la varianza. •Adam (Adaptative moment estimation) [25]: la optimizaci´on de Adam es un m´etodo de descenso de gradiente estoc´astico que se basa en la estimaci´on adaptativa de momentos de primer y segundo orden. Seg´un [26], el m´etodo es eficiente desde el punto de vista computacional, tiene pocos requisitos de memoria, es invariable al reescalado diagonal de los gradientes y se adapta bien a los problemas que son grandes en t´erminos de datos/par´ametros. Es interesante destacar que aunque la convergencia con Adam es m´as r´apida, el algoritmo SGD generaliza mejor [27]. 2.2.8 Regularizaci´on La regularizaci´on es un conjunto de t´ecnicas que influye en el aprendizaje para que el algoritmo generalice mejor. Dicho de otro modo, son t´ecnicas ´utiles para evitar o reducir el sobreajuste. De esta manera, nuestro algoritmo se adaptar´a mejor a los datos que nunca se han visto, as´ı mejorar´a la precisi´on cuando se enfrente a datos completamente nuevos del dominio del problema. La regularizaci´on consiste en penalizar los coeficientes de los pesos en los nodos de forma que estos sean m´as suaves. Es decir, les impide o dificulta tomar valores extremos. Se penalizar´a m´as o menos seg´un un valor que ser´a el coeficiente de regularizaci´on, par´ametro que tambi´en se debe optimizar y ser´a un hiperpar´ametro m´as de nuestro modelo. Sobreajuste Uno de los aspectos m´as importantes a la hora de entrenar redes neuronales es evitar el sobreajuste. Como es un aspecto muy conocido aqui haremos una breve recapitulaci´on. El sobreajuste se refiere al fen´omeno en el que una red neuronal modela muy bien los datos de entrenamiento pero falla cuando ve nuevos datos del mismo dominio del problema. El 31
sobreajuste est´a causado por el ruido en los datos de entrenamiento que la red neuronal capta durante el entrenamiento y lo aprende como un concepto subyacente de los datos. Este ruido aprendido, sin embargo, es ´unico para cada conjunto de entrenamiento. En cuanto el modelo ve nuevos datos del mismo dominio del problema, pero que no contienen este ruido, el rendimiento de la red neuronal empeora mucho. Cuando las redes neuronales alcanzan suficiente complejidad se adaptan perfectamente a los datos de entrenamiento aprendiendo asi su ruido en el modelo. Esto significa que la red neuronal en un momento determinado del periodo de entrenamiento no mejora su capacidad de resolver el problema sino que simplemente empieza a aprender alguna regularidad aleatoria contenida en el conjunto de patrones de entrenamiento [28]. Figura 2.17: Ejemplo de sobreajuste. Regularizaci´on L2 Esta t´ecnica se conoce como decaimiento del peso, weight decay, o regresi´on Ridge. funcion de coste =funcion de coste +λ p X j=1 β2 j(2.31) Cabe mencionar que: λes el par´ametro que decide cuanto queremos penalizar y βson los coeficientes de los pesos que anteriormente hemos denominado w. El t´ermino a˜nadido con lambda es la regularizaci´on, en este caso se utiliza la norma eucl´ıdea sobre los pesos de las conexiones entre neuronas. 32
Cuando λ= 0 la regularizaci´on queda sin efecto y se obtiene el mismo resultado que si no se aplicase. Sin embargo, cuando λ→ ∞ la penalizaci´on aumenta y los estimadores de los coeficientes tienden a 0. Regularizaci´on L1 La regularizaci´on L1, tambi´en conocida como regresi´on Lasso, var´ıa un poco respecto de la regularizaci´on L2, cambia el t´ermino de regularizaci´on. funcion de coste =funcion de coste +λ p X j=1 |βj|(2.32) Este realiza selecci´on de variables ya que los coeficientes pueden tomar valor 0, produce modelos dispersos. Figura 2.18: Comparaci´on L1(izquierda) y L2(derecha). Imagen de [29] Intuitivamente, los pesos m´as peque˜nos reducen el impacto de las neuronas ocultas. En ese caso, esas neuronas ocultas se vuelven despreciables y la complejidad general de la red neuronal se reduce en ambas regularizaciones L1 y L2. Los modelos menos complejos suelen evitar modelar el ruido en los datos y, por tanto, evitar el sobreajuste. Hay que tener cuidado a la hora de elegir el t´ermino de regularizaci´on λ. El objetivo es encontrar el equilibrio adecuado entre la baja complejidad del modelo y la precisi´on. Si el valor de λes demasiado alto, el modelo ser´a sencillo, pero corre el riesgo de que infraajuste a los datos. El modelo no aprender´a lo suficiente sobre los datos de entrenamiento para hacer predicciones ´utiles. Si el valor de λes demasiado bajo, el modelo ser´a m´as complejo y se corre el riesgo de sobreajustar los datos. Aprender´a demasiado sobre las particularidades de los datos de entrenamiento y no podr´a generalizar a nuevos datos. En keras podemos realizar la regularizaci´on en distintos puntos [30]. 33
•kernel regularizer: Regularizador para aplicar una penalizaci´on en el n´ucleo de la capa. Esto es, la funci´on con los pesos que realiza la suma ponderada. •bias regularizer: Regularizador para aplicar una penalizaci´on sobre el sesgo de la capa. El t´ermino de sesgo, bias o t´ermino independiente que hemos denominado b. •activity regularizer: Regularizador para aplicar una penalizaci´on a la salida de la capa. En la salida de la neurona tras pasarla por la funci´on de activaci´on. Dropout Es una t´ecnica de regularizaci´on relativamente sencilla. Consiste en no actualizar los pesos de una neurona, con probabilidad p, en una etapa de entrenamiento. De esta manera, se evita que un nodo o neurona tenga demasiada responsabilidad en la tarea de la clasificaci´on [31]. Figura 2.19: A la izquierda: modelo de red neuronal con dos capas ocultas. A la derecha: ejemplo de red regularizada con dropout. Imagen de [32] En cada iteraci´on el conjunto de neuronas que actualizan sus pesos es diferente. Esto permite capturar m´as aleatoriedad y generalizar mejor. La probabilidad de que una neurona de una capa se actualice o no es otro hiperpar´ametro que debemos seleccionar. Esta probabilidad es la misma para todas las neuronas de una misma capa. Parada temprana Conocido como early stop. Entrenar demasiado tiempo puede provocar un sobreajuste, lo que significa que el ´optimo local de la red neuronal s´olo puede ser eficaz en sus datos de entrenamiento. La soluci´on es simplemente entrenar la red neuronal durante menos tiempo. Es dif´ıcil decidir cu´ando es mejor dejar de formarse si s´olo se observa la curva de aprendizaje de entrenamiento por s´ı misma. El inicio del sobreajuste puede detectarse mediante una validaci´on cruzada en la que los datos disponibles se dividen en subconjuntos de entrenamiento validaci´on y prueba. El subconjunto de entrenamiento se utiliza para calcular el gradiente y actualizar los pesos de la red. Se monitoriza el error de validaci´on con este conjunto durante el entrenamiento. 34
El error de validaci´on suele disminuir durante la fase inicial del entrenamiento al igual que el error en el conjunto de entrenamiento. Sin embargo, cuando la red empieza a sobreajustar los datos, el error en el conjunto de validaci´on comienza a aumentar. Cuando el error de validaci´on aumenta durante un n´umero determinado de iteraciones el entrenamiento se detiene, y se devuelven los pesos al m´ınimo del error de validaci´on. Se detiene el entrenamiento cuando el error de validaci´on del modelo no haya mejorado de forma apreciable durante un n´umero determinado de ´epocas. Este es el proceso b´asico para utilizar la parada temprana: 1. Se pone un contador a 0 y se selecciona una paciencia (un valor entero para el n´umero de ´epocas que est´a dispuesto a esperar a que el modelo mejore antes de terminar el entrenamiento). 2. Se entrena la red neuronal durante una ´epoca. 3. Se eval´ua su rendimiento frente a un conjunto de validaci´on reservado. 4. Si el modelo no tiene m´as rendimiento que en la ´ultima ´epoca, incrementa el contador. 5. Si el contador es igual a su paciencia, se detiene el entrenamiento. en caso contrario, se vuelve al paso 1. ´ Esta es una t´ecnica sencilla que puede mejorar el rendimiento de la red neuronal con un m´ınimo esfuerzo[33]. Figura 2.20: Detenci´on temprana del entrenamiento. 35
A continuaci´on veremos las ventajas e inconvenientes que supone usar esta t´ecnica [28]. Ventajas de la detenci´on temprana: •Frente a otras regularizaciones, solo se lleva a cabo una vez el proceso de entrenamiento, no una vez para cada posible λ. •La medida de rendimiento de la validaci´on cruzada se aplica directamente al conjunto real de par´ametros de la red que va a utilizar. Destacamos las siguientes desventajas: •M´etodo ad hoc. •Depende de los detalles del m´etodo de optimizaci´on. •Algunos de los datos de entrenamiento se utilizan s´olo para decidir cu´ando parar y esto parece un desperdicio. 2.2.9 Hiperpar´ametros Ya hemos visto todos los hiperpar´ametros necesarios para ajustar una red neuronal. Tambi´en hemos visto que los pesos de las neuronas se entrenan autom´aticamente con algoritmos de optimizaci´on. Ahora debemos elegir que hiperpar´ametros usar en nuestro modelo. •N´umero de capas ocultas. •N´umero de neuronas en cada capa. •Funci´on de activaci´on en cada capa. •Dropout en cada capa. •Regularizaci´on en cada capa. •Optimizador. •Tama˜no de lote (batch). •N´umero de ´epocas. Con estos hiperpar´ametros definimos un espacio de b´usqueda. Limitamos el n´umero m´aximo de capas ocultas a 3. El n´umero de neuronas en cada capa ser´a un m´ultiplo de 2 entre 2 y 32. Las funciones de activaci´on podr´an ser alguna de las definidas en Funciones de activaci´on. El dropout tomar´a un valor distinto en cada capa entre 0 y 0.3 con saltos de 0.05. El coeficiente de regularizaci´on tomar´a valores entre 0 y 0.05. Se elige uno de los optimizadores comentados en la secci´on 2.2.7. El n´umero de ´epocas se establece entre 30 y 200. El tama˜no de lote se busca en las potencias de 2 hasta un m´aximo de 25. 36
B´usqueda aleatoria La b´usqueda aleatoria o Random Search sustituye la enumeraci´on exhaustiva de todas las combinaciones por una selecci´on aleatoria. Esto puede aplicarse de forma sencilla al entorno discreto descrito anteriormente, pero tambi´en se generaliza a los espacios continuos y mixtos. Puede superar a la b´usqueda en cuadriculada (Grid Search), especialmente cuando s´olo un peque˜no n´umero de hiperpar´ametros afecta al rendimiento final del algoritmo de aprendizaje autom´atico [34]. Este art´ıculo muestra como la b´usqueda aleatoria es m´as eficiente en la elecci´on de hiperpar´ametros que la b´usqueda en cuadr´ıcula. La b´usqueda aleatoria funciona desplaz´andose de forma iterativa a mejores posiciones en el espacio de b´usqueda, que se muestrean a partir de una hiperesfera que rodea la posici´on actual [35]. Los experimentos aleatorios son m´as eficientes, respecto al tiempo, que los experimentos de cuadr´ıcula para la optimizaci´on de hiperpar´ametros porque no todos los hiperpar´ametros son igualmente importantes. Los experimentos de b´usqueda en cuadr´ıcula asignan demasiados ensayos a la b´usqueda en dimensiones que no importan y sufren de una pobre cobertura en las dimensiones que son importantes, aunque, al explorar todo el espacio en una cantidad de tiempo mucho mayor, logran encontrar mejores resultados. Cuando se utiliza un conjunto de validaci´on relativamente peque˜no, la incertidumbre que conlleva la selecci´on del mejor modelo por validaci´on cruzada puede ser mayor que la incertidumbre en la medici´on del rendimiento del conjunto de pruebas de cualquier modelo [34]. Nosotros utilizaremos la b´usqueda aleatoria por limitaciones en la capacidad de c´omputo, sacrificando de esta manera algo de precisi´on. Figura 2.21: B´usqueda aleatoria frente a b´usqueda en cuadr´ıcula. 2.3 Estimaci´on del error Queremos saber como se comporta el modelo que estamos utilizando. Es decir, como de bueno es el modelo en el contexto del problema. Para ello, vamos a evaluar c´omo predice la clase de las observaciones. Nos limitamos al ´ambito de los clasificadores, aprendizaje inductivo basado en ejemplos. 37
2.3.1 Tarea de inducci´on de un clasificador Dados: •Descripci´on de instancias X, conjunto de atributos y valores •Descripci´on de las hip´otesis H, espacio de ´arboles de decisi´on, reglas o dem´as clasificadores. •Concepto objetivo c:X→0,1. •Ejemplos x, ya sean positivos o negativos, que constituyen el subconjunto D, observaciones disponibles del dominio. Pares (x, c(x)) La tarea de inducci´on de un clasificador consiste en determinar la hip´otesis que coincida con el concepto en todo el espacio de instancias perteneciente al dominio. h∈H/h(x) = c(x)∀x∈X(2.33) Hay muchos motivos por los cuales h(x)6=c(x) para alg´un x∈X. •El concepto cno pertenece al espacio de hip´otesis H. •No se puede describir con el lenguaje de hip´otesis. •Originalmente c∈H, pero hay ejemplos con ruido. •Aunque pertenezca, entrenamos con un subconjunto de X, los ejemplos de D, card(D)card(X) (2.34) Especialmente si Des poco representativo del comportamiento de cen X. En la tarea de clasificaci´on de una observaci´on observamos dos posibles escenarios: •Acierto: la clase predicha coincide con la clase real de la observaci´on. •Fallo: la clase predicha no coincide con la clase real de la observaci´on. Entonces, usamos la tasa de error como medida del rendimiento. La tasa de error es la proporci´on de errores cometidos sobre todo el conjunto de ejemplos. Como ya hemos mencionado, no disponemos de todos los datos del dominio del problema, entonces intentaremos estimar esta tasa de error, que llamaremos error verdadero, con diferentes t´ecnicas. El error verdadero eD, de una hip´otesis hrespecto a un concepto objetivo cy distribuci´on de probabilidad del espacio de instancias Des: eD(h) = Prx∈D[c(x)6=h(x)] (2.35) 38
La tasa de error eS(h) de una hip´otesis hrespecto a un concepto objetivo cy una muestra S⊂X,Ses una muestra de observaciones del dominio X, es: eS(h) = 1 card(S)·X x∈S δ(c(x), h(x)) (2.36) Siendo δ: δ(c(x), h(x)) = 1 sii c(x)6=h(x) (2.37) δ(c(x), h(x)) = 0 sii c(x) = h(x) (2.38) El error de resubstituci´on ero error de entrenamiento, como su propio nombre indica este error se consigue calculando la tasa de error sobre el conjunto de datos utilizado para entrenar el algoritmo. Es inevitablemente optimista, ya que adaptarse mucho a esos datos proporcionar´a una buena estimaci´on del rendimiento cuando la generalizaci´on sobre el dominio puede ser escasa, utilizar esta medida puede conducir a sobreajuste. El error de resubstituci´on estima tasas de error menores que el error verdadero. El error en el conjunto de entrenamiento no es un buen indicador del error sobre datos nuevos. Almacenar todos los datos del dominio ser´ıa el clasificador ´optimo. Como no tenemos disponibles todos los datos del dominio solo se podr´ıan almacenar un peque˜no subconjunto de datos, teniendo una ´optima precisi´on sobre este conjunto pero sin garant´ıas de la capacidad de generalizar. Nos interesa el rendimiento futuro sobre nuevos datos. Para ello se utiliza un conjunto de prueba o test independiente del conjunto de entrenamiento. Se obtiene una buena estimaci´on de la tasa de error cuando la hip´otesis hy el conjunto de muestra Sson independientes. Normalmente solo se dispone de un conjunto de datos etiquetados. Para evitar el sobreajuste o la mala estimaci´on del rendimiento del clasificador debemos utilizar un conjunto de datos con ejemplos independientes que no se hayan utilizado de ninguna manera en la construcci´on del clasificador, ni en la etapa de aprendizaje ni en la etapa de procesamiento (estandarizaci´on, escalado...). Por ello, divideremos en dos subconjuntos, entrenamiento Dy prueba S, los datos de los que disponemos. Calcularemos la estimaci´on del error verdadero, eD, con el error de generalizaci´on, eS. Nos basamos en la suposici´on que los datos de entrenamiento y de prueba son muestras representativas de un mismo problema subyacente. En ocasiones se suele hablar de un tercer conjunto de datos, el conjunto de validaci´on. Este debe ser independiente del conjunto de entrenamiento y de prueba. Se utiliza para optimizar los hiperpar´ametros del clasificador antes de calcular una estimaci´on de su error verdadero con el conjunto de prueba. 39
Normalmente, cuanto m´as grande sea el conjunto de entrenamiento mejor ser´a el clasificador dado que, a priori, ver´a m´as diversidad de instancias del problema a resolver. Cuanto m´as grande sea el conjunto de prueba mejor ser´a la estimaci´on del error verdadero. Y cuanto m´as grande sea el conjunto de validaci´on mejores ser´an los hiperpar´ametros elegidos. 2.3.2 Estratificaci´on Consiste en repartir las observaciones entre las particiones de manera que en todas ellas se respete la distribuci´on inicial de clases. 2.3.3 Hold-out Tambi´en conocido como m´etodo de reserva. Consiste en reservar un cierto n´umero de observaciones del conjunto de datos. Estas observaciones no se utilizan para entrenar el algoritmo de clasificaci´on. Se usar´an en una etapa posterior ´unicamente para estimar el error de generalizaci´on. Existe el dilema de que ambos conjuntos deben ser suficientemente grandes para que el entrenamiento del clasificador sea bueno y tambi´en lo sea la estimaci´on del error. La partici´on m´as com´un usada es 2 3para entrenar y 1 3para estimar el error. Figura 2.22: M´etodo de reserva. Su principal inconveniente es la dependencia de los resultados de la partici´on inicial. Estos subconjuntos podr´ıan no ser representativos, podr´ıamos tener la situaci´on de que una clase no tuviese ningun ejemplo en alguno de los dos conjuntos, ya que la partici´on se realiza de manera aleatoria. Entra en juego la variante estratifiacada. Esta t´ecnica se puede repetir para obtener unos resultados m´as consistentes y menos dependientes de la partici´on inicial. 40
Ventajas •Eficiencia. La diferencia entre la predicci´on y la predicci´on promedio est´a distribuida de manera justa entre los valores de caracter´ıstica de la instancia. •Permite explicaciones comparativas. En lugar de comparar una predicci´on con la predicci´on promedio de todo el conjunto de datos, puede compararse con un subconjunto. •Base te´orica. El valor de Shapley es un m´etodo de explicaci´on con una teor´ıa s´olida como es la teor´ıa de juegos. Las propiedades dan a la explicaci´on una base razonable. Desventajas •Requiere mucho tiempo de c´omputo. Un c´alculo exacto del valor de Shapley es computacionalmente costoso porque hay 2kposibles coaliciones de los valores de la variable y la ausencia de una variable. Se reduce limitando el n´umero de iteraciones M. La disminuci´on de Mreduce el tiempo de c´alculo, pero aumenta la varianza del valor de Shapley. •Puede malinterpretarse. La interpretaci´on debe ser: Dado el conjunto actual de valores de caracter´ısticas, la contribuci´on de un valor de una caracter´ıstica a la diferencia entre la predicci´on real y la predicci´on media es el valor estimado de Shapley. •Siempre se usan todas las variables. El valor de Shapley es el m´etodo de explicaci´on incorrecto si se busca explicaciones dispersas. •Acceso a los datos. No es suficiente acceder a la funci´on de predicci´on. Se necesitan los datos para reemplazar partes de la instancia de inter´es con valores de instancias de datos extra´ıdos al azar. 2.4.3 SHAP SHapley Additive exPlanations [44] es un m´etodo para explicar las predicciones individuales, con base en los Shapley Values. El objetivo de SHAP es explicar la predicci´on de una instancia xcalculando la contribuci´on de cada caracter´ıstica a la predicci´on. El m´etodo de explicaci´on SHAP calcula los valores de Shapley a partir de la teor´ıa de juegos de coalici´on. Una innovaci´on que SHAP trae es que la explicaci´on del valor de Shapley se representa como un m´etodo de atribuci´on de caracter´ısticas aditivas, un modelo lineal. Esta visi´on conecta los valores Shapley y LIME. SHAP especifica la explicaci´on como: g(z0) = φ0+ M X j=1 φjz0 j.(2.49) 47
Donde ges el modelo de explicaci´on, z0∈ {0,1}Mes el vector de coalici´on, Mes el tama˜no m´aximo de la coalici´on y φj∈Res la atribuci´on de caracter´ısticas para una caracter´ıstica j. Como tiene su base en los valores Shapley tambi´en satisface las propiedades de eficiencia, simetr´ıa, simulaci´on y aditividad. SHAP describe las siguientes tres propiedades deseables: •Precisi´on local. f(x) = g(x0) = φ0+ M X j=1 φjx0 j(2.50) Se define φ0=EX(ˆ f(x)) y establece todos x0 jen 1, esta es la propiedad de eficiencia Shapley. Solo con un nombre diferente y usando el vector de coalici´on. f(x) = φ0+ M X j=1 φjx0 j=EX(ˆ f(X)) + M X j=1 φj(2.51) •Ausencia. x0 j= 0 ⇒φj= 0 (2.52) Una caracter´ıstica ausente obtiene una atribuci´on de cero, la denominada “propiedad menor de contabilidad” [44]. •Consistencia. Cuando un modelo cambie de modo que la contribuci´on marginal de un valor de un atributo cambia, el valor de Shapley tambi´en debe cambiar de la misma manera, ya sea aumentar o disminuir. Con esta propiedad se explica el cumplimiento de las propiedades de linealidad y simetr´ıa de Shapley. Podemos diferenciar entre los siguientes m´etodos de SHAP: •KernelSHAP: utiliza la esperanza marginal. Tiene un coste computacional exponencial. Se puede utilizar para aproximar cualquier funci´on. •TreeSHAP: proporciona una estimaci´on eficiente para modelos basados en ´arboles. Utiliza la esperanza condicional. Es una alternativa eficiente a kernel ya que tiene un coste computacional polin´omico. Se beneficia de la propiedad aditiva de SHAP en el ensemble de ´arboles. es un enfoque basado en la retropropagaci´on que atribuye un cambio a las entradas basado en las diferencias entre las entradas y las referencias correspondientes para las activaciones no lineales. •DeepSHAP [45]: es un m´etodo que utiliza un enfoque basado en la retropropagaci´on para aproximar el valor Shapley a las entradas frente a sus respectivas referencias para funciones de activaci´on no lineales. 48
Ventajas •Base te´orica s´olida. Se basa en la teor´ıa de juegos al igual que los valores Shapley y aplica todas sus ventajas. La predicci´on est´a bastante distribuida entre los valores de las caracter´ısticas. Se obtiene explicaciones que comparan la predicci´on con la predicci´on promedio. •Conecta los valores Shapley y LIME. •Implementaci´on r´apida para modelos basados en ´arboles. •Reutilizaci´on de c´alculos. Permite calcular los muchos valores de Shapley necesarios para las interpretaciones del modelo global. Las interpretaciones globales son consistentes con las explicaciones locales, ya que los valores de Shapley son la “unidad at´omica” de las interpretaciones globales. Desventajas •SHAP es lento. Los m´etodos SHAP globales requieren calcular los valores de Shapley para muchas instancias. •KernelSHAP asume independencia de las variables. •TreeSHAP puede dar valores err´oneos de Shapleys para variables no relevantes por el muestreo condicional. •Las desventajas de los valores de Shapley tambi´en se aplican a SHAP: Los valores de Shapley pueden malinterpretarse y se necesita acceso a los datos para calcularlos para nuevos datos. 49
Cap´ıtulo 3 Datos 3.1 Descripci´on Los datos utilizados para la realizaci´on de este trabajo de fin de grado los proporciona el Departamento de Ingenier´ıa El´ectrica de la Universidad de Valladolid, el GIR Adire coordinado por el profesor D. ´ Oscar Duque P´erez, quienes los obtuvieron de forma experimental. A partir de un motor de inducci´on y las herramientas de medida necesarias, se recogi´o el resultado de los datos da˜nando manualmente el motor de forma controlada y obteniendo as´ı su clasificaci´on. Para obtener los datos [46] se utiliza un banco de pruebas que consist´ıa en un motor de inducci´on con las siguientes caracter´ısticas: potencia nominal de 1,1 kW, conexi´on en estrella, tensi´on nominal de 400 V, velocidad nominal de 1410 rpm y corriente nominal de 2,6 A. La carga del motor era un freno magn´etico electromagn´etico Lucas N¨ulle que tambi´en incorpora un medidor de torsi´on mec´anica y velocidad magn´etica. Utilizando diferentes modelos de inversor y a dos niveles de carga se toman las medidas para el motor. Este se va da˜nando manualmente lo que nos proporciona la etiqueta. Contamos con 900 observaciones, siendo la mitad de ellas correspondientes al arm´onico de la banda lateral superior (USH) y la otra mitad al arm´onico de la banda lateral inferior (LSH). Por tanto, tendr´ıamos 450 observaciones dobles. Se toman 2048 valores para cada frecuencia correspondientes a los instantes de tiempo de 2048 milisegundos medidos cuando el motor est´a en funcionamiento. Los datos corresponden a mediciones en las que el motor estaba alimentado por un inversor. Un inversor es un dispositivo que permite controlar la velocidad rotacional de un motor de inducci´on de corriente alterna, el cual es alimentado por una frecuencia y voltaje constante, y entrega al motor una frecuencia y voltaje variable, por lo que tambi´en son denominados variadores de frecuencia [47]. En este experimento se utilizan tres tipos de variadores de frecuencia: •Inversor PowerFlex 40 de Allen Bradley (AB). 50
•Inversor ACS355 de ABB (ABB). •Inversor Altivar 66 de T´el´em´ecanique (TM). Se lleva a cabo el experimento en cada motor con dos niveles de carga para comprobar su posible influencia en el resultado final. •Nivel de carga bajo (NC1) 35 % del par nominal. •Nivel de carga alto (NC2) 60 % del par nominal. Se utiliza un nivel de carga bajo para evitar que sea vac´ıo y un nivel de carga elevado para comprobar su posible influencia. Por ´ultimo, se determina la clasificaci´on del estado del motor seg´un el nivel de perforaci´on del mismo provocado antes de realizar el experimento. El motor se prob´o en cinco condiciones diferentes: desde el motor sano hasta con la barra del rotor totalmente rota. La rotura de la barra se logr´o perforando un agujero en la uni´on de una barra y el anillo final de la jaula. Las clases son ordinales en el rango de R1 a R5: •Motor en buen estado (R1). •Motor levemente da˜nado (R2). •Motor medianamente da˜nado (R3). •Motor severamente da˜nado (R4). •Motor gravemente da˜nado (R5). En la Tabla 3.1 podemos ver la asignaci´on a clases seg´un la perforaci´on, en di´ametro y profundidad, de la barra del rotor. Perforaci´on Estado Profundidad Di´ametro R1 0 0 R2 4.2mm 2.5mm R3 9.4mm 2.5mm R4 17mm 2.5mm R5 17mm 3.5mm Tabla 3.1: Medidas de la perforaci´on seg´un el estado de deterioro. Tanto la distribuci´on de clases como de los distintos grupos de inversores y niveles de carga es homog´enea, es decir, las muestras est´an balanceadas. Contamos con el mismo n´umero de observaciones para cada combinaci´on de variador de frecuencia, nivel de carga y estado del motor. 51
3.2 Procesamiento En este trabajo, el inversor fue configurado para que el transitorio de arranque tuviera una duraci´on de 10 s con una frecuencia final de 50 Hz. La tensi´on se captur´o para calcular la frecuencia fundamental y la velocidad s´ıncrona a lo largo tiempo. Tambi´en se midi´o la velocidad para calcular el deslizamiento del motor. Esta informaci´on permite calcular la componente fundamental y las trayectorias de los arm´onicos relacionados con la falla en el plano tiempo-frecuencia. Aunque el arranque dura 10 s, s´olo se analizan los seis segundos centrales para evitar los efectos de borde al al principio y al final del transitorio que afectan a la cuantificaci´on de la gravedad del fallo. Como ya hemos comentado las observaciones son dobles, para un motor, con un cierto inversor y nivel de carga se obtienen dos observaciones, una para los valores del arm´onico de la banda lateral superior y otra para los valores de la banda lateral inferior. La primera acci´on sobre el conjunto de datos ser´a redimensionarlo de tal manera que cada motor tenga una sola observaci´on con todos los valores y los nombres de las columnas apropiados. Observamos de forma gr´afica la estructura de los datos iniciales. Figura 3.1: Datos iniciales. Observamos de forma gr´afica la estructura de los datos tras la modificaci´on mencionada. Figura 3.2: Modificaci´on de la estructura de los datos. 52
Con esto hacemos un cambio de dimensionalidad. Pasamos a tener la mitad de observaciones, 450, pero el doble de atributos. Por tanto, cada observaci´on tendr´a 4099 variables: •Medida de USH: 2048 variables una por cada valor de esta frecuencia en 2048 milisegundos. •Medida de LSH: 2048 variables una por cada valor de esta frecuencia en 2048 milisegundos. •Modelo de inversor. •Nivel de carga. •Clase. Esta gran cantidad de variables hace que sea relevante uno de los puntos de inter´es planteados entre los objetivos de este trabajo. Concretamente el saber si las redes neuronales son capaces de adaptarse, en este problema concreto, a este problema en concreto, a este elevado n´umero de variables seleccionando aquellas m´as apropiadas para la clasificaci´on o si por el contrario una reducci´on de la dimensi´on en el n´umero de variables ayuda a mejorar la tarea de clasificaci´on. Figura 3.3: Representaci´on de arm´onicos LSH y USH con los datos originales. Podemos ver en la Figura 3.3 la representaci´on del conjunto de datos original, las medidas de ambos arm´onicos frente al tiempo. 53
3.2.1 Conjunto de datos alternativo Se utiliza el conjunto de datos propuesto en [3], como conjunto de datos alternativo con reducci´on de dimensionalidad. En el nuevo conjunto de datos reducido se toman, para cada una de las curvas originales, 7 mediciones construidas a partir de los datos originales para instantes de tiempo entre 250 y 2000 con un intervalo de 250 (250,500, 750, ..., 2000). Cada una de esas 7 mediciones se genera promediando el correspondiente instante con las 5 observaciones inmediatamente anteriores y posteriores al mismo. Constituye una importante reducci´on de dimensionalidad, es decir, del n´umero de variables empleadas, pasando de 2048 variables por banda a tan solo 7. Lo cual proporcionar´a una notable reducci´on del tiempo de ejecuci´on de los algoritmos. Uno de los objetivos de este trabajo es valorar tambi´en si se produce una mejora o empeoramiento en los resultados de clasificaci´on con este conjunto de datos reducido. Figura 3.4: Representaci´on de arm´onicos LSH y USH con los datos alternativos o reducidos. Podemos ver en la Figura 3.4 la representaci´on del conjunto de datos reducido, las medidas de ambos arm´onicos frente al tiempo. Respecto al original se observan lineas mucho m´as rectas porque hay tan solo 7 mediciones por individuo en cada banda. 3.2.2 Subconjuntos de datos A continuaci´on, describimos los subconjuntos con los que vamos a trabajar en la elaboraci´on de modelos. •Conjunto agrupado: se utilizar´an todas las variables. Las variables categ´oricas correspondientes al modelo de inversor y nivel de carga se desdoblar´an como variables dummy 54
por especificaciones del algoritmo. •Subconjuntos modelo de inversor-nivel de carga: se crearan seis subconjuntos, uno para cada combinaci´on modelo de inversor-nivel de carga. 3.3 Terminolog´ıa Para aclarar las ideas a continuaci´on explicaremos como nos vamos a referir a los siguientes t´erminos de forma sin´onima: •Conjunto de datos original, datos completos o conjunto completo: hace referencia al conjunto con todos los instantes de tiempo midiendo los arm´onicos de la onda en los 2048 milisegundos segundos. •Conjunto de datos alternativo, datos reducidos o conjunto reducido: es el conjunto de datos con solo 10 medidas del arm´onico por banda. Los datos con la reducci´on de dimensi´on propuestos en [3]. •Conjunto con variables dummy o conjunto agrupado: representar´a el conjunto con todos los modelos de inversor y niveles de carga representando esta informaci´on en 5 variables dicot´omicas. Podr´a ser tanto de los datos completos como los reducidos. Tambi´en nos referimos a este conjunto como datos agrupados ya que solo hay un grupo. •Por grupos: nos referiremos a los conjuntos separados por grupos seg´un la combinaci´on modelo de inversor y nivel de carga, ya sea para los datos completos o reducidos. Esto dar´a lugar a seis conjuntos. 55
Cap´ıtulo 4 An´alisis y resultados 4.1 Experimento realizado En primer lugar, se realiza una partici´on en dos subconjuntos, uno para entrenamiento, train, y otro para prueba, test. El conjunto de entrenamiento se utilizar´a para entrenar el algoritmo y seleccionar los hiperpar´ametros al mismo tiempo. Como queremos evitar el sobreajuste tambi´en en la selecci´on de hiperpar´ametros, realizaremos otra divisi´on dentro del conjunto de entrenamiento en 5 carpetas, de esta manera, buscaremos el mejor modelo con validaci´on cruzada 5-fold. Posteriormente, se entrena el modelo seleccionado con todos los datos de entrenamiento. Los datos de test se utilizan para estimar el error solo con el mejor modelo ya entrenado. Figura 4.1: Esquema de desarrollo con validaci´on cruzada 5-fold. Imagen de [3] Con el objetivo de obtener un estimador de error m´as preciso y con menos variabilidad, 56
Df Dum Sq Mean Sq F value Pr(>F) Banda 2 0.010499 0.005249 4.142 0.106 Nivel de carga 1 0.005101 0.005101 4.025 0.115 Inversor 2 0.003768 0.001884 1.487 0.329 Banda:Nivel de carga 2 0.000312 0.000156 0.123 0.887 Banda:Inversor 4 0.006653 0.001663 1.312 0.399 Nivel de carga:Inversor 2 0.005321 0.002661 2.099 0.238 Residuales 4 0.005069 0.001267 Tabla 4.7: Anova de tres factores con interacciones de orden 2. Banda-Nivel de carga-Modelo de inversor. Con el objetivo de observar con m´as claridad este an´alisis de la varianza se realiza de nuevo obviando las interacciones e incluyendo solo los efectos principales. Se realiza un proceso de eliminaci´on de las interacciones de forma individual, una a una, comprobando que siguen sin ser significativas, se omiten las tablas intermedias. Df Dum Sq Mean Sq F value Pr(>F) Banda 2 0.010499 0.005249 3.629 0.0585 . Nivel de carga 1 0.005101 0.005101 3.527 0.0849 . Inversor 2 0.003768 0.001884 1.303 0.3076 Residuales 12 0.017356 0.001446 Tabla 4.8: Anova de tres factores sin interacci´on. Banda-Nivel de carga-Modelo de inversor. Podemos concluir de forma clara que la banda y el nivel de carga resultan significativos a nivel 10 %. Sin embargo, el modelo de inversor no lo es, esto nos dice que para cualquier modelo de inversor podemos ajustar un modelo estad´ıstico que proporcione resultados sin diferencias significativas, en cuanto a la tasa de error se refiere, frente a otros inversores. Esto es un gran avance ya que el diferente comportamiento de estos dificultaba en algunos casos su clasificaci´on. Hemos podido comprobar como con esta t´ecnica los resultados podr´ıan ser igual de buenos sea cual sea el modelo de inversor que tenga la m´aquina. Al no resultar significativa la interacci´on entre banda y nivel de carga, podemos elegir ambos individualmente de forma que se optimice el resultado. A la vista de los errores obtenidos lo ´optimo parece ser elegir el nivel de carga alto (NC2) y la banda inferior (LSH) o ambas bandas en alg´un otro caso. Esto se especifica en la siguiente secci´on 4.4 y con los contrastes mostrados a continuaci´on. Para comprobar el buen funcionamiento del an´alisis realizado en la Tabla 4.8 mostramos el gr´afico de residuales. 63
Figura 4.3: Gr´afico de residuales del ANOVA de la Tabla 4.8. Podemos observar como se distribuyen en torno al cero y no parecen tener ning´un patr´on. Podemos decir que siguen una distribuci´on aleatoria luego el an´alisis realizado es correcto. Tambi´en, realizamos los test post-hoc para comprobar en favor de qu´e covariable encontramos las diferencias significativas aunque podemos sospecharlo a la vista de los resultados obtenidos en la secci´on 4.2. Se realiza el test de Duncan a nivel 0.05. Banda Count LS Mean LS Sigma Homogeneous Groups lsh 6 0.218 0.0155259 X ambas 6 0.2225 0.0155259 X ush 6 0.271333 0.0155259 X Tabla 4.9: Parte 1 test post-hoc para la banda. Vemos dos grupos separados: las bandas LSH y ambas se encuentran en el mismo grupo, la banda USH se encuentra en otro grupo diferente. 64
Contrast Sig. Difference ambas - lsh 0.0045 ambas - ush * -0.0488333 lsh - ush * -0.0533333 Tabla 4.10: Parte 2 test post-hoc para la banda. Podemos observar en la Tabla 4.10 como existen diferencias significativas entre las bandas USH y LSH, entre USH y ambas, pero no hay diferencia significativa entre LSH y ambas bandas. Siendo USH la que peores resultados arroja. Sin diferencia significativa entre LSH y ambas bandas intentaremos priorizar LSH por simplicidad. Se realiza el test de Duncan a nivel 0.1. Nivel.de.Carga Count LS Mean LS Sigma Homogeneous Groups NC2 9 0.220444 0.0126769 X NC1 9 0.254111 0.0126769 X Tabla 4.11: Parte 1 test post-hoc para el nivel de carga. Vemos como los dos niveles de carga se encuntran en grupos separados. Contrast Sig. Difference NC1 - NC2 * 0.0336667 Tabla 4.12: Parte 2 test post-hoc para el nivel de carga. Podemos observar en la Tabla 4.12 que existen diferencias significativas entre el nivel de carga alto y el bajo. Siendo el nivel de carga alto, NC2, el que mejores resultados obtiene, lo cual es coherente con la literatura existente [48]. 4.4 Clasificaci´on seg´un modelo de inversor El objetivo ahora es seleccionar que modelo y que datos incorporaremos a este para que el error cometido sea el menor posible. Buscamos para cada inversor cual ser´ıan sus condiciones ´optimas tras el an´alisis realizado en la secci´on anterior. 4.4.1 Inversor AB Encontramos su ´optimo utilizando los datos reducidos por grupos con nivel de carga alto y utilizando ´unicamente la banda LSH. ˆeg= 0.166 65
Mostramos la matriz de confusi´on para el entrenamiento de este modelo promediada en las 20 repeticiones. Se utiliza en cada repetici´on 60 observaciones para entrenamiento y 15 para prueba. Predicho Real R1 R2 R3 R4 R5 R1 2.0 0.05 0.9 0.05 0.0 R2 0.05 2.65 0.3 0.0 0.0 R3 0.6 0.4 2.0 0.0 0.0 R4 0.15 0.0 0.0 2.85 0.0 R5 0.0 0.0 0.0 0.0 3.0 Tabla 4.13: Matriz de confusi´on Inversor AB. Resulta m´as interesante observar la matriz anterior condicionada por filas para estimar la probabilidad de clasificar una instancia de clase X en la columna Y. Se muestra la matriz redondeada a 3 decimales, luego puede que sus valores no sumen exactamente 1. Predicho Real R1 R2 R3 R4 R5 R1 0.667 0.017 0.3 0.017 0.0 R2 0.017 0.883 0.1 0.0 0.0 R3 0.2 0.133 0.667 0.0 0.0 R4 0.05 0.0 0.0 0.95 0.0 R5 0.0 0.0 0.0 0.0 1.0 Tabla 4.14: Matriz de confusi´on Inversor AB condicionada por filas. Podemos apreciar una diferenciaci´on casi perfecta el las clases R4 y R5, las que constituyen m´as da˜nos en el motor. Si agrupamos las clases R1,R2 y R3 frente a estas, haciendo diferenciaci´on entre m´as da˜nados y menos da˜nados obtendr´ıamos una tasa de error muy baja, a pesar de haber utilizado una funci´on de p´erdida categ´orica sin efectos proporcionales. El caso de motor con m´as da˜nos R5 predice correctamente el 100 % de las observaciones de prueba. Por la contra, el motor con menos da˜nos R1 y con da˜nos intermedios R3 tan solo aciertan el 66.7 % de las observaciones de prueba. 66
Interpretabilidad del modelo En primer lugar mostramos el gr´afico obtenido de la importancia de variables. Figura 4.4: Importancia de variables en SHAP para AB. Podemos ver como para las clases de motores con menos da˜nos importan m´as los instantes iniciales x500 y x250. Para los motores con m´as da˜nos, importan m´as los instantes intermedios x1000, x750 y sobre todo respecto de las clases R1 y R2 importan mas los instantes finales x1500. A continuaci´on, mostramos para cada clase los gr´aficos obtenidos con SHAP. Debemos fijarnos en los siguientes puntos para interpretar estos gr´aficos: •Situaci´on del punto: Si el punto est´a a la derecha, SHAP positivo, favorece la clasificaci´on a la clase que estamos estudiando en dicho gr´afico. Si est´a a la izquierda, SHAP negativo, favorece a la no clasificaci´on en esa clase, podr´ıa ser cualquiera de las otras. •Color: Un punto de color rojizo indica que toma valores altos en esa variable. Por el contrario, un color azulado indica que toma valores bajos en esa variable. •Distribuci´on: se trata de ver si los puntos se agrupan segun colores a izquierda o derecha en cada variable. •Orden: las variables se muestran en orden de ¨ımportancia”seg´un los valores SHAP que toman los puntos. La interpretaci´on de estos gr´aficos se realiza de forma local variando levemente los valores de la observaci´on y viendo como cambia as´ı su predicci´on. En la mayor´ıa de estos gr´aficos se toman valores bastante extremos que podr´ıamos considerar incluso outliers, se toma la decisi´on de recortar en torno al valor SHAP -0.5, +0.5 para poder apreciar mejor la influencia de las variables y extraer conclusiones a cerca de su comportamiento. Se muestran a continuaci´on los gr´aficos sobre valores SHAP recortados, sin outiers. 67
•R1. Figura 4.5: SHAP recortado para clase R1 en inversores AB. La variable x500, que es la que m´as peso tiene proporciona una informaci´on de que los valores altos en esa variable, puntos rojos, se sit´uan en los valores SHAP positivos luego favorecen la clasificaci´on de esa variable. Tambi´en se observan m´as puntos azules a la izquierda, SHAP negativos, por lo que valores bajos en esa variable favorecen la no clasificaci´on en R1. Valores bajos de la variable x1000 favorecen a la clasificaci´on como no R1. •R2. Figura 4.6: SHAP recortado para clase R2 en inversores AB. Una vez m´as la variable x500 toma gran importancia. En este caso se observa perfectamente como valores bajos de esta variable motivan a la clasificaci´on como no R2 y valores altos de estas motivan a la clasificaci´on R2. 68
La variable x250 tiene una distribuci´on similar a la anterior. Podr´ıamos concluir que valores bajos en los primeros instantes de tiempo motivar´ıan a una clasificaci´on distinta a R2 y valores altos a una clasificaci´on R2. •R3. Figura 4.7: SHAP recortado para clase R3 en inversores AB. Para la clasificaci´on R3, toman m´as importancia los instantes centrales. Valores bajos en las variables x1250 y x1000 favorecen a la clasificaci´on de la instancia como R3. •R4. Figura 4.8: SHAP recortado para clase R4 en inversores AB. En este caso empiezan a tomar importancia los instantes de x1750 pero tambi´en los de x500, tiene en cuenta momentos del tiempo mayores para la clasificaci´on de esta clase. 69
•R5. Figura 4.9: SHAP recortado para clase R5 en inversores AB. Podemos ver en el gr´afico como es la que mejor se diferencia, que la interpretaci´on sea tan clara puede estar relacionado con que la clase es f´acilmente separable en esa variable. Valores bajos en las variables x1500 y x 1750 favorecen a la clasificaci´on no R5, por el contrario valores altos en estas variables favorecen a la clasificaci´on R5. Con un enfoque m´as general vemos como los individuos de esta clase deben tener valores generalmente m´as altos en todas la variables. Podemos observar como valores bajos en los primeros instantes de tiempo favorecen clasificaciones de motores con menos da˜nos. Valores altos en los primeros instantes de tiempo favorecen clasificaciones de motores con m´as da˜nos. Con los instantes superiores sucede justo al contrario, valores altos indican clasificaciones en motores con m´as da˜nos y valores m´as bajos favorecer´ıan clasificar en motores con menos da˜nos. En los motores R1 y R2 importaron m´as los instantes iniciales, en motores R3 los instantes intermedios y en motores R4 y R5 los instantes finales. 4.4.2 Inversor ABB Encontramos su ´optimo utilizando los datos reducidos por grupos con nivel de carga alto y utilizando la informaci´on de ambas bandas. ˆeg= 0.17 Mostramos la matriz de confusi´on para este modelo promediada en las 20 repeticiones. Se utiliza en cada repetici´on 60 observaciones para entrenamiento y 15 para prueba. 70
Predicho Real R1 R2 R3 R4 R5 R1 1.9 1.1 0.0 0.0 0.0 R2 1.15 1.65 0.15 0.05 0.0 R3 0.0 0.0 3.0 0.0 0.0 R4 0.0 0.05 0.0 2.9 0.05 R5 0.0 0.0 0.0 0.0 3.0 Tabla 4.15: Matriz de confusi´on Inversor ABB. Resulta m´as interesante observar la matriz anterior condicionada por filas para estimar la probabilidad de clasificar una instancia de clase X en la columna Y. Se muestra la matriz redondeada a 3 decimales, luego puede que sus valores no sumen exactamente 1. Predicho Real R1 R2 R3 R4 R5 R1 0.633 0.367 0.0 0.0 0.0 R2 0.383 0.55 0.05 0.017 0.0 R3 0.0 0.0 1.0 0.0 0.0 R4 0.0 0.0 0.017 0.967 0.017 R5 0.0 0.0 0.0 0.0 1.0 Tabla 4.16: Matriz de confusi´on Inversor ABB condicionada por filas. Al igual que en el caso anterior vemos como el modelo es capaz de separar entre los motores con menos da˜nos, R1 a R3, y los motores con m´as da˜nos, R4 y R5, casi a la perfecci´on. Sin embargo, a penas es capaz de decidir entre las clases R1 y R2 ya que las confunde en muchas ocasiones. Tambi´en destacamos el buen comportamiento en R3, R4 y R5 en las cuales predice correctamente casi todas las observaciones. Sin duda nos encontramos ante un modelo muy preciso y que ha funcionado bien ante los datos del laboratorio. Habr´ıa que comprobar su funcionamiento en la pr´actica con datos reales donde hay mucho ruido. Interpretabilidad del modelo A continuaci´on se mostrar´an e interpretar´an los gr´aficos SHAP. Al igual que en la secci´on anterior, encontramos outliers o puntos con valores SHAP muy extremos que nos dificultan la interpretaci´on del gr´afico. Por ello, mostramos directamente los gr´aficos recortados en el entorno -0.5,+0.5 para mejorar la interpretabilidad. 71
•R1. Figura 4.10: SHAP recortado para clase R1 en inversores ABB. Podemos ver como valores bajos en las variables m´as importantes favorecen a la clasificaci´on como R1. Por el contrario, valores altos de estas favorecen a la clasificaci´on como no R1. Resultan m´as importantes los instantes iniciales. 72
•R2. Figura 4.17: SHAP para clase R2 en inversores TM. Vemos claramente como la interpretaci´on de la variable x250 cambia frente a la clase R1, esto podr´ıa ser un factor diferencial entre estos grupos. Valores altos favorecen R2, valores bajos no R2. En el resto de variables se comporta similar a R1, valores bajos favorecen la clasificaci´on R2 y valores altos no R2. •R3. Figura 4.18: SHAP para clase R3 en inversores TM. Con un gr´afico muy similar al de R2, en este caso vemos como valores bajos de la variable x1750 favorecen la clasificaci´on como no R3. Esta variable podr´ıa ser realmente importante para diferenciar este grupo. 79
•R4. Figura 4.19: SHAP para clase R4 en inversores TM. Valores bajos en la variable x1750 favorece la clasificaci´on como R4. Valores altos en ella favorecen no R4. Tambi´en valores altos en las variables x250 y x750 favorecen la clasificaci´on R4. •R5. Figura 4.20: SHAP para clase R5 en inversores TM. Valores bajos en la variable x1750 favorecen la clasificaci´on como no R5, por el contrario, valores altos favorecen R5. Esta variable podr´ıa ser la diferenciadora entre los grupos R4 y R5. Vemos como valores altos en todas las variables favorecen la clasificaci´on R5 puesto que los puntos rojos se sit´uan siempre en la zona de valores positivos. Podemos apreciar como a medida que aumentan los da˜nos en el motor las variables que toman m´as importancia tambi´en van siendo mayores en el tiempo. 80
Para las clases R1, R2, R3 importan m´as los instantes iniciales e intermedios, en R4 y R5 toman m´as importancia los instantes finales. Pese a ser uno de los inversores m´as conflictivos seg´un la literatura existente, hemos logrado obtener una interpretabilidad muy clara, incluso mejor que para los otros inversores. Tambi´en hemos logrado ver que variables son las que diferencian clasificar en un grupo u otro m´as claramente. 4.5 Comparaci´on de m´etodos En este cap´ıtulo compararemos los resultados obtenidos con los mejores resultados obtenidos sobre estos datos con t´ecninas boosting [3]. Resultados de este trabajo AB-NC2 ABB-NC2 TM-NC2 R1 R2 R3 R4 R5 R1 R2 R3 R4 R5 R1 R2 R3 R4 R5 R1 0.667 0.017 0.3 0.017 0.0 R1 0.663 0.367 0.0 0.0 0.0 R1 0.95 0.0 0.05 0.0 0.0 R2 0.017 0.883 0.1 0.0 0.0 R2 0.383 0.55 0.05 0.017 0.0 R2 0.017 0.733 0.25 0.0 0.0 R3 0.2 0.133 0.667 0.0 0.0 R3 0.00 0 1.0 0.0 0.0 R3 0.033 0.517 0.433 0.017 0.0 R4 0.05 0.14 0.0 0.95 0.0 R4 0.0 0.0 0.017 0.967 0.017 R4 0.0 0.0 0.017 0.933 0.05 R5 0.0 0.0 0.0 0.0 1.0 R5 0.0 0.0 0.0 0.0 1.0 R5 0.0 0.0 0.0 0.0 1.0 Tasa Acierto = 0.834 Tasa Acierto = 0.83 Tasa Acierto = 0.811 Resultados de [3] AB-NC2 ABB-NC2 TM-NC2 R1 R2 R3 R4 R5 R1 R2 R3 R4 R5 R1 R2 R3 R4 R5 R1 0.78 0.04 0.17 0.0 0.01 R1 0.54 0.11 0.27 0.0 0.08 R1 0.94 0.0 0.06 0.0 0.0 R2 0.09 0.76 0.15 0.0 0.0 R2 0.13 0.79 0.01 0.07 0.0 R2 0.01 0.66 0.33 0.0 0.0 R3 0.36 0.1 0.54 0.0 0.0 R3 0.16 0.01 0.83 0.0 0.0 R3 0.0 0.29 0.61 0.04 0.06 R4 0.01 0.01 0.0 0.98 0.0 R4 0.0 0.0 0.0 1.0 0.0 R4 0.0 0.0 0.0 1.0 0.0 R5 0.0 0.0 0.0 0.0 1.0 R5 0.0 0.0 0.0 0.0 1.0 R5 0.0 0.0 0.0 0.02 0.98 Tasa Acierto = 0.812 Tasa Acierto = 0.832 Tasa Acierto = 0.838 Tabla 4.19: Comparaci´on de resultados Boosting vs. Redes Neuronales. Podemos observar, en la Tabla 4.19, unos resultados muy similares para ambas t´ecnicas. Ambas logran una tasa de acierto en torno a 0.82. Siendo el peor el modelo de inversor TM para las redes neuronales y el mejor para boosting y por el contrario el modelo de inversor AB el mejor para redes neuronales y el peor para boosting, siempre refiri´endonos a mejor y peor como mejor y peor tasa de acierto. Si que podemos destacar que para todos los modelos aqu´ı presentados y con unos resultados relativamente buenos se utiliza el nivel de carga alto o 2 ya que obten´ıa mejores tasas de acierto. Podemos observar como las redes neuronales son capaces de diferenciar mejor las clases R1 y R2 frente a las t´ecnicas boosting. 81
Cap´ıtulo 5 Conclusiones y trabajo futuro Aqu´ı se exponen las conclusiones obtenidas con la realizaci´on de este trabajo, as´ı como posibles mejoras y trabajos futuros. 5.1 Conclusiones Se ha comprobado en esta memoria que para un mejor funcionamiento de las redes neuronales en este tipo de problema es conveniente un preprocesado de los datos, igual que ocurr´ıa con el boosting. En otras palabras, las redes neuronales no han sido capaces por s´ı mismas de efectuar una selecci´on de variables suficientemente eficiente. Los resultados que se obtuvieron con las redes neuronales considreadas en este trabajo son similares a los obtenidos con boosting. Se mejora en la clasificaci´on para el inversor AB y se pierde ligeramente para el inversor TM. Los resultados obtenidos para el inversor ABB son pr´acticamente id´enticos. Las redes neuronales aqu´ı consideradas obtienen siempre la mejor clasificaci´on a nivel de carga alto (NC2), lo que es coherente con la literatura existente ya que los fallos suelen aparecer de forma m´as clara a niveles de carga altos [48], este hecho tambi´en suced´ıa con las t´ecnicas boosting. Se vuelve a confirmar la hip´otesis de que el nivel de carga alto facilita la detecci´on de fallos y la diagnosis del estado del motor. Tambi´en es destacable la superioridad de los modelos que utilizan la banda LSH sobre los de USH que aparecen tambi´en en boosting. Esto refuerza el mayor inter´es que parece tener LSH. Adem´as, se descarta la hip´otesis de la simetr´ıa entre las bandas. Los resultados parecen decir que LSH es significativamente mejor que USH a la hora de detectar los fallos lo que es una cuesti´on por investigar m´as detenidamente desde el punto de vista el´ectrico. Con la interpretabilidad extra´ıda a trav´es de SHAP se puede ver que variables y de que manera influyen en la predicci´on de nuevas observaciones. En t´erminos generales, podr´ıamos 82
decir que los instantes iniciales del estado transitorio son m´as determinantes en las clasificaciones de motores sanos y los instantes finales son m´as determinantes en la clasificaci´on de los motores m´as da˜nados. Respecto a la interpretabilidad de modelos, podemos destacar la dispersi´on de puntos en las redes neuronales ya que los resultados son est´an m´as mezclados y son menos claros frente a [3] donde los ´arboles proporcionaban unos resultados muy claros. Es interesante destacar que con las t´ecnicas boosting apenas aparec´ıan dos o tres variables con amplia dispersi´on (lo que invita a pensar en una mayor relevancia de la variable) en los gr´aficos SHAP, que eran m´as mucho claros. Con las redes neuronales parece que encuentra patrones en varias variables. Esto puede ser una ventaja porque las redes infieren m´as patrones o una desventaja porque los ´arboles consiguen separan mejor. Para el inversor TM en boosting se utilizan los dos arm´onicos (USH y LSH) mientras que en redes neuronales solo se utiliza LSH. Seg´un la informaci´on proporcionada por los expertos del Departamento de Ingenier´ıa El´ectrica, este inversor es el m´as complejo de estudiar. Hay que notar que en este inversor las redes neuronales obtuvieron una tasa de error de 0.147 utilizando los datos originales frente a 0.189 utilizando los datos reducidos y frente a 0.162 que obtuvo boosting. El modelo basado en los datos originales no se ha detallado en este trabajo sobre todo por problemas en la interpretabilidad, ya que aparec´ıan demasiadas variables en escena. Estos resultados refuerzan la opini´on de los expertos del Departamento de Ingenier´ıa El´ectrica sobre la complejidad del inversor y ponen de manifiesto la necesidad de un estudio m´as exhaustivo sobre el inversor TM. 5.2 Trabajo futuro y posibles mejoras A continuaci´on, se comentan algunas posibles mejoras para este trabajo y trabajos futuros que nacen de este TFG. Los modelos entrenados en este trabajo son relativamente sencillos por la escasez de datos los cuales no pod´ıan satisfacer las necesidades de modelos complejos. Modelos con m´as neuronas y m´as capas podr´ıan funcionar mejor para los datos actuales y conseguir mejores predicciones en cuanto al error se refiere, no obstante, siempre a riesgo de sobreajustar. Tambi´en pueden utilizarse otro tipo de redes neuronales con memoria temporal, como las recurrentes o convolucionales 1D, para tratar los datos como una secuencia temporal. Por un motivo de restricciones temporales y computacionales la b´usqueda de hiperpar´ametros se realiza con una b´usqueda aleatoria en ver de una b´usqueda en cuadr´ıcula. Esta primera es mucho m´as r´apida al fijar el n´umero de iteraciones y encuentra conjuntos de hiperpar´ametros correspondientes con un ´optimo local. La b´usqueda en cuadr´ıcula es capaz de encontrar el ´optimo global a costa de un algoritmo m´as parecido a la fuerza bruta y probar todas las posibilidades, algo que se sal´ıa completamente del tiempo disponible. Puede ser conveniente realizar una b´usqueda m´as exhaustiva de hiper par´ametros. 83
Tras analizar la estructura de las matrices de confusi´on se tiene la intuici´on que penalizar la etapa de entrenamiento seg´un la distancia a la clase podr´ıa tener beneficios en cuanto a la clasificaci´on. En este sentido podr´ıa ser ´util el uso de otras funciones de p´erdida como el MAE en el desarrollo de las redes. Tambi´en puede ser interesante el uso de nuevas t´ecnicas espec´ıficas para datos tabulares, como la que se describe en TabNet [49], que es una red neuronal interpretable, y comprobar si la interpretaci´on que se obtenga mediante esta metodolog´ıa concuerda con la obtenida aqu´ı. Finalmente, otra labor de inter´es a completar, de acuerdo a lo sugerido por los expertos del Departamento de Ingenier´ıa El´ectrica, es la publicaci´on de los resultados aqu´ı obtenidos en revistas internacionales de impacto. 84
Bibliograf´ıa [1] R. Rosa. Estudio sobre la viabilidad de los estad´ısticos de orden superior de la corriente de alimentaci´on como indicadores para determinar el estado de un motor de inducci´on. Trabajo de Fin de Grado, Universidad de Valladolid, 2015. [2] Miljkovi´c Dubravko and Z Hep. Brief review of motor current signature analysis. HDKBR Info-CrSNDT Journal, 15:15–26, 2015. [3] Alejandro Bar´on. Detecci´on y Clasificaci´on de Fallos en Motores mediante Procedimientos Boosting. Trabajo de Fin de Grado, Universidad de Valladolid, 2020. [4] Sahar Zolfaghari, Samsul Bahari Mohd Noor, Mohammad Rezazadeh Mehrjou, Mohammad Hamiruce Marhaban, and Norman Mariun. Broken rotor bar fault detection and classification using wavelet packet signature analysis based on fourier transform and multi-layer perceptron neural network. Applied Sciences, 8(1), 2018. [5] R. A. Fisher. The use of multiple measurements in taxonomic problems. Annals of Eugenics, 7(7):179–188, 1936. [6] D. R. Cox. The regression analysis of binary sequences. Journal of the Royal Statistical Society: Series B (Methodological), 20(2):215–232, 1958. [7] James N. Morgan and John A. Sonquist. Problems in the analysis of survey data, and a proposal. Journal of the American Statistical Association, 58(302):415–434, 1963. [8] Dimitris Bertsimas and Jack Dunn. Optimal classification trees. Machine Learning, 106(7):1039–1082, April 2017. [9] Zhou Yong. Knowledge discovery of interesting classification rules based on adaptive genetic algorithm. International Journal of Computational Intelligence Systems, 10 2007. [10] Warren McCulloch and Walter Pitts. A logical calculus of ideas immanent in nervous activity. Bulletin of Mathematical Biophysics, 5:127–147, 1943. [11] F. Rosenblatt. The perceptron: A probabilistic model for information storage and organization in the brain. Psychological Review, pages 65–386, 1958. [12] Michael A. Nielsen. Neural networks and deep learning, 2018. 85
[13] Prajit Ramachandran, Barret Zoph, and Quoc V. Le. Searching for activation functions. CoRR, abs/1710.05941, 2017. [14] Lu Lu. Dying relu and initialization: Theory and numerical examples. Communications in Computational Physics, 28(5):1671–1706, Jun 2020. [15] The softmax function. http://akashgit.github.io/2017/03/13/unsaturating_ softmax.html.´ Ultimo acceso: 2021-03-25. [16] Ian Goodfellow, Yoshua Bengio, and Aaron Courville. Deep Learning. MIT Press, 2016. http://www.deeplearningbook.org. [17] Kurt Hornik. Approximation capabilities of multilayer feedforward networks. Neural Networks, 4(2):251–257, 1991. [18] A. Cauchy. Methode generale pour la resolution des systemes d’equations simultanees. C.R. Acad. Sci. Paris, 25:536–538, 1847. [19] Sebastian Ruder. An overview of gradient descent optimization algorithms. CoRR, abs/1609.04747, 2016. [20] David E. Rumelhart, Geoffrey E. Hinton, and Ronald J. Williams. Learning Representations by Back-propagating Errors. Nature, 323(6088):533–536, 1986. [21] Garrett B. Goh, Nathan O. Hodas, and Abhinav Vishnu. Deep learning for computational chemistry, 2017. [22] Keras API reference optimizers sgd. https://keras.io/api/optimizers/sgd/.´ Ultimo acceso: 2021-03-25. [23] Ilya Sutskever, James Martens, George Dahl, and Geoffrey Hinton. On the importance of initialization and momentum in deep learning. In Proceedings of the 30th International Conference on International Conference on Machine Learning - Volume 28, ICML’13, page III–1139–III–1147. JMLR.org, 2013. [24] Keras API reference optimizers rmsprop. https://keras.io/api/optimizers/ rmsprop/.´ Ultimo acceso: 2021-03-25. [25] Keras API reference optimizers adam. https://keras.io/api/optimizers/adam/. ´ Ultimo acceso: 2021-03-25. [26] Diederik P. Kingma and Jimmy Ba. Adam: A method for stochastic optimization, 2017. [27] Pan Zhou, Jiashi Feng, Chao Ma, Caiming Xiong, Steven HOI, et al. Towards theoretically understanding why sgd generalizes better than adam in deep learning. arXiv preprint arXiv:2010.05627, 2020. 86
[28] H Jabbar and Rafiqul Zaman Khan. Methods to avoid over-fitting and under-fitting in supervised machine learning (comparative study). Computer Science, Communication and Instrumentation Devices, pages 163–172, 2015. [29] Trevor Hastie, Robert Tibshirani, and Jerome Friedman. The elements of statistical learning – data mining, inference, and prediction. [30] Keras API reference layer weight regularizers. https://keras.io/api/layers/ regularizers/.´ Ultimo acceso: 2021-04-7. [31] Nitish Srivastava, Geoffrey Hinton, Alex Krizhevsky, Ilya Sutskever, and Ruslan Salakhutdinov. Dropout: A simple way to prevent neural networks from overfitting. J. Mach. Learn. Res., 15(1):1929–1958, January 2014. [32] Nitish Srivastava, Geoffrey Hinton, Alex Krizhevsky, Ilya Sutskever, and Ruslan Salakhutdinov. Dropout: A simple way to prevent neural networks from overfitting. Journal of Machine Learning Research, 15(56):1929–1958, 2014. [33] Lutz Prechelt. Early Stopping — But When?, pages 53–67. Springer Berlin Heidelberg, Berlin, Heidelberg, 2012. [34] James Bergstra and Yoshua Bengio. Random search for hyper-parameter optimization. J. Mach. Learn. Res., 13:281–305, February 2012. [35] L. A. Rastrigin. The convergence of the random search method in the extremal control of a many parameter system. Automaton Remote Control, 24:1337–1342, 1963. [36] Tom M. Mitchell. Machine Learning. McGraw-Hill, New York, 1997. [37] Francisco Azuaje, Ian Witten, and Frank E. Witten ih, frank e: Data mining: Practical machine learning tools and techniques. Biomedical Engineering Online - BIOMED ENG ONLINE, 5:1–2, 01 2006. [38] Christoph Molnar. Interpretable Machine Learning. 2019. https://christophm. github.io/interpretable-ml-book/. [39] Tim Miller. Explanation in artificial intelligence: Insights from the social sciences. CoRR, abs/1706.07269, 2017. [40] Marco Tulio Ribeiro, Sameer Singh, and Carlos Guestrin. ”why should i trust you?”: Explaining the predictions of any classifier. In Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, KDD ’16, page 1135–1144, New York, NY, USA, 2016. Association for Computing Machinery. [41] David Alvarez-Melis and Tommi S. Jaakkola. On the robustness of interpretability methods. CoRR, abs/1806.08049, 2018. 87
[42] Lloyd S. Shapley. A Value for n-Person Games. RAND Corporation, Santa Monica, CA, 1952. [43] Erik ˇ Strumbelj and Igor Kononenko. Explaining prediction models and individual predictions with feature contributions. Knowledge and Information Systems, 41:647–665, 12 2013. [44] Scott M Lundberg and Su-In Lee. A unified approach to interpreting model predictions. In I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017. [45] Zeon Trevor Fernando, Jaspreet Singh, and Avishek Anand. A study on the interpretability of neural retrieval models using DeepSHAP. In Proceedings of the 42nd International ACM SIGIR Conference on Research and Development in Information Retrieval. ACM, July 2019. [46] Vanesa Fernandez-Cavero, Luis A. Garc´ıa-Escudero, Joan Pons-Llinares, Miguel A. Fern´andez-Temprano, Oscar Duque-Perez, and Daniel Morinigo-Sotelo. Diagnosis of broken rotor bars during the startup of inverter-fed induction motors using the dragon transform and functional anova. Applied Sciences, 11(9), 2021. [47] Ra´ul Granados. Diagn´ostico de fallos en el rotor de motores el´ectricos en estado transitorio mediante t´ecnicas estad´ısticas. Trabajo de Fin de Grado, Universidad de Valladolid, 2017. [48] Sakthivel Ganesan, Prince Winston David, Praveen Kumar Balachandran, and Devakirubakaran Samithas. Intelligent starting current-based fault identification of an induction motor operating under various power quality issues. Energies, 14(2), 2021. [49] Sercan O Arik and Tomas Pfister. Tabnet: Attentive interpretable tabular learning. arXiv preprint arXiv:1908.07442, 2019. 88
268 m1 = m1 + test1 ["confusion_matrix"] 269 m2 = m2 + test2 ["confusion_matrix"] 270 test = [ test0 , test1 , test2 ] 271 m = [m0 ,m1 ,m2] 272 if verbose: 273 print ("ACC AB = ", test0 [" accuracy "]) 274 print ("ACC ABB = ", test1 [" accuracy "]) 275 print ("ACC TM = ", test2 [" accuracy " ]) 276 elif group ==" both": 277 X_test . iloc [: ,: -5] = standarize ( X_test . iloc [: ,: -5]) 278 test0 = test_model ( model , X_test . loc [( X_test [" model_AB " ]==1) & ( X_test[" level_NC1 "]==1) ], y_test . loc [( X_test [" model_AB " ]==1) & ( X_test [" level_NC1 " ]==1) ], net ) 279 test1 = test_model ( model , X_test . loc [( X_test [" model_ABB " ]==1) & ( X_test[" level_NC1 "]==1) ], y_test . loc [( X_test [" model_ABB " ]==1) & ( X_test [ " level_NC1 " ]==1) ], net ) 280 test2 = test_model ( model , X_test . loc [( X_test [" model_TM " ]==1) & ( X_test[" level_NC1 "]==1) ], y_test . loc [( X_test [" model_TM " ]==1) & ( X_test [" level_NC1 " ]==1) ], net ) 281 test3 = test_model ( model , X_test . loc [( X_test [" model_AB " ]==1) & ( X_test[" level_NC2 "]==1) ], y_test . loc [( X_test [" model_AB " ]==1) & ( X_test [" level_NC2 " ]==1) ], net ) 282 test4 = test_model ( model , X_test . loc [( X_test [" model_ABB " ]==1) & ( X_test[" level_NC2 "]==1) ], y_test . loc [( X_test [" model_ABB " ]==1) & ( X_test [ " level_NC2 " ]==1) ], net ) 283 test5 = test_model ( model , X_test . loc [( X_test [" model_TM " ]==1) & ( X_test[" level_NC2 "]==1) ], y_test . loc [( X_test [" model_TM " ]==1) & ( X_test [" level_NC2 " ]==1) ], net ) 284 errors[0,iter] = 1- test0 [" accuracy "] 285 errors[1,iter] = 1- test1 [" accuracy "] 286 errors[2,iter] = 1- test2 [" accuracy "] 287 errors[3,iter] = 1- test3 [" accuracy "] 288 errors[4,iter] = 1- test4 [" accuracy "] 289 errors[5,iter] = 1- test5 [" accuracy "] 290 n_obs [0, iter] = len( test0 [" yhat "]) 291 n_obs [1, iter] = len( test1 [" yhat "]) 292 n_obs [2, iter] = len( test2 [" yhat "]) 293 n_obs [3, iter] = len( test3 [" yhat "]) 294 n_obs [4, iter] = len( test4 [" yhat "]) 295 n_obs [5, iter] = len( test5 [" yhat "]) 296 m0 = m0 + test0 ["confusion_matrix"] 297 m1 = m1 + test1 ["confusion_matrix"] 298 m2 = m2 + test2 ["confusion_matrix"] 299 m3 = m3 + test3 ["confusion_matrix"] 300 m4 = m4 + test4 ["confusion_matrix"] 301 m5 = m5 + test5 ["confusion_matrix"] 302 test = [test0 ,test1 , test2 , test3 ,test4 , test5] 303 m = [m0 ,m1 ,m2 ,m3 ,m4 ,m5] 304 if verbose: 305 print ("ACC AB NC1 = ", test0 [" accuracy "]) 306 print ("ACC ABB NC1 = ", test1 [" accuracy "]) 95
307 print ("ACC TM NC1 = ", test2 [" accuracy "]) 308 print ("ACC AB NC2 = ", test3 [" accuracy "]) 309 print ("ACC ABB NC2 = ", test4 [" accuracy "]) 310 print ("ACC TM NC2 = ", test5 [" accuracy "]) 311 312 #si no hay que hacer la media ponderada cuando no separamos en grupos: 313 #if group == "no" : 314 315 m_e = np.sum( errors *( n_obs /np . sum( n_obs ))) 316 return {"errors":errors ,"mean_e":m_e ," model ":model, " test ":test ," confusion_matrix":m, " n_obs ": n_obs } 317 318 """ Interpretabilidad de modelos complejos ( SHAP ) """ 319 320 import shap 321 def get_shap ( model , X_train , X_test = None , type="deep ",s=True): 322 #type: " deep " o " kernel " 323 ’’’ 324 if X_test == None : 325 X_test = X_train 326 ’’’ 327 if type =="deep ": 328 explainer = shap . DeepExplainer ( model , np . array ( X_train )) 329 shap_values = explainer . shap_values ( np. array ( X_train )) 330 331 return shap . summary_plot ( shap_values , np . array ( X_train ) , plot_type =" bar", feature_names = X_train .columns , show =s) 332 else: 333 explainer = shap . KernelExplainer ( model . predict_proba , X_train ) 334 shap_values = explainer . shap_values ( X_train ) 335 336 return shap. summary_plot ( shap_values , X_train , plot_type ="bar", feature_names = X_train .columns , show=s) 337 338 def get_shap_class ( model ,X_train , X_test = None , clase = "R1" ,type=" deep" , s= True ): 339 #type: " deep " o " kernel " 340 c= int( clase [ -1]) -1 341 ’’’ 342 if X_test == None : 343 X_test = X_train 344 ’’’ 345 if type =="deep ": 346 explainer = shap . DeepExplainer ( model , np . array ( X_train )) 347 shap_values = explainer . shap_values ( np. array ( X_train )) 348 349 np. save(’./ resultados / shaps / shap_R ’+str (c +1)+’. npy ’, shap_values [c]) 350 aux = shap_values [c] 351 aux_train = X_train . iloc [np . where (np .all(aux <0.5 , axis =1) )] 352 aux = aux[ np . where (np . all(aux <0.5 , axis =1))] 353 aux_train = X_train . iloc [np . where (np .all( aux > -0.5 , axis =1) )] 96
354 aux = aux[ np . where (np . all(aux > -0.5 , axis =1))] 355 np. save(’./ resultados / shaps / shap_clipped_R ’+str(c+1) +’.npy ’, shap_values [c]) 356 357 return shap. summary_plot (aux , aux_train , feature_names = X_train . columns , show=s) 358 else: 359 explainer = shap . KernelExplainer ( model . predict_proba , X_train ) 360 shap_values = explainer . shap_values ( X_train ) 361 # shap . summary_plot ( shap_values [c], X_train , feature_names = X_train . columns) 362 np. save(’./ resultados / shaps / shap_R ’+str (c)+’.npy ’, shap_values [c]) 363 return shap. summary_plot ( shap_values [c], X_train , feature_names = X_train .columns , show=s) 364 365 """ 366 367 368 369 370 Busqueda de hiperparametros """ 371 372 from kerastuner . tuners import RandomSearch 373 import tensorflow as tf 374 from kerastuner import HyperModel 375 from tensorflow import keras 376 from tensorflow . keras import layers 377 from kerastuner . tuners import RandomSearch 378 from kerastuner import HyperModel 379 from keras.regularizers import l1_l2 380 381 class MyHyperModel ( HyperModel ): 382 383 def __init__ ( self , num_classes , input_shape ): 384 self. num_classes = num_classes 385 self. input_shape = input_shape 386 387 def build (self , hp): 388 model = keras. Sequential () 389 reg =hp. Choice (’regularization’,values=[0.0,0.025,0.05]) 390 391 model . add ( layers . Dense ( units = hp . Choice ( ’units_l_0 ’,values = [2 ,4 ,8 ,16 ,32]), 392 activation = hp . Choice (’dense_activation_l_0’, values=[’relu ’,’sigmoid’]) ,# Funcion de activacion 393 input_dim = self . input_shape , activity_regularizer = l1_l2 ( l1 =reg , l2 = reg ) 394 ) 395 ) 396 model . add ( layers . Dropout ( rate =hp . Float (’dropout_0 ’, min_value =0.0 , max_value =0.3 , default =0.05 , step =0.05) )) 97
397 lay =hp. Choice (’layers’, values = [1 ,2 ,3]) 398 for iin range (1 , lay ):# Entre 1 y 3 capas Densas 399 # incluye 0, 1 o 2 capas mas 400 401 model . add ( layers . Dense ( units = hp . Choice ( ’units_l_ ’+str(i),values = [2 ,4 ,8 ,16 ,32]), 402 activation = hp . Choice (’dense_activation_l_’+str (i) , values =[ ’relu ’,’sigmoid’]) ,# Funcion de activacion 403 input_dim = self . input_shape , activity_regularizer = l1_l2 ( l1 =reg , l2 = reg ) 404 ) 405 ) 406 model . add ( layers . Dropout ( rate =hp . Float (’dropout_ ’+str(i), min_value =0.0 , max_value =0.3 , default =0.05 , step =0.05) )) 407 #Capa de salida 408 model . add ( layers . Dense ( self . num_classes , activation = ’softmax’)) 409 410 model . compile( hp . Choice (’ optimizer ’ ,[’adam ’,’rmsprop’]) , ’ categorical_crossentropy’, metrics =[ ’accuracy ’])# Metodo de optimizacion 411 return model 412 413 414 class MyTuner ( RandomSearch ):# esto lo pruebo yo pero no se si ira con RS 415 def run_trial (self , trial , *args , ** kwargs ): 416 # You can add additional HyperParameters for preprocessing and custom training loops 417 # via overriding ‘run_trial ‘ 418 kwargs[’ batch_size ’] = trial . hyperparameters . Choice ( ’batch_size ’, values =[2 ,4 ,8 ,16 ,32]) 419 kwargs[’epochs’] = trial . hyperparameters . Int (’epochs’, 30, 200) 420 super ( MyTuner , self ). run_trial (trial , *args , ** kwargs ) 421 422 def buscar_model (X,y, trials =400 , executions =5 , val_split =0.15) : 423 # maximo de intentos : trials 424 # ejecuciones por intento ( modelo ): executions 425 # Porcentaje de datos que se usan para la validacion del modelo: val_split 426 427 #Se crea el hipermodelo con el espacio de b s q u e d a 428 hypermodel = MyHyperModel ( num_classes =5 , input_shape =X. shape [1]) 429 #se establecen las caracteristicas de la busqueda y el tipo de busqueda 430 #tuner = RandomSearch( 431 tuner = MyTuner (# Para probar tunear epocas y batch_size 432 hypermodel , 433 objective =’val_accuracy’, 434 max_trials = trials , 435 executions_per_trial=executions, 436 directory =’my_dir’, 437 project_name=’helloworld ’, 438 overwrite = True ) 439 #Se entrenan los modelos buscando el mejor 98
440 tuner . search (X, y , 441 #si tuneamos epoch y batch_size no hay que pasarlo como argumento 442 # epochs =100 , 443 # batch_size = 8, 444 # validation_data =( X_test , d_test )) 445 validation_split = val_split ,# validamos con un porcentaje de los mismos datos 446 callbacks =[ tf . keras . callbacks . EarlyStopping ( ’val_loss ’, patience =3)],# Se anade para una parada temprana en las epocas si val_loss no decrece 447 verbose = 0 448 ) 449 model = tuner . get_best_models ( num_models =1) [0] 450 params = tuner . get_best_hyperparameters () [0] 451 return params , tuner 452 453 from sklearn . model_selection import train_test_split 454 def hold_out (X,y , test_size =0.2 , net ="normal",t=400 , e=5 , reps =20 , group = False ) : 455 #reps: repeticiones 456 #t: trials , e: executions 457 # Separacion train - test 458 if group : 459 # Matrices de confusion por grupos 460 m0 = np. zeros ((5 ,5) , dtype = float )#AB -NC1 461 m1 = np. zeros ((5 ,5) , dtype = float )#ABB -NC1 462 m2 = np. zeros ((5 ,5) , dtype = float )#TM -NC1 463 m3 = np. zeros ((5 ,5) , dtype = float )#AB -NC2 464 m4 = np. zeros ((5 ,5) , dtype = float )#ABB -NC2 465 m5 = np. zeros ((5 ,5) , dtype = float )#TM -NC2 466 else: 467 m = np . zeros ((5 ,5) , dtype = float ) 468 469 for rin range (0 ,20): 470 X_train , X_test , y_train , y_test = train_test_split (X, y, test_size = test_size , stratify =y) 471 # normalizar los datos 472 if net ==" dummy ": 473 X_train . iloc [: ,: -5] = standarize ( X_train . iloc [: ,: -5]) 474 X_test . iloc [: ,: -5] = standarize ( X_test . iloc [: ,: -5]) 475 if net =="normal": 476 X_train = standarize ( X_train ) 477 X_test = standarize ( X_test ) 478 #Para LSTM no normalizamos 479 480 # Construimos el modelo con los parametros obtenidos automaticamente 481 params , tuner = buscar_model ( X_train , y_train , trials =t, executions =e) 482 483 model = tuner . hypermodel . build ( params ) 484 model . fit ( X_train , y_train , epochs = params . get (’epochs’), batch_size = params . get( ’batch_size ’),verbose =0) 485 99
486 if group : 487 test0 = test_model ( model , X_test . loc [( X_test [" model_AB " ]==1) & ( X_test[" level_NC1 "]==1) ], y_test . loc [( X_test [" model_AB " ]==1) & ( X_test [" level_NC1 " ]==1) ], net ) 488 test1 = test_model ( model , X_test . loc [( X_test [" model_ABB " ]==1) & ( X_test[" level_NC1 "]==1) ], y_test . loc [( X_test [" model_ABB " ]==1) & ( X_test [ " level_NC1 " ]==1) ], net ) 489 test2 = test_model ( model , X_test . loc [( X_test [" model_TM " ]==1) & ( X_test[" level_NC1 "]==1) ], y_test . loc [( X_test [" model_TM " ]==1) & ( X_test [" level_NC1 " ]==1) ], net ) 490 test3 = test_model ( model , X_test . loc [( X_test [" model_AB " ]==1) & ( X_test[" level_NC2 "]==1) ], y_test . loc [( X_test [" model_AB " ]==1) & ( X_test [" level_NC2 " ]==1) ], net ) 491 test4 = test_model ( model , X_test . loc [( X_test [" model_ABB " ]==1) & ( X_test[" level_NC2 "]==1) ], y_test . loc [( X_test [" model_ABB " ]==1) & ( X_test [ " level_NC2 " ]==1) ], net ) 492 test5 = test_model ( model , X_test . loc [( X_test [" model_TM " ]==1) & ( X_test[" level_NC2 "]==1) ], y_test . loc [( X_test [" model_TM " ]==1) & ( X_test [" level_NC2 " ]==1) ], net ) 493 m0 = m0 + test0 ["confusion_matrix"] 494 m1 = m1 + test1 ["confusion_matrix"] 495 m2 = m2 + test2 ["confusion_matrix"] 496 m3 = m3 + test3 ["confusion_matrix"] 497 m4 = m4 + test4 ["confusion_matrix"] 498 m5 = m5 + test5 ["confusion_matrix"] 499 m = [m0 ,m1 ,m2 ,m3 ,m4 ,m5] 500 else: 501 test = test_model (model , X_test , y_test , net) 502 m = m + test ["confusion_matrix"] 503 504 return m 505 506 """# Building model funtions """ 507 508 from keras . models import Sequential 509 from keras . layers import Dense , LeakyReLU ,Dropout , LSTM , InputLayer 510 from keras . metrics import Accuracy 511 from keras.regularizers import l1_l2 512 from sklearn.metrics import accuracy_score , confusion_matrix 513 514 515 516 # Builds Keras model given a set of hyperparameters (passed as a dictionary ) 517 def build_model (X ,y , params ): 518 519 520 myreg = l1_l2 (l1 = params ["regularizer"][0] , l2 = params ["regularizer"][1]) 521 522 model = Sequential () 523 # Input layer 100
524 model . add ( Dense ( params [" nneurons "][0] , activation = params ["activations" ][0] , input_shape =( X. shape [1] ,) , activity_regularizer = myreg )) 525 model .add ( Dropout ( params [" dropouts "][0]) ) 526 527 #Hidden layers 528 for iin range (1 , len (params[" nneurons " ])): 529 model . add ( Dense ( params [" nneurons "][i], activation = params ["activations" ][i],activity_regularizer=myreg)) 530 model .add ( Dropout ( params [" dropouts "][i])) 531 532 # Output layer 533 model . add ( Dense (5 , activation ="softmax",activity_regularizer=myreg)) 534 535 model . compile( optimizer = params [" optimizer "], loss =" categorical_crossentropy",metrics =[" accuracy "]) 536 model . fit (X,y, epochs = params ["epochs"], batch_size = params [" batch_size "], validation_split = params [" val_split "], verbose = params ["verbose"]) 537 return model 538 539 # Builds Keras model given a set of hyperparameters (passed as a dictionary ) with dummy variables at 6 least columns 540 def build_model_dummy (X,y, params ): 541 542 myreg = l1_l2 (l1 = params ["regularizer"][0] , l2 = params ["regularizer"][1]) 543 544 model = Sequential () 545 # Input layer 546 model . add ( Dense ( params [" nneurons "][0] , activation = params ["activations" ][0] , input_shape =( X. shape [1] ,) , activity_regularizer = myreg )) 547 548 model .add ( Dropout ( params [" dropouts "][0]) ) 549 #Hidden layers 550 for iin range (1 , len (params[" nneurons " ])): 551 model . add ( Dense ( params [" nneurons "][i], activation = params ["activations" ][i],activity_regularizer=myreg)) 552 model .add ( Dropout ( params [" dropouts "][i])) 553 554 # Output layer 555 model . add ( Dense (5 , activation ="softmax")) 556 557 model . compile( optimizer = params [" optimizer "], loss =" categorical_crossentropy",metrics =[" accuracy "]) 558 model . fit (X,y, epochs = params ["epochs"], batch_size = params [" batch_size "], validation_split = params [" val_split "], verbose = params ["verbose"]) 559 return model 560 561 # Builds Keras RNN model given a set of hyperparameters (passed as a dictionary ) 562 def build_LSTM_model (X ,y , params ): 563 564 if " model " in X. columns : 101
565 X = X. drop ([" model "], axis =1) 566 if " level " in X. columns : 567 X = X. drop ([" level "], axis =1) 568 569 X = preprocess_LSTM (X) # reshape data 570 571 lay = len(params[" nneurons " ]) 572 ret = [ True ]*( lay -1) +[ False ] 573 574 model = Sequential () 575 576 # Input layer 577 model . add ( LSTM ( params [" nneurons " ][0] , return_sequences = ret [0] , input_shape =( X. shape [1] , X. shape [2]) )) 578 model .add ( Dropout ( params [" dropouts "][0]) ) 579 580 #Hidden layers 581 for iin range (1 , len (params[" nneurons " ])): 582 model . add ( LSTM ( params [" nneurons "][i], return_sequences = ret [i]) ) 583 model .add ( Dropout ( params [" dropouts "][i])) 584 585 # output layer 586 model . add ( Dense (5 , activation ="softmax")) 587 588 model . compile( optimizer = params [" optimizer "], loss =" categorical_crossentropy",metrics =[" accuracy "]) 589 model . fit (X,y, epochs = params ["epochs"], batch_size = params [" batch_size "], validation_split = params [" val_split "], verbose = params ["verbose"]) 590 return model 591 592 def test_model (model ,X ,y ,net ="normal"): 593 if net == " LSTM ": 594 if " model " in X. columns : 595 X = X. drop ([" model "], axis =1) 596 if " level " in X. columns : 597 X = X. drop ([" level "], axis =1) 598 X = preprocess_LSTM (X) 599 yhat = np . argmax ( model . predict (X) ,axis =1) 600 ytest = np . argmax (y. values , axis =1) 601 #La posicion Cij de la matriz de confusion corresponde con los de la clase real i asignados a la clase j ( real filas , predicho columnas ) 602 return {" accuracy ": accuracy_score ( yhat , ytest ) ," yhat ":yhat ," ytest ":ytest , "confusion_matrix": confusion_matrix (ytest ,yhat , labels =[0 ,1 ,2 ,3 ,4]) } 603 604 """# Saving results """ 605 606 #Crea un fichero el la ubicacion proporcionada 607 def create_file ( arg ): 608 f = open(arg ,"w+") 609 f. close () 610 102
611 # Anade los resultados al final del fichero 612 #Debe ser un fichero ya existente 613 def write_R ( title , m , file, metric=" Error "): 614 file. write (" %s" %title ) 615 try: 616 l = m. shape 617 except: 618 l= len(m) 619 if l == (5 ,5) : 620 file. write ("Matriz global:\n") 621 np. savetxt (file, m, delimiter =",",fmt=" %2.0f") 622 acc_global = prec_mat (m) 623 else: 624 file. write ("Matriz global:\n") 625 np. savetxt (file, np.sum(m, axis =0) , delimiter =",",fmt=" %2.0f") 626 acc_global = prec_mat (np. sum(m, axis =0)) 627 628 if l == 6: 629 file. write (" Matriz AB - NC1 :\n") 630 np. savetxt (file, m [0] , delimiter =",",fmt=" %2.0f") 631 acc_AB_NC1 = prec_mat (m [0]) 632 file. write (" Matriz ABB - NC1 :\ n") 633 np. savetxt (file, m [1] , delimiter =",",fmt=" %2.0f") 634 acc_ABB_NC1 = prec_mat (m [1]) 635 636 file. write (" Matriz TM - NC1 :\n") 637 np. savetxt (file, m [2] , delimiter =",",fmt=" %2.0f") 638 acc_TM_NC1 = prec_mat (m [2]) 639 file. write (" Matriz AB - NC2 :\n") 640 np. savetxt (file, m [3] , delimiter =",",fmt=" %2.0f") 641 acc_AB_NC2 = prec_mat (m [3]) 642 file. write (" Matriz ABB - NC2 :\ n") 643 np. savetxt (file, m [4] , delimiter =",",fmt=" %2.0f") 644 acc_ABB_NC2 = prec_mat (m [4]) 645 file. write (" Matriz TM - NC2 :\n") 646 np. savetxt (file, m [5] , delimiter =",",fmt=" %2.0f") 647 acc_TM_NC2 = prec_mat (m [5]) 648 649 file. write (" Matriz AB :\ n") 650 np. savetxt (file, m [0]+ m[3] , delimiter =",",fmt =" %2.0f") 651 acc_AB = prec_mat (m [0]+ m[3]) 652 file. write (" Matriz ABB :\n") 653 np. savetxt (file, m [1]+ m[4] , delimiter =",",fmt =" %2.0f") 654 acc_ABB = prec_mat (m [1]+m [4]) 655 file. write (" Matriz TM :\ n") 656 np. savetxt (file, m [2]+ m[5] , delimiter =",",fmt =" %2.0f") 657 acc_TM = prec_mat (m [2]+ m[5]) 658 659 file. write (" Matriz NC1 :\n") 660 np. savetxt (file, m [0]+ m [1]+ m [2] , delimiter =",",fmt=" %2.0f") 661 acc_NC1 = prec_mat (m [0]+m [1]+ m[2]) 103
662 file. write (" Matriz NC1 :\n") 663 np. savetxt (file, m [3]+ m [4]+ m [5] , delimiter =",",fmt=" %2.0f") 664 acc_NC2 = prec_mat (m [3]+m [4]+ m[5]) 665 if metric == " Error " or metric == " error ": 666 metrica = 1 667 else: 668 metric = " Acc " 669 metrica = 0 670 file. write (" %s global = %s \n" %(metric , str (abs( metrica - acc_global )))) 671 if l == 2 or l == 6: 672 file. write (" %s NC1 = %s \n" %(metric , str(abs ( metrica - acc_NC1 )))) 673 file. write (" %s NC2 = %s \n" %(metric , str(abs ( metrica - acc_NC2 )))) 674 if l == 3 or l == 6: 675 file. write (" %s AB = %s \n" %(metric , str(abs ( metrica - acc_AB )))) 676 file. write (" %s ABB = %s \n" % (metric , str(abs( metrica - acc_ABB )))) 677 file. write (" %s TM = %s \n" %(metric , str(abs ( metrica - acc_TM )))) 678 if l == 6: 679 file. write (" %s AB - NC1 = %s \n" %(metric , str(abs ( metrica - acc_AB_NC1 ) ))) 680 file. write (" %s ABB -NC1 = %s \n" %(metric , str(abs (metrica - acc_ABB_NC1 )))) 681 file. write (" %s TM -NC1 = %s \n" %(metric , str(abs ( metrica - acc_TM_NC1 )) )) 682 file. write (" %s AB - NC2 = %s \n" %(metric , str(abs ( metrica - acc_AB_NC2 ) ))) 683 file. write (" %s ABB -NC2 = %s \n" %(metric , str(abs (metrica - acc_ABB_NC2 )))) 684 file. write (" %s TM - NC2 = %s \n" %(metric , str(abs ( metrica - acc_TM_NC2 ) ))) 685 686 def guarda_foto ( archivo , imagen , titulo ): 687 plt . title ( titulo ) 688 fig = imagen 689 plt . savefig ( archivo ) 690 plt. clf () 691 692 """# Funciones para la ejecucion """ 693 694 def crea_guarda_shaps (X,y, nombre =""): 695 # Nombre seria por ejemplo "AB - NC1" 696 #Ya esta estandarizado 697 #X = standarize (X) 698 699 700 # Construimos el modelo con los parametros obtenidos automaticamente 701 params , tuner = buscar_model (X,y ,trials =400 , executions =5) 702 703 model = tuner . hypermodel . build ( params ) 704 model . fit (X,y, epochs = params . get (’epochs’) ,batch_size = params .get (’ batch_size ’),verbose =0) 705 104