data-raw/playground.R

In roeder/codaqol: Quality of Life Functions for Compositional Data

library(tidyverse)
library(codaqol)

### Mammalsmilk --------
d.clr <- t(apply(mammalsmilk, 1, function(x){log(x) - mean(log(x))}))
d.pcx <- prcomp(d.clr)

# Sum the total variance
milk_totvar <- sum(d.pcx$sdev^2)
# Calculate the PC1 and PC2 variance
milk_PC1 <- paste("PC1: ", round(sum(d.pcx$sdev[1]^2)/milk_totvar, 3))
milk_PC2 <- paste("PC2: ", round(sum(d.pcx$sdev[2]^2)/milk_totvar, 3))

# The basic biplot function:
biplot(d.pcx, var.axes=T, scale=0, xlab=milk_PC1, ylab=milk_PC2)
# biplot(d.pcx)

# GGPLot PCA
d.pcx$x %>%
  data.frame() %>%
  rownames_to_column() %>%
  ggplot(., aes(x = PC1, y = PC2)) +
  geom_label(aes(label = rowname)) +
  # labs(x = paste0("PC1: ", round(pcvariance[1] , 1), "% var."),
  #      y = paste0("PC2: ", round(pcvariance[2] , 1), "% var.")) +
  guides(colour = FALSE)

### GS count data ------------
gs_pca <- prcomp(sewage_clr)

gs_pcvar <- (gs_pca$sdev^2) / sum(gs_pca$sdev^2)
gs_ax_labels <- paste0("PC", 1:length(gs_pcvar), ": ", round(gs_pcvar, 4) * 100, "%")

biplot(gs_pca, var.axes=T, scale=0, xlab=gs_ax_labels[1], ylab=gs_ax_labels[2])

gs_pca$x %>%
  data.frame() %>%
  rownames_to_column() %>%
  ggplot(., aes(x = PC1, y = PC2)) +
  geom_label(aes(label = rowname)) +
  # labs(x = paste0("PC1: ", round(pcvariance[1] , 1), "% var."),
  #      y = paste0("PC2: ", round(pcvariance[2] , 1), "% var.")) +
  guides(colour = FALSE)

var_sewage <- clr_variance(sewage_clr)