frescImport()
:
Import Frescalo outputs to the workspacefrescS_it()
:
Calculate recorder effort per grid square and time periodfrescTrends()
:
Estimate occurrence trends for each species over timefrescPescPlot()
:
Create plots of trends in relative occupancy for a speciesfrescP_ijt()
:
Estimate probability of occurrence for each species per grid square and
time period
In this document, we will go through a number of functions that we have written in order to extract some useful information and metrics from the output of a Frescalo analysis (Hill, 2012).
These functions use outputs from the
sparta
package, and this tutorial assumes that you have
run the frescalo()
function, that in itself requires some
data wrangling and other pre-functions to be run. See the sparta
vignettes for details.
If you have done that, then you should have a folder on your machine
containing the output from the frescalo()
function. The
folder will be named frescalo_yymmdd, according to the date on
which the analysis was carried out, and the location of the folder
should have been specified in the sinkdir
argument of the
fresalo()
function. Something like this.
The data used for this tutorial are from the island of Öland in southern Sweden. A delightful place with an interesting flora. For 295 vascular plant species, Historical observations (note: the link is to a scanned book) mapped in 1938 were digitised, and these were compared to observations from the modern flora that was published in 2024. Each observation was assigned to a 5 × 5 km grid cell using the old Swedish national grid RT90. As well as being used to assess changes in occurrences and distribution in the modern flora, comparisons of species’ trends using the Frescalo method have also been used in my research (for example Auffret & Svenning 2022). The data are available on the Swedish observation system ArtPortalen.
source("frescFun.R")
Here will follow descriptions of five functions.
frescImport()
: Import Frescalo outputs to the
workspaceTakes a folder containing frescalo outputs and makes list that can be used with rest of the functions.
frescImport(folder)
folder
The location of the folder
A list containing eight objects. The list resembles the output list
from the frescalo()
function, with many items assigned the
same names. However, not all of the items are actually used by
subsequent functions, and we also bring in some additional files that
are needed, but are otherwise not returned automatically by
frescalo()
. The aim is to find a balance between being able
to run the functions below, but also being able to explore the frescalo
output if needed.
out.path <- "data/frescalo_210805"
fresc.res <- frescImport(out.path)
lapply(fresc.res, head) # take a peek
## $paths
## [1] "data/frescalo_210805"
##
## $trend
## Species Time TFactor StDev X Xspt Xest N.0.00 N.0.98
## 1 Acer platanoides 1938 0.411 0.081 30 30 30 91 0
## 2 Acer platanoides 2010 0.568 0.091 76 76 76 91 60
## 3 Alchemilla subcrenata 1938 0.500 0.139 14 14 14 91 0
## 4 Alchemilla subcrenata 2010 0.254 0.073 13 13 13 91 0
## 5 Drosera rotundifolia 1938 0.940 0.284 12 12 12 89 0
## 6 Drosera rotundifolia 2010 0.178 0.091 4 4 4 89 0
##
## $stat
## Location Loc_no No_spp Phi_in Alpha Wgt_n2 Phi_out Spnum_in Spnum_out Iter
## 1 3G6h 1 170 0.693 3.57 13.78 0.85 139.4 215.1 5
## 2 3G6i 2 107 0.699 3.63 13.05 0.85 138.6 213.7 3
## 3 3G7h 3 161 0.748 2.39 17.99 0.85 165.0 216.4 5
## 4 3G7i 4 121 0.723 3.21 15.57 0.85 142.9 207.4 3
## 5 3G8h 5 177 0.764 2.10 15.93 0.85 174.5 218.3 6
## 6 3G8i 6 129 0.734 2.86 16.84 0.85 146.7 204.0 3
##
## $freq
## Location Species Pres Freq Freq1 SDFrq1 Rank Rank1
## 1 3G6h Dryopteris filix-mas 1 1 1 0 1 0.005
## 2 3G6h Epilobium hirsutum 1 1 1 0 2 0.009
## 3 3G6h Geranium dissectum 1 1 1 0 3 0.014
## 4 3G6h Neotinea ustulata 1 1 1 0 4 0.019
## 5 3G6h Orchis militaris 1 1 1 0 5 0.023
## 6 3G6h Pinus sylvestris 1 1 1 0 6 0.028
##
## $lm_stats
## SPECIES NAME b a b_std_err b_tval b_pval
## 1 S1 Acer platanoides 0.002180556 -3.814917 NA NA NA
## 2 S2 Adonis vernalis 0.001972222 -3.456167 NA NA NA
## 3 S3 Aira praecox -0.006736111 13.930583 NA NA NA
## 4 S4 Ajuga pyramidalis -0.002125000 4.640250 NA NA NA
## 5 S5 Alchemilla acutiloba -0.008930556 18.159417 NA NA NA
## 6 S6 Alchemilla filicaulis -0.007069444 14.257583 NA NA NA
## a_std_err a_tval a_pval adj_r2 r2 F_val F_num_df F_den_df Ymin Ymax
## 1 NA NA NA NA 1 NA 1 0 1938 2010
## 2 NA NA NA NA 1 NA 1 0 1938 2010
## 3 NA NA NA NA 1 NA 1 0 1938 2010
## 4 NA NA NA NA 1 NA 1 0 1938 2010
## 5 NA NA NA NA 1 NA 1 0 1938 2010
## 6 NA NA NA NA 1 NA 1 0 1938 2010
## Z_VAL SIG_95
## 1 1.2887048 FALSE
## 2 0.7465314 FALSE
## 3 -3.0297710 TRUE
## 4 -0.7672186 FALSE
## 5 -3.3278560 TRUE
## 6 -3.0932677 TRUE
##
## $log_file
## [1] " Log file for FRESCALO"
## [2] ""
## [3] " Input limits:"
## [4] " Number of samples = 10000"
## [5] " Number of species = 10000"
## [6] " Number of time periods = 100"
##
## $in_data
## V1 V2 V3
## 1 6H0e S1 1938
## 2 5H9e S1 1938
## 3 5H7d S1 1938
## 4 5H6d S1 1938
## 5 5H6e S1 1938
## 6 5H5d S1 1938
##
## $spe_codes
## SPECIES NAME
## 1 S1 Acer platanoides
## 2 S2 Adonis vernalis
## 3 S3 Aira praecox
## 4 S4 Ajuga pyramidalis
## 5 S5 Alchemilla acutiloba
## 6 S6 Alchemilla filicaulis
frescS_it()
: Calculate recorder effort per grid square
and time periodTakes an output list from frescImport()
and returns
recorder effort per grid square and time period. Recorder effort is
defined as the proportion of each grid square’s neighbourhood’s
benchmark species that were observed in the focal grid square during a
particular time period. This is known as Sit in Hill
(2012).
frescS_it(frescalo.results)
frescalo.results
A list output from the
frescImport()
function.
This function first uses the log file from the Frescalo output to
identify the proportion of benchmark species used for the Frescalo
analysis, known as R* in Hill (2012) and specified as the
alpha
argument in sparta::frescalo
. It then
for each grid square identifies the benchmark species for the
neighbourhood (which are common across all time periods) and calculates
the proportion that were observed in each time period.
A data frame with a column with grid square names, and an additional column for each time period in the analysis containing the proportion of benchmark species.
fresc.samp <- frescS_it(fresc.res)
fresc.samp # take a peek
By plotting these outputs, we see that the modern flora has a much higher sampling effort than the historical inventory. In the case of this data set, this is completely expected.
hist(fresc.samp[,2], main=names(fresc.samp)[2], xlim=c(0,1), xlab="Proportion benchmark species")
hist(fresc.samp[,3], main=names(fresc.samp)[3], xlim=c(0,1), xlab="Proportion benchmark species")
frescTrends()
: Estimate occurrence trends for each
species over timeTakes an output list from frescImport()
and returns
metrics of occurrence trends used in the literature.
frescTrends(frescalo.results, return.all=TRUE)
frescalo.results
A list output from the
frescImport()
function.
return.all
Logical. Specify if you not only want trend data for each species, but
also random draws of species trends for nice plotting (see details, as
well as frescPescPlot()
below)
This function returns metrics of relative occupancy and change in
relative occupancy for each species in each time period, as well as
changes over time and a measure of uncertainty, following some
publications that have used the Frescalo method. The term relative
occupancy is taken from Pescott et
al. 2022. It is originally known as Time factor in Hill, 2012, and
called Relative Reporting Rate in Fox et al. 2014. The
function extracts the relative occupancy values from the Frescalo
output, and applies a Z-test on the earliest and most recent time
periods of the data, according to the equation given in Fox et al. 2014.
This value of change over time is (comfortingly) identical to the value
‘b’ in the sparta()
output’s lm_stats object. The function
also calculates linear models based on 100 random draws from the
standard deviations of the relative occupancy values for each species
over time. The distributions of the slopes of these linear models are
split according to their strength into the following categories: “Strong
decline”, “Moderate decline”, “Stable”, “Moderate increase”, “Strong
increase”, according to Pescott et al. 2022.
If return.all=FALSE
a data frame with a column for
species names, columns relating to relative occupancy and z-tests, as
well as the proportion of linear model slopes that fall into the above
categories. If return.all=TRUE
, a list is returned, with
the first object trends
being the above-described data
frame. The second object is a list of the linear model outputs for each
species, and the third gives the time period names - both of these are
required for running the frescPescPlot()
function.
fresc.trends <- frescTrends(fresc.res, return.all=TRUE) # can take some seconds
length(fresc.trends) # just to show that the list contains three items, printing any of it would get a bit messy
## [1] 3
Now to look at the main trends data frame.
fresc.trends$trends
frescPescPlot()
: Create plots of trends in relative
occupancy for a speciesTakes a species name and an output from frescTrends()
and returns a nice plot in the style of Pescott et
al. 2022 and the wonderful Plant Atlas 2020 of Britain and
Ireland.
frescPescPlot(species, trends, point.col="black", line.col="forestgreen")
species
Character. The name of a species that has been
analysed using Frescalo, and for which trends have been estimated using
FrescTrends()
.
trends
An output list from
the FrescTrends()
function.
point.col
Character. The colour to be used for the plot’s points.
line.col
Character. The colour to be used for the plot’s
trend lines.
A plot with a point and whiskers showing the relative occupancy for
the species at each time period, and trend lines for the 100 linear
models run in FrescTrends()
Here we can look at trends in the grassland specialist Ajuga pyramidalis, or pyramidal bugle.
frescPescPlot("Ajuga pyramidalis", fresc.trends, point.col="black", line.col="forestgreen")
We can try another species, and play around with the colours a little. These plots can also be combined with a histogram of the strength of the trend lines calculated in the previous function, and also used in Pescott et al. 2022 and the plant atlas.
my.species <- "Pyrola minor"
frescPescPlot(my.species, fresc.trends, point.col="blue", line.col="brown")
# extract values from relevant columns in the trends data frame
strength.cols <- c("Strong decline", "Moderate decline", "Stable", "Moderate increase", "Strong increase")
my.species.num <- unlist(fresc.trends$trends[fresc.trends$trends$species==my.species ,strength.cols])
barplot(my.species.num, cex.names=0.75) # and plot
frescP_ijt()
: Estimate probability of occurrence for
each species per grid square and time periodTakes an output list from frescImport()
and returns
estimated probability of occurrence for each species per grid square and
time period, as implemented by Eichenberg et al. 2021.
This is known as Pijt in Hill (2012).
frescP_ijt(frescalo.results)
frescalo.results
A list output from the
frescImport()
function.
Equation 2 (and 3) in Hill (2012) estimates probability of occurrence of a given species in a given grid cell at a given time period, based on the species’ time-independent frequency in its neighbourhood (i.e. all observations across all time periods), recorder effort at the grid cell at the specific time period (Sit), and the relative occurrence across the whole study region at the specific time period. By assuming full recorder effort in a specific grid cell at each time period, Eichenberg et al. estimate each species’ probability of occurrence. While Hill uses grid-square probability of occurrence to estimate regional-level trends, Eichenberg uses these estimated regional trends to re-estimate grid-square occurrence probability, therefore assuming that a species’ grid-square occurrence reflects its regional-level relative occupancy. The combination of time-independent frequency and regional-level relative occurrence also means that a species will receive a non-zero probability of occurrence in a grid cell at a given time period if it has occurred somewhere in its neighbourhood during another time period. In an extreme example, a species that has disappeared from a target grid square from one time period to the next, but has a high time-independent neighbourhood frequency (for instance, it was ubiquitous in the first time period), and has increased regionally (relative occurrence is higher in the second time period), will still have a high probability of occurrence in the target grid square in the second time period according to this method.
A data frame with a column with grid square names, species names, and an additional column for each time period in the analysis containing the estimated probability of occurrence.
fresc.prob <- frescP_ijt(fresc.res) # can take some seconds
fresc.prob # take a peek
fresc.prob[fresc.prob$grid.square=="3G6h",] # look at a specific grid square
In Eichenberg et al. (2021), these probabilities of occurrences were summed per grid square to produce an estimate of species richness, but please see ‘Details’ in the box above.
est.rich <- aggregate(fresc.prob[,3:ncol(fresc.prob)], list(fresc.prob$grid.square), FUN="sum")
est.rich # take a peek
End of document