This report provides examples of how to summarize and display data from Motus stations in the National Estuarine Research Reserve System. The examples can also be adapted to any Motus station.
In these examples, we will explore data from the station named NCNERR - Masonboro Island. The report will show you how to prepare the necessary R packages, download Motus data, analyze and plot Motus data, and analyze Motus data in relation to the NERRS’ System Wide Monitoring Program (SWMP).
The goal of these examples is to encourage users to perform similar analyses. To get started, you can copy the code from the .html or .Rmd version of this file and paste into a new R script, and then edit the script to analyze your own data and modify the analyses and figures.
These examples are meant to complement the results displayed in the Motus dashboard. The methods are conducted with R software and require you to have an introductory understanding of R. Advanced users may also wish to explore additional analyses using the motus package in R.
First, you will need to install and attach the following R packages. If you have not yet installed these packages, remove the hashtag that is before the install.packages() function. After you install the packages the first time, you will not need to install them again. However, you must attach the packages using the library() function each time you start a new R session.
#install.packages(c("dplyr", "ggplot2", "lubridate", "SWMPr"))
library(dplyr)
library(ggplot2)
library(lubridate)
library(SWMPr)
Note that to download specific data, you must be logged in with a Birds Canada user account. This account is free. To register for or log into an account, click the login option in the top right corner of the Motus home page.
After you have a Birds Canada account, you need to download Motus data for a specific station. Go to the Motus dashboard. Click “Find a station.” Our example is for station #12461. After you find your station, click “View Station.” Click the “Download” button near the top of the screen and select “Download detections.” Save the .csv file on your computer.
Next, you will load the data file (the .csv file that you downloaded from the Motus dashboard). The following line of code will open a dialogue box that asks you to select the .csv file. It will also save the data set to an object called ‘d’ (d for data). We will refer to d throughout the rest of this example report.
d <- read.csv(file.choose())
Examine the data. The function str() will show you the column names in the data set. NROW() will show you the number of rows in the data set.
str(d)
NROW(d)
## [1] 2506
Note that some stations have commas in the station name, which leads to incorrect formatting when the data are treated as comma separated values, as is done in this report. Therefore, it is necessary for the user to examine their data and manually correct data misalignments if they are present and pose a problem for analysis.
Before analyzing detections with respect to time, we need to convert data to a date-time format (POSIXct), set times to UTC (the default time used by Motus data), and then convert to the local timezone of the Motus station of interest. Then add variables for year, month, and hour so that we can summarize based on those timeframes.
tsStartdt <- force_tz(dmy_hms(d$tsStart), tzone = "UTC")
tsEnddt <- force_tz(dmy_hms(d$tsEnd), tzone = "UTC")
#Convert to local time zone as needed
OlsonNames() #To see all possible timezones
d$tsStartdt <- with_tz(tsStartdt, "US/Eastern")
d$tsEnddt <- with_tz(tsEnddt, "US/Eastern")
d$Year <- year(tsEnddt)
d$Month <- month(tsEnddt)
d$Hour <- hour(tsEnddt)
This section provides a summary of several important measures of the data.
The following code will show the number of species and a list of the species detected by the station.
number.of.species <- d |> select(species) |> unique() |> NROW()
number.of.species
list.of.species <- d |> select(species) |> unique()
knitr::kable(list.of.species, col.names = c("Species"), row.names = FALSE)
Optionally, you might want to remove certain entries in species list, for example, if misaligned data lead to unneeded rows (e.g., a blank, ‘166.38’ or ‘false’):
number.of.species <- d |> select(species) |> filter(!species %in% c('','166.38', 'false')) |> unique() |> NROW()
number.of.species
## [1] 25
list.of.species <- d |> select(species) |> filter(!species %in% c('','166.38', 'false')) |> unique()
knitr::kable(list.of.species, col.names = c("Species"), row.names = FALSE)
Species |
---|
Gray Catbird |
Semipalmated Plover |
Lesser Yellowlegs |
Nelson’s Sparrow |
Savannah Sparrow |
Kirtland’s Warbler |
Bobolink |
Bank Swallow |
Short-billed Dowitcher |
American Robin |
Veery |
Purple Martin |
American Kestrel (Northern) |
Dunlin |
Yellow-rumped Warbler (Myrtle) |
Sora |
Painted Bunting |
Ovenbird |
Swainson’s Thrush |
Barn Swallow |
Piping Plover |
Red Knot |
Northern Waterthrush |
American Kestrel |
Bicknell’s Thrush |
Note that this list only includes individuals that were detected by the Motus station. It does not include individuals that were tagged nearby but not detected by the station. Therefore, the number of species might be lower than the number of species shown on the Motus dashboard.
Save the species list as a .csv file. You can then open the .csv file in a spreadsheet software such as Excel and easily copy and paste the species list into a table in another document.
write.csv(list.of.species, "SpeciesList.csv", row.names = FALSE)
getwd() #Find the location of the saved species list file.
The following code will display the number of unique individuals per species. Note that some individuals might be detected in more than one year, so the row totals might be lower than the sum of values from each year.
number.of.individuals.per.species.per.year <- d |> filter(!species %in% c('', '166.38', 'false')) |> select(species, tagDeployID, Year) |> unique() |> select(species, Year) |> table() |> as.data.frame.matrix()
number.of.individuals.per.species.per.year$Total <- d |> filter(!species %in% c('', '166.38', 'false')) |> select(species, tagDeployID) |> unique() |> group_by(species) |> summarize(Total = n()) |> data.frame() |> pull(Total)
nums <- rbind(number.of.individuals.per.species.per.year, colSums(number.of.individuals.per.species.per.year))
nums <- cbind(Species = c(rownames(number.of.individuals.per.species.per.year), "Total"), nums)
rownames(nums) <- NULL
knitr::kable(nums)
Species | 2022 | 2023 | 2024 | Total |
---|---|---|---|---|
American Kestrel | 0 | 4 | 4 | 4 |
American Kestrel (Northern) | 3 | 2 | 1 | 3 |
American Robin | 1 | 1 | 0 | 1 |
Bank Swallow | 1 | 0 | 0 | 1 |
Barn Swallow | 0 | 4 | 0 | 4 |
Bicknell’s Thrush | 0 | 0 | 1 | 1 |
Bobolink | 2 | 4 | 1 | 6 |
Dunlin | 1 | 1 | 1 | 2 |
Gray Catbird | 4 | 11 | 3 | 13 |
Kirtland’s Warbler | 3 | 4 | 0 | 7 |
Lesser Yellowlegs | 2 | 1 | 0 | 2 |
Nelson’s Sparrow | 1 | 0 | 0 | 1 |
Northern Waterthrush | 0 | 6 | 0 | 6 |
Ovenbird | 0 | 1 | 0 | 1 |
Painted Bunting | 5 | 3 | 0 | 5 |
Piping Plover | 0 | 1 | 0 | 1 |
Purple Martin | 9 | 6 | 1 | 11 |
Red Knot | 0 | 5 | 2 | 5 |
Savannah Sparrow | 2 | 1 | 0 | 3 |
Semipalmated Plover | 2 | 0 | 0 | 2 |
Short-billed Dowitcher | 3 | 3 | 0 | 6 |
Sora | 1 | 2 | 2 | 2 |
Swainson’s Thrush | 0 | 1 | 0 | 1 |
Veery | 3 | 3 | 0 | 3 |
Yellow-rumped Warbler (Myrtle) | 1 | 4 | 1 | 5 |
Total | 44 | 68 | 17 | 96 |
Save the table as a .csv file.
write.csv(nums, "Speciesbyyear.csv")
getwd() #Find the location of the saved species list file.
It can be desirable to report the number of projects that a given station supports and show the connections with other stations.
The following code allows us to count the number of projects that are registered with Motus that had an animal detected by the station.
number.of.projects <- d |> select(tagProjectID) |> unique() |> NROW()
number.of.projects
## [1] 56
Impressively, animals that were detected by the station were part of dozens of registered projects, showing that this station serves a broad research community.
We will use circular plots to examine when animals were detected. We provide examples of monthly plots and hourly plots.
Summarizing data by month can provide important information on the seasons during which animals are present near the station. This can help answer questions like whether they are year-round residents or if they pass by during fall or spring migration.
We must select only one instance of each individual in a given month. This ensures that we don’t double-count the number of individuals (this would incorrectly inflate our counts).
plot_dfm <- d |> filter(!is.na(Year)) |>
select(Year, Month, species, tagDeployID) |> group_by(Year, Month, species, tagDeployID) |> unique()
The next three plots show monthly counts of unique individual animals. The first plot is for all species across all years. The plot shows the number of unique individuals detected per month. The length of the bar indicates the number of individuals, and corresponds to the scale bar on the upper left of the image–in our first plot, the inner ring indicates 10 individuals per month and the outer ring indicates 60 individuals per month. The scale bar will change depending on your data set. The width of the bar does not convey any information.
#Plot by month
ggplot(plot_dfm, aes(x = Month)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
scale_x_continuous(limits = c(1,12),
breaks = seq(1, 12, by = 1),
labels = c("Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec")) +
labs(title = "Monthly Count of Unique Animals", x = NULL, y = "Count")
The second plot shows the count of unique animals per month per year.
#Plot by month and year
ggplot(plot_dfm, aes(x = Month)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
facet_wrap(vars(Year)) +
scale_x_continuous(limits = c(1,12),
breaks = seq(1, 12, by = 1),
labels = c("Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec")) +
labs(title = "Monthly Count of Unique Animals", x = NULL, y = "Count")
The next plot is for a particular species. You can modify the aesthetics of circular plots in many ways.
#Plot by month by species (adjust species name as desired)
plot_dfm |> filter(species == "Kirtland's Warbler") |> ggplot(aes(x = Month)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
facet_wrap(vars(Year)) +
scale_x_continuous(limits = c(1,12),
breaks = seq(1, 12, by = 1),
labels = c("Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec")) +
labs(title = "Monthly Count of Unique Kirtland's Warblers", x = NULL, y = "Count")
Summarizing data by hour can provide important information on the hours during which animals are present near the station. This can help answer questions like whether they migrate at nighttime or daytime.
We must select only one instance of each individual in a given hour in a given month and year. This ensures that we don’t double-count the number of individuals (this would incorrectly inflate our counts). Filter out any detections that don’t have hour data.
plot_dfh <- d |> filter(!is.na(Hour)) |>
select(Year, Month, Hour, species, tagDeployID) |> group_by(Year, Month, Hour, species, tagDeployID) |> unique()
The next three plots show hourly counts of unique animals. The first plot is for all species across all years and months, the second is broken down by month, and the third is broken down by year. The length of the bar indicates the number of species detected in that hour. These plots use a 24-hour clock, so the top is for 00:00 hours and 24:00 hours, and the bottom is for 12:00 hours (i.e., noon). You can modify the aesthetics of the circular plots in many ways.
#Plot by hour
ggplot(plot_dfh, aes(x = Hour)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
scale_x_continuous(limits = c(0,24),
breaks = c(0, 3, 6, 9, 12, 15, 18, 21)) +
labs(title = "Hourly Count of Unique Animals", x = NULL, y = "Count")
#Plot by hour and month
ggplot(plot_dfh, aes(x = Hour)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
facet_wrap(vars(Month)) +
scale_x_continuous(limits = c(0,24),
breaks = c(0, 3, 6, 9, 12, 15, 18, 21)) +
labs(title = "Hourly Count of Unique Animals by Month", x = NULL, y = "Count")
#Plot by hour and year
ggplot(plot_dfh, aes(x = Hour)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
facet_wrap(vars(Year)) +
scale_x_continuous(limits = c(0,24),
breaks = c(0, 3, 6, 9, 12, 15, 18, 21)) +
labs(title = "Hourly Count of Unique Animals by Year", x = NULL, y = "Count")
Sometimes we want to know about individual species. The following plot shows the detections for a single species. You can adjust species name as desired.
plot_dfh |> filter(species == "Kirtland's Warbler") |> ggplot(aes(x = Hour)) +
geom_histogram(fill = "lightgreen", col = "black") +
coord_polar() +
facet_wrap(vars(Year)) +
scale_x_continuous(limits = c(0,24),
breaks = c(0, 3, 6, 9, 12, 15, 18, 21)) +
labs(title = "Hourly Count of Unique Kirtland's Warbler", x = NULL, y = "Count")
In this example, we see that Kirtland’s Warblers were detected during most hours, with the highest counts between 5 and 8am.
We might want to analyze our Motus data in relation to other NERRS data. The NERRS maintains an impressive data monitoring network called the System Wide Monitoring Program (SWMP). In this example, we will download data from Research Creek at the NCNERR (which is near the Masonboro Motus tower) in the range 1 June 2022 - 1 December 2023. We will then plot the Motus detections of Kirtland’s Warblers at the Masonboro Island Motus tower in relation to barometric pressure measured at Research Creek. We might expect a relationship between these variables if we hypothesize that migratory birds prefer to migrate under certain weather conditions. This example shows how to combine Motus and SWMP data in an analysis.
First, we must obtain the SWMP data, as follows. 1. Go to the System Wide Monitoring Program (SWMP) website 2. On the website, go to the Data Graphing and Export option 3. Choose the station of interest: NOCRCMET (Research Creek) 4. Choose range of data: 1 June 2022 - 1 December 2023. Select best available data set. 5. Provide your contact information as requested, and the CDMO will then send you an email with a link to a zipped folder with our requested data. 6. Save the zipped folder on your computer (don’t extract the files–the code will do that for you).
Next, import the SWMP data that is saved in a zipped folder on your computer. The code will open a dialog box that will ask you to select the zipped folder.
swmpd <- import_local(file.choose(), station_code = 'nocrcmet', trace = FALSE, collMethd = c("1", "2"))
Again, we are interested to learn if Kirtland’s Warblers are detected during certain pressure conditions. We will plot barometric pressure in red and detections of Kirtland’s Warblers as vertical black lines.
plot(bp ~ datetimestamp, swmpd, type = "n", xlab = "Date", ylab = "Barometric Pressure (millibars)")
lines(bp ~ datetimestamp, swmpd, col = "red")
kiwas <- d |> filter(species == "Kirtland's Warbler") |> filter(stationName == "NCNERR - Masonboro Island")
abline(v = kiwas$tsEnddt, col = "black")
Let’s zoom in to the fall of 2023 by specifying the xlim (limits of the x-axis):
overplot(swmpd, select = c('bp'), xlab = "Date", xlim = c(ymd_hms("2023-10-11 00:00:00", tz="US/Eastern"), ymd_hms("2023-11-5 00:00:00", tz="US/Eastern")))
kiwas <- d |> filter(species == "Kirtland's Warbler") |> filter(stationName == "NCNERR - Masonboro Island")
abline(v = kiwas$tsEnddt, col = "black")
We see that, in the fall of 2023, Kirtland’s Warblers passed the station on days with rising barometric pressure. This is consistent with Cooper et al. (2023), who found that Kirtland’s Warblers leave the breeding grounds for fall migration when barometric pressure has been increasing over the past 24 hours, which is indicative of calm weather over the next few days.
Let’s add wind speed to the plot to test if the weather was indeed calm when Kirtland’s Warblers passed the station.
overplot(swmpd, select = c('bp', 'wspd'), xlab = "Date", xlim = c(ymd_hms("2023-10-11 00:00:00", tz="US/Eastern"), ymd_hms("2023-11-5 00:00:00", tz="US/Eastern")))
kiwas <- d |> filter(species == "Kirtland's Warbler") |> filter(stationName == "NCNERR - Masonboro Island")
abline(v = kiwas$tsEnddt, col = "black")
It was! The wind speeds were relatively low on the days that Kirtland’s Warblers passed the station. Now see what other questions you can answer!
This report provides examples of ways to summairze and analyze Motus data.
The report was created as part of a NERRS Science Collaborative Grant and included efforts from the Cape Fear Bird Observatory, Althouse and Meade, Inc., and the NCNERR. Ray Danner authored the first version of the report.