flowchart LR
A("Format</br>sound files") --> B("Detect</br>sounds")
B --> C("Measure</br>acoustic features")
C --> D("Mitigate</br>false positives")
D --> E("Summarize</br>USVs counts")
style A fill:#382A5433
style B fill:#395D9C33
style C fill:#3497A933
style D fill:#60CEAC33
style E fill:#DEF5E533
Detecting rat ultrasonic vocalizations
Example code
Objetive
- Example code on how to use ohun to automatically detect rat ultrasonic vocalizations
Analysis flowchart
Load packages
Code
# install sketchy if not installed
if (!requireNamespace("sketchy", quietly = TRUE)) {
install.packages("sketchy")
}
## add 'developer/' to packages to install from github
x <- c(
"maRce10/warbleR",
"ranger",
"maRce10/ohun"
)
sketchy::load_packages(x)1 Custom function
Code
# function to convert seconds to min:s
seg_2_minseg <- function(seg) {
minutos <- seg %/% 60
segundos <- seg %% 60
sprintf("%d:%02d", minutos, segundos) # 2 digits
}2 Prepare sound files
2.1 Set directory with original sound files
- Modify with your own directory paths
Code
# where the original sound files are
sound_file_path <- "./data/raw/sound_files/"
# where to save the consolidated sound files and clips
consolidate_path <- "./data/raw/sound_files/"
# where to save the output data
output_data_path <- "./data/raw"
warbleR_options(wav.path = sound_file_path)These are the sound files used in this example:
Code
files <- list.files(
path = sound_file_path,
full.names = TRUE,
recursive = TRUE,
pattern = ".wav$|.wac$|.mp3$|.flac$"
)[1] "./data/raw/sound_files//A/example_A.flac"
[2] "./data/raw/sound_files//B/example_B.flac"
[3] "./data/raw/sound_files//D/example_D.flac"Sound files can be found at: https://github.com/maRce10/Automatic-detection-of-rat-ultrasonic-vocalizations/tree/master/data/raw/sound_files
2.2 Consolidate sound files
- Put all original sound files in a single folder
Code
dir.create(file.path(consolidate_path, "consolidated_sound_files"))
cns <- consolidate(
path = sound_file_path,
dest.path = file.path(consolidate_path, "consolidated_sound_files"),
parallel = 1,
file.ext = ".wav$|.flac$"
)2.3 Homogenize sound file format
- Convert flac to wav (if needed)
- Downsample files
- Update sound file directory afterwards
Code
# convert flac to wav
if (any(grepl("\\.flac$", files)))
warbleR::wav_2_flac(path = file.path(consolidate_path, "consolidated_sound_files"), reverse = TRUE)
fix_wavs(samp.rate = 200, bit.depth = 16, path = file.path(consolidate_path, "consolidated_sound_files"))
warbleR_options(
wav.path = file.path(consolidate_path, "consolidated_sound_files", "converted_sound_files")
)2.4 Split sound files into 5 min clips
Update sound file directory afterwards
Code
# check files
feature_acoustic_data(path = .Options$warbleR$path)
# segment duration 2 min
clip_duration <- 2 * 60
ssf <- split_acoustic_data(sgmt.dur = clip_duration, cores = 1, path = .Options$warbleR$path)
write.csv(ssf, file.path(consolidate_path, "consolidated_sound_files", "converted_sound_files", "5min_clip_info.csv"), row.names = FALSE)
clips_path <- file.path(consolidate_path, "consolidated_sound_files", "converted_sound_files", "clips")
warbleR_options(
wav.path = clips_path
)3 Automatic detection
3.0.1 Energy detector
- Uses optimized detection parameters
- Saves results in a RDS file
Code
detection <- energy_detector(
path = .Options$warbleR$path,
thinning = 0.5,
bp = c(35, 90),
smooth = 1,
threshold = 2.5,
hold.time = 3,
min.duration = 1,
max.duration = 200,
cores = 1
)
saveRDS(
detection,
file.path(
output_data_path,
"detection.RDS"
)
)3.0.2 Random forest classification
3.0.2.1 Measure acoustic features
- Save features at the end of the analysis
Code
detection <- readRDS(file.path(output_data_path, "detection.RDS"))
# measure spectrographic parameters
spectral_features <- spectro_analysis(
detection,
bp = c(35, 85),
fast = TRUE,
ovlp = 70,
parallel = 1
)
# check if any NA
sapply(spectral_features, function(x)
sum(is.na(x)))
# remove NAs
detection <- detection[complete.cases(spectral_features), ]
spectral_features <- spectral_features[complete.cases(spectral_features), ]
# save acoustic parameters just in case
write.csv(
spectral_features,
file.path(output_data_path, "spectral_features.csv"),
row.names = FALSE
)3.0.2.2 Predict based on pre-defined RF model
- Download model from figshare
Code
model_rds <- tempfile()
# this downloads the 55 kHz USV with bedding detection model
download.file(url = "https://figshare.com/ndownloader/files/57261971", destfile = model_rds)3.0.3 Energy detector
- Uses optimized detection parameters
- Saves results in a RDS file
Code
detection <- energy_detector(
path = .Options$warbleR$path,
thinning = 0.5,
bp = c(35, 90),
smooth = 5,
threshold = 1,
hold.time = 5,
min.duration = 1,
max.duration = 150,
cores = 1
)
saveRDS(
detection,
file.path(
output_data_path,
"detection.RDS"
)
)3.0.4 Random forest classification
3.0.4.1 Measure acoustic features
- Save features at the end of the analysis
Code
detection <- readRDS(file.path(output_data_path, "detection.RDS"))
# measure spectrographic parameters
spectral_features <- spectro_analysis(
detection,
bp = c(35, 85),
fast = TRUE,
ovlp = 70,
parallel = 1
)
# check if any NA
sapply(spectral_features, function(x)
sum(is.na(x)))
# remove NAs
detection <- detection[complete.cases(spectral_features), ]
spectral_features <- spectral_features[!complete.cases(spectral_features), ]
# save acoustic parameters just in case
write.csv(
spectral_features,
file.path(output_data_path, "spectral_features.csv"),
row.names = FALSE
)3.0.4.2 Predict based on pre-defined RF model
- Download model from figshare
Code
model_rds <- tempfile()
# this downloads the 55 kHz USV without bedding detection model
download.file(url = "https://figshare.com/ndownloader/files/57261968", destfile = model_rds)3.0.5 Energy detector
- Uses optimized detection parameters
- Saves results in a RDS file
Code
detection <- energy_detector(
path = .Options$warbleR$path,
thinning = 0.5,
bp = c(20, 30),
smooth = 17,
threshold = 2,
hold.time = 25,
min.duration = 2,
max.duration = 300,
cores = 1
)
saveRDS(
detection,
file.path(
output_data_path,
"detection.RDS"
)
)3.0.6 Random forest classification
3.0.6.1 Measure acoustic features
- Save features at the end of the analysis
Code
detection <- readRDS(file.path(output_data_path, "detection.RDS"))
# measure spectrographic parameters
spectral_features <- spectro_analysis(
detection,
bp = c(20, 30),
fast = TRUE,
ovlp = 70,
parallel = 1
)
# check if any NA
sapply(spectral_features, function(x)
sum(is.na(x)))
# remove NAs
detection <- detection[complete.cases(spectral_features), ]
spectral_features <- spectral_features[!complete.cases(spectral_features), ]
# save acoustic parameters just in case
write.csv(
spectral_features,
file.path(output_data_path, "spectral_features.csv"),
row.names = FALSE
)
# resave detections
saveRDS(
detection,
file.path(
output_data_path,
"detection.RDS"
)
)3.0.6.2 Predict based on pre-defined RF model
- Download model from figshare
Code
model_rds <- tempfile()
# this downloads the 22 kHz USV without bedding detection model
download.file(url = "https://figshare.com/ndownloader/files/57261959", destfile = model_rds)- (from here the code is based on the 55 kHz detection procedure)
Apply the model to the measured acoustic parameters
Code
# read model
rf_model <- readRDS(model_rds)
# read acoustic features
spectral_features <- read.csv(
file.path(output_data_path, "spectral_features.csv"),
stringsAsFactors = FALSE
)
# apply model
detection$class <- predict(object = rf_model, data = spectral_features)$predictions
# keep only true positives
filtered_detection <- detection[detection$class == "true.positive", ]
saveRDS(
filtered_detection,
file.path(
output_data_path,
"random_forest_filtered_detection.RDS"
)
)4 Summarized
4.1 Reassemble detections to original sound files
Code
# read detections
filtered_detection <- readRDS(file.path(output_data_path, "random_forest_filtered_detection.RDS"))
# read clip information
clip_info <- read.csv(
file.path(
consolidate_path,
"consolidated_sound_files",
"converted_sound_files",
"5min_clip_info.csv"
)
)
# reassemble to original (long) sound files
reass_detec <- reassemble_detection(detection = filtered_detection, Y = clip_info, pb = FALSE)4.2 USV counts per minute
- To get counts for each minute the code uses
hop.size = 60
Code
# counts per minute
## include files so it includes those with no detections
count_min <- acoustic_activity(
X = reass_detec,
time.window = 60,
hop.size = 60,
path = file.path(consolidate_path, "/consolidated_sound_files/"),
files = list.files(
path = file.path(consolidate_path, "/consolidated_sound_files/"),
pattern = "\\.wav$"
)
)
# convert to minute rate
count_min$minute <- count_min$start / 60 + 1
wide_count_min <- reshape(count_min[, c("counts", "minute", "sound.files")],
direction = "wide",
idvar = "sound.files",
timevar = "minute")
names(wide_count_min) <- c("sound.files", paste("min", 1:max(count_min$minute)))
wide_count_min$total <- apply(wide_count_min[, -1], 1, sum, na.rm = TRUE)
# print results
wide_count_min| sound.files | min 1 | min 2 | min 3 | min 4 | min 5 | total | |
|---|---|---|---|---|---|---|---|
| 1 | example_A.wav | 166 | 148 | 205 | 155 | 134 | 808 |
| 6 | example_B.wav | 155 | 187 | 210 | 173 | 49 | 774 |
| 11 | example_D.wav | 196 | 248 | 250 | 221 | 70 | 985 |
4.2.1 Minute with highest USV count
- Finds the minute with the highest USV count per sound file
- The function scans across the entire sound file counting sounds every 1 s (
hop.size = 1)
Code
# counting
count_high_min <- acoustic_activity(
X = reass_detec,
time.window = 60,
hop.size = 1,
path = file.path(consolidate_path, "/consolidated_sound_files/")
)
# rate per minute
count_high_min$rate <- count_high_min$rate * 60
# get highest per sound file
sub_count_high_min <- count_high_min[count_high_min$duration == 60, ]
# get the row with the highest rate per sound file
highes_min <- do.call(rbind, lapply(split(sub_count_high_min, sub_count_high_min$sound.files), function(x)
x[which.max(x$rate), ]))
highes_min$`start (min:s)` <- seg_2_minseg(highes_min$start)
highes_min$`end (min:s)` <- seg_2_minseg(highes_min$end)
# order columns
highes_min <- highes_min[, c(
"sound.files",
"start",
"end",
"start (min:s)",
"end (min:s)",
"duration",
"counts",
"rate"
)]
highes_min| sound.files | start | end | start (min:s) | end (min:s) | duration | counts | rate | |
|---|---|---|---|---|---|---|---|---|
| example_A.wav | example_A.wav | 23 | 83 | 0:23 | 1:23 | 60 | 220 | 220 |
| example_B.wav | example_B.wav | 98 | 158 | 1:38 | 2:38 | 60 | 236 | 236 |
| example_D.wav | example_D.wav | 39 | 99 | 0:39 | 1:39 | 60 | 290 | 290 |
4.2.2 Counts before and after a given time
- For the sake of the example, lest’s say we introduced a mate into the cage in the middle of the experiment. The following data simulates the mate entry time for each sound file:
Code
# simulated mate entry times
mate_entry <- data.frame(sound.files = unique(highes_min$sound.files),
time = c(56, 64, 68))
mate_entry| sound.files | time |
|---|---|
| example_A.wav | 56 |
| example_B.wav | 64 |
| example_D.wav | 68 |
- We can now count the number of USVs before and after the mate entry time
- The code also counts the number of USVs in the minute right before and right after the mate entry
Code
# count before and after
out <- lapply(mate_entry$sound.files, function(x) {
# get the sound file
sf <- reass_detec[reass_detec$sound.files == x, ]
sf$mid.point <- (sf$end + sf$start) / 2
# get the time of the mate entry
t <- mate_entry$time[mate_entry$sound.files == x]
# get the calls before and after the mate entry
before <- sum(sf$mid.point < t)
after <- sum(sf$mid.point >= t)
min.before <- sum(sf$mid.point < t & sf$mid.point >= t - 60)
min.after <- sum(sf$mid.point >= t & sf$mid.point < t + 60)
return(
data.frame(
sound.files = x,
mate.entry = t,
count.before = before,
count.after = after,
count.min.before = min.before,
count.min.after = min.after
)
)
})
counts <- do.call(rbind, out)
counts| sound.files | mate.entry | count.before | count.after | count.min.before | count.min.after |
|---|---|---|---|---|---|
| example_A.wav | 56 | 159 | 649 | 159 | 152 |
| example_B.wav | 64 | 172 | 602 | 171 | 187 |
| example_D.wav | 68 | 228 | 757 | 223 | 255 |
Session information
Click to see
─ Session info ───────────────────────────────────────────────────────────────
setting value
version R version 4.3.2 (2023-10-31)
os Ubuntu 22.04.2 LTS
system x86_64, linux-gnu
ui X11
language es_ES:en
collate en_US.UTF-8
ctype en_US.UTF-8
tz America/Costa_Rica
date 2025-08-19
pandoc 3.1.11 @ /usr/lib/rstudio/resources/app/bin/quarto/bin/tools/x86_64/ (via rmarkdown)
─ Packages ───────────────────────────────────────────────────────────────────
package * version date (UTC) lib source
backports 1.5.0 2024-05-23 [1] CRAN (R 4.3.2)
bitops 1.0-9 2024-10-03 [1] CRAN (R 4.3.2)
brio 1.1.5 2024-04-24 [1] CRAN (R 4.3.2)
cachem 1.1.0 2024-05-16 [1] CRAN (R 4.3.2)
checkmate 2.3.3 2025-08-18 [1] CRAN (R 4.3.2)
class 7.3-20 2022-01-13 [4] CRAN (R 4.1.2)
classInt 0.4-11 2025-01-08 [1] CRAN (R 4.3.2)
cli 3.6.5 2025-04-23 [1] CRAN (R 4.3.2)
crayon 1.5.3 2024-06-20 [1] CRAN (R 4.3.2)
curl 6.4.0 2025-06-22 [1] CRAN (R 4.3.2)
DBI 1.2.3 2024-06-02 [1] CRAN (R 4.3.2)
devtools 2.4.5 2022-10-11 [1] CRAN (R 4.3.2)
digest 0.6.37 2024-08-19 [1] CRAN (R 4.3.2)
dplyr 1.1.4 2023-11-17 [1] CRAN (R 4.3.2)
DT * 0.31 2023-12-09 [1] CRAN (R 4.3.2)
dtw 1.23-1 2022-09-19 [1] CRAN (R 4.3.2)
e1071 1.7-16 2024-09-16 [1] CRAN (R 4.3.2)
ellipsis 0.3.2 2021-04-29 [3] CRAN (R 4.1.1)
evaluate 1.0.4 2025-06-18 [1] CRAN (R 4.3.2)
farver 2.1.2 2024-05-13 [1] CRAN (R 4.3.2)
fastmap 1.2.0 2024-05-15 [1] CRAN (R 4.3.2)
fftw 1.0-9 2024-09-20 [1] CRAN (R 4.3.2)
fs 1.6.6 2025-04-12 [1] CRAN (R 4.3.2)
generics 0.1.3 2022-07-05 [1] CRAN (R 4.3.2)
ggplot2 3.5.2 2025-04-09 [1] CRAN (R 4.3.2)
glue 1.8.0 2024-09-30 [1] CRAN (R 4.3.2)
gtable 0.3.6 2024-10-25 [1] CRAN (R 4.3.2)
htmltools 0.5.8.1 2024-04-04 [1] CRAN (R 4.3.2)
htmlwidgets 1.6.4 2023-12-06 [1] CRAN (R 4.3.2)
httpuv 1.6.13 2023-12-06 [1] CRAN (R 4.3.2)
httr 1.4.7 2023-08-15 [1] CRAN (R 4.3.2)
igraph 2.1.4 2025-01-23 [1] CRAN (R 4.3.2)
jsonlite 2.0.0 2025-03-27 [1] CRAN (R 4.3.2)
kableExtra * 1.4.0 2024-01-24 [1] CRAN (R 4.3.2)
KernSmooth 2.23-20 2021-05-03 [4] CRAN (R 4.0.4)
knitr * 1.50 2025-03-16 [1] CRAN (R 4.3.2)
later 1.3.2 2023-12-06 [1] CRAN (R 4.3.2)
lattice 0.20-45 2021-09-22 [4] CRAN (R 4.1.1)
lifecycle 1.0.4 2023-11-07 [1] CRAN (R 4.3.2)
magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.3.2)
MASS 7.3-55 2022-01-13 [4] CRAN (R 4.1.2)
Matrix 1.6-5 2024-01-11 [1] CRAN (R 4.3.2)
memoise 2.0.1 2021-11-26 [3] CRAN (R 4.1.2)
mime 0.13 2025-03-17 [1] CRAN (R 4.3.2)
miniUI 0.1.1.1 2018-05-18 [3] CRAN (R 4.0.1)
NatureSounds * 1.0.5 2025-01-17 [1] CRAN (R 4.3.2)
ohun * 1.0.3 2025-08-18 [1] local
packrat 0.9.2 2023-09-05 [1] CRAN (R 4.3.2)
pbapply 1.7-4 2025-07-20 [1] CRAN (R 4.3.2)
pillar 1.11.0 2025-07-04 [1] CRAN (R 4.3.2)
pkgbuild 1.4.8 2025-05-26 [1] CRAN (R 4.3.2)
pkgconfig 2.0.3 2019-09-22 [3] CRAN (R 4.0.1)
pkgload 1.4.0 2024-06-28 [1] CRAN (R 4.3.2)
profvis 0.3.8 2023-05-02 [1] CRAN (R 4.3.2)
promises 1.2.1 2023-08-10 [1] CRAN (R 4.3.2)
proxy 0.4-27 2022-06-09 [1] CRAN (R 4.3.2)
purrr 1.1.0 2025-07-10 [1] CRAN (R 4.3.2)
R6 2.6.1 2025-02-15 [1] CRAN (R 4.3.2)
ranger * 0.16.0 2023-11-12 [1] CRAN (R 4.3.2)
RColorBrewer 1.1-3 2022-04-03 [1] CRAN (R 4.3.2)
Rcpp 1.1.0 2025-07-02 [1] CRAN (R 4.3.2)
RCurl 1.98-1.17 2025-03-22 [1] CRAN (R 4.3.2)
remotes 2.5.0 2024-03-17 [1] CRAN (R 4.3.2)
rjson 0.2.23 2024-09-16 [1] CRAN (R 4.3.2)
rlang 1.1.6 2025-04-11 [1] CRAN (R 4.3.2)
rmarkdown 2.28 2024-08-17 [1] CRAN (R 4.3.2)
Rraven * 1.0.14 2017-11-17 [1] CRAN (R 4.3.2)
rstudioapi 0.17.1 2024-10-22 [1] CRAN (R 4.3.2)
scales 1.4.0 2025-04-24 [1] CRAN (R 4.3.2)
seewave * 2.2.3 2023-10-19 [1] CRAN (R 4.3.2)
sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.3.2)
sf 1.0-21 2025-05-15 [1] CRAN (R 4.3.2)
shiny 1.8.0 2023-11-17 [1] CRAN (R 4.3.2)
signal 1.8-1 2024-06-26 [1] CRAN (R 4.3.2)
sketchy 1.0.5 2025-06-13 [1] local (/home/marce/Dropbox/R_package_testing/sketchy)
stringi 1.8.7 2025-03-27 [1] CRAN (R 4.3.2)
stringr 1.5.1 2023-11-14 [1] CRAN (R 4.3.2)
svglite 2.1.3 2023-12-08 [1] CRAN (R 4.3.2)
systemfonts 1.1.0 2024-05-15 [1] CRAN (R 4.3.2)
testthat 3.2.3 2025-01-13 [1] CRAN (R 4.3.2)
tibble 3.3.0 2025-06-08 [1] CRAN (R 4.3.2)
tidyselect 1.2.1 2024-03-11 [1] CRAN (R 4.3.2)
tuneR * 1.4.7 2024-04-17 [1] CRAN (R 4.3.2)
units 0.8-7 2025-03-11 [1] CRAN (R 4.3.2)
urlchecker 1.0.1 2021-11-30 [1] CRAN (R 4.3.2)
usethis 3.1.0 2024-11-26 [1] CRAN (R 4.3.2)
vctrs 0.6.5 2023-12-01 [1] CRAN (R 4.3.2)
viridisLite 0.4.2 2023-05-02 [1] CRAN (R 4.3.2)
warbleR * 1.1.36 2016-04-19 [1] Github (maRce10/warbleR@a938839)
xaringanExtra 0.8.0 2024-05-19 [1] CRAN (R 4.3.2)
xfun 0.52 2025-04-02 [1] CRAN (R 4.3.2)
xml2 1.3.6 2023-12-04 [1] CRAN (R 4.3.2)
xtable 1.8-4 2019-04-21 [3] CRAN (R 4.0.1)
yaml 2.3.10 2024-07-26 [1] CRAN (R 4.3.2)
[1] /home/marce/R/x86_64-pc-linux-gnu-library/4.3
[2] /usr/local/lib/R/site-library
[3] /usr/lib/R/site-library
[4] /usr/lib/R/library
──────────────────────────────────────────────────────────────────────────────