Support vector regression : propiedades y aplicaciones

TRABAJO FIN DE GRADO SUPPORT VECTOR REGRESSION: PROPIEDADES Y APLICACIONES Realizado por: Juan José Martín Guareño Supervisado por: Dr. Rafael Blanquero Bravo y Dr. Emilio Carrizosa Priego FACULTAD DE MATEMÁTICAS DEPARTAMENTO DE ESTADÍSTICA E INVESTIGACIÓN OPERATIVA Support Vector Regression: Propiedades y Aplicaciones 3 Índice general 1. Introducción 7 2. Support Vector Machines 9 2.1. SVM Caso Lineal Separable . . . . . . . . . . . . . . . . . . . . . . 9 2.2. SVM Caso Lineal No Separable . . . . . . . . . . . . . . . . . . . . 11 2.3. SVMCasoNoLineal.......................... 12 3. Support Vector Regression 15 3.1. Funciones de Pérdida . . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.2. FuncionaldeRiesgo........................... 16 3.3. CasoLineal ............................... 18 3.4. CasoNoLineal ............................. 21 3.4.1. Truco del Kernel . . . . . . . . . . . . . . . . . . . . . . . . 22 3.4.2. MatrizdeGram......................... 23 3.4.3. FunciónKernel ......................... 24 3.4.4. Núcleo de Mercer . . . . . . . . . . . . . . . . . . . . . . . . 24 3.4.5. MatrizKernel.......................... 25 3.4.6. Construcción del Núcleo . . . . . . . . . . . . . . . . . . . . 26 3.4.7. Tipos de Kernel . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.4.8. Kernels Asimétricos . . . . . . . . . . . . . . . . . . . . . . . 33 3.5. Estimación de Parámetros del SVR . . . . . . . . . . . . . . . . . . 35 3.5.1. Validación Cruzada . . . . . . . . . . . . . . . . . . . . . . . 36 4. Support Vector Regression en R 37 4.1. Funciones Específicas SVR . . . . . . . . . . . . . . . . . . . . . . . 37 4.2. EjemploPráctico ............................ 39 4.2.1. Análisis Descriptivo . . . . . . . . . . . . . . . . . . . . . . . 40 4.3. Regresión mediante SVR . . . . . . . . . . . . . . . . . . . . . . . . 44 4.3.1. KernelLineal .......................... 45 4.3.2. Kernel Polinomial . . . . . . . . . . . . . . . . . . . . . . . . 47 4.3.3. KernelRadial.......................... 50 4.3.4. Kernel Sigmoidal . . . . . . . . . . . . . . . . . . . . . . . . 52 4.4. Comparación de Modelos . . . . . . . . . . . . . . . . . . . . . . . . 56 Bibliografía 60 Support Vector Regression: Propiedades y Aplicaciones 5 Abstract Statistical learning plays a key role in many areas of science, finance and industry. In particular, supervised learning plays a key role in the fields of statistics, data mining and artificial intelligence, intersecting with areas of engineering and other disciplines. Mathematical optimization has played a crucial role in supervised learning. Techniques from very diverse fields within mathematical optimization have been shown to be useful. Support Vector Machine (SVM) and Support Vector Regression (SVR) are ones of the main exponents as application of the mathematical optimization to supervised learning. SVM and SVR are state of the art methods for supervised learning and regression. These geometrical optimization problems can be written as convex quadratic optimization problems with linear constraints, in principle solvable by any nonlinear optimization procedure. In this work we analyze SVMs and SVRs: how the problems are obtained and expressed in a manageable way. On the one hand, we describe the techniques used by the algorithms of supports vectors dedicated to the classification, in linear and nonlinear cases. On the other hand, we focus on the theoretical development of the techniques in the field of support vector regression. We pay more attention to the nonlinear case, where the algorithm of support vector shows its full potential, using a kernel function to calculate a nonlinear approximation function. Finally we bring these theoretical procedures into practice with the help of the statistical language and environment R. Support Vector Regression: Propiedades y Aplicaciones 7 Capítulo 1 Introducción El aprendizaje estadístico desempeña un papel clave en muchas áreas de la ciencia, las finanzas y la industria. Algunos ejemplos de problemas de aprendizaje son: predecir si un individuo que ha tenido un accidente de tráfico volverá a tener otro accidente, en base a su conducción; predecir el valor de una acción en seis meses desde el momento del estudio, sobre la base de las medidas de rendimiento de la empresa y los datos económicos; o identificar los factores de riesgo para una enfermedad, basándose en variables clínicas y demográficas. En particular, el aprendizaje supervisado es una parte importante del aprendizaje estadístico, con gran influencia en los campos de la minería de datos y la inteligencia artificial. Partiendo de un conjunto C , constituido por ejemplos de una población específica de los cuales se conocen una serie de variables, una vez determinada la variable respuesta y las variables explicativas, se buscan procedimientos para predecir el valor de la variable respuesta, conocido el valor de las variables explicativas de un futuro ejemplo de dicha población. En el caso de que la variable respuesta fuese cualitativa estaríamos ante un problema de clasificación, y si la variable respuesta fuese cuantitativa tendríamos un problema de regresión. La optimización matemática es crucial en el aprendizaje supervisado, donde han demostrado ser útiles técnicas precedentes de muy diversos campos dentro de la misma. Las máquinas de vectores soporte (SVM) y la regresión de vectores soporte (SVR) son unos de los principales exponentes de la aplicación de la optimización matemática al aprendizaje supervisado. donde han demostrado ser útiles técnicas precedentes de muy diversos campos. Tanto para el caso de las (SVM) como (SVR), el objetivo es realizar la predicción a partir de un problema de optimización geométrica que se puede escribir como un problema de optimización cuadrático convexo con restricciones lineales, en principio resoluble mediante cualquier procedimiento de optimización no lineal. En este documento se muestra cómo podemos obtener dichos problemas y las transformaciones necesarias para facilitar su resolución. En el capítulo 2 describiremos las técnicas utilizadas por los algoritmos de vectores soportes dedicados a la clasificación, en los casos Lineal y No Lineal. 8 CAPÍTULO 1. Introducción El capítulo 3 se centra en el desarrollo teórico de las técnicas de los vectores soporte en el ámbito de la regresión. Se presenta una primera parte donde se realiza el estudio de la función de regresión a partir de una función lineal y una segunda parte más interesante donde el algoritmo de los vectores soporte muestra todo su potencial, con el uso de una función kernel para el cálculo de una función de aproximación no lineal. En el capítulo 4 llevaremos estos procedimientos teóricos a la práctica con la ayuda del programa estadístico R. Support Vector Regression: Propiedades y Aplicaciones 9 Capítulo 2 Support Vector Machines Las máquinas de vectores soporte (SVM, del inglés Support Vector Machines) tienen su origen en los trabajos sobre la teoría del aprendizaje estadístico y fueron introducidas en los años 90 por Vapnik y sus colaboradores, [ 1 ]. Aunque originariamente las SVMs fueron pensadas para resolver problemas de clasificación binaria, actualmente se utilizan para resolver diversos tipos de problemas, por ejemplo, la regresión, en el cual nos centraremos más adelante. De hecho, desde su introducción, han ido ganando un merecido reconocimiento gracias a sus sólidos fundamentos teóricos. La idea es seleccionar un hiperplano de separación que equidiste de los ejemplos más cercanos de cada clase para, de esta forma, conseguir lo que se denomina un margen máximo a cada lado del hiperplano. Además, a la hora de definir el hiperplano, sólo se consideran los ejemplos de entrenamiento que distan del hiperplano la distancia margen. Estos ejemplos reciben el nombre de vectores soporte. 2.1. SVM Caso Lineal Separable Dado un conjunto separable de ejemplos C = {(x1, y1)..., (xn, yn)} , donde xi∈ Rd e yi∈ { +1 ,− 1 } se puede definir un hiperplano de separación, como una función lineal que es capaz de separar dicho conjunto: D(x) = (w1x1+... +wdxd) + b=hw, xi+b , donde wi∈R∀i= 1, ..., d yb∈R. Dicho hiperplano deberá cumplir las siguientes desigualdades: hw, xii+b≥0si yi= +1 ∀i= 1, ..., n hw, xii+b≤0si yi=−1∀i= 1, ..., n O escrito de forma más compacta de la siguiente forma: yi(hw, xii+b)≥0∀i= 1, ..., n 16 CAPÍTULO 3. Support Vector Regression Buscamos una función de regresión, en principio lineal, es decir, de la siguiente forma: f(x)=(w1x1+... +wdxd) + b=hw, xi+b, donde wi∈R∀i= 1, ..., d yb∈R. Lo que pretendemos es definir una función de pérdida que ignore los errores asociados con los puntos que caen dentro de una banda que está a una cierta distancia de la función de regresión lineal. Es decir, si el punto verifica |y−f ( x ) | ≤ ε la función pérdida deberá anularse. Los ejemplos que disten exactamente ε de nuestra función se denominarán vectores soporte. Siguiendo esta estrategia nosotros podemos definir las siguientes funciones de pérdida: Lε 1(y, f(x)) = máx {0,|y−f(x)|−ε} Lε 2(y, f(x)) = máx {0,(y−f(x))2−ε} Dichas funciones son denominadas funciones de pérdida ε -insensible lineal y cuadrática respectivamente. En lo que sigue centraremos el estudio de la SVR tomando como función perdida la función de pérdida ε -insensible lineal, se puede hacer un estudio similar para la otra función de pérdida. 3.2. Funcional de Riesgo Hasta ahora el algoritmo SV para la regresión puede parecer bastante extraño y casi no relacionado con otros métodos existentes de estimación de la función de regresión. Sin embargo una vez expresado con notación matemática más estándar se observa la conexión. Para más simplicidad vamos a tener solo en cuenta solo al caso lineal, sabiendo que podría extenderse al caso no lineal utilizando el método del kernel. Supongamos que disponemos de un conjunto C = {(x1, y1), ..., (xn, yn)} , donde xi∈Rd e yi∈R . Asumimos que este conjunto de datos se ha obtenido a partir de cierta función de probabilidad P ( x, y ). Nuestro objetivo será encontrar una función fque minimice el funcional: R[f] = Zc(x, y, f(x))dP(x, y) Vamos a denotar como c ( x, y, f ( x )) a la función que penaliza los errores de estimación. Cuando nos refiramos a esta función lo haremos como función coste. 3.2. Funcional de Riesgo 17 Basándonos en los datos empíricos y dado que no conocemos la medida de probabilidad, solo podemos usar los datos de S para la estimación de una función que minimiza R [ f ]. Una posible aproximación consiste en la sustitución de la integración por la estimación empírica para conseguir el llamado funcional de riesgo empírico: Remp[f] := 1 n n X i=1 c(xi, yi, f(xi)) Un primer intento sería encontrar la función f := argmin f∈H ( Remp [ f ]) para alguna hipotética clase H . Sin embargo si la clase H es muy amplia, hay una gran cantidad de ejemplos, cuando se trata de unos datos en espacios de muy alta dimensión esto puede no ser una buena idea, ya que se dará sobreajuste y por tanto malas propiedades de generalización. Por lo tanto le añadimos un termino de control de capacidad, que en el caso del algoritmo SV es kwk2 , lo que conduce al funcional de riesgo generalizado: Rreg[f] := λ 2kwk2+Remp[f] donde λ > 0es llamada la constante de regularización. La pregunta que nos surge ahora es, ¿ Qué función de coste debemos usar en Rreg ? La función coste utilizada para el caso descrito de SVR es: c(x, y, f(x)) = Lε 1(y, f(x)) = |y−f(x)|ε, con |z|ε= m´ax{0,|z|−ε}. Las funciones de coste del tipo |y−f ( x ) |p ε con p > 1pueden no ser deseable, ya que el aumento superlineal conduce a una pérdida de las propiedades de robustez del estimador (véase ejemplos en [ 5 ]). En estos casos la derivada de la función coste puede crecer sin límite. Para p < 1tendríamos una función de coste no convexa. Para el caso c ( x, y, f ( x )) = Lε 2 ( y, f ( x )), es decir, para la función ε− insensible cuadrática, tendríamos una aproximación en forma de mínimos cuadráticos que, a diferencia de nuestra función de coste estándar del SV, conduce a una inversión de la matriz en lugar de un problema de programación cuadrática. Por un lado vamos a pretender evitar el uso de funciones complejas ya que esto puede conducir a problemas de optimización difíciles, por otro lado podríamos usar una función de coste que se adapte mejor a los datos. Bajo el supuesto de que las muestran fuesen generadas por una dependencia funcional subyacente, más un ruido aditivo, es decir: yi = fv ( xi ) + ξi con densidad p ( ξ ), la función de coste óptima un sentido de probabilidad máxima sería: c(x, y, f(x)) = −log p(y−f(x)) 18 CAPÍTULO 3. Support Vector Regression Algunos ejemplos de funciones de coste asociadas a modelos de densidad conocidos son: Funciones de Coste ε-insensitive c(x, y, f(x)) = |y−f(x)|ε Laplace c(x, y, f(x)) = |y−f(x)| Gauss c(x, y, f(x)) = (y−f(x))2 Huber c(x, y, f(x)) =          1 2σ(y−f(x))2si |y−f(x)| ≤ σ |y−f(x)|− σ 2si |y−f(x)|> σ Polinomial 1 p|y−f(x)|p Polinomial a trozos c(x, y, f(x)) =            1 pσp−1(y−f(x))psi |y−f(x)| ≤ σ |y−f(x)|−σp−1 psi |y−f(x)|> σ Cuadro 3.1: Tabla correspondiente a algunas funciones de coste asociadas a funciones de densidad comunes. Las funciones coste deben ser funciones convexas. En caso de que una función no sea convexa habría que encontrar una aproximación convexa con el fin de hacer frente a la situación de manera eficiente. El requisito de convexidad se pide porque queremos asegurar la existencia y unicidad de un mínimo en los problemas de optimización, y que los óptimos locales sean globales [2]. 3.3. Caso Lineal Vamos a considerar en este caso un conjunto de igual estructura que en el caso de SVM, modificando el posible valor de la variable respuesta. Sea el conjunto C = {(x1, y1), ..., (xn, yn)} , donde xi∈Rd e yi∈R . Así pues podríamos plantear el problema de forma análoga al caso de SVM pero con las mencionadas diferencias como: mín 1 2kwk2+C n X i=1 (ξi+ξ∗ i) s.a yi−(hw, xii+b)≤ε+ξii= 1, ..., n (hw, xii+b)−yi≤ε+ξ∗ ii= 1, ..., n ξi≥0, ξ∗ i≥0i= 1, ..., n. 3.3. SVR Caso Lineal 19 La constante C > 0determina el equilibrio entre la regularidad de f y la cuantía hasta la cual toleramos desviaciones mayores que ε . Consideraremos ξi y ξ∗ i las variables que controlan el error cometido por la función de regresión al aproximar el i-ésimo ejemplo. Figura 3.1: Ahora la idea es englobar a todos los ejemplos en una banda en torno a nuestra función de predicción f. Así pues un valor muy grande de la constante C , en el caso límite ( C→ ∞ ) estaríamos considerando que el conjunto está perfectamente representado por nuestro hiperplano predictor ( ξi→ 0). Por contra, un número demasiado pequeño para C permitiría valores de ξi elevados, es decir, estaríamos admitiendo un número muy elevado de ejemplos mal representados. Tras la obtención del problema primal pasamos a plantear el problema dual asociado. Para ello determinamos la función lagrangiana. La idea clave es construir una función de Lagrange con la función objetivo y las restricciones correspondientes, mediante la introducción de un conjunto de variables duales. Nuestro objetivo es obtener el problema dual a partir de las condiciones de dualidad fuerte [ 6 ]. Así definimos L := L ( w, b, ξ, ξ∗, α, α∗, η, η∗ ), donde w , b , ξ , ξ∗ son las variables originales del problema y α , α∗ , η , η∗ son las variables duales asociadas a las restricciones. Esta función tiene un punto de silla respecto a las variables del primal. Una demostración de esto se puede encontrar en [6]. Con esto, L:= 1 2kwk2+C n X i=1 (ξi+ξ∗ i)− n X i=1 αi(ε+ξi−yi+hw, xii+b) − n X i=1 α∗ i(ε+ξ∗ i+yi−hw, xii−b)− n X i=1 (ηiξi+η∗ iξ∗ i) lo que buscamos es hallar el: m´ax α,α∗,η,η∗{m´ın w,b,ξ,ξ∗L(w, b, ξ, ξ∗, α, α∗, η, η∗)}. Se entiende que las variables del dual son positivas. Ademas se deduce de la 20 CAPÍTULO 3. Support Vector Regression condición de punto silla que las derivadas parciales respecto de las variables del primal deberían anularse para el óptimo. δbL= n X i=1 (α∗ i−αi) = 0 (3) δwL=w− n X i=1 (α∗ i−αi)xi= 0 ⇒w= n X i=1 (α∗ i−αi)xi δξiL=C−αi−ηi= 0 ⇒ηi=C−αi≥0 δξ∗ iL=C−α∗ i−η∗ i= 0 ⇒η∗ i=C−α∗ i≥0 Como la única restricción tanto para η como para η∗ es ser positivas, podemos expresarlas en función de αi y α∗ i . Sustituyendo las ecuaciones anteriores en la función de Lagrange obtenemos: m´ax αi,α∗ i{1 2 n X i,j=1 (αi−α∗ i)(αj−α∗ j)hxi, xji+C n X i=1 (ξi+ξ∗ i) − n X i=1 αi(ε+ξi−yi+h n X i=1 (α∗ i−αi)xi!, xii+b) − n X i=1 α∗ i(ε+ξ∗ i+yi+h n X i=1 (α∗ i−αi)xi!, xii−b) − n X i=1 ((C−αi)ξi+ (C−α∗ i)ξ∗ i)} Simplificando la expresión anterior e imponiendo la restricción (3) aun no utilizada, obtenemos el siguiente problema de optimización dual: máx −1 2 n X i,j=1 (αi−α∗ i)(αj−α∗ j)hxi, xji−ε n X i=1 (αi+α∗ i) + n X i=1 yi(αi−α∗ i) s.a n X i=1 (αi−α∗ i)=0 αi, α∗ i∈[0, C] A partir de estas expresiones, también estaríamos extrayendo una expresión de nuestra función de predicción: f(x) = n X i=1 (αi−α∗ i)hxi, xi+b De esta forma estaríamos obteniendo la función buscada sin depender la resolución del problema de la dimensión en la que se encuentra nuestros ejemplos de entrada y pasaría a depender únicamente de los vectores soporte. Para completar nuestra función de regresión deberíamos calcular b para ello usamos las condiciones de complementariedad de holgura de Karush-Kuhn-Tucker 3.4. SVR Caso No Lineal 21 (KKT) [ 9 ]. Estas establecen que en la solución óptima el producto entre las variables de holgura y las restricciones duales deben anularse, es decir: αi(ε+ξi−yi+hw, xii+b)=0 α∗ i(ε+ξ∗ i+yi−hw, xii−b)=0 (C−αi)ξi= 0 (C−α∗ i)ξ∗ i= 0, esto nos permite deducir varias conclusiones. En primer lugar, solo los ejemplos ( xi, yi ) tales que αi = 0 ó α∗ i = 0 estarían fuera de nuestro tubo o banda construido. Por otra parte αiα∗ i = 0, es decir, no se pueden activar a la vez las dos variables duales asociadas a un mismo ejemplo. Por último en los casos que αi, α∗ i∈ (0 , C )tendríamos que la variable ξi, ξ∗ i correspondiente se debe anular, por lo que podríamos despejar el valor de bde las primeras restricciones. b=yi−hw, xii−εsi αi∈(0, C) b=yi−hw, xii+εsi α∗ i∈(0, C) 3.4. Caso No Lineal De igual forma que nos planteamos la búsqueda de la función de regresión como una función lineal, podíamos plantearnos el caso que quisiésemos una aplicación más general. Consideramos el mismo conjunto de antes: C = {(x1, y1), ..., (xn, yn)} , donde xi∈Rd e yi∈R . Sea Φ: X→ F la función que hace corresponder a cada punto de entrada x un punto en el espacio de características F , donde F es un espacio de Hilbert. Este espacio de características puede ser de dimensión elevada o incluso infinita. Nuestro objetivo va a ser trasladar nuestros ejemplos a este nuevo espacio de características y aquí encontrar la función que mejor aproxime las imágenes de nuestro conjunto. En este caso la función que estamos buscando será de la siguiente forma: f(x) = hw, Φ (x)i+b Plantearíamos el problema primal al igual que en el caso lineal a diferencias que en este caso no nos depende directamente de los ejemplos del conjunto sino de sus imágenes por una cierta función φ. mín 1 2kwk2+C n X i=1 (ξi+ξ∗ i) s.a yi−hw, φ(xi)i ≤ ε+ξii= 1, ..., n hw, φ(xi)i−yi≤ε+ξ∗ ii= 1, ..., n ξi≥0, ξ∗ i≥0i= 1, ..., n 22 CAPÍTULO 3. Support Vector Regression Como ya hemos dicho anteriormente, la complejidad de este problema va a depender de la dimensión en la que se encuentren nuestros ejemplos, y en este caso, tras ser transformados por la función φ , estos ejemplos podrían tener una dimension extremadamente alta, lo que complicaría en exceso la resolución de este problema primal. Así pues, procederemos a construir el problema dual asociado. Como ya lo transformamos para el caso lineal, por analogía podríamos deducir el dual buscado. máx −1 2 n X i,j=1 (αi−α∗ i)(αi−α∗ i)hφ(xi), φ(xj)i−ε n X i=1 (αi+α∗ i) + n X i=1 yi(αi−α∗ i) s.a n X i=1 yi(αi−α∗ i)=0 αi, α∗ i∈[0, C] 3.4.1. Truco del Kernel ¿Podríamos resolver este problema sin llegar a conocer explícitamente la función φ? La respuesta es que sí. Tras plantear el problema dual, observamos que la función objetivo solo depende del producto interno de las imágenes de nuestros ejemplos. El algoritmo del truco del kernel es ampliamente utilizado en los algoritmos de cálculo de productos internos de la forma hΦ(x),Φ(x0)ien el espacio de características F. El truco consiste en que, en lugar de calcular estos productos internos en F , lo que realmente utilizamos computacionalmente, debido a su posible alta dimensionalidad, es definir una función kernel , K : X×X→R que asigna a cada par de elementos del espacio de entrada X , un valor real correspondiente al producto escalar de las imágenes de dichos elementos en el nuevo espacio F, es decir, K(x, x0) = hΦ(x),Φ(x0)i donde Φ : X→ F . Los primeros en aplicarlo para este tipo de problemas fueron Cortes y Vapnik en [1]. Ejemplo: Consideramos un espacio de entrada bidimensional, X⊆R2 y sea Φla función definida de la siguiente forma: Φ : x= (x1, x2)7−→ Φ(x) = (x2 1, x2 2,√2x1x2)∈ F =R3 A la hora de evaluar el producto interno en el espacio Fnos queda: hφ(x), φ(z)i=h(x2 1, x2 2,√2x1x2),(z2 1, z2 2,√2z1z2)i =x2 1z2 1+x2 2z2 2+ 2x1x2z1z2 = (x1z1+x2z2)2 =hx, zi2 3.4. SVR Caso No Lineal 23 Luego en este caso la función kernel sería: K(x, z) = hx, zi2 Problema a resolver aplicando el kernel: Una vez hayamos fijado el kernel a utilizar pasaríamos a resolver el problema de la forma: máx −1 2 n X i,j=1 (αi−α∗ i)(αi−α∗ i)K(xi, xj)−ε n X i=1 (αi+α∗ i) + n X i=1 yi(αi−α∗ i) s.a n X i=1 yi(αi−α∗ i) = 0 αi, α∗ i∈[0, C], y la función de predicción: f(x) = n X i=1 (αi−α∗ i)K(xi, x) + b 3.4.2. Matriz de Gram Dado un conjunto de puntos C = {x1, ..., xl} la matriz de Gram se define como G cuyas entradas son: gij = hxi, xji . Si usamos una función kernel K para evaluar los productos internos en el espacio de características con la función característica φ . Las entradas asociadas a la matriz de Gram son: kij = hφ ( xi ) , φ ( xj ) i. En estos casos la matriz también es conocida como matriz kernel y presenta la forma característica: K:=       K(x1, x1)K(x1, x2)··· K(x1, xl) K(x2, x1)K(x2, x2)··· K(x2, xl) . . .. . ..... . . K(xl, x1)K(xl, x2)··· K(xl, xl)       Para el caso que usemos una función kernel, esta matriz es simétrica por construcción y además es semidefinida positiva, es decir, si consideramos el caso general de la matriz kernel tenemos: kij =K(xi, xj) = hφ(xi), φ(xj)i,con i, j = 1, ..., l. Dado α∈Rlcualquiera se tiene: v0K v = l X i,j=1 vivjKij = l X i,j=1 vivjhφ(xi), φ(xj)i =*l X i,j=1 viφ(xi), l X i,j=1 vjφ(xj)+ =      l X i,j=1 viφ(xi)      2 ≥0 24 CAPÍTULO 3. Support Vector Regression 3.4.3. Función Kernel Un producto interno en el espacio de características tiene asociado un kernel equivalente en el espacio de entradas: K(x, z) = hφ(x), φ(z)i Hemos visto cómo construir una matriz de evaluaciones, aplicando la función núcleo a todos los pares un conjunto de entradas. La matriz resultante es semidefinida positiva. Ahora vamos a introducir un método alternativo para demostrar que una función candidata a ser kernel, de verdad lo es. Esto representa una de las herramientas teóricas necesarias para la creación de nuevos y más complejos núcleos a partir de otros más simples. Una de las observaciones más importantes es la relación con las matrices semidefinidas positivas. Como hemos visto, la matriz kernel evaluada para cualquier conjunto de datos es semidefinida positiva. Dada esta peculiaridad de las funciones kernel nos lleva a definir lo siguiente: una función K : X×X→R se dice finitamente semidefinida positiva, si es simétrica y las matrices formadas con cualquier conjunto finito del espacio X son semidefinidas positivas. Esta definición no requiere que Xsea siquiera un espacio vectorial. Teorema 3.1 Sea una función K : X×X→R , ya sea continua o con dominio finito. Se puede descomponer en: K(x, z) = hφ(x), φ(z)i con cierta función característica φ en un espacio de Hilbert F aplicado a los argumentos y seguido de la evaluación del producto interno en F si y solo si K es finitamente semidefinida positiva. Una prueba de este resultado la podemos encontrar en [14]. Dada una función K que satisfaga la propiedad anterior, es decir, que sea finitamente semidefinida positiva, llamamos al correspondiente espacio Fk como el espacio de Hilbert reproducido por el kernel, del inglés (RKHS). 3.4.4. Núcleo de Mercer Ahora sí estamos en las condiciones óptimas para pasar al teorema de Mercer como consecuencia del análisis anterior. El teorema de Mercer se utiliza generalmente para la construcción del espacio de características para un kernel válido. Dado que ya hemos logrado esto con la construcción (RKHS), nosotros no necesitamos el teorema de Mercer en sí mismo. Lo incluimos ya que define el espacio de características en términos de un vector de características explícito, en lugar de utilizar el espacio de funciones de nuestra construcción (RKHS). 3.4. SVR Caso No Lineal 25 Teorema 3.2 Sea X un conjunto compacto en Rn . Supongamos que K es una función continua y simétrica tal que el operador Tk:L2(X)−→ L2(X) (Tkf)(·) = ZX K(·, x)f(x)dx es positivo. Esto es: ZX×X K(x, z)f(x)f(z)dx dz ≥0∀f∈L2(X), entonces nosotros podemos desarrollar K ( x, z )en una serie uniformemente convergente, en X×X, en términos de unas funciones φj, satisfaciendo hφj, φii=δij K(x, z) = ∞ X j=1 φj(x)φj(z) Además la serie P∞ j=1 kφjk2 L2(X)es convergente. La prueba de este resultado la podemos encontrar en [14]. 3.4.5. Matriz Kernel Dado un conjunto de datos de entrenamiento Ce = {x1, ..., xl} y una función kernel K ( ·,· )nosotros introducimos facilmente la matriz de Gram o en este caso también llamada matriz kernel con las entradas Kij = H ( xi, xj ). En la sección anterior hemos visto que si la función K es un kernel válido la matriz kernel asociada a cada conjunto de entrenamiento Ce tiene que ser semidefinida positiva. Lo que llamabamos ser la función finitamente semidefinida positiva (FPSF). Este hecho nos permite manipular los kernels sin tener en cuenta necesariamente el espacio de características correspondiente. A condición de que mantengamos la propiedad, tenemos garantizado un kernel válido, es decir, que existe un espacio de características para la cual es la función núcleo correspondiente. Este razonamiento sobre la medida de similitud implica que la función del núcleo puede ser más natural que la realización de una construcción explícita de su espacio de características. La modularidad intrinseca de máquinas kernel también significa que cualquier función kernel se puede utilizar siempre que produzca matrices simétricas y semidefinidas positivas, y podemos aplicar cualquier algoritmo de núcleo, siempre y cuando podamos aceptar como entrada una matriz adecuada a nuestro conjunto de datos. En vista de nuestra caracterización de los núcleos, en cuanto a la propiedad (FPSF), queda claro que la matriz del núcleo es el ingrediente básico en la teoría de los métodos con kernels. Esta matriz contiene toda la información disponible para llevar a cabo la etapa de aprendizaje, con la única excepción de las etiquetas de salida en el caso de aprendizaje supervisado. Vale la pena tener en cuenta que es solo a través de la matriz kernel que el algoritmo de aprendizaje obtiene la 32 CAPÍTULO 3. Support Vector Regression En definitiva el all −subsets kernel es definido por la función φ: φ:x7−→ (φA(x))A⊆{1,...,n} con el correspondiente kernel K⊆(x, z)dado por: K⊆(x, z) =hφ(x), φ(z)i= =X A⊆{1..n}Y i∈A xizi ANOVA Kernels Tanto el kernel polinomial, como all-subsets kernel, tienen limitado el control de qué características usar y la forma de ponderar sus parámetros. Nosotros ahora presentamos un método que permite una mayor libertad en la especificación del conjunto de monomios. El ANOVA kernel Kd es como el all-subsets, excepto que se restringe a subconjuntos cardinalidad d . Podemos usar la notación xi para la expresión xi1 1xi2 2···xin n , donde~ i= (i1, ..., in)∈ {0,1}ny además la restricción: n X j=1 ij=d para el caso d= 0 hay una característica con valor constante 1 correspondiente al conjunto vacío. La diferencia entre el ANOVA y el polinomial es la exclusión de las coordenadas repetidas. Para más información acerca de la historia de este método ver la sección 9.10 de [14]. La función característica del núcleo ANOVA de grado dviene dado por: φd:x7−→ (φA(x))|A|=d donde para cada subconjunto Aes definido por: φA(x) = Y i∈A xi=xiA donde iAes la función indicador del conjunto A. La dimensión del resultante φA es claramente n d , ya que es el número de tales subconjuntos, mientras que el producto interno resultante viene dado por: Kd(x, z) = hφd(x), φd(z)i=X |A|=d φA(x)φA(z) =X 1≤i1<i2<···<id≤n (xi1zi1)(xi2zi2)···(xidzid) =X 1≤i1<i2<···<id≤n d Y j=1 xijzij 3.4. SVR Caso No Lineal 33 Perceptrón Multicapa ó Sigmoidal El perceptrón más establecido, con una capa oculta, también tiene una representación válida en los kernels, K(x, x0) = tanh(γhx, x0i+R) para ciertos valores escalares γ y R que aquí serían los parámetros de nuestro kernel. Splines Los Splines son una opción muy usada en el modelado debido a su flexibilidad. Un Spline finito, de orden κcon Nnodos situados en τses de la forma: K(x, x0) = κ X r=0 xrx0r+ N X s=1 (x−τs)κ +(x0−τs)κ +dτ Se define un Spline infinito en el intervalo [0,1) como: K(x, x0) = κ X r=0 xrx0r+Z1 0(x−τs)κ +(x0−τs)κ +dτ En el caso particular de que κ= 1,(S∞ 1)el kernel está dado por: K(x, x0) = 1 + hx, x0i+1 2hx, x0imin(x, x0)−1 6min(x, x0)3, donde la solucion es cúbica a trozos. En principio los kernels splines son kernels univariantes, pero gracias a la propiedad de que a partir del producto tensorial de kernels se obtiene un nuevo kernel, podemos construir kernels splines multidimensionales. Series de Fourier Las series de Fourier consideran una expansión en un espacio de caracteristicas 2N+ 1 dimensional. Este kernel esta definido en el intervalo [-π 2,π 2]. K(x, x0) = sin(N+1 2)(x−x0) sin(1 2(x−x0)) Al igual que los kernels splines, pueden usarse de base para construir kernels multidimensionales. 3.4.8. Kernels Asimétricos Hay veces que, dado el conjunto C = {x1, ..., xn} que representan l objetos en Rd , solo disponemos de una medida de similitud entre los objetos no simétrica, que podemos acumular en una matriz S∈ Mn×n . Como hemos dicho que la medida de similitud no es simétrica, sij 6=sji para algún par i,j. Algunos ejemplos reales en los que se obtiene una medida de similitud no simétrica pueden ser: 34 CAPÍTULO 3. Support Vector Regression Página web que enlaza a otra; una puede enlazar a la otra sin que esta enlace sobre la misma. Niños en la clase que quieran estar sentados juntos; dos niños pueden no tener la misma preferencia en estar sentados juntos. Aparición en diferentes textos, la palabra i puede aparecer en textos donde aparece la palabra j, pero no al contrario. Tanto para la clasificación como para la regresión las funciones de decisión son respectivamente de la forma: f1(x) = n X i=1 αiyiK(xi, x) + b f2(x) = n X i=1 (αi−α∗ i)K(xi, x) + b En caso de que la matriz K propuesta a ser matriz kernel sea asimétrica procederemos mediante una serie de técnicas para la simetrización de la correspondiente matriz. Interpretación de Asimetría Hay una particular forma de elegir sij , que denotaremos con ∧ , y definimos como: sij =|xi∧xj| |xi|=Pn k=1 |m´ın(xik, xjk)| Pn k=1 |xik| En particular sij ≥0ysij ≤1 Asumimos que acumulamos nuestro conjunto de datos C en una matriz de datos X. Esta medida es útil para diversos tipos de ejemplos. Para el caso de que la matriz X sea una matriz de términos por documentos, la medida |xi| mediría el número de documentos indexados por el término i , y |xi∧xj| el número de documentos indexados por ambos. Por lo tanto sij puede ser interpretado como el grado en el que el tema representado por el término i es un subconjunto del tema representado por el término j. Esta medida numérica de subsethood es debida a Kosko [7]. En el caso de ser una matriz de cocitación |xi| es el número de citas recibidas por el autor o página web i u |xi∧xj| son los autores o webs que se citan mutuamente. Este tipo de problemas tienen en común que las normas de los individuos calculadas por |xi| siguen la ley de Zipf [ 12 ]: hay algunos individuos con unas normas muy grandes (muy citados) y hay una gran cantidad de individuos con normas muy pequeñas. Esta asimetría puede interpretarse como una especie de jerarquía de los individuos que se organizan en una especie de árbol. En lo alto se encuentran los individuos con grandes normas que corresponden a individuos muy famosos. En la base se encuentran individuos con norma pequeña que son individuos raros. Con este tipo de norma las matrices de similitud que se generan son asimétricas y frente a este tipo de matrices se usan unos métodos para poder obtener unas 3.5. Estimación de Parámetros del SVR 35 matrices kernels con todas las propiedades necesarias. Describimos a continuación algunos métodos conocidos: 1) Una forma sencilla para lograr la simetrización de nuestra matriz es considerar la matriz de valores medios. De la descomposición: S=1 2(S+ST) + 1 2(S−ST) nos quedamos con la parte simétrica, quedandonos con la matriz K=1 2(S+ST). Sin embargo de esta forma despreciamos información importante, además nosotros ignoramos la parte asimétrica. 2) El método sugerido por Schölkopf [ 13 ] consiste en tomar como matriz kernel K = STS . Este método gana sentido cuando tenemos matrices de cocitación porque Kij = 1 cuando hay un único k tal que Ski = Skj = 1, es decir, existe un autor que cita simultaneamente a i y a j . Sin embargo vamos a perder un caso de similitud que se produce cuando dos autores citan a su vez a un tercero. Esta información la recoge SST. 3) Para el caso de la clasificación podemos utilizar la etiqueta del conjunto al que corresponden para la construcción de nuestra matriz kernel. La construcción de la matriz Kla hace a partir de la matriz de similitud dada de la forma: K={kij}n i,j=1, kij =(m´ax(sij, sji)si xi, xj∈C1óxi, xj∈C2 m´ın(sij, sji)c.c siendo C1 y C2 las dos clases diferenciadas a clasificar. Este método es denominado pick-out. De esta forma conseguimos que la matriz K sea simétrica. Para una mayor información de porque consideramos la matriz kernel de esta forma podemos ayudarnos de [11]. Sin embargo con estas transformaciones nos aseguraríamos la simetría de la matriz, pero no necesariamente cumple la propiedad semidefinida positiva. En este caso sustituiríamos la matriz K por K + λI , con λ > 0suficientemente grande, para que los autovalores se hagan no negativos. 3.5. Estimación de Parámetros del SVR Los problemas primal y dual descritos anteriormente, dependen de una serie de parámetros que pueden hacer variar la solución. Además, cuando introducimos el kernel, éste depende a su vez de unos parámetros de los cuales depende nuestra solución. Es necesario dar un procedimiento de selección de dichos parámetros. El procedimiento más simple y utilizado consiste en definir una malla sobre el espacio de los parámetros y seleccionar aquellos valores de los parámetros que proporcionan mejores estimadores del error esperado de predicción. Para ello hace falta tener un procedimiento de estimación de dicho error esperado, como se propone a continuación. 36 CAPÍTULO 3. Support Vector Regression 3.5.1. Validación Cruzada Probablemente el método más simple y más ampliamente utilizado para la estimación del error de predicción, sea la validación cruzada. Este método estima directamente la generalización del error cuando se aplica el método, de función de predicción ˆ f ( x ), a una muestra de test. Se denota por Err = E [ L ( y, ˆ f ( x ))], siendo L(y, ˆ f(x)) una función de perdida adecuada. Idealmente, si tuviésemos un conjunto suficientemente grande, podríamos seleccionar un conjunto de validación para evaluar nuestro modelo de predicción. Puesto que a menudo el conjunto no es lo suficientemente grande como nosotros requeriríamos esto no suele ser posible. En casos en los que no dispongamos un conjunto amplio, nosotros podemos usar la validación cruzada para evaluar nuestro predictor. Cuando aplicamos la validación cruzada con k -pliegues, subdividimos nuestro conjunto en k partes de aproximadamente el mismo tamaño. Utilizamos k− 1partes para ajustar el modelo y la restante para validarlo. Realizamos este mismo procedimiento k veces, seleccionando cada vez un conjunto de validación distinto. En cada caso nosotros calculamos el error de predicción del modelo y por último promediamos estos errores. En otras palabras, sea κ : { 1 , ..., n} 7−→ { 1 , ..., k} una función que indica la partición a la que se le asigna i -ésima observación por la aleatorización. Para j∈ { 1 , ..., n} denotamos por ˆ f−κ(j) ( x )la función ajustada con el conjunto de datos a excepción del subconjunto κ ( j )-ésimo. La correspondiente estimación de error mediante la validación cruzada con k-pliegues es: CV (ˆ f) = 1 n n X i=1 L(yi,ˆ f−κ(i)(xi)). Las elecciones más comunes son k = 5 y k = 10. El caso k = n es conocido como leave - one - out . En este caso κ ( i ) = i , y la función de predicción para el caso i -ésimo se realiza con la totalidad del conjunto excepto el i -ésimo elemento. Esta validación cruzada es muy costosa, y por eso se reserva para los casos donde el conjunto de datos es muy pequeño. Hay algunos casos especiales en los que es beneficioso utilizar otros k -pliegues. Un estudio sobre la elección de estos k lo podemos encontrar en [4]. En nuestro caso, como utilizamos la validación cruzada para la estimación de los parámetros de nuestro predictor, la situación que nos queda es la siguiente. Sea θnuestro conjunto de parámetros a estimar. La estimación de error sería: CV (ˆ f, θ) = 1 n n X i=1 L(yi,ˆ f−κ(i)(xi, θ)). La función CV ( θ )nos proporciona una estimación de la curva de error test, y buscamos el parámetro ˆ θ que minimice dicha función. Nosotros elegiríamos como función predictora aquella formada con estos parámetros, es decir, f ( x, ˆ θ )ajustada con todo el conjunto de datos. Support Vector Regression: Propiedades y Aplicaciones 37 Capítulo 4 Support Vector Regression en R En la actualidad existe una gran diversidad de repositorios web y de paquetes software de libre distribución dedicados a la implementación de SVM, SVR y muchas de sus variantes. Para apreciar el funcionamiento práctico de estos algoritmos, practicaremos la regresión mediante SVR a ciertos datos muestrales obtenidos del repositorio UCI Machine learning repository [8]. Para abordar este tipo de regresión con este conjunto de datos, utilizaremos el programa R , en el cual podremos realizar tanto la regresión, como un estudio descriptivo. Usaremos la librería e 1071 [ 10 ], paquete software pensado para resolver problemas de clasificación y regresión mediante máquinas de vectores soporte, la cual podemos instalar fácilmente en R. R implementa diversas formulaciones SVR con la posibilidad de usar tres importantes tipos de kernels y la posibilidad de usar tecnicas como la validación cruzada descrita en la sección 3.5.1, para la elección de los parámetros del kernel. 4.1. Funciones Específicas SVR Comenzaremos describiendo las funciones que utilizaremos de la librería e1071 dedicadas a la resolución de problemas de regresión mediante máquinas de vectores soporte. • El comando llamado svm es el utilizado por esta librería a la hora de obtener tanto la clasificación como la regresión mediante máquinas de vectores soporte. ## S3 method for class 'formula' svm(formula, data = NULL, ..., subset, na.action = na.omit, scale = TRUE) ## Default S3 method: svm(x, y = NULL, scale = TRUE, type = NULL, kernel = "radial", degree = 3, gamma = if (is.vector(x)) 1 else 1 / ncol(x), coef0 = 0, cost = 1, nu = 0.5, 38 CAPÍTULO 4. Support Vector Regression en R class.weights = NULL, cachesize = 40, tolerance = 0.001, epsilon = 0.1, shrinking = TRUE, cross = 0, probability = FALSE, fitted = TRUE, ..., subset, na.action = na.omit) Se proporciona una interfaz de fórmula. Si la variable de predictor ” y ”incluye factores, la interfaz de fórmula debe ser utilizada para obtener una matriz de modelo correcto. El modelo de probabilidad para la clasificación se ajusta a una distribución logística utilizando máxima verosimilitud para los valores de decisión de todos los clasificadores binarios, y calcula las probabilidades de clase a posteriori para el problema de optimización cuadrática utilizando multiclase. El modelo de regresión probabilística asume errores de Laplace para los predictores, y estima el parámetro de escala utilizando máxima verosimilitud. Un objeto de la clase svm contiene el modelo ajustado, incluyendo: .SV →los vectores de soporte resultantes. .index → índice de los vectores de soporte resultantes en la matriz de datos. Este índice se refiere a los datos preprocesados. .coefs → los correspondientes coeficientes veces las etiquetas de formación. .rho →el termino independiente. .sigma → en caso de un modelo de regresión probabilístico, el parámetro de escala de la distribución de Laplace mediante estimación de máxima verosimilitud. .probA, probB → vectores numéricos de longitud k(k−1) 2 , siendo k , el número de las clases que contienen los parámetros de la distribución logística. • El comando llamado tune.svm se utiliza para obtener una estimación del error cometido por los modelos según los parámetros específicos. tune.svm (method, x, y = NULL, data = list(), validation.x = NULL, validation.y = NULL, ranges = NULL, predict.func = predict, tunecontrol = tune.control(), ...) best.tune(...) Se le proporciona el conjunto de datos junto al kernel con los posibles valores de los coeficientes específicos del modelo. Esta función calcula el error como medida de rendimiento, el error de clasificación se utiliza para la clasificación, y el error cuadrático medio para la regresión. Este error es calculado mediante validación cruzada. La validación cruzada se realiza aleatoriamente a partir del conjunto de datos, una vez construidas estos pliegues, se mantienen constantes durante el proceso de formación. Un objeto de la clase tune.svm contiene: 4.2. Ejemplo Práctico 39 .best.parameters →tabla de datos 1×k,knúmero de parámetros. .best.performance →mejor rendimiento del modelo. .performances → si se solicita, un conjunto de datos con todas las combinaciones de parámetros, junto con los resultados de rendimiento correspondientes. .train.ind → lista de índices de vectores utilizados para la división en conjuntos de entrenamiento y validación. .best.model → modelo completo asociado al mejor conjunto de parámetros. Además de éstas, usaremos otras funciones que no son específicas de máquinas de vectores soporte, como predict , utilizada por R para predecir la variable respuesta a partir de un modelo. Esta función se utiliza tanto para la clasificación como para la regresión en nuestro caso. Para la representación de las gráficas de este trabajo utilizaremos un mapa de calor obtenido con la función filled.contour. 4.2. Ejemplo Práctico Trabajaremos con el fichero de título Boston Housing Data obtenido de [ 8 ]. Estos estudios tienen su origen en la biblioteca StatLib mantenido en la Universidad Carnegie Mellon. Su fecha de creación data del 7 de julio de 1993 y los creadores son David Harrison y Daniel L Rubinfeld [ 3 ]. Son unos estudios realizados en Boston que recoge una serie de características de las zonas de Boston y el valor medio de la vivienda en cada zona. Para ello dividieron el territorio de Boston en 506 zonas. Se recogían un total de catorce variables de cada zona: 1. CRIM Tasa de criminalidad. 2. ZN Proporción de zona residencial por cada 25000 pies cuadrados. 3. INDUS Proporción de hectáreas de comercio al pormenor. 4. CHAS (1) Si delimita con el río Charles, (0) si no. 5. NOX Concentración de óxido nitroso. 6. RM Promedio de habitaciones por vivienda. 7. AGE Proporción de viviendas ocupadas por sus propietarios construidas antes de 1940. 8. DIS Distancia ponderadas a 5 centros de empleo en Boston. 9. RAD Índice de accesibilidad a carreteras centrales. 10. TAX Valor total de la tasa de impuestos a la propiedad por 10000 $. 40 CAPÍTULO 4. Support Vector Regression en R 11. PTRATIO Relación alumno-profesor en la zona. 12. B1000( Bk − 0 , 63) 2 donde Bk es la proporción de personas de color. 13. LSTAT Estado inferior de la población. 14. MEDV Valor medio en miles de dolares, de las viviendas ocupadas. Del conjunto de variables que disponemos trece son variables continuas entre ellas la variable a la que le queremos aplicar la regresión, la variable número 14, MEDV. Y una variable binaria que toma los valores 0 o 1. En este conjunto de datos no disponemos de valores perdidos. 4.2.1. Análisis Descriptivo Comenzaremos con un análisis descriptivo de las variables para tener una idea de como son las variables con la cual vamos a trabajar. datosHousing<- read.table(file="housing.data.txt",header=FALSE) colnames(datosHousing)<- c("CRIM","ZN","IND","CHAS","NOX","RM", "AGE","DIS","RAD","TAX","PTRA","B", "LST","MEDV") summary(datosHousing) ## CRIM ZN IND CHAS ## Min. : 0.00632 Min. : 0.00 Min. : 0.46 Min. :0.000 ## 1st Qu.: 0.08204 1st Qu.: 0.00 1st Qu.: 5.19 1st Qu.:0.000 ## Median : 0.25651 Median : 0.00 Median : 9.69 Median :0.000 ## Mean : 3.61352 Mean : 11.36 Mean :11.14 Mean :0.069 ## 3rd Qu.: 3.67708 3rd Qu.: 12.50 3rd Qu.:18.10 3rd Qu.:0.000 ## Max. :88.97620 Max. :100.00 Max. :27.74 Max. :1.000 ## NOX RM AGE DIS ## Min. :0.3850 Min. :3.561 Min. : 2.90 Min. : 1.130 ## 1st Qu.:0.4490 1st Qu.:5.886 1st Qu.: 45.02 1st Qu.: 2.100 ## Median :0.5380 Median :6.208 Median : 77.50 Median : 3.207 ## Mean :0.5547 Mean :6.285 Mean : 68.57 Mean : 3.795 ## 3rd Qu.:0.6240 3rd Qu.:6.623 3rd Qu.: 94.08 3rd Qu.: 5.188 ## Max. :0.8710 Max. :8.780 Max. :100.00 Max. :12.127 ## RAD TAX PTRA B ## Min. : 1.000 Min. :187.0 Min. :12.60 Min. : 0.32 ## 1st Qu.: 4.000 1st Qu.:279.0 1st Qu.:17.40 1st Qu.:375.38 ## Median : 5.000 Median :330.0 Median :19.05 Median :391.44 ## Mean : 9.549 Mean :408.2 Mean :18.46 Mean :356.67 ## 3rd Qu.:24.000 3rd Qu.:666.0 3rd Qu.:20.20 3rd Qu.:396.23 ## Max. :24.000 Max. :711.0 Max. :22.00 Max. :396.90 4.2. Análisis descriptivo 41 ## LST MEDV ## Min. : 1.73 Min. : 5.00 ## 1st Qu.: 6.95 1st Qu.:17.02 ## Median :11.36 Median :21.20 ## Mean :12.65 Mean :22.53 ## 3rd Qu.:16.95 3rd Qu.:25.00 ## Max. :37.97 Max. :50.00 La matriz de correlaciones: ## CRIM ZN IND CHAS NOX RM AGE DIS RAD TAX PTRA B LST MEDV ## CRIM 1.0 -0.2 0.4 -0.1 0.4 -0.2 0.4 -0.4 0.6 0.6 0.3 -0.4 0.5 -0.4 ## ZN -0.2 1.0 -0.5 0.0 -0.5 0.3 -0.6 0.7 -0.3 -0.3 -0.4 0.2 -0.4 0.4 ## IND 0.4 -0.5 1.0 0.1 0.8 -0.4 0.6 -0.7 0.6 0.7 0.4 -0.4 0.6 -0.5 ## CHAS -0.1 0.0 0.1 1.0 0.1 0.1 0.1 -0.1 0.0 0.0 -0.1 0.0 -0.1 0.2 ## NOX 0.4 -0.5 0.8 0.1 1.0 -0.3 0.7 -0.8 0.6 0.7 0.2 -0.4 0.6 -0.4 ## RM -0.2 0.3 -0.4 0.1 -0.3 1.0 -0.2 0.2 -0.2 -0.3 -0.4 0.1 -0.6 0.7 ## AGE 0.4 -0.6 0.6 0.1 0.7 -0.2 1.0 -0.7 0.5 0.5 0.3 -0.3 0.6 -0.4 ## DIS -0.4 0.7 -0.7 -0.1 -0.8 0.2 -0.7 1.0 -0.5 -0.5 -0.2 0.3 -0.5 0.2 ## RAD 0.6 -0.3 0.6 0.0 0.6 -0.2 0.5 -0.5 1.0 0.9 0.5 -0.4 0.5 -0.4 ## TAX 0.6 -0.3 0.7 0.0 0.7 -0.3 0.5 -0.5 0.9 1.0 0.5 -0.4 0.5 -0.5 ## PTRA 0.3 -0.4 0.4 -0.1 0.2 -0.4 0.3 -0.2 0.5 0.5 1.0 -0.2 0.4 -0.5 ## B -0.4 0.2 -0.4 0.0 -0.4 0.1 -0.3 0.3 -0.4 -0.4 -0.2 1.0 -0.4 0.3 ## LST 0.5 -0.4 0.6 -0.1 0.6 -0.6 0.6 -0.5 0.5 0.5 0.4 -0.4 1.0 -0.7 ## MEDV -0.4 0.4 -0.5 0.2 -0.4 0.7 -0.4 0.2 -0.4 -0.5 -0.5 0.3 -0.7 1.0 Podemos representar un histograma de la variable MEDV: Histogram of MEDV MEDV Frequency 10 20 30 40 50 0 50 100 150 48 CAPÍTULO 4. Support Vector Regression en R ## ## - sampling method: 10-fold cross validation ## ## - best parameters: ## degree cost epsilon ## 3 1 0.1 ## ## - best performance: 0.2627534 ## ## - Detailed performance results: ## degree cost epsilon error dispersion ## 1 3 0.1 0.1 0.3953308 0.1817405 ## 2 4 0.1 0.1 0.4495784 0.2071792 ## 3 5 0.1 0.1 0.6334360 0.6379059 ## 4 3 1.0 0.1 0.2627534 0.1754837 ## 5 4 1.0 0.1 0.4032249 0.4316404 ## 6 5 1.0 0.1 2.3424848 6.4478723 ## 7 3 10.0 0.1 0.5558749 1.0074161 ## 8 4 10.0 0.1 2.3552565 6.3333634 ## 9 5 10.0 0.1 3.9880200 10.2653011 ## 10 3 100.0 0.1 1.7685265 3.9239705 ## 11 4 100.0 0.1 6.1934556 16.3795169 ## 12 5 100.0 0.1 6.2235547 12.2730292 ## 13 3 0.1 1.0 0.5642490 0.1301769 ## 14 4 0.1 1.0 0.6259962 0.1446387 ## 15 5 0.1 1.0 0.7316682 0.4238583 ## 16 3 1.0 1.0 0.4450963 0.1004615 ## 17 4 1.0 1.0 1.3929350 2.7675810 ## 18 5 1.0 1.0 1.4406510 2.9467190 ## 19 3 10.0 1.0 0.3915978 0.1273010 ## 20 4 10.0 1.0 3.1440790 8.5417739 ## 21 5 10.0 1.0 1.9549975 4.4892543 ## 22 3 100.0 1.0 0.3744535 0.1233382 ## 23 4 100.0 1.0 2.7446738 7.2553780 ## 24 5 100.0 1.0 1.9533168 4.6003029 En este caso el mejor modelo seria: mejorajuste1_Norm=tuneado1_Norm$best.model print(mejorajuste1_Norm) ## ## Call: ## best.svm(x = MEDV_N ~ ., data = datosHousingNorm, degree = (3:5), ## cost = 10^(-1:2), epsilon = 10^(-1:0), kernel = "polynomial", ## scale = TRUE) 4.3. Regresión mediante SVR 49 ## ## ## Parameters: ## SVM-Type: eps-regression ## SVM-Kernel: polynomial ## cost: 1 ## degree: 3 ## gamma: 0.07692308 ## coef.0: 0 ## epsilon: 0.1 ## ## ## Number of Support Vectors: 372 Generando así un mapa de calor y un error medio: d1 <- interp(puntCP1Norm,puntCP2Norm,ypredNorm1,duplicate="mean") filled.contour(d1, main = "Mapa de Calor Kernel Polinomial", xlab="C.P.1",ylab="C.P.2") -2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Polinomial C.P.1 C.P.2 50 CAPÍTULO 4. Support Vector Regression en R errorModeloL1Medio(mejorajuste1_Norm,datosHousingNorm,MEDV_N) ## 0.1988483 Se tiene también que el error cuadrático medio cometido con esta función de predicción es: ## 0.1098226 4.3.3. Kernel Radial Ahora realizamos el estudio con el kernel radial. Tenemos en este caso un nuevo parámetro γjunto al coste y epsilon. tuneado2_Norm <- tune.svm(MEDV_N~.,data=datosHousingNorm, kernel= "radial",gamma = 10^(-4:-1), cost = 10^(-1:1), epsilon = 10^(-1:0), scale=TRUE) summary(tuneado2_Norm) ## ## Parameter tuning of 'svm': ## ## - sampling method: 10-fold cross validation ## ## - best parameters: ## gamma cost epsilon ## 0.1 10 0.1 ## ## - best performance: 0.1259211 ## ## - Detailed performance results: ## gamma cost epsilon error dispersion ## 1 1e-04 0.1 0.1 0.9798190 0.31440489 ## 2 1e-03 0.1 0.1 0.7591207 0.30945188 ## 3 1e-02 0.1 0.1 0.4154568 0.23749867 ## 4 1e-01 0.1 0.1 0.3686023 0.24588329 ## 5 1e-04 1.0 0.1 0.7537387 0.30924793 ## 6 1e-03 1.0 0.1 0.3964701 0.22075054 ## 7 1e-02 1.0 0.1 0.2219599 0.15499423 ## 8 1e-01 1.0 0.1 0.1768023 0.13373057 ## 9 1e-04 10.0 0.1 0.3967259 0.22037508 ## 10 1e-03 10.0 0.1 0.2799216 0.17562604 ## 11 1e-02 10.0 0.1 0.1692663 0.10806448 ## 12 1e-01 10.0 0.1 0.1259211 0.07722572 ## 13 1e-04 0.1 1.0 0.9798601 0.26939380 4.3. Regresión mediante SVR 51 ## 14 1e-03 0.1 1.0 0.8495537 0.25519538 ## 15 1e-02 0.1 1.0 0.6377322 0.20342812 ## 16 1e-01 0.1 1.0 0.6232611 0.22389760 ## 17 1e-04 1.0 1.0 0.8451833 0.25397370 ## 18 1e-03 1.0 1.0 0.6000316 0.17935163 ## 19 1e-02 1.0 1.0 0.4555450 0.12728025 ## 20 1e-01 1.0 1.0 0.4036266 0.11584329 ## 21 1e-04 10.0 1.0 0.5969894 0.17833236 ## 22 1e-03 10.0 1.0 0.4684074 0.12180408 ## 23 1e-02 10.0 1.0 0.3233276 0.06902243 ## 24 1e-01 10.0 1.0 0.3581760 0.07587053 En este caso el mejor modelo seria: mejorajuste2_Norm=tuneado2_Norm$best.model print(mejorajuste2_Norm) ## ## Call: ## best.svm(x = MEDV_N ~ ., data = datosHousingNorm, gamma = 10^(-4:-1), ## cost = 10^(-1:1), epsilon = 10^(-1:0), kernel = "radial", ## scale = TRUE) ## ## ## Parameters: ## SVM-Type: eps-regression ## SVM-Kernel: radial ## cost: 10 ## gamma: 0.1 ## epsilon: 0.1 ## ## ## Number of Support Vectors: 340 Generando asi un mapa de calor y un error medio: d2 <- interp(puntCP1Norm,puntCP2Norm,ypredNorm2,duplicate="mean") filled.contour(d2, main = "Mapa de Calor Kernel Radial", xlab="C.P.1",ylab="C.P.2") 52 CAPÍTULO 4. Support Vector Regression en R -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Radial C.P.1 C.P.2 errorModeloL1Medio(mejorajuste2_Norm,datosHousingNorm,MEDV_N) ## 0.08626746 Con error cuadrático medio: ## 0.03465188 Apreciamos que ambos errores de aproximación son cercanos a 0. También se observa que el mapa de calor aproximado es similar al de los datos de entrada. 4.3.4. Kernel Sigmoidal Por último realizamos el estudio con el kernel Sigmoidal. Tenemos el parámetro γyRjunto al coste y epsilon. tuneado3_Norm <- tune.svm(MEDV_N~.,data=datosHousingNorm, kernel= "sigmoid",coef0 = c(0,5), gamma = 10^(-3:0), cost = 10^(-1:1),scale=TRUE) summary(tuneado3_Norm) 4.3. Regresión mediante SVR 53 ## ## Parameter tuning of 'svm': ## ## - sampling method: 10-fold cross validation ## ## - best parameters: ## gamma coef0 cost ## 0.001 0 10 ## ## - best performance: 0.3112602 ## ## - Detailed performance results: ## gamma coef0 cost error dispersion ## 1 0.001 0 0.1 8.521782e-01 2.324651e-01 ## 2 0.010 0 0.1 4.936465e-01 1.539728e-01 ## 3 0.100 0 0.1 7.704087e-01 2.336415e-01 ## 4 1.000 0 0.1 1.923297e+01 2.951157e+00 ## 5 0.001 5 0.1 1.021740e+00 2.555519e-01 ## 6 0.010 5 0.1 1.021382e+00 2.554847e-01 ## 7 0.100 5 0.1 1.013524e+00 2.538353e-01 ## 8 1.000 5 0.1 3.557984e+01 2.988301e+00 ## 9 0.001 0 1.0 4.931048e-01 1.539584e-01 ## 10 0.010 0 1.0 3.126837e-01 1.106426e-01 ## 11 0.100 0 1.0 5.270987e+01 1.417341e+01 ## 12 1.000 0 1.0 1.882462e+03 3.321538e+02 ## 13 0.001 5 1.0 1.021378e+00 2.554857e-01 ## 14 0.010 5 1.0 1.017759e+00 2.547409e-01 ## 15 0.100 5 1.0 9.471725e-01 2.421115e-01 ## 16 1.000 5 1.0 3.583226e+03 2.627001e+02 ## 17 0.001 0 10.0 3.112602e-01 1.096166e-01 ## 18 0.010 0 10.0 3.270285e-01 1.165002e-01 ## 19 0.100 0 10.0 5.234178e+03 1.425616e+03 ## 20 1.000 0 10.0 1.884112e+05 3.431191e+04 ## 21 0.001 5 10.0 1.017718e+00 2.547481e-01 ## 22 0.010 5 10.0 9.841693e-01 2.501437e-01 ## 23 0.100 5 10.0 7.419826e-01 2.053973e-01 ## 24 1.000 5 10.0 3.566550e+05 2.853332e+04 En este caso el mejor modelo seria: mejorajuste3_Norm=tuneado3_Norm$best.model print(mejorajuste3_Norm) ## ## Call: ## best.svm(x = MEDV_N ~ ., data = datosHousingNorm, gamma = 10^(-3:0), 54 CAPÍTULO 4. Support Vector Regression en R ## coef0 = c(0, 5), cost = 10^(-1:1), kernel = "sigmoid", ## scale = TRUE) ## ## ## Parameters: ## SVM-Type: eps-regression ## SVM-Kernel: sigmoid ## cost: 10 ## gamma: 0.001 ## coef.0: 0 ## epsilon: 0.1 ## ## ## Number of Support Vectors: 390 Generando asi un mapa de calor y un error medio: d3 <- interp(puntCP1Norm,puntCP2Norm,ypredNorm3,duplicate="mean") filled.contour(d3, main = "Mapa de Calor Kernel Sigmoidal", xlab="C.P.1",ylab="C.P.2") -2 -1 0 1 2 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Sigmoidal C.P.1 C.P.2 4.3. Regresión mediante SVR 55 errorModeloL1Medio(mejorajuste3_Norm,datosHousingNorm,MEDV_N) ## 0.3298638 Y error cuadrático medio asociado: ## 0.3033291 Finalizaremos el estudio con un resumen con la comparación de los modelos utilizados. 56 CAPÍTULO 4. Support Vector Regression en R 4.4. Comparación de Modelos Tras realizar el mismo estudio con cada uno de los kernels predefinidos en la librería e1071 podemos comparar su eficacia frente a nuestro conjunto de datos. Original Aproximación -2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor Normalizado C.P.1 C.P.2 -3 -2 -1 0 1 2 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Lineal C.P.1 C.P.2 E.L1 = 0.3235691, ECM=0.2869733 -2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor Normalizado C.P.1 C.P.2 -2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Polinomial C.P.1 C.P.2 E.L1 = 0.1988483, ECM=0.1098226 -2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor Normalizado C.P.1 C.P.2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Radial C.P.1 C.P.2 E.L1 = 0.0862674, ECM=0.0346518 -2 -1 0 1 2 3 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor Normalizado C.P.1 C.P.2 -2 -1 0 1 2 -4 -2 0 2 4 6 -3 -2 -1 0 1 2 3 Mapa de Calor kernel Sigmoidal C.P.1 C.P.2 E.L1 = 0.3298638, ECM=0.3033291 4.4. Comparación de Modelos 57 Apreciamos que ambos errores alcanzan su mínimo para el caso del kernel Radial. Por lo que cabe esperar que nuestra función de aproximación sera mejor si utilizamos el kernel Radial con los mejores parámetros estimados mediante validación cruzada.