pca_plots.qmd

# PCA Plots

```{r}
#| label: setup
#| echo: false
#| message: false
source("_common.R")
library(knitr)

```

## Introduction

This tutorial is based on a nice blog post by Claus Wilke: <https://clauswilke.com/blog/2020/09/07/pca-tidyverse-style/>

Functions from the Tidyverse packages (`dplyr` , `tidyr`, `ggplot2`) and `broom` will be used in this tutorial. As example data a lipidomics dataset is imported that contains the lipid concentrations as columns and with one additional column indicating the experimental group

```{r}
#| label: load-lib
#| echo: true
#| message: false
library(here)
library(tidyverse)
library(broom)
library(cowplot)

```

```{r}
#| label: read-data
#| echo: true
#| message: false

d_long <- read_csv(here("data/Metabolites-1614644-supplementary.csv")) 
```

```{r}
#| label: tbl-data
#| echo: false
#| tbl-cap: "Test dataset in tidy format"
kable(head(d_long, n = 5),digits = 2)
```

## Transform and Scale the Data

```{r}
#| label: prcomp
#| echo: true
#| message: false

pca_res <- d_long %>% 
  select(where(is.numeric)) %>% # retain only numeric columns
  log2() |> 
  prcomp(scale = TRUE)

```

Compare this to plotting non-normalized concentrations

```{r}
#| label: pcaplot-1
#| echo: true
#| message: false

# Use broom to make the data ggplot ready and include annotations
pca_annot <- pca_res |> broom::augment(d_long)

pca_contrib <- pca_res |> tidy(matrix = "eigenvalues")

ggplot(data = pca_annot, aes(.fittedPC1, .fittedPC2, color = Group)) + 
geom_point(size = 1.5) + 
stat_ellipse(level = 0.95) +
theme_light() +
theme(aspect.ratio=1)

```

```{r}
#| label: pcaplot-2
#| echo: true
#| message: false
#| fig-width: 6
#| fig-height: 3

# Use broom to make the data ggplot ready and include annotations
pca_annot <- pca_res |> broom::augment(d_long)

pca_12 <- ggplot(data = pca_annot, aes(.fittedPC1, .fittedPC2, color = Group)) + 
  geom_point(size = 1.5) + 
  stat_ellipse(level = 0.95)+
  theme_light() +
  theme(legend.position = "none")


pca_34 <- ggplot(data = pca_annot, aes(.fittedPC3, .fittedPC4, color = Group)) + 
  geom_point(size = 1.5) + 
  stat_ellipse(level = 0.95) +
  theme_light() +
theme(legend.position = "none")

cowplot::plot_grid(pca_12, pca_34, nrow = 1, 
                   rel_widths = c(1,1), 
                   labels = c("A","B"))

```
```{r}
#| label: pcaplot-3
#| echo: true
#| message: false
#| fig-width: 6
#| fig-height: 3


# PCA plot of merged


 
m_log <- d_long  |> 
  dplyr::select(-Group) |> 
  dplyr::select(where(~!any(is.na(.)))) |> 
  column_to_rownames("SubjectID") |> 
  log2() 

d_annot <- d_long |> 
  select(SubjectID, Group)|> 
  column_to_rownames("SubjectID")

d_annot$Group <- factor(d_annot$Group)


# get pca result with annotation
pca_res <- PCAtools::pca(mat = m_log, center = TRUE, 
                         metadata = d_annot , 
                         scale = TRUE, 
                         transposed = TRUE)
```

```{r}
#| label: pcaplot-4
#| echo: true
#| message: false
#| fig-width: 6
#| fig-height: 6
#| 
  p_PCA <- PCAtools::biplot(
    pca_res,
    x = 1,
    y = 2 ,
    colby = "Group",
    #colkey = colby_color,
    legendPosition = 'right',
    pointSize = 3,
    lengthLoadingsArrowsFactor = 1.5,
    lab = NULL,
    sizeLoadingsNames = 2.5,
    alphaLoadingsArrow = 0.1,
    colLoadingsNames = "grey35",
    showLoadings = TRUE,
    ntopLoadings = 20,
    drawConnectorsLoadings = T,
    boxedLoadingsNames = T,
    widthLoadingsArrows = 0.5,
    widthConnectorsLoadings = 0.001,
    max.overlaps = 33,
    ellipseType = "t",
    # xlim = c(-22, 25),
    # ylim = c(-21, 20),
    ellipse = TRUE,
    ellipseLevel = 0.95,
    ellipseLineSize = 1,
    ellipseFill = TRUE,
    ellipseAlpha = .05,
    gridlines.minor = FALSE,
  )

p_PCA
    # theme(
    #   aspect.ratio = 1,
    #   panel.grid.major = element_line(size = .5),
    #   legend.position = c(.91, .93),
    #   legend.title = element_blank(),
    #   legend.background = element_rect(fill = alpha("white", 0.0), inherit.blank = TRUE)
    # )

```