Membangun plot skor menggunakan komponen utama

Aug 18 2020

Saya mencoba membuat plot skor dari dua komponen utama pertama. Saya mulai dengan membagi data menjadi tiga bingkai data berdasarkan class. Saya kemudian mengubah data dan melakukan PCA.

Data saya adalah sebagai berikut:

14      1   82.0 12.80   7.60   1070   105   400
14      1   82.0 11.00   9.00    830   145   402
14      1  223.6 17.90  10.35   2200   135   500
15      1  164.0 14.50   9.80   1946   138   500
15      1  119.0 12.90   7.90   1190   140   400
15      1   74.5  7.50   6.30    653   177   350
15      1   74.5 11.13   8.28    930   113   402
16      1  279.5 14.30   9.40   1575   230   700
16      1   82.0  7.80   6.70    676   175   525
16      1   67.0 11.00   8.30    920   106   300
16      2  112.0 11.70   8.00   1353   140   560
16      2  149.0 12.80   8.70   1550   170   550
16      2  119.0  8.50   7.40    888   175   250
16      2  119.0 13.30   9.60   1275   157   450
16      2  238.5 14.90   8.90   1537   183   700
16      2  205.0 12.00   7.90   1292   201   600
16      2   82.0  9.40   6.20    611   209   175
16      2  119.0 15.95  10.25   1350   145   450
16      2  194.0 16.74  10.77   1700   120   450
17      2  336.0 22.20  10.90   3312   135   450
17      3  558.9 23.40  12.60   4920   152   600
17      3  287.0 14.30   9.40   1510   176   800
17      3  388.0 23.72  11.86   3625   140   500
17      3  164.0 11.90   9.80    900   190   600
17      3  194.0 14.40   9.20   1665   175   600
17      3  194.0 14.40   8.90   1640   175   600
17      3  186.3  9.70   8.00   1081   205   600
17      3  119.0  8.00   6.50    625   196   400
17      3  119.0  9.40   6.95    932   165   250
17      3   89.4 14.55   9.83   1378   146   400

Kolom 1 type:, Kolom 2 class:, Kolom 3 v1:, Kolom 4 v2:, Kolom 5 v3:, Kolom 6 v4:, Kolom 7 v5:, Kolom-8:v6

Kode saya adalah sebagai berikut:

data <- read.csv("data.csv")
result <- split(data, data$class);

data1 <- result[[1]][,3:8];
data1Logged <- log10(data1)
pca.data1Logged = prcomp( ~ v1 + 
                         v2 + 
                         v3 + 
                         v4 + 
                         v5 + 
                         v6, 
                       data = data1Logged, scale. = FALSE );

data2 <- result[[2]][,3:8];
data2Logged <- log10(data2)
pca.data2Logged = prcomp( ~ v1 + 
                         v2 + 
                         v3 + 
                         v4 + 
                         v5 + 
                         v6, 
                       data = data2Logged, scale. = FALSE );

data3 <- result[[3]][,3:8];
data3Logged <- log10(data3)
pca.data3Logged = prcomp( ~ v1 + 
                         v2 + 
                         v3 + 
                         v4 + 
                         v5 + 
                         v6, 
                       data = data3Logged, scale. = FALSE );

Untuk masing-masing dari ketiganya class, saya ingin memiliki plot skor untuk PC1 dan PC2:

pca.data1Logged$x[,1:2]
pca.data2Logged$x[,1:2] pca.data3Logged$x[,1:2]

Ini yang terbaik yang bisa saya ketahui:

opar <- par(mfrow = c(1,3))
plot(pca.data1Logged$x[,1:2]) plot(pca.data2Logged$x[,1:2])
plot(pca.data3Logged$x[,1:2])
par(opar)

Tapi saya ingin plot ini diskalakan, diwarnai, ditumpangkan, dll. Saya sudah mulai membaca tentang ggplot, tetapi saya tidak memiliki pengalaman untuk melakukan ini. Saya ingin yang berikut ini:

https://cran.r-project.org/web/packages/ggfortify/vignettes/plot_pca.html

Masalah di atas adalah bahwa saya telah memecah data menjadi 3 bingkai data terpisah, jadi tidak ada judul untuk "class1", "class2," class3 ".

Jawaban

3 BappaDas Aug 18 2020 at 12:51

Anda dapat menggunakan factoextradan FactoMineRmenyukai

library("factoextra")
library("FactoMineR")

