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
<- austraits::get_versions()[["doi"]][1]
most_recent
most_recent
[1] "10.5281/zenodo.11188867"
<- austraits::load_austraits(doi = most_recent) austraits
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 $traits %>%
austraits::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) dplyr
Second, extract the numeric trait data you want to analyse. Only field-collected data on adult plants with reported coordinates is used:
<-
trait_data %>% join_locations())$traits %>%
(austraits rename(latitude = `latitude (deg)`, longitude = `longitude (deg)`) %>%
::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(
dplyrID = row_number(),
latitude = as.numeric(latitude),
longitude = as.numeric(longitude)
%>%
) distinct(value, dataset_id, location_name, taxon_name, .keep_all = TRUE) %>%
::filter(!is.na(latitude)) %>%
dplyr::filter(!is.na(longitude)) %>%
dplyr::left_join(woody_taxa, by = "taxon_name") dplyr
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:
<- terra::rast("export/wc2.1_2.5m_bio_12.tif") target_raster
Set the target CRS:
<- terra::crs(target_raster) target_crs
Load australia as a vector file to bound observations, using the target crs from above:
<- terra::vect("export/AUS_2021_AUST_SHP_GDA2020/", crs=target_crs) australia
Crop the extent of the target raster with the vector file above and confirm the data are as expected:
<- target_raster %>%
precip_data_cropped ::crop(australia)
terra
plot(precip_data_cropped)
Project the AusTraits data into a spatial format and confirm the data plot as expected:
<- terra::vect(trait_data,
occurrence_data_species_family_vect 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 ::extract(
terrax = 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 ::values(extracted_tmp) %>%
terramutate(
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)")