install.
packages("ggplot2")
library(ggplot2)
set.seed(435)
X_1 <- sample(x = 10:40, size = 100, replace = TRUE)
X_1
X_2 <- 2.5 * X_1 + 10
X_2
X_2 <- X_2 + rnorm(n = 100, mean = 10, sd = 30)
datos <- data.frame(X_1, X_2)
ggplot(data = datos, aes(x = X_1, y = X_2)) +
geom_point() +
geom_segment(aes(x = 10, y = 2.5 * 10 + 10, xend = 40, yend = 2.5 * 40
+ 10),
colour = "firebrick", arrow = arrow(ends = "both")) +
theme_bw()
######## Ejemplo cálculo eigenvectors y eigenvalues
datos <- data.frame(X1 = c(2.5, 0.5, 2.2, 1.9, 3.1, 2.3, 2, 1, 1.5,
1.1),
X2 = c(2.4, 0.7, 2.9, 2.2, 3, 2.7, 1.6, 1.1, 1.6,
0.9))
datos
###En primer lugar se resta a cada valor la media de la variable a la
que
##pertenece. Con esto se consigue centralizar las variables y que su
media sea 0.
datos_centrados <- datos
datos_centrados$X1 <- datos$X1 - mean(datos$X1)
datos_centrados$X2 <- datos$X2 - mean(datos$X2)
datos_centrados
##se calcula la matriz de correlaciones
matriz_cov <- cov(datos_centrados)
matriz_cov
## Dado que la matriz de covarianzas es cuadrada, se pueden obtener
sus
## correspondientes eigenvectors y eigenvalues.
eigen <- eigen(matriz_cov)
eigen$values
eigen$vectors
##Los eigenvectors ordenados de mayor a menor eigenvalues se
corresponden con
##las componentes principales.
#Una vez obtenidos los eigenvectors (componentes principales) se
calcula el valor
#que toma cada componente para cada observación en función de las
variables
#originales (principal component scores). Para ello, simplemente se
tienen
#que multiplicar los eigenvectors transpuestos por los datos
originales centrados
#y también transpuestos.
t_eigenvectors <- t(eigen$vectors)
t_eigenvectors
t_datos_centrados <- t(datos_centrados)
t_datos_centrados
#Producto matricial
pc_scores <- t_eigenvectors %*% t_datos_centrados
rownames(pc_scores) <- c("PC1", "PC2")
# Se vuelve a transponer para que los datos estén en modo tabla
t(pc_scores)
## realizar comentariop12
datos_recuperados <- t(eigen$vectors %*% pc_scores)
datos_recuperados[, 1] <- datos_recuperados[, 1] + mean(datos$X1)
datos_recuperados[, 2] <- datos_recuperados[, 2] + mean(datos$X2)
datos_recuperados
##Ejemplo cálculo directo de PCA con R
##comentar sobre los datos-p13
data("USArrests")
data()# para ver base de datos del pquete basico de R
head(USArrests) ## realizar comentarios
apply(X = USArrests, MARGIN = 2, FUN = mean) ##realizar comentarios
apply(X = USArrests, MARGIN = 2, FUN = var) ##realizar comentarios
##La función prcomp() es una de las múltiples funciones en R que
##realizan PCA.
pca <- prcomp(USArrests, scale = TRUE)
names(pca) ##realizar comentarios p14
pca$center
pca$scale
##rotation contiene el valor de los loadings phi para cada componente
(eigenvector).
pca$rotation
##primera componente
##pc1 = -0.5358995
Mueder-0.5831836Assault-0.2781909UrbanPop-0.5434321Rape
#La función prcomp() calcula automáticamente el valor de las
componentes principales
#para cada observación (principal component scores) multiplicando los
datos por los
#vectores de loadings. El resultado se almacena en la matriz x.
head(pca$x)
dim(pca$x)
#Mediante la función biplot() se puede obtener una representación
bidimensional de las
#dos primeras componentes.
biplot(x = pca, scale = 0, cex = 0.8, col = c("blue4", "brown3"))
#La imagen especular, cuya interpretación es equivalente, se puede
obtener invirtiendo
#el signo de los loadings y de los principal component scores.
pca$rotation <- -pca$rotation
pca$x <- -pca$x
biplot(x = pca, scale = 0, cex = 0.8, col = c("blue4", "brown3"))
## Una vez calculadas las componentes principales, se puede conocer la
varianza explicada
## por cada una de ellas, la proporción respecto al total y la
proporción de varianza acumulada.
library(ggplot2)
pca$sdev^2
prop_varianza <- pca$sdev^2/sum(pca$sdev^2)
prop_varianza
ggplot(data = data.frame(prop_varianza, pc = 1:4),aes(x = pc, y =
prop_varianza)) +
geom_col(width = 0.3) +
scale_y_continuous(limits = c(0, 1)) +
theme_bw() +
labs(x = "Componente principal", y = "Proporción de varianza
explicada")
prop_varianza_acum <- cumsum(prop_varianza)
prop_varianza_acum
##p18
ggplot(data = data.frame(prop_varianza_acum, pc = 1:4),
aes(x = pc, y = prop_varianza_acum, group = 1)) +
geom_point() +
geom_line() + theme_bw() +
labs(x = "Componente principal", y = "Proporción de varianza explicada
acumulada")
######################################################################
########################
#p35 Algoritmo de t-SNE
#Ejemplo con tsne
library(readr)
library(dplyr)
# Carga de datos
datos <- read_csv(paste0("http://archive.ics.uci.edu/ml/machine-
learning-",
"databases/optdigits/optdigits.tra"))
##datos
dim(datos)
#str(datos)
# La última columna contiene el número real al que se corresponde la
observación. # Se renombra como "numero"
datos <- datos %>% rename(numero = `0_26`)
# La función tsne() recibe como argumento una matriz, no un
data.frames
datos <- data.matrix(datos)
# Debido a los requerimientos computacionales del t-SNE, se limita
este ejemplo # únicamente a 1000 observaciones.
datos <- datos[1:1000,]
##p39
install.packages("tsne")
library(tsne)
library(ggplot2)
# También se limita el número de iteraciones (epoch) a 100, aunque los
# resultados podrían mejorar si se aumentara
set.seed(321)
tsne_reduction <- tsne(datos, k = 2, perplexity = 30, epoch = 100)
# Para poder representar el verdadero número al que corresponde cada
imagen,
# se adjunta la variable "numero" del set de datos
resultados <- as.data.frame(tsne_reduction)
colnames(resultados) <- c("dim_1", "dim_2")
resultados$numero <- as.character(datos[ ,"numero"])
ggplot(data = resultados, aes(x = dim_1, y = dim_2)) +
geom_point(aes(color = numero)) +
theme_bw()
##p Ejemplo con R tsne
library(readr)
library(dplyr)
# Carga de datos
datos <- read_csv(paste0("http://archive.ics.uci.edu/ml/machine-
learning-", "databases/optdigits/
optdigits.tra"))
dim(datos)
# La última columna contiene el número real al que se corresponde la
observación.
# Se renombra como "numero"
datos <- datos %>% rename(numero = `0_26`)
install.packages("Rtsne")
library(Rtsne)
library(ggplot2)
tsne <- Rtsne(X = datos[, -65], is_distance = FALSE, dims = 2,
perplexity = 30,
theta = 0.5, max_iter = 500)
# El objeto devuelto por Rtsne() almacena los valores de las
dimensiones en el
# elemento Y. Como en este caso se ha especificado que la reducción se
haga
# a dos dimensiones (k=2), Y tiene solo dos columnas.
head(tsne$Y)
# Para poder representar el verdadero número al que corresponde cada
imagen,
# se adjunta la variable "numero" del set de datos
resultados <- as.data.frame(tsne$Y)
colnames(resultados) <- c("dim_1", "dim_2")
resultados$numero <- as.character(datos$numero)
ggplot(data = resultados, aes(x = dim_1, y = dim_2)) +
geom_point(aes(color = numero)) +
theme_bw()
pca <- prcomp(x = datos[,-65])
resultados <- as.data.frame(pca$x[, 1:2])
resultados$numero <- as.character(datos$numero)
ggplot(data = resultados, aes(x = PC1, y = PC2)) +
geom_point(aes(color = numero)) +
theme_bw()
#Se reduce de nuevo la dimensionalidad, pero esta vez empleando PCA y
representando
#las dos primeras componentes
pca <- prcomp(x = datos[,-65])
resultados <- as.data.frame(pca$x[, 1:2])
resultados$numero <- as.character(datos$numero)
ggplot(data = resultados, aes(x = PC1, y = PC2)) +
geom_point(aes(color = numero)) +
theme_bw()
#Véase ahora la reducción a un espacio de 3 dimensiones
library(scatterplot3d)
library(RColorBrewer)
tsne <- Rtsne(X = datos[, -65], is_distance = FALSE, dims = 3,
perplexity = 30,
theta = 0.5, max_iter = 500)
resultados <- as.data.frame(tsne$Y)
colnames(resultados) <- c("dim_1", "dim_2", "dim_3")
resultados$numero <- as.factor(datos$numero)
colores <- brewer.pal(n = 10, name = "Set3")
colores <- colores[as.numeric(resultados$numero)]
scatterplot3d(x = resultados$dim_1,
y = resultados$dim_2,
z = resultados$dim_3,
pch = 20, color = colores, cex.lab = 0.8,
grid = TRUE, box = FALSE)
legend("bottom", legend = levels(resultados$numero),
col = colores, pch = 16,
inset = -0.23, xpd = TRUE, horiz = TRUE)
