Illustrate a 3D visualisation of the three main PCs using plot3d() package in R?

Question

I have a dataset (view/download) here. I need to run PCA on this data set and need to illustrate a 3D visualisation of the three main PCs using plot3d() package.

I have looked for more information before coming here on which package to use, some people use princomp and some us prcomp (in university we have personally been using prcomp). The confusing part is were to move forward with this?

The full question is "Illustrate a 3D visualisation of the three main PCs using plot3d() package in R. Use the following colours for samples (data points) belonging to the various immune subgroups/subtypes (C1: red, C2: yellow, C3: green, C4: cyan, C5: blue, C6: purple)"

After reading in the CSV and loading rgl library

In the current state I have removed the rows with the subgroups: see code below

# Remove last row from dataset
dataMinusLastRow = data[1:(nrow(data)-1),]

# Convert dataset to numeric and apply
coltemp <- colnames(dataMinusLastRow)
pcaData = t(apply(dataMinusLastRow, 1, as.numeric))
colnames(pcaData) <- coltemp

# PCA
pc = prcomp(t(pcaData), scale=T, center=T)
summary(pc) # From this we can see that PC184 accounts for 95% of cumulative proportion

# Getting the cumulative proportion in numeric dataset
vars <- apply(pc$x, 2, var)  
props <- vars / sum(vars)
pc.cs <- cumsum(props) # Cumulative Proportion
pc.index<-min(which(pc.cs>0.95)) # Index of which cumulative proportion > 95%

I think the main part which confuses me is that it asks for the three main PC's which I am not sure what they are and have had not much support from university with information regarding this.

Thanks to anyone who can help in any shape/form!

Answer 1

Regarding your question on how to decide which PCs to choose for plotting, prcomp is ordering PCs by their proportions in the variance. For example, using the dataset iris:

iris_pc <- prcomp(iris[,1:4], scale  = TRUE, center = TRUE)

and you get

> summary(iris_pc)
Importance of components:
                          PC1    PC2     PC3     PC4
Standard deviation     1.7084 0.9560 0.38309 0.14393
Proportion of Variance 0.7296 0.2285 0.03669 0.00518
Cumulative Proportion  0.7296 0.9581 0.99482 1.00000

BTW, using the princomp instead of prcomp gives similar results:

iris_pc2 <- princomp(iris[,1:4], cor=TRUE, scores=TRUE)

> summary(iris_pc)
Importance of components:
                          Comp.1    Comp.2     Comp.3      Comp.4
Standard deviation     1.7083611 0.9560494 0.38308860 0.143926497
Proportion of Variance 0.7296245 0.2285076 0.03668922 0.005178709
Cumulative Proportion  0.7296245 0.9581321 0.99482129 1.000000000

So, this means if you want to represent the main 3 PCs, you have to pick the first three PCs.

However, based on your dataset, if you take a look at the output of prcomp, you will see:

pc = prcomp(t(pcaData), scale=T, center=T)

> summary(pc)
Importance of components:
                          PC1    PC2     PC3     PC4     PC5     PC6     PC7    PC8     PC9    PC10   PC11    PC12
Standard deviation     9.3653 7.9806 5.46230 5.02379 3.95681 3.36495 3.21719 3.1042 2.98449 2.73437 2.6529 2.56841
Proportion of Variance 0.1993 0.1447 0.06781 0.05736 0.03558 0.02573 0.02352 0.0219 0.02024 0.01699 0.0160 0.01499

You can see that the first three PCs count for only 41% of the total variance. This can be explain by the diversity of your samples (coming from different cancers for example). So it means that your graphic representation won't be perfect to separate clusters.

Here, the graphic representation of the cumulative proportion of variances accross PC:

plot(pc, type = "l")

Next for the plotting and color attribution, I get a lot of trouble to extract the values from your Subgroups row.... for some reasons, I was not able to convert them into a character vector... weird ... anyway, I re-defined groups based on the code you put in patients ids:

patients = as.character(colnames(data))
grp = c(rep("C1",length(grep(".C1_",patients))),
        rep("C2",length(grep(".C2_",patients))),
        rep("C3",length(grep(".C3_",patients))),
        rep("C4",length(grep(".C4_",patients))),
        rep("C5",length(grep(".C5_",patients))),
        rep("C6",length(grep(".C6_",patients))))

For plotting results of PCA, you can use plot3d from rgl package for example. It gives you an interactive plot.

library(rgl)
colors = c("red","yellow", "green","cyan", "blue","purple")
plot3d(pc$x[,1:3], col =colors[as.factor(grp)], main = "pc / prcomp")
legend3d("topright",legend = unique(grp),col = colors, pch = 16, cex = 2)

And here the snapshot of your plot:

Hope that it helps you to figure it out the solution of your problem