How to perform PCA on proteomic data set

Question

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!