library(dplyr)
library(stringr)
library(terra)
library(purrr)
library(ggplot2)31 Analysis example: Using AusTraits for trait-climate analyses
31.1 Introduction
A more involved analysis is merging trait data with climatic parameters.
Please note the data summary and analysis presented here is similar to that in a manuscript currently in review, Towers et al, “Revisiting the role of mean annual precipitation in shaping functional trait distributions at a continental scale”, doi: 10.1101/2023.06.29.546983
The code presented here must be run offline on your machine, as the climate data files required are not hosted on the traits.build-book GitHub repository.
This example uses the most recent AusTraits release, version 6.0.0
most_recent <- austraits::get_versions()[["doi"]][1]
most_recent[1] "10.5281/zenodo.11188867"
austraits <- austraits::load_austraits(doi = most_recent)Loading data from 'data/austraits/austraits-6.0.0.rds'
31.2 Extract data from AusTraits
For the analysis here, two subsets of AusTraits data are required.
First, extracting data on woodiness. There are three traits for which AusTraits has near-complete taxon coverage. Woodiness & plant growth form are in Wenk_2022 and life history is in Wenk_2023_1:
woody_taxa <-
austraits$traits %>%
dplyr::filter(trait_name == "woodiness_detailed") %>%
dplyr::filter(dataset_id == "Wenk_2022") %>%
dplyr::select(taxon_name, value) %>%
dplyr::mutate(woodiness = ifelse(stringr::str_detect(value, "herb"), "herbaceous", "woody")) %>%
dplyr::select(-value) %>%
dplyr::distinct(taxon_name, .keep_all = TRUE)Second, extract the numeric trait data you want to analyse. Only field-collected data on adult plants with reported coordinates is used:
trait_data <-
(austraits %>% join_locations())$traits %>%
rename(latitude = `latitude (deg)`, longitude = `longitude (deg)`) %>%
dplyr::filter(!is.na(latitude)) %>%
dplyr::filter(trait_name == "leaf_mass_per_area") %>%
dplyr::filter(basis_of_record == "field") %>%
dplyr::filter(life_stage == "adult") %>%
dplyr::group_by(dataset_id, location_name, taxon_name) %>%
dplyr::mutate(value = mean(as.numeric(value))) %>%
dplyr::ungroup() %>%
dplyr::mutate(
ID = row_number(),
latitude = as.numeric(latitude),
longitude = as.numeric(longitude)
) %>%
distinct(value, dataset_id, location_name, taxon_name, .keep_all = TRUE) %>%
dplyr::filter(!is.na(latitude)) %>%
dplyr::filter(!is.na(longitude)) %>%
dplyr::left_join(woody_taxa, by = "taxon_name")31.3 Download and open spatial climatic data
Download the climatic layers of interest from WorldClim:
https://www.worldclim.org/data/worldclim21.html
Download a shapefile of Australia:
https://www.abs.gov.au/statistics/standards/australian-statistical-geography-standard-asgs-edition-3/jul2021-jun2026/access-and-downloads/digital-boundary-files
Load a target raster, to set the baseline Coordinate Reference System (CRS) and resolution:
target_raster <- terra::rast("export/wc2.1_2.5m_bio_12.tif")Set the target CRS:
target_crs <- terra::crs(target_raster)Load australia as a vector file to bound observations, using the target crs from above:
australia <- terra::vect("export/AUS_2021_AUST_SHP_GDA2020/", crs=target_crs)Crop the extent of the target raster with the vector file above and confirm the data are as expected:
precip_data_cropped <- target_raster %>%
terra::crop(australia)
plot(precip_data_cropped)Project the AusTraits data into a spatial format and confirm the data plot as expected:
occurrence_data_species_family_vect <- terra::vect(trait_data,
geom = c("longitude", "latitude"),
crs = terra::crs(precip_data_cropped))
plot(occurrence_data_species_family_vect)Use the function terra::extract to extract climate variables for each AusTraits data point:
extracted_tmp <-
terra::extract(
x = precip_data_cropped, # climatic data
y = occurrence_data_species_family_vect, # AusTraits trait data
method = "simple",
cells = TRUE, bind = TRUE, ID = FALSE, xy = TRUE)Turn extracted values into a dataframe and separate data for herbaceous versus woody taxa:
extracted <-
terra::values(extracted_tmp) %>%
mutate(
value_herbs = ifelse(woodiness == "herbaceous", value, NA),
value_woody = ifelse(woodiness == "woody", value, NA)
)31.4 Plot your data
ggplot() +
geom_point(
data = extracted,
aes(x = wc2.1_2.5m_bio_12,
y = value_herbs),
shape = 16,
alpha = 0.4,
color = "darkseagreen3"
) +
geom_point(
data = extracted,
aes(x = wc2.1_2.5m_bio_12,
y = value_woody),
shape = 16,
alpha = 0.4,
color = "coral"
) +
scale_x_continuous(
trans = "log10"
) +
scale_y_continuous(
trans = "log10"
) +
labs(
x = "Annual precipitation (mm)",
y = "Leaf mass per area (g/m2)")