monogenea

# Tue Feb  4 19:43:33 2020 ------------------------------
setwd("~/Documents/Tutorials/birdsong")
library(parallel)
library(tidyverse)
library(abind)
library(caret)
library(tuneR)
library(warbleR)
source("funs.R")

monogenea

#### Download HQ male song recordings > 30s long from Europe ####
query <- querxc("type:song type:male len_gt:30 q_gt:C area:europe")
query$Species <- with(query, paste(Genus, Specific_epithet))
# Select top 50 most abundant bird species
speciesCount <- sort(table(query$Species), decreasing = T)
topSpecies <- names(speciesCount)[1:50]
query <- query[query$Species %in% topSpecies, ]
# Downsample to min size among the 50 classes
balancedClasses <- lapply(topSpecies, function(x){
   set.seed(100)

monogenea

#### Pre-processing ####
# Read files
fnames <- list.files("mp3/", full.names = T, patt = "*.mp3")

# Write metadata for Kaggle dataset
ids <- str_extract(fnames, pattern = "[0-9]{4,}")
query$Path <- fnames[match(query$Recording_ID, ids)]
write.csv(query, "metadata.csv")

# Play random file - setWavPlayer in macOS if "permission denied"

monogenea

# Encode species from fnames regex
species <- str_extract(fnames, patt = "[A-Za-z]+-[a-z]+") %>%
      gsub(patt = "-", rep = " ") %>% factor()

# Stratified sampling: train (80%), val (10%) and test (10%)
set.seed(100)
idx <- createFolds(species, k = 10)
valIdx <- idx$Fold01
testIdx <- idx$Fold02
# Define samples for train, val and test

monogenea

# Define targets and augment data
target <- model.matrix(~0+species)

targetTrain <- do.call("rbind", lapply(1:(dim(Xtrain)[1]/length(fnamesTrain)),
                                       function(x) target[-c(valIdx, testIdx),]))
targetVal <- do.call("rbind", lapply(1:(dim(Xval)[1]/length(fnamesVal)),
                                     function(x) target[valIdx,]))
targetTest <- do.call("rbind", lapply(1:(dim(Xtest)[1]/length(fnamesTest)),
                                      function(x) target[testIdx,]))
# Assemble Xs and Ys

monogenea

# Plot spectrogram from random training sample - range is 0-22.05 kHz
image(train$X[sample(dim(train$X)[1], 1),,,],
      xlab = "Time (s)",
      ylab = "Frequency (kHz)",
      axes = F)
# Generate mel sequence from Hz points, standardize to plot
freqs <- c(0, 1, 5, 15, 22.05)
mels <- 2595 * log10(1 + (freqs*1e3) / 700) # https://en.wikipedia.org/wiki/Mel_scale
mels <- mels - min(mels)
mels <- mels / max(mels)

monogenea

# Fri Feb  7 15:49:46 2020 ------------------------------
setwd("~/Documents/Tutorials/birdsong")
library(keras)
use_condaenv("plaidml")
use_backend("plaidml")
k_backend() # plaidml
library(tidyverse)
library(caret)
library(e1071)
library(pheatmap)

monogenea

# Build model
model <- keras_model_sequential() %>% 
        layer_conv_2d(input_shape = dim(train$X)[2:4], 
                      filters = 16, kernel_size = c(3, 3),
                      activation = "relu") %>% 
        layer_max_pooling_2d(pool_size = c(2, 2)) %>% 
        layer_dropout(rate = .2) %>% 
        
        layer_conv_2d(filters = 32, kernel_size = c(3, 3),
                      activation = "relu") %>% 

monogenea

# Print summary
summary(model)
model %>% compile(optimizer = optimizer_adam(decay = 1e-5),
                  loss = "categorical_crossentropy",
                  metrics = "accuracy")

history <- fit(model, x = train$X, y = train$Y,
               batch_size = 16, epochs = 50,
               validation_data = list(val$X, val$Y))

monogenea

# Grep species, set colors for heatmap
speciesClass <- gsub(colnames(train$Y), pat = "species", rep = "")
cols <- colorRampPalette(rev(brewer.pal(n = 7, name = "RdGy")))
        
# Validation predictions
predProb <- predict(model, val$X)
predClass <- speciesClass[apply(predProb, 1, which.max)]
trueClass <- speciesClass[apply(val$Y, 1, which.max)]

# Plot confusion matrix