Estadística multivariante: algunos problemas con la alta dimensionalidad

TRABAJO FIN DE GRADO Grado en Matemáticas Estadística multivariante: algunos problemas con la alta dimensionalidad. Andrés Romero Pavón Dirigido por: Dr. Juan Manuel Muñoz Pichardo FACULTAD DE MATEMÁTICAS Departamento de Estadística e Investigación Operativa Sevilla, Julio 2024 A mi familia, en especial a mis padres y mi hermano por haber estado en los momentos más duros. A mis amigos, por siempre acompañarme. Y, por último a D.Juan Manuel Muñoz, por haberme prestado su ayuda siempre que lo he necesitado. Muchísimas gracias. Resumen El presente trabajo aborda la problemática asociada a la alta dimensionalidad en el ámbito de la estadística multivariante. La investigación se centra en describir y presentar algunas soluciones y estrategias enfocadas a resolver este problema, bien desde el enfoque de la reducción de la dimensión, bien directamente para abordar técnicas de clasificación o predicción, planteando desafíos sustanciales en la inferencia estadística y el análisis de datos. Un problema muy presente en la actualidad al que tiene que enfrentarse la estadística es el análisis de conjuntos de datos de altas dimensiones, especialmente cuando “p≫N”, es decir, cuando la dimensión o número de variables es superior al tamaño muestral o número de casos. La alta dimensionalidad presenta obstáculos críticos en la identificación de patrones y la interpretación de relaciones entre variables. En este contexto, se exploran en profundidad las limitaciones de los métodos clásicos de estadística multivariante y se proponen enfoques innovadores para hacer frente a la falta de observaciones en relación con la abundancia de variables. Se examinan técnicas de reducción de dimensionalidad, como el análisis de componentes principales, así como métodos más avanzados como la regularización y la selección de características. Además, se aborda la problemática de la selección de modelos y la mejora de la eficiencia computacional en entornos de alta dimensionalidad. Los resultados de este trabajo proporcionan ideas valiosas para investigadores y profesionales en campos donde la alta dimensionalidad es una característica inherente, como la genómica, la minería de datos y la inteligencia artificial. Se espera que estas contribuciones fomenten el avance en la comprensión y el tratamiento de la alta dimensionalidad en el análisis estadístico multivariante. Abstract This paper addresses the problems associated with high dimensionality in the field of multivariate statistics. Entitled “Multivariate statistics: Problems with high dimensionality’’, the research focuses on situations where the number of explanatory variables psignificantly exceeds the number of observations N, posing substantial challenges in statistical inference and data analysis. High dimensionality presents critical obstacles in identifying patterns and interpreting relationships between variables. In this context, the limitations of classical multivariate statistical methods are explored in depth and innovative approaches are proposed to cope with the lack of observations in relation to the abundance of variables. Dimensionality reduction techniques, such as principal component analysis, as well as more advanced methods such as regularisation and feature selection are discussed. In addition, the issue of model selection and improving computational efficiency in high dimensionality environments is addressed. The results of this paper provide valuable insights for researchers in fields where high dimensionality is an inherent feature, such as genomics, data mining and artificial intelligence. It is hoped that these contributions will further advance the understanding and treatment of high dimensionality in multivariate statistical analysis. Índice general 1. Introducción - presentación 1 1.1. Nociones y técnicas estadísticas previas. . . . . . . . . . . . . . . . . . . . 2 1.1.1. Análisis Discriminante Lineal. . . . . . . . . . . . . . . . . . . . . . 2 1.1.2. Regresión Logística. . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 1.1.3. Regularización en regresión. . . . . . . . . . . . . . . . . . . . . . . 5 1.1.3.1. Regularización Ridge. . . . . . . . . . . . . . . . . . . . . 5 1.1.3.2. Regularización LASSO. . . . . . . . . . . . . . . . . . . . 6 1.1.3.3. Regularización Elastic Net. . . . . . . . . . . . . . . . . . 7 1.2. Problemas con la alta dimensionalidad y posibles formas de abordarlo. . . 9 1.3. Introducción a las técnicas básicas para afrontar el problema p≫N. . . . 10 1.3.1. Análisis de Componentes Principales. . . . . . . . . . . . . . . . . . 10 1.3.2. Análisis Discriminante Lineal Diagonal y Nearest Shrunken Centroids. 10 1.3.3. Clasificadores lineales con regularización cuadrática y L1. . . . . . . 11 2. Análisis teórico de algunas técnicas para combatir el problema de la alta dimensionalidad. 13 2.1. Reducción de la dimensión: Análisis de Componentes Principales. . . . . . 13 2.1.1. Planteamiento del problema . . . . . . . . . . . . . . . . . . . . . . 15 2.1.1.1. Descomposición en valores singulares SVD . . . . . . . . . 15 2.1.2. Análisis de Componentes Principales mediante descomposición SVD. 16 2.1.3. Determinación del número de componentes a retener . . . . . . . . 19 2.1.4. Interpretación de las componentes . . . . . . . . . . . . . . . . . . . 20 2.2. Técnicas de clasificación . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.2.1. Análisis Discriminante Lineal Diagonal (DLDA). . . . . . . . . . . . 21 2.2.2. Clasificadores Lineales con Regularización Cuadrática. . . . . . . . 25 2.2.2.1. Análisis Discriminante Regularizado. . . . . . . . . . . . . 25 2.2.2.2. Regresión Logística con Regularización Cuadrática. . . . . 27 2.2.2.3. Clasificador de vectores de soporte. . . . . . . . . . . . . . 29 i 2.2.2.4. Selección de características. . . . . . . . . . . . . . . . . . 32 2.2.2.5. Estrategias computacionales cuando p≫N......... 32 2.2.3. Clasificadores Lineales con Regularización L1............. 34 2.2.3.1. Regularización LASSO y Elastic Net . . . . . . . . . . . . 34 3. IMPLEMENTACIÓN EN R. 39 3.1. Introducción................................... 39 3.2. Análisis de Componentes Principales. . . . . . . . . . . . . . . . . . . . . . 40 3.2.1. Cálculos alternativos. . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3.2.1.1. Selección del número de componentes. . . . . . . . . . . . 45 3.3. Técnicas de clasificación . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 3.3.1. Análisis Discriminante Lineal Diagonal . . . . . . . . . . . . . . . . 53 3.4. Clasificadores Lineales con Regularización Cuadrática . . . . . . . . . . . . 55 3.4.1. Análisis Discriminante Regularizado . . . . . . . . . . . . . . . . . . 55 3.4.2. Regresión Logística con Regularización Cuadrática . . . . . . . . . 58 3.5. Clasificadores Lineales con Regularización L1................. 60 3.5.1. Regularización LASSO y Elastic Net . . . . . . . . . . . . . . . . . 60 Bibliografía 63 ii Capítulo 1 Introducción - presentación Un problema muy presente en la actualidad al que tiene que enfrentarse la Estadística es el análisis de conjuntos de datos de altas dimensiones, especialmente cuando p≫N, es decir, cuando la dimensión o el número de variables es superior al tamaño muestral o número de observaciones. El principal objetivo de este proyecto será el de abordar la problemática existente con la alta dimensionalidad cuando se tiene un conjunto de datos con Nobservaciones ypvariables, con p≫N. Concretamente, se llevará a cabo el análisis teórico de este problema, así como la descripción y la propuesta de soluciones y estrategias enfocadas a su resolución, bien con el objetivo de reducir la dimensionalidad para obtener una mayor interpretabilidad de los datos, o bien para llevar a cabo técnicas de clasificación o predicción. Estos problemas han adquirido una importancia creciente, especialmente en el área de la genómica y otras áreas de la biología computacional. Veremos que, en este contexto, la alta varianza y el “sobreajuste” (overfitting) son problemas importantes a tratar. Supongamos, por ejemplo, que queremos realizar un modelo de regresión lineal múltiple sobre un conjunto de datos donde el número de variables es mucho mayor que el número de observaciones. Sabemos que el número de parámetros βja estimar es excesivamente alto, por lo que puede resultar bastante costoso computacionalmente. Además, a priori, veremos también que se incluyen en el modelo variables que no aportan apenas información a la hora de explicar la variable respuesta Y. En estos casos, estaríamos considerando un modelo que se sobreajusta a los datos, disminuyendo así el sesgo y con una varianza excesiva (overfitting). Por otra parte, si el modelo tiene un número insuficiente de variables, es decir, si las variables seleccionadas no son capaces de explicar lo suficientemente bien el modelo, la varianza será menor y se tendrá un mayor sesgo (underfitting), siendo su capacidad de predicción mucho menor. El objetivo es encontrar un equilibrio, reduciendo el número de variables hasta seleccionar aquellas que, aun sin perder apenas información, mantengan un alto porcentaje de la varianza original y cuya ausencia no tenga un efecto significativo sobre la capacidad predictiva del modelo. Para resolver esta cuestión, se estudiarán técnicas de reducción de dimensionalidad, como pueden ser el Análisis Discriminante Lineal (LDA) o el Análisis de Componentes Principales (PCA), entre otros . Se estudiará sobre todo este último para el caso en el que el número de variables pes mucho mayor que el número de observaciones N. 1 1.1. NOCIONES Y TÉCNICAS ESTADÍSTICAS PREVIAS. siendoλ1yλ2dos parámetros de ajuste que controlan la cantidad de regularización. La parte “L1” permite la “selección” de variables, mientras que la parte “L2” estabiliza la regularización “L1” y elimina las deficiencias de la regresión LASSO. Esto es, elimina la limitación sobre el número de variables seleccionadas e incentiva el efecto agrupación. Sobre este método, cabe destacar la capacidad de seleccionar variables agrupadas. Recordemos que el procedimiento LASSO no posee esta propiedad. Si b βens(λ1, λ2)=[b β1··· b βp]Tes una solución del criterio elastic net simple con parámetros de regularización (λ1, λ2), y para dos componentes βiyβjse verifica b βi(λ1, λ2)×b βj(λ1, λ2)>0(idéntico signo), entonces se tiene que D[λ1, λ2](i, j) = 1 ||y||1|b βi(λ1, λ2)−b βj(λ1, λ2)| ≤ 1 λ2q2(1 −ρij) siendo ρij =ZT iZjel coeficiente de correlación muestral entre los regresores iyj. Así, si se tiene que ambos regresores están altamente correlados (hay variables agrupadas), es decir, ρij ≃1, entonces D[λ1, λ2](i, j)≃0⇒b βi(λ1, λ2)≃b βj(λ1, λ2), es decir, concluimos que los coeficientes son muy similares. Estimador Elastic Net Sea Zuna matriz de datos e yun vector. Vamos a definir en primer lugar el conjunto de datos aumentado (Z∗, y∗). Z∗= (1 + λ2)−1/2"Z λ1/2 2Ip#;y∗="y 0p# Dada una matriz Zestandarizada y un vector de salida ycentrado, dados los parámetros (λ1, λ2)y el conjunto de datos aumentado (Z∗, y∗), el estimador elastic net simple se obtiene resolviendo un problema tipo LASSO c β∗: m´ın β∗((y∗−Z∗β∗)T(y∗−Z∗β∗) + λ1 √1 + λ2||β∗||1) Sin embargo, se puede comprobar empíricamente que la solución simple no proporciona muy buenos resultados. Zou y Hastie (2005) consideran el estimador elastic net corregido dado por b βenet =q1 + λ2c β∗. Como podemos observar, se trata de un coeficiente re-escalado del coeficiente elastic net simple. Por tanto, preserva la propiedad de selección de variables. b βens =1 √1 + λ2c β∗⇒b βenet = (1 + λ2)b βens 8CAPÍTULO 1. INTRODUCCIÓN - PRESENTACIÓN 1.2. PROBLEMAS CON LA ALTA DIMENSIONALIDAD Y POSIBLES FORMAS DE ABORDARLO. Otro resultado de interés es: b βenet : m´ın ββT1 1 + λ2[ZTZ+λ2Ip]β−2yTZβ+λ1||β||1 mientras que, como se vió anteriormente, b βLASSO : m´ın β(y−Zβ)T(y−Zβ) + λ1||β||1= m´ın ββT(ZTZ)β−2yTZβ+λ1||β||1 donde ZTZes la matriz de correlación muesstral de los regresores. Por tanto, b βenet se puede obtener como un problema LASSO. 1.2. Problemas con la alta dimensionalidad y posibles formas de abordarlo. Como se ha mencionado en la breve introducción, en este trabajo se estudiarán aquellos problemas existentes con la alta dimensionalidad que surgen, especialmente, cuando se tiene un conjunto de datos con un número de variables pmucho mayor que el número de observaciones N. Se estudiarán, desde el punto de vista del álgebra lineal y la estadística, métodos como la Descomposición en Valores Singulares (ó SV D, por sus siglas en inglés, "Singular Value Descomposition") para el Análisis de Componentes Principales, o el cálculo de clasificadores para un enfoque predictivo. Veamos un pequeño ejemplo donde nos enfrentemos ante este problema. Supongamos que tenemos un conjunto de datos expresado en una matriz Xque consta de Nobservaciones (filas) y pvariables (columnas). En ámbitos como la genómica, o la biología computacional, es muy usual encontrarnos con el problema p≫N, es decir, tenemos un conjunto de observaciones relativamente pequeño, en el que se evalúan un gran número de variables. A menudo disponemos de datos de salida que caracterizan a cada una de las observaciones, como pueden ser algún tipo de cáncer, especie biológica, tiempo de supervivencia u otros. Denominaremos Yia la variable de salida de la observación i-ésima. Un procedimiento estadístico habitual consiste en construir un modelo de predicción que utilice la matriz de datos de entrada Xpara predecir el vector de salida Y. Queremos construir un modelo de regresión lineal que sea capaz de predecir la variable de salida Yen función de los datos de entrada X. Supongamos que X∈RN×p, con p≫N. Y=Xβ+ε(1.1) con ε= [ε1,··· , εN]′, donde εies el error asociado a la i-ésima observación y además εi∼N(0, σ2). Para estimar los parámetros del modelo, se utiliza el método de mínimos cuadrados. En este caso, la función de error cuadrático está dada por S(β) = n X i=1 ε2 i=εTε= (Y−Xβ)T(Y−Xβ)(1.2) Los estimadores por mínimos cuadrados deben satisfacer la ecuación ∂S ∂β = 0 CAPÍTULO 1. INTRODUCCIÓN - PRESENTACIÓN 9 1.3. INTRODUCCIÓN A LAS TÉCNICAS BÁSICAS PARA AFRONTAR EL PROBLEMA p≫N. Se tiene entonces que la condición anterior se reduce a XTXˆ β=XTY Si la matriz (XTX)−1existe, el estimador ˆ βviene dado por ˆ β= (XTX)−1XTY. Pero aquí surge un problema, y es que el rango de la matriz XTXes como máximo N. Como la dimensión de esta matriz es p×p, y se tiene que p≫Npor hipótesis, obtenemos que la matriz XTXes singular, y la inversa no existe. La mayoría de procedimientos simples, como pueden ser la Regresión Lineal, el Análisis Discriminante Lineal o la Regresión Logística, entre otros, requieren la hipótesis de que p<N. Es por ello que habría que buscar otros procedimientos alternativos para abordar este problema, y es lo que se tratará a lo largo de este proyecto. 1.3. Introducción a las técnicas básicas para afrontar el problema p≫N. En esta subsección, se va a introducir un resumen teórico de las técnicas básicas que se van a estudiar para afrontar el problema de la alta dimensionalidad. En primer lugar, se verá el Análisis de Componentes Principales. 1.3.1. Análisis de Componentes Principales. El análisis de componentes principales (PCA, del inglés Principal Component Analysis) es una técnica utilizada para extraer información relevante de un conjunto de datos que “a priori” puede ser confuso. Además, se emplea para reducir la dimensión de un conjunto complejo de datos y así revelar posibles patrones o estructuras en ellos. Su objetivo principal es encontrar un conjunto más pequeño de nuevas variables no correlacionadas (denominadas componentes principales) que capturen la mayor parte de la variabilidad presente en los datos originales . La técnica del PCA puede ser llevada a cabo utilizando la descomposición en valores singulares SVD sobre la matriz de datos, que veremos en el próximo capítulo. 1.3.2. Análisis Discriminante Lineal Diagonal y Nearest Shrunken Centroids. Otras de las técnicas que se estudiarán serán el Análisis discriminante lineal diagonal yNearest Shrunken Centroids. El Análisis Discriminante Lineal (LDA, del inglés Linear Discriminant Analysis) es una técnica utilizada en estadística y aprendizaje automático para encontrar una combinación lineal de rasgos o variables predictoras que caracterizan dos o más clases de objetos. La combinación resultante puede ser utilizada como un clasificador lineal, o, más comúnmente, para la reducción de dimensiones antes de la posterior clasificación. En el caso del LDA diagonal, se asume que la matriz de covarianza de las variables predictoras dentro de cada clase es diagonal, lo que significa que no hay correlación entre ellas. 10 CAPÍTULO 1. INTRODUCCIÓN - PRESENTACIÓN 1.3. INTRODUCCIÓN A LAS TÉCNICAS BÁSICAS PARA AFRONTAR EL PROBLEMA p≫N. Esta simplificación reduce el número de parámetros a estimar y puede ser útil cuando se trabaja con conjuntos de datos de alta dimensionalidad. El análisis discriminante lineal diagonal se utiliza para encontrar una combinación lineal de las variables predictoras que maximice la separación entre los distintos grupos (véase Análisis Discriminante Lineal, Wikipedia [23]). Por su parte, Nearest Shrunken Centroids oNSC es un algoritmo de clasificación utilizado principalmente en el análisis de datos de alta dimensionalidad, como puede ser el caso de la clasificación de microarrays de expresión génica en biología molecular. El objetivo del NSC es encontrar los centroides representativos de cada clase y luego “encoger” estos centroides hacia el origen, lo que reduce el impacto de las variables menos relevantes en la clasificación. Las observaciones futuras se asignan a aquella clase cuya distancia entre la observación y el centroide (“encogido”) es mínima . Se utilizará como referencia el artículo de Hastie y otros [20]. 1.3.3. Clasificadores lineales con regularización cuadrática y L1. En último lugar, veremos algunos clasificadores lineales con regularización tanto cuadrática como L1. Para estos últimos, estudiaremos diversas aplicaciones. Clasificadores lineales con regularización cuadrática. Análisis Discriminante Regularizado : Una extensión del análisis discriminante clásico que utiliza técnicas de regularización para controlar la complejidad del modelo y evitar el sobreajuste. Ayuda a mejorar la capacidad de generalización del modelo al penalizar coeficientes grandes. Regresión Logística con Regularización Cuadrática : Una variante de la regresión logística que incluye términos de penalización cuadrática para controlar la complejidad del modelo y prevenir el sobreajuste. Esto se logra agregando términos de regularización a la función de costo. Clasificador de vectores de soporte : Del inglés Support vector classifier. Es un algoritmo de clasificación que busca encontrar el hiperplano óptimo que maximiza el margen entre las clases en un espacio de características de alta dimensión. El SVC puede manejar datos no linealmente separables mediante el uso de funciones kernel. Veremos, además, métodos como la selección de características (Feature Selection) o los llamados “atajos computacionales cuando p≫N” (Computational Shortcuts when p≫N). Clasificadores lineales con regularización L1 Regularización Lasso : Un método de análisis de regresión que realiza selección de variables y regularización para mejorar la exactitud e interpretabilidad del modelo estadístico producido por este. Elastic Net : Un método de regresión regularizado que combina las regularizaciones L1yL2de los modelos Lasso y Ridge. Esto permite beneficiarse de las propiedades de selección de variables de la regresión Lasso mientras controla la complejidad del modelo de manera más suave. CAPÍTULO 1. INTRODUCCIÓN - PRESENTACIÓN 11 Capítulo 2 Análisis teórico de algunas técnicas para combatir el problema de la alta dimensionalidad. A lo largo de todo este capítulo, se estudiarán diferentes técnicas para combatir el problema de la alta dimensionalidad. En particular, se abordarán el Análisis de Componentes Principales mediante la descomposición de la matriz en valores singulares, el Análisis Discriminante Lineal Diagonal, así como algunos clasificadores lineales con reegularización cuadrática y L1. Este estudio se enfocará desde el punto de vista teórico en primer lugar, y luego se verán algunos ejemplos para ayudar a su comprensión. 2.1. Reducción de la dimensión: Análisis de Componentes Principales. La técnica principal que se abordará en este capítulo es la técnica de Análisis de Componentes Principales oPCA. El Análisis de Componentes Principales (PCA) es una técnica fundamental en el análisis de datos utilizada para reducir la dimensionalidad de un conjunto de datos con un número elevado de variables. Su objetivo principal es encontrar un conjunto de nuevas variables no correlacionadas (a las que se les llamará componentes principales) que capturen la mayor parte de la variabilidad presente en los datos originales. Técnicamente, el PCA busca la proyección según la cual los datos queden mejor representados en términos de mínimos cuadrados. 13 2.1. REDUCCIÓN DE LA DIMENSIÓN: ANÁLISIS DE COMPONENTES PRINCIPALES. Figura 2.1: (Fuente: Análisis de Componentes Principales, Wikipedia) En la Figura 2.1 vemos un ejemplo del PCA de una distribución normal multivariante. Los vectores muestran los autovectores de la matriz de correlación escalados mediante la raíz cuadrada del correspondiente autovalor, y desplazados para que su origen coincidan con la media estadística. En el Análisis de Componentes Principales encontramos las direcciones en los datos con mayor variación, es decir, los vectores propios correspondientes a los autovalores más grandes de la matriz de covarianza y se proyectan los datos en estas direcciones. La motivación para hacer esto es que la mayor parte de la información de segundo orden va en estas direcciones. La elección del número de las direcciones se hace, a menudo, por prueba y error, pero también existen métodos basados en algunos principios o criterios objetivos. Como se recoge en Uriel y Aldás ([21]), cuando se examina la información proporcionada por un conjunto de pvariables cuantitativas, X1, X2, . . . , Xp, que pueden estar correlacionadas en mayor o menor grado, surge la duda de si sería suficiente analizar únicamente las menos correlacionadas entre ellas o alguna combinación de estas. Además, cuando el número de variables pes elevado, analizar todas en conjunto puede ser complicado. Por ello, es útil reducir el número de variables a r(donde r < p), que contengan la mayor parte de la información de las pvariables originales. El Análisis de Componentes Principales (PCA) es una técnica utilizada para convertir las variables originales en un conjunto nuevo de variables llamadas componentes principales, que no están correlacionadas entre sí. Estas componentes son el resultado de combinar linealmente las variables originales y se pueden ordenar según la cantidad de información que capturan del conjunto de variables original. La primera componente principal es la que recoge una mayor cantidad de información, luego le sigue la segunda y así sucesivamente. En términos generales, el objetivo principal del PCA es reducir la dimensionalidad del problema, reduciendo el número de variables de parcomponentes, siendo r < p, asegurando que se pierda la menor cantidad de información posible. Esta reducción puede ser útil para facilitar la aplicación de otras técnicas estadísticas como la regresión o el análisis de conglomerados, que se llevarían a cabo sobre el conjunto de componentes 14CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.1. REDUCCIÓN DE LA DIMENSIÓN: ANÁLISIS DE COMPONENTES PRINCIPALES. principales que representan a las variables originales. 2.1.1. Planteamiento del problema En esta sección se llevará a cabo el planteamiento del problema al que nos enfrentamos. Para ello, nos apoyaremos en los trabajos de Uriel y Aldás ([21]) y en los de Madsen y otros ([13]), así como en el trabajo realizado por la Universidad de Murcia ([4]). Supongamos que disponemos de un conjunto de datos con nobservaciones (individuos) y con pvariables. Podemos expresar esto mediante una matriz X, donde X∈Rn×p. Cabe destacar que el PCA se aplica sobre variables numéricas. Para datos categóricos se utiliza otro tipo de técnicas como el análisis de correspondencias. 2.1.1.1. Descomposición en valores singulares SVD En álgebra lineal la factorización de una matriz es la descomposición de la misma como producto de dos o más matrices según una forma canónica. En esta sección nos centraremos en la descomposición en valores singulares SVD apoyándonos en el artículo de Madsen y otros ([13]). Una matriz real X, donde X∈Rn×p, posee una descomposición de la forma X=UΓVT,(2.1) donde Ues una matriz de dimensión n×pcon columnas ortonormales (UTU=I), mientras que Ves una matriz ortonormal (VTV=I) de dimensión p×p, y Γes una matriz diagonal de dimensión p×pcuyos elementos diagonales son positivos o cero, los llamados valores singulares. A partir de esta matriz Xpodemos construir dos matrices simétricas, definidas positivas, XXTyXTX, y aplicando la descomposición vista anteriormente en (2.1), se tendría que XXT=UΓVTVΓUT=UΓ2UT(2.2) XTX=VΓ2VT(2.3) Supongamos que se tiene n≥p. Vamos a ver que las matrices XXT, cuya dimensión es n×nyXTX, cuya dimensión es p×pcomparten pautovalores y los restantes n−p autovalores de XXTserán nulos. Usando la descomposición vista anteriormente en (2.3), podemos identificar los autovectores y autovalores de XTXcomo las columnas de Vy los elementos diagonales de Γ al cuadrado, respectivamente. (Esto último muestra que los autovalores de XTXdeben ser no negativos). Denotando un autovector, por ejemplo, por vy el elemento diagonal correspondiente por γ, se tiene XTXv =γ2v luego multiplicando a ambos lados por Xobtenemos XXTXv =γ2Xv CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 15 2.1. REDUCCIÓN DE LA DIMENSIÓN: ANÁLISIS DE COMPONENTES PRINCIPALES. Pero esto significa que tenemos un autovector u=Xv cuyo autovalor asociado es γ2para la matriz XXT, pues XXTXv =γ2Xv Acabamos de ver que XXTyXTXcomparten pautovalores. Lo siguiente que se probará será que los restantes n−pautovalores de XXTson nulos. Para ver esto, consideremos un autovector de XXT,u⊥:XXTu⊥=β⊥u⊥, que es ortogonal a los pautovectores ui ya vistos anteriormente, es decir, UTu⊥= 0. Usando la descomposición XXT=UΓ2UT, inmediatamente se ve que los autovalores restantes β⊥deben ser todos nulos. XXTu⊥=UΓ2UTu⊥= 0. El rango Rde XXTviene determinado por la dimensión más pequeña de X,(R ≤ p). Esto asegura que XXTtiene, como mucho, pautovalores mayores que cero. Concuerda con lo que acabamos de ver, donde hemos concluido que los autovalores restantes (n−p) de XXTson cero. Es evidente que a la hora de calcular los autovalores no nulos, estas dos matrices cuadradas pueden ser intercambiadas. Esta es una propiedad de la que podemos aprovecharnos cuando estamos trabajando con matrices que tienen muchas más columnas(variables) que filas (observaciones), para así facilitar los cálculos y sea menos costoso computacionalmente. 2.1.2. Análisis de Componentes Principales mediante descomposición SVD. Sea X(con componentes Xj, j = 1, . . . , p ) un vector aleatorio con distribución de probabilidad PX. Sea también {Xα|α= 1, . . . , n}una muestra de X. Denotaremos por convención para la matriz Xlas columnas como las variables j= 1, . . . , p y las filas como las observaciones α= 1, . . . , n : en otras palabras, los elementos de nuestra matriz vendrán dados por xα,j. Esto es, la matriz X∈Rn×p. Antes de aplicar el análisis de componentes principales, en primer lugar, hemos de centrar los datos. Para ello, se resta a cada columna la media de la variable a la que pertenece. Con esto se consigue que todas las variables tengan media cero y evitar el sesgo. Cuando las varianzas de las variables son muy distintas, también se suelen escalar los datos a una varianza unitaria para eliminar el efecto de las distintas unidades en las que puedan estar medidas estos datos, y que ninguna tenga más representación que otra por estar en unidades más altas. Para cada variable, definimos una nueva variable con media 0 que viene dada por lo siguiente : Xj:= Xj−1 n n X α=1 Xα,j =Xj−Xjj= 1,··· , p y reformulando, se tiene entonces que la matriz de covarianzas muestral es : C≡1 nXTX de donde hemos eliminado la media de los datos: Xα,j := Xα,j −Xj. 16CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.1. REDUCCIÓN DE LA DIMENSIÓN: ANÁLISIS DE COMPONENTES PRINCIPALES. Podemos usar SV D para llevar a cabo el PCA. Para ello, descomponemos Xusando SV D, es decir, X=UΓVT y vemos que podemos expresar la matriz de covarianzas como C=1 nXTX=1 nVΓ2VT Con estas operaciones, podríamos calcular U,VyΓcomo los autovalores y autovectores respectivos de las matrices XTXyXXT, para así obtener la descomposición de Xen valores singulares. Para simplificar estas operaciones, también es importante señalar el hecho de que XV =UΓ. Si se denota la matriz de vectores propios ordenados según los autovalores por ˜ U, entonces podemos utilizar el PCA para hacer una transformación de los datos como sigue: Y=X˜ U. Los autovectores son precisamente las llamadas “componentes principales”. Si seleccionamos únicamente, por ejemplo, las primeras dcolumnas de Y, habríamos proyectado los datos de pdimensiones en ddimensiones. En este caso, se tiene que Ves una matriz p×pyUes de dimensión n×p. Sabiendo que esta descomposición ordena los valores singulares en orden descendente, tenemos que si p<n, las pcolumnas de Vcorresponden a los autovectores asociados a los autovalores ordenados de C. Por otra parte, si p≥n, únicamente nos interesan las primeras ncolumnas de V, pues corresponden a los autovectores asociados a los primeros nautovalores de C, todos ellos ordenados de mayor a menor, pues los p−nautovalores restantes son nulos. El resto de las p−ncolumnas no interesan para el cálculo y se corresponden con los autovectores asociados a los p−nautovalores nulos de la matriz XTX. Estos autovectores o columnas se pueden ignorar a la hora de hacer la transformación de nuestros datos. Esta transformación de los datos puede escribrise como sigue : Y=X˜ U=UΓVT˜ U, donde VT˜ Ues una matriz p×p, si n>póp×nsi p≥n, con unos en la diagonal y cero en los elementos fuera de la diagonal. Esto se debe a que todos los autovectores son ortogonales entre sí y además están normalizados, es decir, ambas matrices poseen columnas ortonormales. Para concluir, hay que destacar que esta transformación Yde los datos en términos de las componentes principales tiene dimensión, a priori, n×psi n>p ón×nsi p≥n. Esto quiere decir que si tenemos un conjunto de datos con el problema que queríamos abordar de p≫N, al aplicar PCA obtenemos una simplificación de los datos donde, a lo sumo, se tienen Ncomponentes principales. Es decir, el número máximo de componentes principales es min(n, p). Es un problema posterior el elegir un número k < N conveniente de componentes principales que expliquen un gran porcentaje de la varianza original. El Análisis de Componentes Principales, como se ha mencionado anteriormente, busca combinaciones lineales de las variables originales de manera que extraigan la máxima varianza del conjunto de variables original. Atendiendo nuevamente a los trabajos de Uriel y Aldás ([21]), concluimos que para obtener las componentes principales, se realiza una descomposición en valores singulares de la matriz de covarianzas o de correlaciones. Si las variables originales tienen variabilidades CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 17 2.2. TÉCNICAS DE CLASIFICACIÓN Figura 2.5: Gráficas de error para un conjunto de datos. (Fuente: Elements of Statistical Learning, cap.18) La Figura 2.5 muestra las distintas gráficas de error para un conjunto de datos. Estas gráficas muestran cómo varía el error de clasificación para el conjunto de entrenamiento, validación cruzada de 10 iteraciones, y para un conjunto de test para los distintos valores del parámetro ∆. El valor ∆ = 4.34 es el valor óptimo elegido después de realizar la validación cruzada, obteniéndose una reducción de dimensionalidad pasando a tener 43 variables (genes). Figura 2.6: (Fuente: Elements of Statistical Learning, cap.18) 24CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.2. TÉCNICAS DE CLASIFICACIÓN La Figura 2.5 demuestra la efectividad de la reducción. Sin reducción, cometemos 5/20 errores en los datos de test y algunos errores en los datos de entrenamiento y validación cruzada. Los centroides reducidos logran cero errores en los datos de test para una abanico bastante amplio de valores de ∆. La Figura 2.6 muestra los cuatro centroides dkj de cada gen para las distintas clases en los datos (gris). Las barras azules son las versiones reducidas de estos centroides, d′ kj obtenidas mediante la umbralización suave de las barras grises, usando ∆=4.3. Las puntuaciones discriminantes (2.4) se pueden usar para construir estimaciones de probabilidad de clase: ˆpk(x∗) = e1 2δk(x∗) PK l=1 e1 2δl(x∗).(2.10) Estas probabilidades permiten clasificar observaciones, o en su lugar, decidir no clasificar una observación si la probabilidad es baja. Nótese que se pueden usar otras formas de selección de características, incluyendo también la umbralización dura vista en esta sección. Esto muestra la importancia de llevar a cabo algún tipo de selección de características con análisis discriminante lineal diagonal en problemas de alta dimensionalidad. 2.2.2. Clasificadores Lineales con Regularización Cuadrática. En esta sección del capítulo 2, nos apoyaremos en el libro “Elements of Statistical Learning” por Hastie, Tibshirani y Friedman (véase [6]) para analizar y estudiar diversos métodos de clasificación lineal con regularización cuadrática. Estos métodos, como pueden ser el Análisis Discriminante Regularizado o la Regresión Logística con Regularización Cuadrática son útiles en problemas de alta dimensionalidad donde el número de variables pes mucho mayor que el número de observaciones N. 2.2.2.1. Análisis Discriminante Regularizado. El Análisis Discriminante Lineal (LDA, por sus siglas en inglés) es una técnica clásica de clasificación supervisada utilizada para encontrar la combinación lineal de características que mejor separa dos o más clases de objetos. El análisis discriminante lineal conlleva la inversión de una matriz de covarianza de dimensión p×p. Sin embargo, cuando p≫N, esta matriz puede ser enorme, de hecho tiene rango máximo N < p, y por lo tanto no es invertible. En estos casos, el LDA tradicional se vuelve ineficaz debido a problemas de sobreajuste y singularidad en la matriz de covarianza. Para abordar estos problemas, Jerome Friedman, en 1989 (véase Friedman [3]) propuso el Análisis Discriminante Regularizado (RDA, del inglés Regularized Discriminant Analysis), una combinación entre el análisis discriminante lineal LDA y el análisis discriminante cuadrático QDA capaz de superar estos problemas de singularidad regularizando la estimación de la matriz de covarianzas ˆ Σ. La matriz de covarianza regularizada tiene la forma ˆ Σk(α) = αb Σk+ (1 −α)ˆ Σ,(2.11) CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 25 2.2. TÉCNICAS DE CLASIFICACIÓN donde ˆ Σes la matriz de covarianzas muestral tal y como se usa en LDA, y b Σkes la matriz de covarianzas muestral restringida a la clase k. Aquí, el parámetro a determinar α∈[0,1] permite una combinación de modelos entre LDA y QDA. En la práctica, αse puede elegir basando el desempeño del modelo en datos de validación, o mediante validación cruzada. Figura 2.7: Error de entrenamiento y de test para un conjunto de datos (Fuente: Elements of Statistical Learning, cap.4) La Figura 2.7 muestra los resultados del RDA aplicado a un cierto conjunto de datos. Tanto el error de entrenamiento como el de test mejoran con el aumento de α, aunque el error de test aumenta bruscamente después de α= 0.9. La gran diferencia entre el error de entrenamiento y de test se debe en parte al hecho de que hay muchas observaciones repetidas sobre un pequeño número de individuos, diferentes en los conjuntos de entrenamiento y test. Algunas modificaciones similares permiten que ˆ Σpueda ser reducida hacia la covarianza escalar, ˆ Σ(γ) = γˆ Σ + (1 −γ)σ2I(2.12) para γ∈[0,1]. Aquí σ2Ies una matriz escalar que asegura que la matriz de covarianza no sea singular. Reemplazando ˆ Σen (2.11) por ˆ Σ(γ), se obtiene una familia más general de covarianzas ˆ Σk(α, γ)dependientes de un par de parámetros. Cuando nos encontramos en el caso en el que estamos trabajando con conjuntos de datos de dimensión muy alta, por ejemplo un array de expresión génica, los métodos que se utilizarán se centran en el caso γ= 0 en (2.12), y otros métodos regularizados del análisis discriminante lineal. La idea principal es encontrar un balance entre ambas mediante la introducción de un parámetro γde regularización. Aquí, por ejemplo, se utiliza una versión de RDA que regulariza la matriz de covarianzas ˆ Σhacia su diagonal: ˆ Σ(γ) = γˆ Σ + (1 −γ)diag(ˆ Σ),con γ∈[0,1].(2.13) 26CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.2. TÉCNICAS DE CLASIFICACIÓN Nótese que para γ= 0, el método se corresponde con el LDA diagonal, que es la versión “sin reducción” de los centroides reducidos más cercanos. La forma de la reducción en (2.13) es similar a la regresión ridge, que reduce la matriz de covarianza total de las variables hacia una matriz diagonal (escalar). Utilizando el resumen teórico visto en la sección 1.1.1, y suponiendo que el vector aleatorio Xsigue una distribución normal multivariante para cada clase i, se tiene que la regla óptima de clasificación viene dada, en este caso, por : x0se asigna aGisiTi(x0) = m´ax j=1,...,k Tj(x0), siendo x0una nueva observación que queremos clasificar, y siendo Tj(x0)una función discriminante tal y como se tiene definida en la sección 1.1.1. En el contexto de alta dimensionalidad, los problemas se agravan debido a la maldición de la dimensionalidad. La regularización se vuelve esencial para obtener modelos que generalicen bien, pues estos sufren de sobreajuste. Esto es, se ajusta demasiado a los datos de entrenamiento pero falla en generalizar a nuevos datos. Cuando p≫N, el coste computacional de invertir la matriz de tamaño p×pse puede combatir utilizando varios métodos, entre ellos la descomposición en valores singulares vista en la sección anterior, y así trabajar en una dimensión de tamaño N, menor que p. En estos casos, la selección de los parámetros de regularización es crucial. El valor óptimo de estos parámetros normalmente se suele elegir mediante validación cruzada. En conclusión, el Análisis Discriminante Regularizado es una técnica poderosa para la clasificación en situaciones de alta dimensionalidad. Al introducir términos de regularización, se mejora la estabilidad y la precisión del modelo, permitiendo su aplicación en contextos donde el número de dimensiones es mucho mayor que el número de observaciones. 2.2.2.2. Regresión Logística con Regularización Cuadrática. La regresión logística es un método de clasificación muy utilizado. Sin embargo, en contextos de alta dimensionalidad, donde el número de características (p) es mucho mayor que el número de observaciones (N), los modelos de regresión logística pueden sufrir de sobreajuste y otros problemas debido a la singularidad de las matrices involucradas. La regularización cuadrática, también conocida como Regresión Ridge ó L2, es una técnica que introduce una penalización en los coeficientes del modelo para mejorar su estabilidad y capacidad de generalización. El modelo de regresión logística surge cuando se quieren modelar las probabilidades a posteriori de las Kclases mediante funciones lineales en x, teniéndose que la suma de estas sea igual a uno y permanezcan en [0,1]. El modelo tiene la forma CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 27 2.2. TÉCNICAS DE CLASIFICACIÓN log Pr(G= 1|X=x) Pr(G=K|X=x)=β10 +βT 1x(2.14) log Pr(G= 2|X=x) Pr(G=K|X=x)=β20 +βT 2x(2.15) . . . (2.16) log Pr(G=K−1|X=x) Pr(G=K|X=x)=β(K−1)0 +βT K−1x(2.17) El modelo se especifica en términos de K−1log-odds o transformaciones logit (reflejando la condición de que las probabilidades sumen uno). Aunque el modelo utiliza la última clase Kcomo denominador en las razones de probabilidades, la elección del denominador es arbitraria en el sentido de que las estimaciones no varían bajo esta elección. Un cálculo simple muestra que Pr(G=k|X=x) = exp(βk0+βT kx) 1 + PK−1 ℓ=1 exp(βℓ0+βT ℓx), k = 1, . . . , K −1(2.18) Pr(G=K|X=x) = 1 1 + PK−1 ℓ=1 exp(βℓ0+βT ℓx)(2.19) y claramente suman a uno. Agrupamos el conjunto de parámetros a estimar θ= {β10, βT 1, . . . , β(K−1)0, βT K−1}y denotamos las probabilidades Pr(G=k|X=x) = pk(x;θ). Cuando K= 2, estamos ante una regresión logística binaria. Este modelo es bastante simple, ya que hay una sola función lineal. Es muy utilizado en aplicaciones bioestadísticas donde la variable respuesta es una variable binaria, esto es, tiene dos clases. Por ejemplo, el que los pacientes sobrevivan o mueran, tengan una enfermedad cardíaca o no, o alguna condición está presente o ausente. Los modelos de regresión logística suelen ajustarse mediante máxima verosimilitud, utilizando la verosimilitud condicional de Gdado X. Dado que Pr(G|X)especifica completamente la distribución condicional, la distribución multinomial es la apropiada. La log-verosimilitud para Nobservaciones es ℓ(θ) = N X i=1 log pgi(xi;θ),(2.20) donde pk(xi;θ) = Pr(G=k|X=xi;θ). Discutimos detalladamente el caso de dos clases, ya que los algoritmos se simplifican bastante. Considerando el caso de dos clases, codificamos gimediante una variable respuesta yibinaria(0/1), donde yi= 1 cuando gi= 1, y yi= 0 cuando gi= 2. Sea p1(x;θ) = p(x;θ) yp2(x;θ) = 1 −p(x;θ). La log-verosimilitud puede escribirse como ℓ(β) = N X i=1 {yilog p(xi;β) + (1 −yi)log(1 −p(xi;β))}(2.21) 28CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.2. TÉCNICAS DE CLASIFICACIÓN = N X i=1 nyiβTxi−log(1 + eβTxi)o.(2.22) Aquí β={β10, β1}, y asumimos que el vector de entradas xiincluye el término constante 1para simplificar los cálculos incluyendo el término independiente. Para maximizar la log-verosimilitud y obtener los parámetros del modelo, igualamos sus derivadas a cero. Estas ecuaciones son entonces ∂ℓ(β) ∂β = N X i=1 xi(yi−p(xi;β)) = 0, La regresión logística, tal y como la acabamos de ver, puede modificarse de manera similar para tratar el caso p≫N. Con Kclases, utilizamos una versión simétrica del modelo de regresión logística multiclase: Pr(G=k|X=x) = expβk0+xTβk PK ℓ=1 exp(βℓ0+xTβℓ). Se tienen Kvectores de coeficientes de parámetros de log-odds β1, β2, . . . , βK. Regularizamos el ajuste maximizando la log-verosimilitud añadiendo una penalización m´ax {β0k,βk}K 1"N X i=1 log Pr(gi|xi)−λ 2 K X k=1 ∥βk∥2 2#. Esta regularización resuelve automáticamente la redundancia en la parametrización, es decir, el parámetro de regularización penaliza los coeficientes de tal manera que se evita el problema de sobreparametrización y colinealidad. Además, obliga a que PK k=1 ˆ βkj = 0, j = 1, . . . , p (Se puede comprobar fácilmente derivando la expresión). Nótese que los términos constantes βk0no están regularizados y por lo tanto, se debe imponer una restricción sobre ellos, por ejemplo, que PK k=1 ˆ βk0= 0. Esto es porque considerando los βk0=βk0+ α, las probabilidades no varían, teniéndose así infinitas soluciones. Para corregir este problema, se impone la condición anterior. El problema de optimización resultante es convexo y puede resolverse mediante el algoritmo de Newton u otros métodos numéricos. Este método ofrece una alternativa atractiva al clasificador de vectores de soporte, que se verá a continuación, especialmente en el caso multiclase. Bajo estas condiciones, y una vez aplicada la regularización, si se tiene una nueva observación x0, la regla de clasificación vendrá dada por : x0se asigna a Gisi c Pr[Gi|X=x0] = m´ax jc Pr[Gj|X=x0], j = 1, . . . , k. 2.2.2.3. Clasificador de vectores de soporte. El clasificador de vectores de soporte (SVM, del inglés Support Vector Classifier) es una técnica de aprendizaje supervisado usada principalmente para problemas de clasificación. La idea básica detrás de SVM es encontrar el hiperplano óptimo que separa las diferentes clases en el espacio de características. CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 29 2.2. TÉCNICAS DE CLASIFICACIÓN Cuando se tiene un problema de alta dimensionalidad (donde p≫N), nos encontramos ante una situación especial. Veremos el funcionamiento de SVM y cómo se adapta este modelo a situaciones de alta dimensionalidad. Supongamos que tenemos un conjunto de entrenamiento {(xi, yi)}N i=1, donde xi∈Rpes un vector de características e yi∈ {−1,1}es la etiqueta de clase. El objetivo de SVM es encontrar el hiperplano que maximiza la separación(margen) entre las dos clases. El hiperplano se puede definir como: βTx+β0= 0 donde β∈Rpes un vector de pesos y β0es el término independiente. La regla de clasificación inducida por el hiperplano se puede expresar de la siguiente forma : G(x) = sign[βTx+β0]. Como las clases son separables, podemos encontrar una función f(x) = βTx+β0con yif(xi)>0,∀i. Figura 2.8: Clasificador de vector de soporte.(Fuente: Elements of Statistical Learning, cap.12) Queremos entonces encontrar el hiperplano que crea el mayor margen entre los puntos de entrenamiento para las clases 1 y -1 (ver el panel izquierdo de la Figura 2.8). El problema de optimización m´ax β,β0,∥β∥=1 M sujeto a yixT iβ+β0≥M, i = 1, . . . , N, (2.23) resuelve esto que buscamos. Los puntos en la figura están alejados M unidades del hiperplano por cada lado, y por tanto esta distancia es de 2M unidades de longitud. Se le llama el margen. Los vectores de soporte son aquellos puntos que están en el borde del margen. El margen γse puede escribir como: M=1 ∥β∥. 30CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.2. TÉCNICAS DE CLASIFICACIÓN Para maximizar el margen, necesitamos minimizar ∥β∥. El problema de optimización se puede reformular como m´ın β,β0∥β∥ sujeto a yixT iβ+β0≥1, i = 1, . . . , N, (2.24) de donde hemos eliminado la restricción ∥β∥= 1. Este es un problema de optimización cuadrática con restricciones lineales. El caso en el que las clases no sean separables (véase el panel derecho de la figura 2.8) no se estudiará. En el caso de alta dimensionalidad, donde el número de características pes mucho mayor que el número de muestras N, el SVM sigue siendo aplicable, pero se pueden considerar ciertas condiciones especiales de regularización. Para evitar el sobreajuste debido a la alta dimensionalidad, se introduce un término de regularización. El problema de optimización se convierte en: m´ın β,β0"1 2∥β∥2+C N X i=1 ξi#sujeto a yi(βTxi+β0)≥1−ξi, ξi≥0,∀i(2.25) donde ξison variables de holgura que permiten errores en la clasificación, y Ces un hiperparámetro que controla el equilibrio entre maximizar el margen y minimizar los errores de clasificación. El caso separable se corresponde con C=∞. Cuando p≫N, este método es especialmente interesante porque, en general, las clases son perfectamente separables por un hiperplano a menos que hayan vectores de características idénticos en diferentes clases. Sin ninguna regularización, el clasificador de vectores de soporte encuentra el hiperplano de separación con el margen más grande. Esto es, el hiperplano que produce la mayor distancia entre las clases en los datos de entrenamiento. Un hecho sorprendente es que en la mayoría de ocasiones, cuando p≫Nel clasificador de vectores de soporte no regularizado funciona tan bien como el mejor clasificador regularizado. El sobreajuste aquí no parece ser un problema, en parte debido a la insensibilidad de la pérdida por clasificación errónea. Existen muchos métodos diferentes para generalizar el clasificador de vectores de soporte de dos clases a K > 2clases. En el enfoque “one versus one” (ovo), calculamos todos los K 2clasificadores por cada par de clases. Para cada observación del conjunto de test, la clase asignada es aquella que más veces se repite al utilizar todos los clasificadores. En el “one versus all” (ova), cada clase se compara con todas las demás. Se tienen K clasificadores. Para clasificar una observación del conjunto de test, calculamos la distancia signada al hiperplano para cada uno de los Kclasificadores. La clase ganadora es aquella clase con la que se obtiene la mayor distancia. Tibshirani y Hastie (2007) plantearon un clasificador al que llamaron margin tree, en el que el clasificador de vectores de soporte se utiliza en un árbol binario, al igual que en CART (Classification and Regression Trees). Las clases se organizan de forma jerárquica, lo que puede ser útil para clasificar a pacientes con diferentes tipos de tumores, por ejemplo. Las tasas de error son insensibles a la elección de C (el parámetro de regularización de (2.25), para valores de C > 0.001. Cuando p > N, el hiperplano de vectores de soporte puede separar perfectamente los datos de entrenamiento eligiendo C=∞ CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 31 2.2. TÉCNICAS DE CLASIFICACIÓN 2.2.2.4. Selección de características. La selección de características es un requisito importante para un clasificador cuando p es grande. Ni el análisis discriminante, ni la regresión logística, ni el clasificador de vectores de soporte realizan la selección de características automáticamente, porque todos usan regularización cuadrática. Todas las variables o características tienen pesos no nulos en ambos modelos. Se han propuesto métodos ad-hoc para la selección de características, por ejemplo, eliminar genes (variables) con coeficientes pequeños y reajustar el clasificador. Esto se hace de manera regresiva, comenzando con los pesos más pequeños y avanzando hacia los pesos más grandes. Esto se conoce como eliminación recursiva de características (Guyon et al., 2002). Los tres métodos discutidos en esta sección (RDA, LR y SVM) se pueden adaptar para crear límites de decisión no lineales mediante el uso de lo que se llama “kernels”. Los kernels permiten a estos métodos manejar límites de decisión que no son lineales en la estructura de los datos. Esto significa que pueden capturar relaciones más complejas entre las características de entrada y las salidas deseadas en un problema. Generalmente, la principal motivación para utilizar kernels es aumentar la complejidad del modelo. Con p≫Nlos modelos ya son suficientemente complejos y existe un riesgo elevado de sobreajuste. Sin embargo, a pesar de la alta dimensionalidad, los kernels radiales, también conocidos como kernels RBF (Radial Basis Function), son particularmente útiles en problemas de alta dimensionalidad. Estos kernels ayudan a reducir la influencia de valores atípicos al disminuir el impacto de los productos internos entre puntos de datos que están muy separados en el espacio de características. Esto proporciona robustez al modelo frente a datos atípicos, lo cual es común en problemas con muchas dimensiones. En resumen, aunque los modelos en problemas de alta dimensionalidad ya son complejos, los kernels radiales pueden ofrecer resultados superiores. Esto se debe a su capacidad para manejar relaciones no lineales de manera efectiva y para mitigar los efectos negativos de valores atípicos en los datos. Estos métodos ayudan a identificar y eliminar características irrelevantes, mejorando así la interpretabilidad y el rendimiento del modelo. 2.2.2.5. Estrategias computacionales cuando p≫N. Las técnicas computacionales discutidas en esta sección se aplican a cualquier método que ajuste un modelo lineal con regularización cuadrática en los coeficientes. Cuando p≫N, los cálculos se pueden realizar en un espacio de dimensión Nen lugar de p, mediante la descomposición en valores singulares introducida en la Sección 14.5. Aquí está la intuición geométrica: así como dos puntos en un espacio tridimensional siempre yacen en una línea, Npuntos en un espacio de dimensión pyacen en un subespacio afín de dimensión N−1. Dada la matriz de datos N×pX, sea X=UDVT=RVT(2.26) la descomposición en valores singulares (SVD) de X; es decir, Ves p×Ncon columnas ortonormales, Ues N×Nortogonal, y Duna matriz diagonal con elementos d1≥d2≥ . . . ≥dN. La matriz Res N×N, con filas rT i. 32CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 2.2. TÉCNICAS DE CLASIFICACIÓN Como ejemplo, primero consideremos las estimaciones de una regresión ridge (resumen teórico incluido en la sección 1.1.3.1): ˆ β= (XTX+λI)−1XTy(2.27) Reemplazando Xpor RVTy después de algunas simplificaciones, esto se puede expresar como ˆ β=V(RTR+λI)−1RTy.(2.28) Así ˆ β=Vˆ θ, donde ˆ θes el estimador de la regresión ridge utilizando las Nobservaciones (ri, yi),i= 1,2, . . . , N. En otras palabras, podemos reducir la matriz de datos Xa una matriz de dimensión menor R, y trabajar con las filas de R. Este “truco” se puede probar que reduce el coste computacional de los cálculos de O(p3)aO(pN2)cuando p≫N. Estos resultados pueden generalizarse a todos los modelos que sean lineales y que tengan regularizaciones cuadráticas. Consideremos cualquier problema de aprendizaje supervisado donde se use una función lineal f(X) = β0+XTβpara modelar un parámetro en la distribución condicional de Y|X. Para ello, ajustamos los parámetros βminimizando alguna función de pérdida PN i=1 L(yi, f(xi)) sobre datos con una regularización cuadrática en β. La regresión logística puede ser un ejemplo útil a tener en cuenta. Entonces, como mencionan Hastie y otros en [6], se tiene el siguiente resultado: Sea f∗(ri) = θ0+rT iθcon ridefinido en (2.26), y consideremos el par de problemas de optimización: (ˆ β0,ˆ β) = arg m´ın β0,β∈Rp N X i=1 L(yi, β0+xT iβ) + λβTβ(2.29) (ˆ θ0,ˆ θ) = arg m´ın θ0,θ∈RN N X i=1 L(yi, θ0+rT iθ) + λθTθ(2.30) Entonces ˆ β0=ˆ θ0, y ˆ β=Vˆ θ. El teorema dice que podemos simplemente reemplazar los pvectores xipor los Nvectores ri, y realizar nuestro ajuste regularizado como antes, pero con menos variables predictivas. El vector ˆ θde dimensión Nse puede transformar luego a un vector de dimensión pusando simplemente una multiplicación matricial. Geométricamente, estamos rotando las características o variables a un sistema de coordenadas en el que todas menos las primeras Ncoordenadas son cero. Estas rotaciones se pueden hacer porque la penalización o regularización cuadrática es invariante bajo rotaciones, y los modelos lineales son equivalentes. Este resultado se puede aplicar a algunos de los métodos de aprendizaje vistos en este trabajo, como la regresión logística regularizada o el análisis discriminante lineal. También puede aplicarse a las redes neuronales con regularización cuadrática. Sin embargo, no se puede aplicar a métodos como la regresión LASSO, ya que usa regularizaciones no cuadráticas ( L1) en los coeficientes. Este tipo de regularización se abordará en la siguiente sección. CAPÍTULO 2. ANÁLISIS TEÓRICO DE ALGUNAS TÉCNICAS PARA COMBATIR EL PROBLEMA DE LA ALTA DIMENSIONALIDAD. 33 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. ## Mean :-0.01906 Mean :-0.02781 Mean :-0.01992 Mean :-0.3287 ## 3rd Qu.: 0.31000 3rd Qu.: 0.35250 3rd Qu.: 0.16249 3rd Qu.: 0.6925 ## Max. : 0.94000 Max. : 2.24000 Max. : 1.15000 Max. : 1.5000 ## 5 6 7 8 ## Min. :-0.82500 Min. :-0.700000 Min. :-0.92000 Min. :-0.70500 ## 1st Qu.:-0.22500 1st Qu.:-0.156250 1st Qu.:-0.24625 1st Qu.:-0.20499 ## Median : 0.00000 Median : 0.000000 Median : 0.00000 Median : 0.00000 ## Mean : 0.02609 Mean : 0.006718 Mean : 0.01969 Mean :-0.02313 ## 3rd Qu.: 0.21000 3rd Qu.: 0.184995 3rd Qu.: 0.24750 3rd Qu.: 0.16001 ## Max. : 1.71500 Max. : 1.160000 Max. : 0.94000 Max. : 0.72496 ## 9 10 ## Min. :-0.9000000 Min. :-1.65498 ## 1st Qu.:-0.3600000 1st Qu.:-0.47747 ## Median : 0.0000000 Median : 0.00000 ## Mean : 0.0007803 Mean : 0.01924 ## 3rd Qu.: 0.2950000 3rd Qu.: 0.46752 ## Max. : 0.9900000 Max. : 1.49002 Aquí se puede ver un pequeño resumen de las primeras 10 variables a considerar en nuestro estudio. Como observamos, la media está muy cercana a cero en casi todas ellas. De todas formas, más adelante centraremos las variables para evitar cometer errores de interpretación. Por un lado, en la variable nci.data se tienen almacenado los datos, y en nci.labs las correspondientes etiquetas para cada observación. Vamos a ver cuántos tipos de cáncer diferentes nos encontramos en nuestro conjunto de datos y el número de observaciones asociadas a cada uno de ellos. table(nci.labs) ## nci.labs ## BREAST CNS COLON K562A-repro K562B-repro LEUKEMIA ## 7 5 7 1 1 6 ## MCF7A-repro MCF7D-repro MELANOMA NSCLC OVARIAN PROSTATE ## 1 1 8 9 6 2 ## RENAL UNKNOWN ## 9 1 Como se puede observar, nos encontramos 14 tipos de cáncer diferentes, aunque uno está etiquetado como “UNKNOW”. Vemos que los más comunes son “NSCLC” y “RENAL”, seguidos de “MELANOMA”. 3.2. Análisis de Componentes Principales. Para tratar de corregir el problema de la alta dimensionalidad en este conjunto de datos, se procede al Análisis de Componentes Principales sobre los mismos. En primer lugar, se hallan las componentes principales y posteriormente se decidirá el número de ellas a retener en el modelo: 40 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. Aplicamos la función “princomp” al dataset para realizar el Análisis de Componentes Principales. # princomp(nci.data) ## Error in princomp.default(nci.data): 'princomp'can only be used with more units than variables Aquí surge un problema y es que “princomp” intenta calcular los autovalores y autovectores de la matriz de covarianza muestral, lo cual requiere que n>p(en nuestro caso no se tiene, y por ello la función da error). 3.2.1. Cálculos alternativos. Dada la imposibilidad de aplicar el método clásico de PCA, se ha de aplicar la descomposición en valores singulares SVD que se ha visto en el capítulo 2. Vamos a utilizar la descomposición SVD sobre la matriz de datos para llevar a cabo el PCA. En primer lugar, se centrarán los datos. Z<- scale(nci.data,center=TRUE,scale=FALSE) summary(colMeans(Z)) #Para comprobar que las variables están centradas. ## Min. 1st Qu. Median Mean 3rd Qu. Max. ## -1.214e-16 -5.746e-18 3.253e-19 2.731e-19 6.336e-18 1.041e-16 Como se observa, la media de cada variable es prácticamente 0, por lo que se comprueba lo que queríamos. El cálculo directo de la matriz de covarianza muestral (ˆ Σ=1 n−1ZTZ) es inviable en este caso, se trabajaría con una matriz de dimensiones 6830 ×6830. Se puede recurrir para combatir este problema a la descomposión en valores singulares de Z, la cual puede ser calculada de forma eficiente a pesar de estas dimensiones. Se tiene entonces que Z=UDVT, donde Ues una matriz de dimensión n×ncuyas columnas son los llamados vectores singulares izquierdos, Ves una matriz de dimensión p×ncuyas columnas son los vectores singulares derechos y la matriz Des una matriz diagonal de dimensión n×nde los valores singulares. Con esta descomposición, se verifica que los autovectores Ljdel PCA asociados a los autovalores no nulos coinciden con las columnas de V, esto es, Lj=Vj. Además, los autovalores λjdel PCA verifican λj=d2 j1 n−1, con djlos elementos diagonales de la matriz D. La función “svd” del paquete base de R Program proporciona esta descomposición. SVD<- svd(t(Z)) str(SVD) ## List of 3 ## $ d: num [1:64] 200 149 133 107 102 ... ## $ u: num [1:6830, 1:64] -0.0051 -0.00164 -0.00251 0.00494 -0.00337 ... ## $ v: num [1:64, 1:64] -0.0991 -0.1079 -0.1255 -0.1873 -0.2514 ... CAPÍTULO 3. IMPLEMENTACIÓN EN R. 41 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. Por tanto, los autovalores del PCA son aval<- (SVD$dˆ2)/(nrow(nci.data)-1) aval ## [1] 6.332156e+02 3.529278e+02 2.799189e+02 1.830830e+02 1.635573e+02 ## [6] 1.490968e+02 1.222882e+02 1.197912e+02 1.121777e+02 9.171077e+01 ## [11] 8.928951e+01 8.512997e+01 7.782686e+01 7.514521e+01 7.109609e+01 ## [16] 6.846864e+01 6.708114e+01 6.182238e+01 6.155363e+01 6.002420e+01 ## [21] 5.839505e+01 5.497665e+01 5.359226e+01 5.034073e+01 5.015011e+01 ## [26] 4.716814e+01 4.513922e+01 4.381074e+01 4.320373e+01 4.226847e+01 ## [31] 4.075683e+01 3.990295e+01 3.724347e+01 3.684533e+01 3.557327e+01 ## [36] 3.520311e+01 3.449044e+01 3.397331e+01 3.225367e+01 3.172713e+01 ## [41] 3.103290e+01 3.038929e+01 2.976259e+01 2.893812e+01 2.831753e+01 ## [46] 2.677309e+01 2.646792e+01 2.517931e+01 2.327441e+01 2.283683e+01 ## [51] 2.201186e+01 2.082772e+01 2.016363e+01 1.946066e+01 1.826767e+01 ## [56] 1.775399e+01 1.669647e+01 1.536295e+01 1.434969e+01 1.241896e+01 ## [61] 1.042527e+01 9.940037e+00 8.913814e+00 1.812288e-28 y los autovectores del PCA son, salvo parámetro multiplicativo, avec=SVD$u str(avec) # Los 64 autovectores de dimensión 6830 ## num [1:6830, 1:64] -0.0051 -0.00164 -0.00251 0.00494 -0.00337 ... En realidad esta idea la explota el procedimiento “prcomp()” que realiza el PCA a partir de la Descomposición en Valores Singulares. Veámoslo y luego se compararán los resultados: PCA=prcomp(nci.data) #Por defecto, se tiene center=TRUE, scale=FALSE str(PCA) ## List of 5 ## $ sdev : num [1:64] 25.2 18.8 16.7 13.5 12.8 ... ## $ rotation: num [1:6830, 1:64] -0.0051 -0.00164 -0.00251 0.00494 -0.00337 ... ## ..- attr(*, "dimnames")=List of 2 ## .. ..$ : chr [1:6830] "1" "2" "3" "4" ... ## .. ..$ : chr [1:64] "PC1" "PC2" "PC3" "PC4" ... ## $ center : Named num [1:6830] -0.0191 -0.0278 -0.0199 -0.3287 0.0261 ... ## ..- attr(*, "names")= chr [1:6830] "1" "2" "3" "4" ... ## $ scale : logi FALSE ## $ x : num [1:64, 1:64] -19.8 -21.5 -25.1 -37.4 -50.2 ... ## ..- attr(*, "dimnames")=List of 2 ## .. ..$ : chr [1:64] "V1" "V2" "V3" "V4" ... ## .. ..$ : chr [1:64] "PC1" "PC2" "PC3" "PC4" ... ## - attr(*, "class")= chr "prcomp" 42 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. names(PCA) ## [1] "sdev" "rotation" "center" "scale" "x" Como vimos en el capítulo anterior, la varianza de la componente principal i-ésima se corresponde con el autovalor i-ésimo λi. Con esta información, obtenemos los autovalores no nulos como la varianza de las componentes, o lo que es lo mismo, la desviación típica al cuadrado de estas. PCA$sdevˆ2 ## [1] 6.332156e+02 3.529278e+02 2.799189e+02 1.830830e+02 1.635573e+02 ## [6] 1.490968e+02 1.222882e+02 1.197912e+02 1.121777e+02 9.171077e+01 ## [11] 8.928951e+01 8.512997e+01 7.782686e+01 7.514521e+01 7.109609e+01 ## [16] 6.846864e+01 6.708114e+01 6.182238e+01 6.155363e+01 6.002420e+01 ## [21] 5.839505e+01 5.497665e+01 5.359226e+01 5.034073e+01 5.015011e+01 ## [26] 4.716814e+01 4.513922e+01 4.381074e+01 4.320373e+01 4.226847e+01 ## [31] 4.075683e+01 3.990295e+01 3.724347e+01 3.684533e+01 3.557327e+01 ## [36] 3.520311e+01 3.449044e+01 3.397331e+01 3.225367e+01 3.172713e+01 ## [41] 3.103290e+01 3.038929e+01 2.976259e+01 2.893812e+01 2.831753e+01 ## [46] 2.677309e+01 2.646792e+01 2.517931e+01 2.327441e+01 2.283683e+01 ## [51] 2.201186e+01 2.082772e+01 2.016363e+01 1.946066e+01 1.826767e+01 ## [56] 1.775399e+01 1.669647e+01 1.536295e+01 1.434969e+01 1.241896e+01 ## [61] 1.042527e+01 9.940037e+00 8.913814e+00 1.798566e-28 Comparación de los autovalores no nulos obtenidos directamente y a través de la orden prcomp(): round(cbind(aval,PCA$sdevˆ2)[1:10,],4) ## aval ## [1,] 633.2156 633.2156 ## [2,] 352.9278 352.9278 ## [3,] 279.9189 279.9189 ## [4,] 183.0830 183.0830 ## [5,] 163.5573 163.5573 ## [6,] 149.0968 149.0968 ## [7,] 122.2882 122.2882 ## [8,] 119.7912 119.7912 ## [9,] 112.1777 112.1777 ## [10,] 91.7108 91.7108 Vemos que no existe literalmente diferencia alguna entre los resultados obtenidos mediante prcomp() y realizando los cálculos directamente. La matriz de cargas o matriz que contiene a los autovectores, de dimensiones p×n, es decir, n= 64 autovectores de dimensión 6830 ×1cada uno, CAPÍTULO 3. IMPLEMENTACIÓN EN R. 43 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. dim(PCA$rotation) ## [1] 6830 64 viene dada por “PCA$rotation”. str(PCA$rotation) ## num [1:6830, 1:64] -0.0051 -0.00164 -0.00251 0.00494 -0.00337 ... ## - attr(*, "dimnames")=List of 2 ## ..$ : chr [1:6830] "1" "2" "3" "4" ... ## ..$ : chr [1:64] "PC1" "PC2" "PC3" "PC4" ... Veamos, por ejemplo, las primeras 15 componentes de PC1 (pues PC1 tiene 6830 componentes) PCA$rotation[1:15,1]# Primera componente principal ##123456 ## -0.005096247 -0.001642354 -0.002509243 0.004940063 -0.003365039 0.001382038 ## 7 8 9 10 11 12 ## -0.002909218 -0.001253347 -0.001191334 -0.003507481 -0.001068655 0.001386408 ## 13 14 15 ## 0.001290305 -0.003362776 -0.011424883 Comparación de los primeros 15 autovectores obtenidos directamente y a través de la orden “prcomp()”: for (i in 1:15) cat("Discrepancia",i,":",sum(PCA$rotation[,i]-avec[,i]),"\n") ## Discrepancia 1 : 1.607995e-15 ## Discrepancia 2 : -8.653527e-15 ## Discrepancia 3 : 9.022615e-15 ## Discrepancia 4 : -2.678409e-14 ## Discrepancia 5 : 5.186963e-14 ## Discrepancia 6 : 5.119745e-14 ## Discrepancia 7 : 2.699179e-14 ## Discrepancia 8 : -9.110396e-14 ## Discrepancia 9 : 2.960471e-14 ## Discrepancia 10 : -4.27435e-14 ## Discrepancia 11 : 3.318613e-14 ## Discrepancia 12 : -1.522626e-14 ## Discrepancia 13 : 3.870392e-15 ## Discrepancia 14 : -1.598048e-15 ## Discrepancia 15 : -2.416645e-14 44 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. Como se puede observar, la suma de las diferencias entre las componentes de cada autovector es prácticamente cero. Por tanto , cumple con lo que estábamos probando. 3.2.1.1. Selección del número de componentes. Vamos a representar en un gráfico de sedimentación los autovalores según la posición que ocupen, ordenados de mayor a menor (eje abscisas) y el valor de cada uno de estos autovalores. Además, el gráfico incluye una recta con altura igual a la media de los autovalores que tomaremos como umbral. Según uno de los criterios vistos en teoría, se seleccionarán únicamente aquellos autovalores cuyo valor esté por encima del umbral establecido. plot(aval,type="h",ylab="Valor", xlab="Autovalores") abline(h=mean(aval),col="red",lty=2,lwd=2) 0 10 20 30 40 50 60 0 100 300 500 Autovalores Valor sum(aval>mean(aval)) ## [1] 17 Como se puede observar, este criterio indica que debemos seleccionar las primeras 17 componentes principales. Vamos a calcular la proporción de varianza total acumulada por estas 17 componentes. CAPÍTULO 3. IMPLEMENTACIÓN EN R. 45 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. cumsum(aval)[17]/sum(aval) ## [1] 0.6448598 Atendiendo al resultado obtenido, las primeras 17 componentes principales explican aproximadamente el 65% de la varianza total, un resultado que, según el objetivo del estudio, podría considerarse como “pobre”. Si se quiere explicar un mayor porcentaje, por ejemplo el 90%, veamos con cuántas componentes debemos quedarnos : plot(cumsum(aval)/sum(aval),type="s",xlab="Autovalores", ylab="Porcentaje de varianza total acumulada") abline(h=0.9,col="red",lty=2) 0 10 20 30 40 50 60 0.2 0.4 0.6 0.8 1.0 Autovalores Porcentaje de varianza total acumulada Vamos a representar en una tabla, además, los correspondientes autovalores con el porcentaje de varianza total acumulado. cbind(aval,100*cumsum(aval)/sum(aval)) ## aval ## [1,] 6.332156e+02 14.89294 ## [2,] 3.529278e+02 23.19364 ## [3,] 2.799189e+02 29.77720 ## [4,] 1.830830e+02 34.08323 46 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. ## [5,] 1.635573e+02 37.93002 ## [6,] 1.490968e+02 41.43671 ## [7,] 1.222882e+02 44.31287 ## [8,] 1.197912e+02 47.13030 ## [9,] 1.121777e+02 49.76867 ## [10,] 9.171077e+01 51.92567 ## [11,] 8.928951e+01 54.02571 ## [12,] 8.512997e+01 56.02793 ## [13,] 7.782686e+01 57.85838 ## [14,] 7.514521e+01 59.62576 ## [15,] 7.109609e+01 61.29791 ## [16,] 6.846864e+01 62.90826 ## [17,] 6.708114e+01 64.48598 ## [18,] 6.182238e+01 65.94001 ## [19,] 6.155363e+01 67.38773 ## [20,] 6.002420e+01 68.79947 ## [21,] 5.839505e+01 70.17289 ## [22,] 5.497665e+01 71.46592 ## [23,] 5.359226e+01 72.72638 ## [24,] 5.034073e+01 73.91037 ## [25,] 5.015011e+01 75.08988 ## [26,] 4.716814e+01 76.19925 ## [27,] 4.513922e+01 77.26091 ## [28,] 4.381074e+01 78.29132 ## [29,] 4.320373e+01 79.30745 ## [30,] 4.226847e+01 80.30158 ## [31,] 4.075683e+01 81.26016 ## [32,] 3.990295e+01 82.19866 ## [33,] 3.724347e+01 83.07461 ## [34,] 3.684533e+01 83.94120 ## [35,] 3.557327e+01 84.77786 ## [36,] 3.520311e+01 85.60582 ## [37,] 3.449044e+01 86.41702 ## [38,] 3.397331e+01 87.21606 ## [39,] 3.225367e+01 87.97465 ## [40,] 3.172713e+01 88.72086 ## [41,] 3.103290e+01 89.45074 ## [42,] 3.038929e+01 90.16548 ## [43,] 2.976259e+01 90.86548 ## [44,] 2.893812e+01 91.54609 ## [45,] 2.831753e+01 92.21211 ## [46,] 2.677309e+01 92.84180 ## [47,] 2.646792e+01 93.46431 ## [48,] 2.517931e+01 94.05652 ## [49,] 2.327441e+01 94.60392 ## [50,] 2.283683e+01 95.14103 ## [51,] 2.201186e+01 95.65874 ## [52,] 2.082772e+01 96.14860 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 47 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. ## [53,] 2.016363e+01 96.62284 ## [54,] 1.946066e+01 97.08055 ## [55,] 1.826767e+01 97.51019 ## [56,] 1.775399e+01 97.92776 ## [57,] 1.669647e+01 98.32045 ## [58,] 1.536295e+01 98.68178 ## [59,] 1.434969e+01 99.01928 ## [60,] 1.241896e+01 99.31137 ## [61,] 1.042527e+01 99.55657 ## [62,] 9.940037e+00 99.79035 ## [63,] 8.913814e+00 100.00000 ## [64,] 1.812288e-28 100.00000 which(100*cumsum(aval)/sum(aval)>90)[1] ## [1] 42 Atendiendo al resultado obtenido anteriormente, para explicar el 90% de la varianza total hacen falta 42 componentes. Si el umbral establecido fuese el 0.8, con 30 componentes seleccionadas sería suficiente. Vamos a analizar en primer lugar el PCA con 30 componentes. ncomp<- 42 L<- PCA$rotation[,1:ncomp] #Podemos utilizar "avec" o bien "PCA$rotation".,→ dim(L) # Para ver el número de autovectores. ## [1] 6830 42 R<- cor(L) corrplot::corrplot(R, method = "number", number.cex = 0.75) 48 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.2. ANÁLISIS DE COMPONENTES PRINCIPALES. 1.00 0.00 0.01 −0.01 −0.01 0.01 0.00 0.00 0.01 0.00 0.00 0.00 0.00 −0.01 −0.01 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.01 0.00 1.00 0.00 0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.01 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.01 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 −0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 −1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1 PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8 PC9 PC10 PC11 PC12 PC13 PC14 PC15 PC16 PC17 PC18 PC19 PC20 PC21 PC22 PC23 PC24 PC25 PC26 PC27 PC28 PC29 PC30 PC31 PC32 PC33 PC34 PC35 PC36 PC37 PC38 PC39 PC40 PC41 PC42 PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8 PC9 PC10 PC11 PC12 PC13 PC14 PC15 PC16 PC17 PC18 PC19 PC20 PC21 PC22 PC23 PC24 PC25 PC26 PC27 PC28 PC29 PC30 PC31 PC32 PC33 PC34 PC35 PC36 PC37 PC38 PC39 PC40 PC41 PC42 Como se puede observar y como ya sabíamos, se comprueba que las correlaciones entre componentes principales dos a dos es nula. Los nuevos datos ya transformados, teniéndose 30 variables en vez de 6830 variables, vendrán dados por : dataPCA <- Z%* % L dim(dataPCA) ## [1] 64 42 Esto es, se tienen las mismas 64 observaciones que en un principio, pero proyectadas sobre las 30 componentes principales seleccionadas. Vamos a representar los datos sobre algunos componentes principales. En la Figura de la izquierda, representaremos los datos sobre las componentes 1 y 2, lo haremos sobre las componentes 1 y 3 en la figura del centro y sobre las componentes 2 y 3 en la de la derecha. Se aplica color a los puntos según la clase indicada en la variable “nci.labs”. par(mfrow=c(1,3)) plot(PCA$x[,1:2], col=Cols(nci.labs), pch=19, xlab="Z1",ylab="Z2")#Z1 vs Z2 plot(PCA$x[,c(1,3)], col=Cols(nci.labs), pch=19, xlab="Z1",ylab="Z3")# Z1 vs Z3 plot(PCA$x[,c(2,3)], col=Cols(nci.labs), pch=19, xlab="Z2",ylab="Z3")#Z2 vs Z3 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 49 3.4. CLASIFICADORES LINEALES CON REGULARIZACIÓN CUADRÁTICA ## 1 1 3 0 0 3 ## MCF7A-repro MCF7D-repro MELANOMA NSCLC OVARIAN PROSTATE ## 1 0 4 1 1 1 ## RENAL UNKNOWN ## 4 0 Como anteriormente en el LDLA, la función “rda” da error con aquellas clases que aparecen en menos de 2 observaciones. Para ello, descarto aquellas clases que tienen menos de 2 observaciones en el conjunto de entrenamiento. class_counts <- table(train_labels) # Filtramos aquellas clases con al menos 2 observaciones valid_classes <- names(class_counts[class_counts >= 2]) # Filtramos los datos para incluir sólo las clases válidas en el RDA. train_data <- train_data[train_labels %in % valid_classes, ] train_labels <- train_labels[train_labels %in % valid_classes] test_data <- test_data[test_labels %in % valid_classes, ] test_labels <- test_labels[test_labels %in % valid_classes] Lo siguiente que se hará será convertir los datos de entrenamiento y las etiquetas a matrices numéricas: train_data_matrix <- t(as.matrix(train_data)) train_labels_vector <- as.matrix(as.numeric(train_labels)) test_data_matrix <- t(as.matrix(test_data)) test_labels_vector <- as.matrix(as.numeric(test_labels)) Consideramos los siguientes posibles valores de alpha y delta para realizar validación cruzada. alpha_values <- seq(0,0.99,by = 0.1) delta_values <- seq(0,3,by = 0.1) Haciendo uso de la función “rda” incluida en el paquete “rda” (Guo y otros, [5]), construimos el modelo en primer lugar con alpha y deta predeterminados: modelo <- rda(train_data_matrix, train_labels_vector) predictions <- predict.rda(modelo,train_data_matrix, train_labels_vector, xnew = test_data_matrix, type="class") Realizamos ahora validación cruzada para encontrar los mejores valores de alpha y delta. 56 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.4. CLASIFICADORES LINEALES CON REGULARIZACIÓN CUADRÁTICA crossvalidation <- rda.cv(modelo, train_data_matrix, train_labels_vector, alpha = alpha_values, delta = delta_values) ## Fold 1 : ## Fold 2 : ## Fold 3 : # Encontrar la combinación de alpha y delta que minimiza el error cv_errors <- crossvalidation$cv.err min_error_idx <- which(cv_errors == min(cv_errors), arr.ind = TRUE) best_alpha <- alpha_values[min_error_idx[1,1]] best_delta <- delta_values[min_error_idx[1,2]] cat("Mejor valor de alpha:", best_alpha, "\n") ## Mejor valor de alpha: 0 cat("Mejor valor de delta:", best_delta, "\n") ## Mejor valor de delta: 0.3 Construimos el modelo con los valores de alpha y delta obtenidos: modelobueno <- rda(train_data_matrix, train_labels_vector, xnew = test_data_matrix, ynew = test_labels_vector, alpha = best_alpha, delta = best_delta) prediccion <- predict.rda(modelobueno, train_data_matrix, train_labels_vector, test_data_matrix, alpha = best_alpha, delta = best_delta, type = "class") accuracyrda <- mean(prediccion == test_labels_vector) cat("El modelo presenta una precisión del",round(accuracyrda*100,2), " %"),→ ## El modelo presenta una precisión del 22.22 % Este clasificador de RDA, atendiendo a lo obtenido anteriormente, presenta un resultado bastante malo, aproximadamente del 22 %. Este porcentaje tan bajo se debe al escaso número de observaciones para estimar la matriz de covarianzas, lo que provoca que las estimaciones no sean las adecuadas, y consecuentemente, las clasificaciones de la predicción tampoco. CAPÍTULO 3. IMPLEMENTACIÓN EN R. 57 3.4. CLASIFICADORES LINEALES CON REGULARIZACIÓN CUADRÁTICA 3.4.2. Regresión Logística con Regularización Cuadrática Instalamos y cargamos, en primer lugar, el paquete “glmnet” (véase Friedman y otros, [17]) para implementar esta técnica: library(glmnet) Realizamos las transformaciones necesarias de los datos y los dividimos en dos conjuntos de entrenamiento y test para llevar a cabo la clasificación y evaluar el rendimiento del modelo. # Convertir los datos a data frame y las etiquetas a formato numérico. nci_data <- data.frame(NCI60$data) nci_labels <- as.factor(NCI60$labs) # Verificamos la distribución de clases table(nci_labels) ## nci_labels ## BREAST CNS COLON K562A-repro K562B-repro LEUKEMIA ## 7 5 7 1 1 6 ## MCF7A-repro MCF7D-repro MELANOMA NSCLC OVARIAN PROSTATE ## 1 1 8 9 6 2 ## RENAL UNKNOWN ## 9 1 # Dividir los datos en conjunto de entrenamiento y prueba set.seed(123) train_indices <- sample(1:nrow(nci_data), size = 0.7 *nrow(nci_data)) train_data <- nci_data[train_indices, ] test_data <- nci_data[-train_indices, ] train_labels <- nci_labels[train_indices] test_labels <- nci_labels[-train_indices] A la hora de llevar a cabo la Regresión Logística con Regularización, la función “glmnet” requiere que cada clase se vea representada en, al menos, 2 observaciones en el conjunto de entrenamiento. Para ello, descarto aquellas clases que tienen menos de 2 observaciones en el conjunto de entrenamiento. # Verificar la distribución de clases en el conjunto de entrenamiento table(train_labels) ## train_labels ## BREAST CNS COLON K562A-repro K562B-repro LEUKEMIA ## 6 4 4 1 1 3 58 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.4. CLASIFICADORES LINEALES CON REGULARIZACIÓN CUADRÁTICA ## MCF7A-repro MCF7D-repro MELANOMA NSCLC OVARIAN PROSTATE ## 0 1 4 8 5 1 ## RENAL UNKNOWN ## 5 1 train_class_counts <- table(train_labels) valid_classes <- names(train_class_counts[train_class_counts >1]) train_data <- train_data[train_labels %in % valid_classes, ] train_labels <- train_labels[train_labels %in % valid_classes] Se convierten las variables etiquetas a formato numérico para poder aplicar la regresión logística multinomial, ya que “glmnet” lo requiere. Esto lo hacemos con “as.numeric”. # Verificar la nueva distribución de clases table(train_labels) ## train_labels ## BREAST CNS COLON K562A-repro K562B-repro LEUKEMIA ## 6 4 4 0 0 3 ## MCF7A-repro MCF7D-repro MELANOMA NSCLC OVARIAN PROSTATE ## 0 0 4 8 5 0 ## RENAL UNKNOWN ## 5 0 table(test_labels) ## test_labels ## BREAST CNS COLON K562A-repro K562B-repro LEUKEMIA ## 1 1 3 0 0 3 ## MCF7A-repro MCF7D-repro MELANOMA NSCLC OVARIAN PROSTATE ## 1 0 4 1 1 1 ## RENAL UNKNOWN ## 4 0 # Convertir las etiquetas en formato numérico para glmnet train_labels_numeric <- as.numeric(as.factor(train_labels)) test_labels_numeric <- as.numeric(as.factor(test_labels)) Por último, se utiliza “cv.glmnet” para realizar una validación cruzada y encontrar el mejor valor de lambda (parámetro de regularización). Se entrena finalmente el modelo con el mejor valor de lambda obtenido. CAPÍTULO 3. IMPLEMENTACIÓN EN R. 59 3.5. CLASIFICADORES LINEALES CON REGULARIZACIÓN L1 #Aplicar la validación cruzada con "cv.glmnet". Aplicamos alpha = 0 #para que la regularización sea Ridge. Si alpha = 1, #la regularización sería LASSO. modelo_cv <- cv.glmnet(as.matrix(train_data), train_labels_numeric, nfolds=10,family = "multinomial",alpha = 0) # Obtener el mejor lambda, o mejor dicho, el valor de lambda que # proporciona el menor error de validación cruzada. best_lambda <- modelo_cv$lambda.min best_lambda ## [1] 2.911386 # Entrenar el modelo final con el mejor lambda final_model <- glmnet(as.matrix(train_data), train_labels_numeric, family = "multinomial", lambda = best_lambda, alpha = 0) Una vez que se tiene el modelo construido, lo evaluamos realizando predicciones sobre él utilizando el conjunto de test y comparándolo con los valores originales. predictions <- predict(final_model, as.matrix(test_data), type= "class") # Evaluar la precisión del modelo accuracy <- mean(predictions == test_labels_numeric) print(paste("Precisión del modelo de Regresión Logística con", "Regularización:",round(accuracy *100,2), " %")) ## [1] "Precisión del modelo de Regresión Logística con Regularización: 70 %" La predicción, al igual que anteriormente en LDLA, no es lo suficientemente buena debido a que hay clases del conjunto de datos que no se ven representadas a la hora de entrenar el modelo, pues se dispone de un número muy bajo de observaciones. Sin embargo, parece que mejora ligeramente la precisión, alcanzándose un 70%. 3.5. Clasificadores Lineales con Regularización L1 3.5.1. Regularización LASSO y Elastic Net Cargamos nuevamente la librería “glmnet” (véase Friedman y otros, [17]), que utilizaremos para llevar a cabo la clasificación LASSO y Elastic Net, así como los datos que utilizaremos para los clasificadores. 60 CAPÍTULO 3. IMPLEMENTACIÓN EN R. 3.5. CLASIFICADORES LINEALES CON REGULARIZACIÓN L1 library(glmnet) # Cargar los datos del paquete ISLR library(ISLR) data("NCI60") Transformamos tanto el conjunto de datos como las etiquetas a un dataframe y a una variable de tipo factor, respectivamente. # Convertir los datos a un data frame nci_data <- data.frame(NCI60$data) nci_labels <- as.factor(NCI60$labs) Dividimos el conjunto de datos en dos conjuntos de entrenamiento y test para construir los clasificadores. Posteriormente, y como hemos hecho con los demás clasificadores, filtramos únicamente aquellas clases con 2 o más observaciones para construir el clasificador. # Dividir los datos en conjunto de entrenamiento y prueba set.seed(123) train_indices <- sample(1:nrow(nci_data), size = 0.7 *nrow(nci_data)) train_data <- nci_data[train_indices, ] test_data <- nci_data[-train_indices, ] train_labels <- nci_labels[train_indices] test_labels <- nci_labels[-train_indices] # Filtrar clases con al menos 2 observaciones en el conjunto # de entrenamiento class_counts <- table(train_labels) valid_classes <- names(class_counts[class_counts >= 2]) train_data <- train_data[train_labels %in % valid_classes, ] train_labels <- train_labels[train_labels %in % valid_classes] train_labels <- as.numeric(as.factor(train_labels)) test_labels <- as.numeric(as.factor(test_labels)) Ajustamos el modelo LASSO haciendo uso de la función “cv.glmnet” y evaluamos la precisión de este realizando predicciones sobre él con el conjunto de test anteriormente creado: lasso_model <- cv.glmnet(as.matrix(train_data), train_labels, family = "multinomial", alpha = 1) # Mejor lambda best_lambda_lasso <- lasso_model$lambda.min CAPÍTULO 3. IMPLEMENTACIÓN EN R. 61 3.5. CLASIFICADORES LINEALES CON REGULARIZACIÓN L1 # Hacer predicciones en el conjunto de test lasso_predictions <- predict(lasso_model, newx = as.matrix(test_data), s = best_lambda_lasso, type = "class") # Evaluar el desempeño del modelo lasso_accuracy <- mean(lasso_predictions == test_labels) cat("Precisión del modelo LASSO en el conjunto de prueba:", lasso_accuracy, "\n") ## Precisión del modelo LASSO en el conjunto de prueba: 0.6 Hacemos lo mismo que para el clasificador LASSO, pero en este caso, para Elastic Net, combinando así la penalización Ridge y la penalización LASSO. Ajustamos el modelo eligiendo, mediante validación cruzada, el lambda que minimiza el error de clasificación. # Ajustar el modelo Elastic Net elastic_net_model <- cv.glmnet(as.matrix(train_data), train_labels, family = "multinomial", alpha = 0.5) # Mejor lambda best_lambda_elastic_net <- elastic_net_model$lambda.min Una vez tenemos el modelo construido, procedemos como anteriormente y evaluamos la precisión de este realizando predicciones sobre el conjunto de test. # Hacer predicciones en el conjunto de test elastic_net_predictions <- predict(elastic_net_model, newx = as.matrix(test_data), s = best_lambda_elastic_net, type = "class") # Evaluar el desempeño del modelo elastic_net_accuracy <- mean(elastic_net_predictions == test_labels) cat("Precisión del modelo Elastic Net en el conjunto de prueba:",+ elastic_net_accuracy, "\n") ## Precisión del modelo Elastic Net en el conjunto de prueba: 0.7 Vemos que la precisión es mejor en el modelo con regularización Elastic Net, pues intenta encontrar la mejor combinación entre LASSO y Ridge, aumentando así la precisión que se conseguiría con LASSO por separado. 62 CAPÍTULO 3. IMPLEMENTACIÓN EN R. Bibliografía [1] Efron,B.;Hastie,T.;Johnstone,I. yTibshirani,R. (2004). Least angle regression. Ann. Statist. [2] Friedman,J.;Hastie,T. yTibshirani,R. (2010). Regularization Paths for Generalized Linear Models via Coordinate Descent. Journal of Statistical Software. [3] Friedman,J.H. (1989). «Regularized Discriminant Analysis». Journal of the American Statistical Association,84(405), pp. 165–175. doi: 10.1080/01621459.1989. 10478752. [4] Grupo Innovación Docente: Estadística en Ciencias Sociales. Universidad de Murcia. (Consultado en Mayo, 2024). «Análisis de Componentes Principales.». https://www.um.es/estadempresa/multivarianteR/ACP/resumen_ACP.html. [5] Guo,Y.;Hastie,T. yTibshirani,R. (2023). rda: Shrunken Centroids Regularized Discriminant Analysis. https://CRAN.R-project.org/package=rda. R package version 1.2-1. [6] Hastie,T.;Tibshirani,R. yFriedman,J. H. (2009). The elements of Statistical Learning. Springer, 2.aedición. [7] Hastie, T., Tibshirani, R., Friedman, J. H. (Página web). «Datasets for "The Elements of Statistical Learning"». Disponible en https://hastie.su.domains/ ElemStatLearn/datasets/nci.info.txt. [8] Hoerl,A.E. yKennard,R. (1970). Ridge regression: Biased estimation for nonorthogonal problems. Technometrics. [9] Izenman,A. J. (2013). Modern Multivariate Statistical Techniques Regression, Classification and Manifold Learning. Springer. [10] James,G.;Witten,D.;Hastie,T. yTibshirani,R. (2021). ISLR: Data for an Introduction to Statistical Learning with Applications in R. https://CRAN.R-project.org/package=ISLR. R package version 1.4. [11] Kuhn,M. (2022). caret: Classification and Regression Training. https://CRAN.R-project.org/package=caret. R package version 6.0-92. [12] Kuhn,M. yRamey,J. (2021). sparsediscrim: Sparse and Regularized Discriminant Analysis. https://CRAN.R-project.org/package=sparsediscrim. R package version 0.3.0. 63 Bibliografía [13] Madsen,R.E.;Hansen,L.K. yWinter,0. (February 2004). Singular Value Decomposition and Principal Component Analysis. [14] Posit team (2024). RStudio: Integrated Development Environment for R. Posit Software, PBC, Boston, MA. http://www.posit.co/. [15] R Core Team (2021). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/. [16] Seber,G. A. (1984). Multivariate Observations. John Wiley and Sons. [17] Simon,N.;Friedman,J.;Hastie,T. yTibshirani,R. (2011). «Regularization Paths for Cox’s Proportional Hazards Model via Coordinate Descent». Journal of Statistical Software,39(5), pp. 1–13. doi: 10.18637/jss.v039.i05. https://www.jstatsoft.org/v39/i05/. [18] Tibshirani,R. (1996). Regression shrinkage and selection via the lasso. J. R. Statist. Soc.. [19] —— (2011). Regression shrinkage and selection via the lasso: a retrospective. J. R. Statist. Society : Series B (Statistical Methodology). [20] Tibshirani,R.;Hastie,T.;Narasimhan,B. yChu,G. (2002). Diagnosis of multiple cancer types by shrunken centroids of gene expression. Proceedings of the National Academy of Sciences - PNAS, 1.aedición. [21] Uriel,J., E. y Aldás (2005). Análisis Multivariante Aplicado. Thomson Editores Spain. [22] Wickham,H. (2016). ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York. ISBN 978-3-319-24277-4. https://ggplot2.tidyverse.org. [23] Wikipedia (Consultado en Enero, 2024). «Análisis Discriminante Lineal.». https://es.wikipedia.org/wiki/AnÃąlisis_discriminante_lineal. [24] Zou,H. yHastie,T. (2005). Regularization and variable selection via the elastic net. Journal of the Royal Statistical Society Series B 67 (2). 64