2011-07-05 13 views
23

Me pregunto si es posible trazar los resultados de pca biplot con ggplot2. Supongamos que si quiero mostrar los siguientes resultados de biplot con ggplot2Trazando pca biplot con ggplot2

fit <- princomp(USArrests, cor=TRUE) 
summary(fit) 
biplot(fit) 

Cualquier ayuda será muy apreciada. Gracias

+2

[Este] (http://groups.google.com/group/ ggplot2/browse_thread/thread/5fea365578c3910f/47a e63e7ff18508e) en la lista de distribución de ggplot2 podría ser un buen lugar para comenzar. – joran

+0

Recomiendo en su lugar aceptar la respuesta de MYaseen208 sobre el paquete 'ggbiplot'. Empecé a ajustar la respuesta de crayola (que es genial, pero innecesario dado el paquete) para hacer cosas ya disponibles en 'ggbiplot' (por ejemplo, eliminar etiquetas). –

Respuesta

41

Tal vez esto help-- que es una adaptación de código que escribí hace algún tiempo. Ahora dibuja flechas también.

PCbiplot <- function(PC, x="PC1", y="PC2") { 
    # PC being a prcomp object 
    data <- data.frame(obsnames=row.names(PC$x), PC$x) 
    plot <- ggplot(data, aes_string(x=x, y=y)) + geom_text(alpha=.4, size=3, aes(label=obsnames)) 
    plot <- plot + geom_hline(aes(0), size=.2) + geom_vline(aes(0), size=.2) 
    datapc <- data.frame(varnames=rownames(PC$rotation), PC$rotation) 
    mult <- min(
     (max(data[,y]) - min(data[,y])/(max(datapc[,y])-min(datapc[,y]))), 
     (max(data[,x]) - min(data[,x])/(max(datapc[,x])-min(datapc[,x]))) 
     ) 
    datapc <- transform(datapc, 
      v1 = .7 * mult * (get(x)), 
      v2 = .7 * mult * (get(y)) 
      ) 
    plot <- plot + coord_equal() + geom_text(data=datapc, aes(x=v1, y=v2, label=varnames), size = 5, vjust=1, color="red") 
    plot <- plot + geom_segment(data=datapc, aes(x=0, y=0, xend=v1, yend=v2), arrow=arrow(length=unit(0.2,"cm")), alpha=0.75, color="red") 
    plot 
} 

fit <- prcomp(USArrests, scale=T) 
PCbiplot(fit) 

Es posible que desee cambiar el tamaño del texto, así como la transparencia y los colores, al gusto; sería fácil convertirlos en parámetros de la función. Nota: se me ocurrió que esto funciona con prcomp pero su ejemplo es con princomp. Puede, de nuevo, necesitar adaptar el código en consecuencia. Nota 2: el código para geom_segment() se toma prestado de la publicación de la lista de correo vinculada desde el comentario a OP.

PC biplot

+0

Me gustaría agregar nombres de observaciones y flechas para las variables. ¿Alguna idea? – MYaseen208

+0

Listo - ¡espero que ayude! – crayola

+1

Pequeña actualización para la versión 0.9 de ggplot2, ahora necesita agregar library ("ggplot2") y library ("grid") para trazar flechas. –

4

Esto hará que los estados trazada, aunque no las variables de

fit.df <- as.data.frame(fit$scores) 
fit.df$state <- rownames(fit.df) 

library(ggplot2) 
ggplot(data=fit.df,aes(x=Comp.1,y=Comp.2))+ 
    geom_text(aes(label=state,size=1,hjust=0,vjust=0)) 

enter image description here

+2

Buen intento. ¿Cómo agregar variables con flechas? – MYaseen208

+0

@Henry ¿Alguna solución similar para pls biplot? http://stackoverflow.com/questions/39137287/plotting-pls-biplot-with-ggplot2 – aelwan

7

Si utiliza la excelente FactoMineR paquete para la PCA, es posible encontrar esto útil para la toma de parcelas con ggplot2

# Plotting the output of FactoMineR's PCA using ggplot2 
# 
# load libraries 
library(FactoMineR) 
library(ggplot2) 
library(scales) 
library(grid) 
library(plyr) 
library(gridExtra) 
# 
# start with a clean slate 
rm(list=ls(all=TRUE)) 
# 
# load example data from the FactoMineR package 
data(decathlon) 
# 
# compute PCA 
res.pca <- PCA(decathlon, quanti.sup = 11:12, quali.sup=13, graph = FALSE) 
# 
# extract some parts for plotting 
PC1 <- res.pca$ind$coord[,1] 
PC2 <- res.pca$ind$coord[,2] 
labs <- rownames(res.pca$ind$coord) 
PCs <- data.frame(cbind(PC1,PC2)) 
rownames(PCs) <- labs 
# 
# Just showing the individual samples... 
ggplot(PCs, aes(PC1,PC2, label=rownames(PCs))) + 
    geom_text() 
# 
# Now get supplementary categorical variables 
cPC1 <- res.pca$quali.sup$coor[,1] 
cPC2 <- res.pca$quali.sup$coor[,2] 
clabs <- rownames(res.pca$quali.sup$coor) 
cPCs <- data.frame(cbind(cPC1,cPC2)) 
rownames(cPCs) <- clabs 
colnames(cPCs) <- colnames(PCs) 
# 
# Put samples and categorical variables (ie. grouping 
# of samples) all together 
p <- ggplot() + opts(aspect.ratio=1) + theme_bw(base_size = 20) 
# no data so there's nothing to plot... 
# add on data 
p <- p + geom_text(data=PCs, aes(x=PC1,y=PC2,label=rownames(PCs)), size=4) 
p <- p + geom_text(data=cPCs, aes(x=cPC1,y=cPC2,label=rownames(cPCs)),size=10) 
p # show plot with both layers 
# 
# clear the plot 
dev.off() 
# 
# Now extract variables 
# 
vPC1 <- res.pca$var$coord[,1] 
vPC2 <- res.pca$var$coord[,2] 
vlabs <- rownames(res.pca$var$coord) 
vPCs <- data.frame(cbind(vPC1,vPC2)) 
rownames(vPCs) <- vlabs 
colnames(vPCs) <- colnames(PCs) 
# 
# and plot them 
# 
pv <- ggplot() + opts(aspect.ratio=1) + theme_bw(base_size = 20) 
# no data so there's nothing to plot 
# put a faint circle there, as is customary 
angle <- seq(-pi, pi, length = 50) 
df <- data.frame(x = sin(angle), y = cos(angle)) 
pv <- pv + geom_path(aes(x, y), data = df, colour="grey70") 
# 
# add on arrows and variable labels 
pv <- pv + geom_text(data=vPCs, aes(x=vPC1,y=vPC2,label=rownames(vPCs)), size=4) + xlab("PC1") + ylab("PC2") 
pv <- pv + geom_segment(data=vPCs, aes(x = 0, y = 0, xend = vPC1*0.9, yend = vPC2*0.9), arrow = arrow(length = unit(1/2, 'picas')), color = "grey30") 
pv # show plot 
# 
# clear the plot 
dev.off() 
# 
# Now put them side by side 
# 
library(gridExtra) 
grid.arrange(p,pv,nrow=1) 
# 
# Now they can be saved or exported... 
# 
# tidy up by deleting the plots 
# 
dev.off() 

Y esto es lo que las parcelas finales parece, tal vez el tamaño del texto en la trama izquierda podría ser un poco más pequeño:

enter image description here

16

Aquí es la forma más sencilla a través ggbiplot:

library(ggbiplot) 
fit <- princomp(USArrests, cor=TRUE) 
biplot(fit) 

enter image description here

ggbiplot(fit, labels = rownames(USArrests)) 

enter image description here

+3

Dado que esto no está en CRAN, así es como se obtiene el paquete: 'library (devtools); install_github ("vqv/ggbiplot") '. Esta es definitivamente la mejor respuesta; Me pregunto si podría ser oscurecido por el feo 'biplot' inicial. Esto es lo que vi por primera vez en una pantalla pequeña, casi lo ignoro antes de desplazarme hacia 'ggbiplot'. –

5

también puede utilizar factoextra que también tiene un backend ggplot2:

library("devtools") 
install_github("kassambara/factoextra") 
fit <- princomp(USArrests, cor=TRUE) 
fviz_pca_biplot(fit) 

enter image description here

O ggord:

install_github('fawda123/ggord') 
library(ggord) 
ggord(fit)+theme_grey() 

enter image description here

O ggfortify:

devtools::install_github("sinhrks/ggfortify") 
library(ggfortify) 
ggplot2::autoplot(fit, label = TRUE, loadings.label = TRUE) 

enter image description here