I edited the plotPCA code a little to do just this! First I added a "nPC" variable to the function:

function (object, intgroup = "condition", ntop = 500, returnData = FALSE, nPC = 2)

It defaults to two, so it behaves just as the original if you don't specify.

Then you need to change the following line:

d <- data.frame(PC1 = pca$x[, 1], PC2 = pca$x[, 2], group = group, intgroup.df, name = colnames(object))

to:

d <- cbind(pca$x[,seq_len(min(nPC, ncol(pca$x))), drop = FALSE], data.frame(group = group, intgroup.df, name = colnames(object)))

It's just pulling more columns from the pca$x data frame in the prcomp object, while keeping everything else in the data frame. I didn't change anything to be able to plot these additional PCs, as I usually use this to get a conveniently formatted data frame using returnData=T. But because the data is there, it wouldn't be too hard to change.

Hope this helps!

Thanks you, Dr. Michael. Kind Regards,