#PCA analysis
df.pca <- PCA(df[,-c(1,2)], graph = T)
# Visualize
# Use habillage to specify groups for coloring
fviz_pca_ind(df.pca,
             label = "none", # hide individual labels
             habillage = as.factor(df$class), # color by groups
             palette = c("#00AFBB", "#E7B800", "#FC4E07"),
             addEllipses = TRUE # Concentration ellipses, legend.title = "Class")

Anda dapat mengubah Dim1 dan 2 ke PC1 dan 2 secara manual. Untuk itu, Anda dapat mencatat nilai "Dim1 (63.9%)" dan "Dim2 (23.3%)" dari plot ini dan gunakan kode berikut untuk mengubah Dim1 dan 2 menjadi PC1 dan 2 sejenisnya

fviz_pca_ind(df.pca,
             label = "none", # hide individual labels
             habillage = as.factor(df$class), # color by groups
             palette = c("#00AFBB", "#E7B800", "#FC4E07"),
             addEllipses = TRUE, # Concentration ellipses
             xlab = "PC1 (63.9%)", ylab = "PC2 (23.3%)", legend.title = "Class")

Jika Anda ingin log mengubah data, maka Anda dapat menggunakan

df[,3:8] <- log10(df[,3:8]) 

df.pca <- PCA(df, graph = T)

fviz_pca_ind(df.pca,
             label = "none", # hide individual labels
             habillage = as.factor(df$class), # color by groups
             palette = c("#00AFBB", "#E7B800", "#FC4E07"),
             addEllipses = TRUE, # Concentration ellipses
legend.title = "Class")

Untuk mengubah Dim1 dan 2 ke PC1 dan 2 secara manual, Anda dapat menggunakan kode berikut

fviz_pca_ind(df.pca,
             label = "none", # hide individual labels
             habillage = as.factor(df$class), # color by groups
             palette = c("#00AFBB", "#E7B800", "#FC4E07"),
             addEllipses = TRUE, # Concentration ellipses
             xlab = "PC1 (64.9%)", ylab = "PC2 (22.6%)", legend.title = "Class")

Data

df =
structure(list(Type = c(14L, 14L, 14L, 15L, 15L, 15L, 15L, 16L, 
16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 17L, 17L, 
17L, 17L, 17L, 17L, 17L, 17L, 17L, 17L, 17L), class = c(1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), v1 = c(82, 82, 
223.6, 164, 119, 74.5, 74.5, 279.5, 82, 67, 112, 149, 119, 119, 
238.5, 205, 82, 119, 194, 336, 558.9, 287, 388, 164, 194, 194, 
186.3, 119, 119, 89.4), v2 = c(12.8, 11, 17.9, 14.5, 12.9, 7.5, 
11.13, 14.3, 7.8, 11, 11.7, 12.8, 8.5, 13.3, 14.9, 12, 9.4, 15.95, 
16.74, 22.2, 23.4, 14.3, 23.72, 11.9, 14.4, 14.4, 9.7, 8, 9.4, 
14.55), v3 = c(7.6, 9, 10.35, 9.8, 7.9, 6.3, 8.28, 9.4, 6.7, 
8.3, 8, 8.7, 7.4, 9.6, 8.9, 7.9, 6.2, 10.25, 10.77, 10.9, 12.6, 
9.4, 11.86, 9.8, 9.2, 8.9, 8, 6.5, 6.95, 9.83), v4 = c(1070L, 
830L, 2200L, 1946L, 1190L, 653L, 930L, 1575L, 676L, 920L, 1353L, 
1550L, 888L, 1275L, 1537L, 1292L, 611L, 1350L, 1700L, 3312L, 
4920L, 1510L, 3625L, 900L, 1665L, 1640L, 1081L, 625L, 932L, 1378L
), v5 = c(105L, 145L, 135L, 138L, 140L, 177L, 113L, 230L, 175L, 
106L, 140L, 170L, 175L, 157L, 183L, 201L, 209L, 145L, 120L, 135L, 
152L, 176L, 140L, 190L, 175L, 175L, 205L, 196L, 165L, 146L), 
    v6 = c(400L, 402L, 500L, 500L, 400L, 350L, 402L, 700L, 525L, 
    300L, 560L, 550L, 250L, 450L, 700L, 600L, 175L, 450L, 450L, 
    450L, 600L, 800L, 500L, 600L, 600L, 600L, 600L, 400L, 250L, 
    400L)), class = "data.frame", row.names = c(NA, -30L))
2 jay.sf Aug 18 2020 at 12:53

Anda dapat mengikat hasil terpisah dan menambahkan kolom warna yang Anda gunakan plot.

rb <- rbind(cbind(pca.data1Logged$x[,1:2], d=2), cbind(pca.data2Logged$x[,1:2], d=3),
            cbind(pca.data3Logged$x[,1:2], d=4))

plot(rb, col=rb[,"d"], pch=20, main="PCA Plot")
legend("bottomleft", paste("data", 1:3), col=2:4, pch=20)


Data:

data <- read.table(header=F, text="14      1   82.0 12.80   7.60   1070   105   400
14      1   82.0 11.00   9.00    830   145   402
14      1  223.6 17.90  10.35   2200   135   500
15      1  164.0 14.50   9.80   1946   138   500
15      1  119.0 12.90   7.90   1190   140   400
15      1   74.5  7.50   6.30    653   177   350
15      1   74.5 11.13   8.28    930   113   402
16      1  279.5 14.30   9.40   1575   230   700
16      1   82.0  7.80   6.70    676   175   525
16      1   67.0 11.00   8.30    920   106   300
16      2  112.0 11.70   8.00   1353   140   560
16      2  149.0 12.80   8.70   1550   170   550
16      2  119.0  8.50   7.40    888   175   250
16      2  119.0 13.30   9.60   1275   157   450
16      2  238.5 14.90   8.90   1537   183   700
16      2  205.0 12.00   7.90   1292   201   600
16      2   82.0  9.40   6.20    611   209   175
16      2  119.0 15.95  10.25   1350   145   450
16      2  194.0 16.74  10.77   1700   120   450
17      2  336.0 22.20  10.90   3312   135   450
17      3  558.9 23.40  12.60   4920   152   600
17      3  287.0 14.30   9.40   1510   176   800
17      3  388.0 23.72  11.86   3625   140   500
17      3  164.0 11.90   9.80    900   190   600
17      3  194.0 14.40   9.20   1665   175   600
17      3  194.0 14.40   8.90   1640   175   600
17      3  186.3  9.70   8.00   1081   205   600
17      3  119.0  8.00   6.50    625   196   400
17      3  119.0  9.40   6.95    932   165   250
17      3   89.4 14.55   9.83   1378   146   400")

names(data) <- c("sth", "class", paste0("v", 1:6))