Define Libraries

library("stringr")
library("dplyr")
library("reshape2")
library("vegan")

Define Path

dir.wrk <- str_replace(getwd(), "/scripts", "")
dir.data <- file.path(dir.wrk, "data/data_permutation_withReplacement")
dir.output <- file.path(dir.wrk, "data/data_processed")

Function: Log Transformation

TransformLog <- function(dat) {
    dml <- vegan::decostand(dat[, 5:10], "log")
    dml <- cbind(dat[, -c(5:10)], dml)
    
    return(dml)
}

Function: Permanova Test

getPermanova <- function(dat){
     # GET DISTANCE MATRIX ---
     dist_dml <- vegan::vegdist(x=as.matrix(dat[,5:10]), method="euclidean", binary=FALSE, diag=TRUE, upper=TRUE, na.rm=FALSE)

     # PERMANOVA TEST ---
     set.seed(12345)
     y_permanova <- vegan::adonis(dist_dml ~ IncomeGroup+Ethnicity+EducationLevel+GeoRegion, 
                                   data=dat, permutations=999, method="euclidean", parallel=8)
     
     return(y_permanova)
}

Perform Permanova

nperm <- 1000

for (i in 1:nperm) {
    # DEFINE FILE ---
    file.dat <- file.path(dir.data, paste("abundance_permmute_", i, ".tsv", sep = ""))
    
    # LOAD DATA ---
    dat <- read.delim(file.dat, header = TRUE, stringsAsFactors = FALSE)
    
    # LOG TRANSFORM ---
    dml <- TransformLog(dat)
    
    # PERFORM PERMANOVA TEST ---
    pnova <- getPermanova(dml)
    dm <- as.data.frame(pnova$aov.tab)
    
    # WRITE OUTPUT ---
    file.output <- file.path(dir.data, paste("stats_permanova_permmute_", i, ".tsv", 
        sep = ""))
    write.table(dm, file.output, sep = "\t", row.names = TRUE, col.names = NA, quote = FALSE)
}