Introduction
The secpRod
package provides a number of methods to
calculate secondary production of populations, both those with clear
cohort structure and those where cohort structure is not possible to
discern. First, we walk through the available methods using simulated
populations that mimic the data structure of real sampling regimes. The
parameters to create these simulated populations is outlined in
companion article, vignette(“Simulating the sampling of populations”).
These simulated populations skip some of the steps of moving from
initial sampling processes (e.g., length to mass conversions, sample
subsetting, etc.) that are often part of the process. More complete
examples that showcase helper functions within secpRod
that
deal with these processes can be found below in More complete examples. Here, we just
showcase the secondary production methods available.
The methods within the package use two main objects,
sampleInfo
the sample-level information (e.g., site, date ,
replicate) and density and length or mass. An additional dataframe,
taxaInfo
that houses species-specific information such as
length to mass conversions, and method specific information such as
cohort production interval estimates (CPIs), production:biomass ratios
(PBs), and growth equations. The columns contain information for the
calculation of production. These include, but are not limited to:
taxonID: a character string that matches the name of taxonID from sampleInfo
massForm: a character string that is coercible to a formula for the conversion from length to mass (e.g.,
afdm_mg~(a*lengthClass^b)
)a: numeric variable for the coefficients used in massForm
b: numeric variable for the coefficients used in massForm
percAsh: numeric integer 0–100
-
method: character string for the method to use. Must be one of the following:
- ‘is’: increment-summation method. This cohort-based method calculates interval and annual production, , from field data as the sum of all interval production, , plus initial biomass, :
‘rs’: removal-summation method is similar to the increment-summation method but instead calculates the prodcution lost during a sampling interval as the product of the decrease in density, and the mean individual mass, plus the increase in biomass, .
‘sf’: size frequency method
growthForm: a character string that is coercible to a formula. This formula is parsed and will reference existing information (e.g., size information) if variables are included in the formula. Further, an additional data
min.cpi: integer of the minimum estimated cohort production interval for adjusting annual production estimates using the size-frequency method
max.cpi: integer of the maximum estimated cohort production interval for adjusting annual production estimates using the size-frequency method
pb: numeric of the production to biomass (PB) ratio for the specific taxa. This can take three forms 1) a single value, 2) a vector of values the same length as
bootNum
, and 3) a string of a distribution to randomly sample, not including then =
(e.g., ‘rnorm(mean = 5, sd = 0.5)’, ‘runif(min = 3, max = 8)’). The function will automatically sample thebootNum
values. Fair warning: We suggest you explicitly name the parameters in the function call string. Possible unwanted and unknown things may happen otherwise. The tests for this feature are still under development, it will not reject any nonsensical values (e.g, negative, NA, or Inf). If your distribution samples negative values you will get negative productions estimates–at least until themin.growth
feature is implemented.min.growth: a minimal value of production to assign when density is >0 and production estimates for a taxon are <0.
notes: notes for researcher use. This column will be maintained in output summaries.
…: additional information or data to pass to the function. This is mainly used to pass environmental data to inform growth rate information by included information such as temperature, food availability, etc.
taxaInfo <- data.frame(
taxonID = c("sppX"),
massForm = c("afdm_mg~(a*lengthClass^b)*percAsh"),
# a is NA for this example because we simulated mass growth and don't need to convert from length to mass
a = c(NA),
# b is NA for this example because we simulated mass growth and don't need to convert from length to mass
b = c(NA),
# percAsh is NA for this example because we simulated mass growth and don't need to convert from length to mass
percAsh = c(NA),
# method can accept one or more values. This allow comparisons among different methods.
method = c("is"),
g.a = c(NA),
growthForm = c("log(g_d) = 1 - 0.25*log(afdm_mg)"),
min.cpi = c(335),
max.cpi = c(365),
pb = c("runif(min = 3, max = 8)"),
min.growth = c(0.001),
notes = c("This is here for your use. No information will be used, but this column will be maintained in some summaries. See *** for more information.")
)
A multi-method comparison with simulated data
First, we load and visualize the simulated cohort data set, which contains a single species and follows the progression of density and distribution of individual masses through it life cycle.
The data set is available within the secpRod
package
with:
data('singleCohortSim', package = 'secpRod')
The sampleInfo
for this species includes:
#> # A tibble: 10 × 5
#> taxonID dateID repID density afdm_mg
#> <chr> <dbl> <int> <dbl> <list>
#> 1 sppX 1 1 463 <dbl [463]>
#> 2 sppX 1 2 458 <dbl [458]>
#> 3 sppX 1 3 479 <dbl [479]>
#> 4 sppX 1 4 422 <dbl [422]>
#> 5 sppX 1 5 639 <dbl [639]>
#> 6 sppX 1 6 575 <dbl [575]>
#> 7 sppX 1 7 518 <dbl [518]>
#> 8 sppX 1 8 677 <dbl [677]>
#> 9 sppX 1 9 566 <dbl [566]>
#> 10 sppX 1 10 468 <dbl [468]>
summary_stats <- singleCohortSim %>%
unnest(afdm_mg, keep_empty = TRUE) %>%
group_by(dateID) %>%
dplyr::summarise(
massMean = mean(afdm_mg, na.rm = TRUE),
massSD= sd(afdm_mg, na.rm = TRUE),
larvalDensityMean = mean(density, na.rm = TRUE)
)
sim_plot =
ggplot(summary_stats, aes(x = dateID)) +
stat_halfeye(data = splitCohortSim, aes(x = dateID, y = density),
color = 'green')+
stat_halfeye(data = splitCohortSim %>% unnest(afdm_mg), aes(x = dateID, y = afdm_mg*100),
color = 'red')+
geom_path(aes(y = larvalDensityMean), color = 'green') +
geom_path(aes(y = massMean * 100), color = 'red') +
scale_y_continuous(
name = "Larval Density",
sec.axis = sec_axis(~./100, name = "Mean Mass (mg)"),
) +
theme_minimal() +
labs(title = "Larval Density and Mean Mass over Time", x = "Day")+
theme(axis.title.y.right = element_text(color = 'red'))
More complete examples (#more-complete-examples)
Currently, the analysis is fully implemented for the size-frequency
and production:biomass ratio methods. This quick tutorial will outline
steps to calculate secondary production for whole communities. First,
load the package with library(secpRod)
.
The package comes with a number of data sets of macroinvertebrate species and community data:
simulated data set of a single univoltine species accessed with
data("univoltine")
a full community data set from Junker and Cross (2014) that can be accessed with
data("wbtData")
The simulated data set is a single data frame of artificially sampled size-frequency data through time.
data("univoltine", package = 'secpRod')
#clean this data set to align
head(univoltine)
#> # A tibble: 6 × 5
#> taxonID dateID repID lengthClass n_m2
#> <chr> <date> <chr> <dbl> <dbl>
#> 1 sppA 2023-01-31 1 0.5 0
#> 2 sppA 2023-01-31 1 1 0
#> 3 sppA 2023-01-31 1 1.5 0
#> 4 sppA 2023-01-31 1 2 0
#> 5 sppA 2023-01-31 1 2.5 0
#> 6 sppA 2023-01-31 1 3 0
For this simple example, we can create the required ‘taxaInfo’ object:
# clean simulated data to workable form
taxaInfo <- data.frame(
taxonID = c("sppA"),
massForm = c("afdm_mg~(a*lengthClass^b)*percAsh"),
a = c(0.0025),
b = c(2.692),
percAsh = c(0.958),
method = c("sf"),
g.a = c(NA),
growthForm = c("log(g_d) = 1 - 0.25*log(afdm_mg) - "),
min.cpi = c(335),
max.cpi = c(365),
pb = c("runif(min = 3, max = 8)"),
min.growth = c(0.001),
notes = c("This is here for your use. No information will be used, but this column will be maintained in some summaries. See *** for more information.")
)
This object represents the taxonomic information for calculating production. Each row contains information for a single species/taxon (in this case the only row because species == 1 ). The columns contain information for calculation of production. These include, but are not limited to:
taxonID: a character string that matches the name of taxonID from sampleInfo
massForm: a character string that is coercible to a formula for the conversion from length to mass (e.g.,
afdm_mg~(a*lengthClass^b)
)a: numeric variable for the coefficients used in massForm
b: numeric variable for the coefficients used in massForm
percAsh: numeric integer 0–100
method: character string for the method to use. Must be one of the following: ‘sf’ (this will be updated as methods are finalized.)
growthForm: description forthcoming…
min.cpi: integer of the minimum estimated cohort production interval for adjusting annual production estimates using the size-frequency method
max.cpi: integer of the maximum estimated cohort production interval for adjusting annual production estimates using the size-frequency method
pb: numeric of the production to biomass (PB) ratio for the specific taxa. This can take three forms 1) a single value, 2) a vector of values the same length as
bootNum
, and 3) a string of a distribution to randomly sample, not including then =
(e.g., ‘rnorm(mean = 5, sd = 0.5)’, ‘runif(min = 3, max = 8)’). The function will automatically sample thebootNum
values. Fair warning: We suggest you explicitly name the parameters in the function call string. Possible unwanted and unknown things may happen otherwise. The tests for this feature are still under development, it will not reject any nonsensical values (e.g, negative, NA, or Inf). If your distribution samples negative values you will get negative productions estimates–at least until themin.growth
feature is implemented.min.growth: description forthcoming…
notes: notes for researcher use. This column will be maintained in output summaries.
This ‘wbtData’ object is a list with two (2) elements:
sampleInfo: a list with 32 data frames. One (1) for each taxonomic entity
taxaInfo: a data frame with 32 rows and 13 columns. Each row represents the taxonomic information for a single taxonomic entity and the columns contain information for calculation of production as described above.
Examples
Single-species walkthrough
A quick walkthrough of the calculation of secondary production for a single species.
## isolate a single species data frame from sampleInfo object
acentrella <- wbtData[['sampleInfo']][['Acentrella sp.']]
## let's take a look at the data set
head(acentrella, 10)
#> # A tibble: 10 × 5
#> taxonID repID dateID lengthClass n_m2
#> <chr> <dbl> <dttm> <dbl> <dbl>
#> 1 Acentrella sp. 1 2008-07-17 00:00:00 0.5 0
#> 2 Acentrella sp. 1 2008-07-17 00:00:00 1 0
#> 3 Acentrella sp. 1 2008-07-17 00:00:00 2 62.5
#> 4 Acentrella sp. 1 2008-07-17 00:00:00 3 72.9
#> 5 Acentrella sp. 1 2008-07-17 00:00:00 4 72.9
#> 6 Acentrella sp. 1 2008-07-17 00:00:00 5 52.1
#> 7 Acentrella sp. 1 2008-07-17 00:00:00 6 0
#> 8 Acentrella sp. 1 2008-07-17 00:00:00 7 0
#> 9 Acentrella sp. 1 2008-07-17 00:00:00 8 0
#> 10 Acentrella sp. 1 2008-07-17 00:00:00 9 0
The data contain all the replicates of density and body length distributions in long format.
The first step is to convert length to mass for estimating biomass
patterns. Here we use the convert_length_to_mass()
function, which adds a column of the individual masses based on the
length-to-mass formula and coefficients in taxaInfo.
We can take a look at what this looks like:
acentrellaInfo = subset(taxaInfo, taxonID == "Acentrella sp.")
acentrellaInfo
#> # A tibble: 1 × 13
#> # Groups: taxonID [1]
#> taxonID massForm a b percAsh method g.a growthForm min.cpi max.cpi
#> <chr> <chr> <dbl> <dbl> <dbl> <chr> <lgl> <chr> <dbl> <dbl>
#> 1 Acentre… afdm_mg… 0.0076 2.69 4.9 sf NA "log(g_d)… 335 365
#> # ℹ 3 more variables: pb <dbl>, min.growth <dbl>, notes <chr>
Information on its use can be viewed with
?convert_length_to_mass
acentrellaMass = convert_length_to_mass(taxaSampleList = acentrella,
taxaInfo = wbtData[['taxaInfo']])
head(acentrellaMass)
#> # A tibble: 6 × 6
#> taxonID repID dateID lengthClass n_m2 afdm_mg
#> <chr> <dbl> <dttm> <dbl> <dbl> <dbl>
#> 1 Acentrella sp. 1 2008-07-17 00:00:00 0.5 0 0.00577
#> 2 Acentrella sp. 1 2008-07-17 00:00:00 1 0 0.0372
#> 3 Acentrella sp. 1 2008-07-17 00:00:00 2 62.5 0.240
#> 4 Acentrella sp. 1 2008-07-17 00:00:00 3 72.9 0.716
#> 5 Acentrella sp. 1 2008-07-17 00:00:00 4 72.9 1.55
#> 6 Acentrella sp. 1 2008-07-17 00:00:00 5 52.1 2.83
From here you can view the size frequency histograms using the
plot_cohorts()
function. Check out the function options
with ?plot_cohorts
.
plot_cohorts(taxaSampleListMass = acentrellaMass,
param = 'length',
massClass = 'afdm_mg')
These figures can be helpful for identifying cohort structures and getting ballpark cohort production intervals (CPI) for species when estimating production using the size frequency method. In the future, I hope to implement an approach to delineate cohorts and estimate growth parameters from size-frequency data, but this is a feature still under development.
The next step is to estimate production. Again, only size-frequency is fully operational but I will work to update the other methods soon!
The function calc_production()
is the workhorse function
that will estimate community production.
To apply it to our single species example, we input the sample information and taxa information along with how many bootstraps we would like, as:
calc_production(
taxaSampleListMass = acentrellaMass,
infoCols = c(1:3),
taxaInfo = acentrellaInfo,
bootNum = 10,
wrap = TRUE,
taxaSummary = 'full'
)
#> $P.boots
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]
#> P.ann.samp 7894.01 7079.587 6613.092 8324.897 7190.444 6737.9 8233.515
#> B.ann.samp 260.319 249.3877 264.5394 281.8144 246.7723 220.5484 290.1389
#> N.ann.samp 665.1926 663.8497 581.4014 701.705 596.9849 596.5369 792.4515
#> [,8] [,9] [,10]
#> P.ann.samp 8605.992 7213.12 6350.352
#> B.ann.samp 262.7356 249.5579 229.219
#> N.ann.samp 754.534 649.3973 508.9685
#>
#> $taxaSummary
#> $taxaSummary$summaryType
#> [1] "full"
#>
#> $taxaSummary$taxonID
#> [1] "Acentrella sp."
#>
#> $taxaSummary$method
#> [1] "sf"
#>
#> $taxaSummary$P.ann.samp
#> [1] 7370.212
#>
#> $taxaSummary$P.uncorr.samp
#> [1] 7067.326
#>
#> $taxaSummary$cpi
#> [1] 350
#>
#> $taxaSummary$pb
#> [1] 29.03778
#>
#> $taxaSummary$meanN
#> [1] 657.723
#>
#> $taxaSummary$meanB
#> [1] 253.8146
#>
#> $taxaSummary$meanIndMass
#> [1] 0.3858989
#>
#> $taxaSummary$Nmean
#> 0.5 1 2 3 4 5
#> 2008-07-17 0.000000 11.45833 48.95833 33.33333 39.583333 30.208333
#> 2008-08-13 0.000000 12.35703 62.91970 105.36194 56.619625 7.291667
#> 2008-09-23 53.028394 399.27515 545.08029 132.82990 7.291667 33.333333
#> 2008-10-26 78.012266 582.35127 666.79595 173.86885 3.882576 0.000000
#> 2008-11-23 1.041667 290.29353 555.27152 123.57480 15.136670 2.083333
#> 2008-12-17 208.701307 593.97944 287.82091 103.54200 14.204545 2.083333
#> 2009-02-03 20.014881 250.22631 142.14720 27.80258 4.166667 2.083333
#> 2009-04-07 90.277778 239.72446 332.67915 102.71991 27.272727 4.924242
#> 2009-05-16 0.000000 33.44907 125.79807 32.69121 18.750000 5.208333
#> 2009-06-23 0.000000 19.50921 73.99330 80.57125 23.639759 1.041667
#> 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
#> 2008-07-17 12.500000 2.083333 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 2.083333 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 30.208333 20.833333 3.125000 3.125 0 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 2.840909 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0.000000 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 2.083333 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 1.041667 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 6.875000 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 5.208333 0.000000 1.041667 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 8.798759 0.000000 0.000000 0.000 0 0 0 0 0 0 0 0 0 0 0
#> 21 22 23 24 25 26 27 28 29 30
#> 2008-07-17 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 0 0 0 0 0 0 0 0 0 0
#>
#> $taxaSummary$Nsd
#> 0.5 1 2 3 4 5
#> 2008-07-17 0.000000 29.64635 45.01854 27.25197 46.791335 31.614388
#> 2008-08-13 0.000000 24.55841 130.78970 151.69187 58.651399 9.882118
#> 2008-09-23 99.904573 640.40507 820.70373 166.25512 7.030714 34.652708
#> 2008-10-26 210.757419 434.06890 508.27227 144.97113 9.218577 0.000000
#> 2008-11-23 3.294039 195.29622 348.76995 107.84317 43.107151 4.392052
#> 2008-12-17 388.902788 748.41641 349.28272 114.55605 16.517014 6.588078
#> 2009-02-03 44.522062 535.29716 316.10626 53.51328 7.283395 6.588078
#> 2009-04-07 222.149872 399.20601 492.05437 170.97389 58.965796 9.319089
#> 2009-05-16 0.000000 102.16768 255.46931 57.94451 32.125093 11.251286
#> 2009-06-23 0.000000 31.63802 93.27213 107.03519 32.020550 3.294039
#> 6 7 8 9 10 11 12 13 14 15 16 17 18 19
#> 2008-07-17 15.372183 4.392052 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 4.392052 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 56.299275 32.200064 5.031728 7.030714 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 8.983743 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0.000000 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 6.588078 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 3.294039 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 21.740659 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 11.251286 0.000000 3.294039 0.000000 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 24.384816 0.000000 0.000000 0.000000 0 0 0 0 0 0 0 0 0 0
#> 20 21 22 23 24 25 26 27 28 29 30
#> 2008-07-17 0 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 0 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 0 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 0 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 0 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 0 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 0 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 0 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 0 0 0 0 0 0 0 0 0 0 0
#>
#> $taxaSummary$Bmean
#> 0.5 1 2 3 4 5
#> 2008-07-17 0.000000000 0.4267083 11.77357 23.86831 61.470320 85.519414
#> 2008-08-13 0.000000000 0.4601757 15.13101 75.44435 87.926563 20.642617
#> 2008-09-23 0.305806299 14.8690067 131.08164 95.11277 11.323480 94.366250
#> 2008-10-26 0.449884307 21.6867611 160.35199 124.49868 6.029385 0.000000
#> 2008-11-23 0.006007126 10.8105310 133.53244 88.48566 23.506256 5.897891
#> 2008-12-17 1.203547188 22.1197943 69.21556 74.14118 22.058727 5.897891
#> 2009-02-03 0.115422630 9.3184280 34.18375 19.90802 6.470560 5.897891
#> 2009-04-07 0.520617562 8.9273390 80.00313 73.55253 42.352756 13.940469
#> 2009-05-16 0.000000000 1.2456435 30.25209 23.40852 29.117520 14.744727
#> 2009-06-23 0.000000000 0.7265230 17.79401 57.69299 36.710995 2.948945
#> 6 7 8 9 10 11 12 13 14 15 16 17 18 19
#> 2008-07-17 57.799582 14.5857 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 9.633264 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 139.682324 145.8570 31.33826 43.02554 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 13.136269 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0.000000 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 9.633264 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 4.816632 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 31.789770 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 24.083159 0.0000 10.44609 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 40.685167 0.0000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 20 21 22 23 24 25 26 27 28 29 30
#> 2008-07-17 0 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 0 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 0 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 0 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 0 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 0 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 0 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 0 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 0 0 0 0 0 0 0 0 0 0 0
#>
#> $taxaSummary$Bsd
#> 0.5 1 2 3 4 5
#> 2008-07-17 0.0000000 1.1040302 10.82612 19.51375 72.66387 89.499937
#> 2008-08-13 0.0000000 0.9145551 31.45248 108.61887 91.08177 27.976152
#> 2008-09-23 0.5761338 23.8486848 197.36394 119.04687 10.91824 98.101384
#> 2008-10-26 1.2154045 16.1647259 122.23000 103.80648 14.31585 0.000000
#> 2008-11-23 0.0189962 7.2728313 83.87267 77.22104 66.94258 12.433845
#> 2008-12-17 2.2427404 27.8710270 83.99598 82.02779 25.64984 18.650768
#> 2009-02-03 0.2567516 19.9344664 76.01766 38.31815 11.31063 18.650768
#> 2009-04-07 1.2811029 14.8664317 118.32990 122.42575 91.57001 26.382223
#> 2009-05-16 0.0000000 3.8047243 61.43560 41.49113 49.88816 31.852250
#> 2009-06-23 0.0000000 1.1782000 22.43021 76.64248 49.72581 9.325384
#> 6 7 8 9 10 11 12 13 14 15 16 17 18 19
#> 2008-07-17 71.08046 30.74935 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 20.30870 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 260.32597 225.43702 50.45939 96.80008 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 41.54053 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0.00000 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 30.46305 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 15.23153 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 100.52808 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 52.02557 0.00000 33.03342 0.00000 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 112.75457 0.00000 0.00000 0.00000 0 0 0 0 0 0 0 0 0 0
#> 20 21 22 23 24 25 26 27 28 29 30
#> 2008-07-17 0 0 0 0 0 0 0 0 0 0 0
#> 2008-08-13 0 0 0 0 0 0 0 0 0 0 0
#> 2008-09-23 0 0 0 0 0 0 0 0 0 0 0
#> 2008-10-26 0 0 0 0 0 0 0 0 0 0 0
#> 2008-11-23 0 0 0 0 0 0 0 0 0 0 0
#> 2008-12-17 0 0 0 0 0 0 0 0 0 0 0
#> 2009-02-03 0 0 0 0 0 0 0 0 0 0 0
#> 2009-04-07 0 0 0 0 0 0 0 0 0 0 0
#> 2009-05-16 0 0 0 0 0 0 0 0 0 0 0
#> 2009-06-23 0 0 0 0 0 0 0 0 0 0 0
#>
#> $taxaSummary$datesInfo
#> dateID N n_m2_mean afdm_mg_m2_mean
#> 1 2008-07-17 10 178.1250 255.4436
#> 2 2008-08-13 10 246.6333 209.2380
#> 3 2008-09-23 10 1228.1304 706.9621
#> 4 2008-10-26 10 1507.7518 326.1530
#> 5 2008-11-23 10 987.4015 262.2388
#> 6 2008-12-17 10 1212.4149 204.2700
#> 7 2009-02-03 10 447.4826 80.7107
#> 8 2009-04-07 10 804.4733 251.0866
#> 9 2009-05-16 10 222.1467 133.2977
#> 10 2009-06-23 10 207.5540 156.5586
#> 11 2009-07-16 NA 192.8395 206.0011
This will output a lot of information though it has a simple structure:
P.boots–vectors of bootstrapped estimates of annual production, annual biomass, and annual abundance.
taxaSummary–this comes in ‘full’ and ‘short’ versions. As you see above the full can be a lot. The ‘short’ is a paired down version of this
calc_production(
taxaSampleListMass = acentrellaMass,
infoCols = c(1:3),
taxaInfo = acentrellaInfo,
bootNum = 10,
wrap = TRUE,
taxaSummary = 'short'
)
#> $P.boots
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]
#> P.ann.samp 6577.237 7689.972 8032.055 7666.049 7863.849 7524.68 7949.663
#> B.ann.samp 240.9445 270.1543 284.5143 232.0916 291.3699 249.4617 274.3334
#> N.ann.samp 588.0506 694.8745 691.2424 639.103 809.3192 647.6252 729.6995
#> [,8] [,9] [,10]
#> P.ann.samp 6265.627 7721.312 7122.196
#> B.ann.samp 214.3935 269.127 252.255
#> N.ann.samp 451.6489 691.9201 685.241
#>
#> $taxaSummary
#> $taxaSummary$summaryType
#> [1] "short"
#>
#> $taxaSummary$taxonID
#> [1] "Acentrella sp."
#>
#> $taxaSummary$method
#> [1] "sf"
#>
#> $taxaSummary$P.ann.samp
#> [1] 7370.212
#>
#> $taxaSummary$cpi
#> [1] 350
#>
#> $taxaSummary$pb
#> [1] 29.03778
#>
#> $taxaSummary$meanN
#> [1] 657.723
#>
#> $taxaSummary$meanB
#> [1] 253.8146
#>
#> $taxaSummary$meanIndMass
#> [1] 0.3858989
#>
#> $taxaSummary$datesInfo
#> dateID N n_m2_mean afdm_mg_m2_mean
#> 1 2008-07-17 10 178.1250 255.4436
#> 2 2008-08-13 10 246.6333 209.2380
#> 3 2008-09-23 10 1228.1304 706.9621
#> 4 2008-10-26 10 1507.7518 326.1530
#> 5 2008-11-23 10 987.4015 262.2388
#> 6 2008-12-17 10 1212.4149 204.2700
#> 7 2009-02-03 10 447.4826 80.7107
#> 8 2009-04-07 10 804.4733 251.0866
#> 9 2009-05-16 10 222.1467 133.2977
#> 10 2009-06-23 10 207.5540 156.5586
#> 11 2009-07-16 NA 192.8395 206.0011
You can also turn this off with ‘none’.
An example using the PB method
The pb method allows for multiple options to set the production to biomass ratio:
A single value. This will use the same value for each bootstrap sample. The variability is derived entirely from resampling the sample units (e.g. surbers).
a character string for a distribution from which to sample excluding the
n=
, such asrnorm(mean = 5, sd = 0.5)
. The function will automatically samplebootNum
values.a numeric vector the same length as
bootNum
. This is a more explicit option to set all the values of pb that will be used. In the future, if the lengthbootNum
an error will be returned for the taxa.
acentrellaInfo$method = "pb"
acentrellaInfo$pb = "runif(min = 3, max = 8)"
calc_production(
taxaSampleListMass = acentrellaMass,
infoCols = c(1:3),
taxaInfo = acentrellaInfo,
bootNum = 10,
wrap = TRUE,
taxaSummary = 'short'
)
#> $P.boots
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]
#> P.ann.samp 1088.466 1712.322 1742.765 1987.615 1809.562 2060.79 1218.325
#> B.ann.samp 270.4973 226.5819 251.4287 259.6744 259.8879 290.8709 283.018
#> N.ann.samp 686.2019 607.9946 633.1745 629.4606 639.5113 736.3917 806.5242
#> [,8] [,9] [,10]
#> P.ann.samp 1831.69 1489.918 1000.54
#> B.ann.samp 269.8947 256.5597 261.3129
#> N.ann.samp 720.5662 697.1245 774.2984
#>
#> $taxaSummary
#> $taxaSummary$summaryType
#> [1] "short"
#>
#> $taxaSummary$taxonID
#> [1] "Acentrella sp."
#>
#> $taxaSummary$method
#> [1] "pb"
#>
#> $taxaSummary$P.ann.samp
#> [1] 1968.429
#>
#> $taxaSummary$pb
#> [1] 7.755384
#>
#> $taxaSummary$meanN
#> [1] 657.723
#>
#> $taxaSummary$meanB
#> [1] 253.8146
#>
#> $taxaSummary$meanIndMass
#> [1] 0.3858989
#>
#> $taxaSummary$datesInfo
#> dateID N n_m2_mean afdm_mg_m2_mean
#> 1 2008-07-17 10 178.1250 255.4436
#> 2 2008-08-13 10 246.6333 209.2380
#> 3 2008-09-23 10 1228.1304 706.9621
#> 4 2008-10-26 10 1507.7518 326.1530
#> 5 2008-11-23 10 987.4015 262.2388
#> 6 2008-12-17 10 1212.4149 204.2700
#> 7 2009-02-03 10 447.4826 80.7107
#> 8 2009-04-07 10 804.4733 251.0866
#> 9 2009-05-16 10 222.1467 133.2977
#> 10 2009-06-23 10 207.5540 156.5586
#> 11 2009-07-16 NA 192.8395 206.0011
Applying this to multiple taxa
The current internal code to do this with just
calc_production
is still buggy. But you can run production
by splitting the full sampleInfo into a list by taxa and using
apply
or functions from the (purrr
package)[https://purrr.tidyverse.org/]
The code below shows how to accomplish this. Again, only size frequency and pb are accepted currently.
# load needed packages
library(magrittr)
library(purrr)
# The sampleInfo is already a list by taxonID
str(wbtData[['sampleInfo']], list.len = 5)
#> list<tibble[,5]> [1:32]
#> $ Acentrella sp. : tibble [3,100 × 5] (S3: tbl_df/tbl/data.frame)
#> ..$ taxonID : chr [1:3100] "Acentrella sp." "Acentrella sp." "Acentrella sp." "Acentrella sp." ...
#> ..$ repID : num [1:3100] 1 1 1 1 1 1 1 1 1 1 ...
#> ..$ dateID : POSIXct[1:3100], format: "2008-07-17" "2008-07-17" ...
#> ..$ lengthClass: num [1:3100] 0.5 1 2 3 4 5 6 7 8 9 ...
#> ..$ n_m2 : num [1:3100] 0 0 62.5 72.9 72.9 ...
#> $ Ameletus spp. : tibble [3,038 × 5] (S3: tbl_df/tbl/data.frame)
#> ..$ taxonID : chr [1:3038] "Ameletus spp." "Ameletus spp." "Ameletus spp." "Ameletus spp." ...
#> ..$ repID : num [1:3038] 1 1 1 1 1 1 1 1 1 1 ...
#> ..$ dateID : POSIXct[1:3038], format: "2008-07-17" "2008-07-17" ...
#> ..$ lengthClass: num [1:3038] 0.5 1 2 3 4 5 6 7 8 9 ...
#> ..$ n_m2 : num [1:3038] 0 0 0 0 0 0 0 0 0 0 ...
#> $ Attenella spp. : tibble [3,038 × 5] (S3: tbl_df/tbl/data.frame)
#> ..$ taxonID : chr [1:3038] "Attenella spp." "Attenella spp." "Attenella spp." "Attenella spp." ...
#> ..$ repID : num [1:3038] 1 1 1 1 1 1 1 1 1 1 ...
#> ..$ dateID : POSIXct[1:3038], format: "2008-07-17" "2008-07-17" ...
#> ..$ lengthClass: num [1:3038] 0.5 1 2 3 4 5 6 7 8 9 ...
#> ..$ n_m2 : num [1:3038] 0 0 41.7 20.8 0 ...
#> $ Baetis spp. : tibble [3,100 × 5] (S3: tbl_df/tbl/data.frame)
#> ..$ taxonID : chr [1:3100] "Baetis spp." "Baetis spp." "Baetis spp." "Baetis spp." ...
#> ..$ repID : num [1:3100] 1 1 1 1 1 1 1 1 1 1 ...
#> ..$ dateID : POSIXct[1:3100], format: "2008-07-17" "2008-07-17" ...
#> ..$ lengthClass: num [1:3100] 0.5 1 2 3 4 5 6 7 8 9 ...
#> ..$ n_m2 : num [1:3100] 0 0 20.8 31.2 93.8 ...
#> $ Barbaetis spp. : tibble [3,100 × 5] (S3: tbl_df/tbl/data.frame)
#> ..$ taxonID : chr [1:3100] "Barbaetis spp." "Barbaetis spp." "Barbaetis spp." "Barbaetis spp." ...
#> ..$ repID : num [1:3100] 1 1 1 1 1 1 1 1 1 1 ...
#> ..$ dateID : POSIXct[1:3100], format: "2008-07-17" "2008-07-17" ...
#> ..$ lengthClass: num [1:3100] 0.5 1 2 3 4 5 6 7 8 9 ...
#> ..$ n_m2 : num [1:3100] 0 0 20.8 0 0 ...
#> [list output truncated]
#> @ ptype: tibble [0 × 5] (S3: tbl_df/tbl/data.frame)
#> ..$ taxonID : chr(0)
#> ..$ repID : num(0)
#> ..$ dateID : 'POSIXct' num(0)
#> - attr(*, "tzone")= chr "UTC"
#> ..$ lengthClass: num(0)
#> ..$ n_m2 : num(0)
# if you have a large data frame with all taxa the code below shows how to do this:
## install a package to split into a named_list
# devtools::install_github('jimjunker1/junkR')
# library(junkR)
## Convert to data frame
# sampleInfoDf = wbtData[['sampleInfo']] %>% bind_rows
## Convert back to named list split by taxonID
# sampleInfoList = sampleInfoDf %>% junkR::named_group_split(taxonID)
sampleInfoList = wbtData[['sampleInfo']]
sampleInfoListMass = sampleInfoList %>%
purrr::map(
~convert_length_to_mass(.x,
taxaInfo = taxaInfo)
)
fullProduction = sampleInfoListMass %>%
purrr::map2(., list(taxaInfo),
~calc_production(.x,
infoCols = c(1:3),
taxaInfo = .y,
bootNum = 10,
wrap = TRUE,
taxaSummary = 'short'
)
)
str(fullProduction, list.len = 5)
#> List of 32
#> $ Acentrella sp. :List of 2
#> ..$ P.boots :List of 30
#> .. ..$ : num 6968
#> .. ..$ : num 244
#> .. ..$ : num 612
#> .. ..$ : num 7217
#> .. ..$ : num 257
#> .. .. [list output truncated]
#> .. ..- attr(*, "dim")= int [1:2] 3 10
#> .. ..- attr(*, "dimnames")=List of 2
#> .. .. ..$ : chr [1:3] "P.ann.samp" "B.ann.samp" "N.ann.samp"
#> .. .. ..$ : NULL
#> ..$ taxaSummary:List of 10
#> .. ..$ summaryType: chr "short"
#> .. ..$ taxonID : chr "Acentrella sp."
#> .. ..$ method : chr "sf"
#> .. ..$ P.ann.samp : num 7370
#> .. ..$ cpi : num 350
#> .. .. [list output truncated]
#> $ Ameletus spp. :List of 2
#> ..$ P.boots :List of 30
#> .. ..$ : num 796
#> .. ..$ : num 41.6
#> .. ..$ : num 31.7
#> .. ..$ : num 736
#> .. ..$ : num 41.9
#> .. .. [list output truncated]
#> .. ..- attr(*, "dim")= int [1:2] 3 10
#> .. ..- attr(*, "dimnames")=List of 2
#> .. .. ..$ : chr [1:3] "P.ann.samp" "B.ann.samp" "N.ann.samp"
#> .. .. ..$ : NULL
#> ..$ taxaSummary:List of 10
#> .. ..$ summaryType: chr "short"
#> .. ..$ taxonID : chr "Ameletus spp."
#> .. ..$ method : chr "sf"
#> .. ..$ P.ann.samp : num 660
#> .. ..$ cpi : num 350
#> .. .. [list output truncated]
#> $ Attenella spp. :List of 2
#> ..$ P.boots :List of 30
#> .. ..$ : num 2617
#> .. ..$ : num 112
#> .. ..$ : num 71.7
#> .. ..$ : num 2814
#> .. ..$ : num 125
#> .. .. [list output truncated]
#> .. ..- attr(*, "dim")= int [1:2] 3 10
#> .. ..- attr(*, "dimnames")=List of 2
#> .. .. ..$ : chr [1:3] "P.ann.samp" "B.ann.samp" "N.ann.samp"
#> .. .. ..$ : NULL
#> ..$ taxaSummary:List of 10
#> .. ..$ summaryType: chr "short"
#> .. ..$ taxonID : chr "Attenella spp."
#> .. ..$ method : chr "sf"
#> .. ..$ P.ann.samp : num 2534
#> .. ..$ cpi : num 350
#> .. .. [list output truncated]
#> $ Baetis spp. :List of 2
#> ..$ P.boots :List of 30
#> .. ..$ : num 3568
#> .. ..$ : num 139
#> .. ..$ : num 221
#> .. ..$ : num 3821
#> .. ..$ : num 133
#> .. .. [list output truncated]
#> .. ..- attr(*, "dim")= int [1:2] 3 10
#> .. ..- attr(*, "dimnames")=List of 2
#> .. .. ..$ : chr [1:3] "P.ann.samp" "B.ann.samp" "N.ann.samp"
#> .. .. ..$ : NULL
#> ..$ taxaSummary:List of 10
#> .. ..$ summaryType: chr "short"
#> .. ..$ taxonID : chr "Baetis spp."
#> .. ..$ method : chr "sf"
#> .. ..$ P.ann.samp : num 3833
#> .. ..$ cpi : num 350
#> .. .. [list output truncated]
#> $ Barbaetis spp. :List of 2
#> ..$ P.boots :List of 30
#> .. ..$ : num 1787
#> .. ..$ : num 66.4
#> .. ..$ : num 92.4
#> .. ..$ : num 1495
#> .. ..$ : num 52.7
#> .. .. [list output truncated]
#> .. ..- attr(*, "dim")= int [1:2] 3 10
#> .. ..- attr(*, "dimnames")=List of 2
#> .. .. ..$ : chr [1:3] "P.ann.samp" "B.ann.samp" "N.ann.samp"
#> .. .. ..$ : NULL
#> ..$ taxaSummary:List of 10
#> .. ..$ summaryType: chr "short"
#> .. ..$ taxonID : chr "Barbaetis spp."
#> .. ..$ method : chr "sf"
#> .. ..$ P.ann.samp : num 1512
#> .. ..$ cpi : num 350
#> .. .. [list output truncated]
#> [list output truncated]