2
$\begingroup$

I'm trying to perform Principal Component Analysis using R on a proteomics dataset.

As the dataset contains a lot missing values I tried different approaches.

  1. I ran PCA using prcomp() and autoplot() after deleting all samples with missing values. This worked out nicely, but the plots are not representative as I had to delete more than 70% of the samples.

  2. I imputed the missing values and ran a PCA (again prcomp() and autoplot()) which worked out, but provided me with a PCA plot showing a strong horseshoe effect (which seems off as the PCA I performed on the complete samples of my data doesn't show such a distribution).

  3. I then came across this article (https://www.sciencedirect.com/science/article/pii/S1574954121000261) which introduces a modified PCA, InDaPCA (PCA of Incomplete Data). I transformed the code provided by the authors according to my data but again got errors/warnings because of missing values.

This is the structure of my dataframe, the code I used to perform InDaPCA as well as the resulting error:

num [1:261, 1:1132] 13.18 13.18 14.14 8.32 9.11 ...
 - attr(*, "dimnames")=List of 2
  ..$ : NULL
  ..$ : chr [1:1132] "protein.A0A075B6H9" "protein.A0A075B6I0" "protein.A0A075B6I1" "protein.A0A075B6I9" ...

InDaPCA <- function(protein_data_filtered) {
  
  #scaling
  X <- scale(protein_data, center = TRUE, scale = TRUE)
  
  #correlation
  C <- cor(X, use = "pairwise.complete.obs")
  
  #Eigenvalue
  Eigenvalues <- eigen(C)$values
  Eigenvalues.pos <- Eigenvalues[Eigenvalues > 0]
  Eigenvalues.pos.as.percent <- 100 * Eigenvalues.pos / sum(Eigenvalues.pos)
  
  #Eigenvectors
  V <- eigen(C)$vectors
  
  #Principal components
  X2 <- X
  X2[is.na(X2)] <- 0
  PC <- as.matrix(X2) %*% V
  
  #object.standardized
  PCstand1 <- PC[, Eigenvalues > 0] / sqrt(Eigenvalues.pos)[col(PC[, Eigenvalues > 0])]
  PCstand2 <- PCstand1 / sqrt(nrow(PC) - 1)
  
  #loadings
  loadings <- cor(X, PC, use = "pairwise.complete.obs")
  
  #arrows for biplot
  arrows <- cor(X, PC, use = "pairwise.complete.obs") * sqrt(nrow(X) - 1)
  
  #output
  PCA <- list()
  PCA$Correlation.matrix <- C
  PCA$Eigenvalues <- Eigenvalues
  PCA$Positive.Eigenvalues <- Eigenvalues.pos
  PCA$Positive.Eigenvalues.as.percent <- Eigenvalues.pos.as.percent
  PCA$Square.root.of.eigenvalues <- sqrt(Eigenvalues.pos)
  PCA$Eigenvectors <- V
  PCA$Component.scores <- PC
  PCA$Variable.scores <- loadings
  PCA$Biplot.objects <- PCstand2
  PCA$Biplot.variables <- arrows
  
  return(PCA)
}

pca_result <- InDaPCA(protein_data_filtered)

Error in eigen(C) : unendliche oder fehlende Werte in 'x'
Called from: eigen(C)  

I would highly appreciate any suggestions on how to perform PCA of datasets with many NA values and/or how to handle the horseshoe effect.

Thanks a lot in advance!

$\endgroup$
1
  • $\begingroup$ Not a direct answer to your problem but how about using a package tailored for proteomics datasets? For example msqrob2 workflow uses limma's plotMDS() after normalizing the data. $\endgroup$
    – haci
    Commented Sep 14, 2023 at 8:53

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.