DEB_euler: Dynamic Energy Budget model

View source: R/DEB_euler.R

DEB_eulerR Documentation

Dynamic Energy Budget model

Description

Implementation of the Standarad Dynamic Energy Budget model of Kooijman Note does not do a proper integration, i.e. rate of change is the difference between two fixed time periods which is ok as long as the timestep is very small relative to the dynamics - hourly time steps or shorter are ok for long-lived (months to years) animals Michael Kearney Dec 2015

Usage

DEB_euler(
  step = 1/24,
  z = 7.997,
  del_M = 0.242,
  p_Xm = 13290 * step,
  kap_X = 0.85,
  v = 0.065 * step,
  kap = 0.886,
  p_M = 32 * step,
  p_T = 0,
  E_G = 7767,
  kap_R = 0.95,
  k_J = 0.002 * step,
  E_Hb = 73590,
  E_Hj = E_Hb,
  E_Hp = 186500,
  h_a = 2.16e-11/(step^2),
  s_G = 0.01,
  T_REF = 20 + 273.15,
  T_A = 8085,
  T_AL = 18721,
  T_AH = 90000,
  T_L = 288,
  T_H = 315,
  T_A2 = T_A,
  T_AL2 = T_AL,
  T_AH2 = T_AH,
  T_L2 = T_L,
  T_H2 = T_H,
  E_0 = 1040000,
  f = 1,
  E_sm = 1116,
  K = 1,
  andens_deb = 1,
  d_V = 0.3,
  d_E = 0.3,
  d_Egg = 0.3,
  stoich_mode = 0,
  mu_X = 525000,
  mu_E = 585000,
  mu_V = 5e+05,
  mu_P = 480000,
  mu_N = 244000/5,
  kap_X_P = 0.1,
  n_X = c(1, 1.8, 0.5, 0.15),
  n_E = c(1, 1.8, 0.5, 0.15),
  n_V = c(1, 1.8, 0.5, 0.15),
  n_P = c(1, 1.8, 0.5, 0.15),
  n_M_nitro = c(1, 4/5, 3/5, 4/5),
  h_N = 384238,
  clutchsize = 2,
  clutch_ab = c(0.085, 0.7),
  minclutch = 0,
  batch = 1,
  lambda = 1/2,
  VTMIN = 26,
  VTMAX = 39,
  ma = 1e-04,
  mi = 0,
  mh = 0.5,
  arrhenius = matrix(data = matrix(data = c(rep(T_A, 8), rep(T_AL, 8), rep(T_AH, 8),
    rep(T_L, 8), rep(T_H, 8)), nrow = 8, ncol = 5), nrow = 8, ncol = 5),
  arrhenius2 = matrix(data = matrix(data = c(rep(T_A2, 8), rep(T_AL2, 8), rep(T_AH2,
    8), rep(T_L2, 8), rep(T_H2, 8)), nrow = 8, ncol = 5), nrow = 8, ncol = 5),
  acthr = 1,
  X = 10,
  E_pres = E_0/3e-09,
  V_pres = 3e-09,
  E_H_pres = 0,
  q_pres = 0,
  hs_pres = 0,
  p_surv = 1,
  E_s_pres = 0,
  p_B_pres = 0,
  E_R = 0,
  E_B = 0,
  stages = 7,
  stage = 0,
  breeding = 0,
  Tb = 33,
  starve_mode = 1,
  fdry = 0.3,
  L_b = 0.42,
  L_j = 1.376,
  S_instar = rep(1.618, stages),
  spawnday = 1,
  day = 1,
  metab_mode = 0,
  age = 0
)

Arguments

step

= 1/24, step size (days)

z

= 7.997, Zoom factor (cm)

del_M

= 0.242, Shape coefficient (-)

p_Xm

= 13290*step, Surface area-specific maximum feeding rate J/cm2/h

kap_X

= 0.85, Digestive efficiency (decimal %)

v

= 0.065*step, Energy conductance (cm/h)

kap

= 0.886, fraction of mobilised reserve allocated to soma

p_M

= 32*step, Volume-specific somatic maintenance (J/cm3/h)

p_T

= 0, (Structural-)Surface-area-specific heating cost (J/cm2/h)

E_G

= 7767, Cost of structure (J/cm3)

kap_R

= 0.95, Fraction of reproduction energy fixed in eggs

k_J

= 0.002*step, Maturity maintenance rate coefficient (1/h)

E_Hb

= 7.359e+04, Maturity at birth (J)

E_Hj

= E_Hb, Maturity at metamorphosis (J)

E_Hp

= 1.865e+05, Maturity at puberty

h_a

= 2.16e-11*(step^2), Weibull ageing acceleration (1/h2)

s_G

= 0.01, Gompertz stress coefficient (-)

T_REF

= 20+273.15, Reference temperature for rate correction (deg C)

T_A

= 8085 Arrhenius temperature

T_AL

= 18721, Arrhenius temperature for decrease below lower boundary of tolerance range T_L

T_AH

= 90000, Arrhenius temperature for decrease above upper boundary of tolerance range T_H

T_L

= 288, Lower boundary (K) of temperature tolerance range for Arrhenius thermal response

T_H

= 315, Upper boundary (K) of temperature tolerance range for Arrhenius thermal response

T_A2

= 8085 Arrhenius temperature for maturity maintenance (causes 'Temperature Size Rule' effect)

T_AL2

= 18721, Arrhenius temperature for decrease below lower boundary of tolerance range T_L for maturity maintenance (causes 'Temperature Size Rule' effect)

T_AH2

= 90000, Arrhenius temperature for decrease above upper boundary of tolerance range T_H for maturity maintenance (causes 'Temperature Size Rule' effect)

T_L2

= 288, Lower boundary (K) of temperature tolerance range for Arrhenius thermal response for maturity maintenance (causes 'Temperature Size Rule' effect)

T_H2

= 315, Upper boundary (K) of temperature tolerance range for Arrhenius thermal response for maturity maintenance (causes 'Temperature Size Rule' effect)

E_0

= 1.04e+06, Energy content of the egg (derived from core parameters) (J)

f

= 1, functional response (-), usually kept at 1 because gut model controls food availability such that f=0 when gut empty

E_sm

= 1116, Maximum volume-specific energy density of stomach (J/cm3)

K

= 500, Half saturation constant (#/cm2)

andens_deb

= 1, Animal density (g/cm3)

d_V

= 0.3, Dry mass fraction of structure

d_E

= 0.3, Dry mass fraction of reserve

d_Egg

= 0.3, Dry mass fraction of egg

stoich_mode

= 0, adjust chemical indices to chemical potentials (0) or vice versa (1), or leave as is (2)

mu_X

= 525000, Molar Gibbs energy (chemical potential) of food (J/mol)

mu_E

= 585000, Molar Gibbs energy (chemical potential) of reserve (J/mol)

mu_V

= 500000, Molar Gibbs energy (chemical potential) of structure (J/mol)

mu_P

= 480000, Molar Gibbs energy (chemical potential) of faeces (J/mol)

mu_N

= 244e3/5, Molar Gibbs energy (chemical potential) of nitrogenous waste (J/mol), synthesis from NH3, Withers page 119

kap_X_P

= 0.1, Faecation efficiency of food to faeces (-)

n_X

= c(1,1.8,0.5,.15), Chem. indices of C, O, H and N in food

n_E

= c(1,1.8,0.5,.15), Chem. indices of C, O, H and N in reserve

n_V

= c(1,1.8,0.5,.15), Chem. indices of C, O, H and N in structure

n_P

= c(1,1.8,0.5,.15), Chem. indices of C, O, H and N in faeces

n_M_nitro

= c(1,4/5,3/5,4/5), Chem. indices of C, O, H and N in nitrogenous waste

h_N

= 384238, molar enthalpy of nitrogenous waste (combustion frame of reference) (J/mol), overridden if n_M_nitro specified as urea, uric acid or ammonia

clutchsize

= 2, Clutch size (#), overridden by clutch_ab

clutch_ab

= c(0,0), paramters for relationship between length (cm) and clutch size: clutch size = a*L_w-b, make a and b zero if fixed clutch size

minclutch

= 0, Minimum clutch size if not enough in reproduction buffer for clutch size predicted by clutch_ab - if zero, will not operate

batch

= 1, Invoke Pequerie et al.'s batch laying model?

lambda

= 1/2

VTMIN

= 26, Voluntary thermal maximum, degrees C, controls whether repro event can occur at a given time

VTMAX

= 39, Voluntary thermal maximum, degrees C, controls whether repro event can occur at a given time

arrhenius

= matrix(data = matrix(data = c(rep(T_A,8),rep(T_AL,8),rep(T_AH,8),rep(T_L,8),rep(T_H,8)), nrow = 8, ncol = 5), nrow = 8, ncol = 5), Stage-specific 5-parameter Arrhenius thermal response for DEB model (T_A,T_AL,T_AH,T_L,T_H)

arrhenius2

= matrix(data = matrix(data = c(rep(T_A2,8),rep(T_AL2,8),rep(T_AH2,8),rep(T_L2,8),rep(T_H2,8)), nrow = 8, ncol = 5), nrow = 8, ncol = 5), Stage-specific 5-parameter Arrhenius thermal response for maturity maintenance (causes 'Temperature Size Rule' effect) DEB model (T_A2,T_AL2,T_AH2,T_L2,T_H2)

acthr

= 1

X

= 11

E_pres

= E_0/3e-9

V_pres

= 3e-9

E_H_pres

= 0

q_pres

= 0

hs_pres

= 0

p_surv

= 1

E_s_pres

= 0

E_R

= 0

E_B

= 0

stages

= 7, how many life stages?

stage

= 0

breeding

= 0

Tb

= 33

starve_mode

= 1, Determines how reproduction buffer is used during starvation, where 0 means it is not used, 1 means it is used before structure is mobilised and 2 means it is used to maximise reserve density

fdry

= 0.3, Dry mass fraction of food

S_instar

= rep(1.6, stages), stress at instar n: L_n^2/ L_n-1^2 (-)

p_B_past

= 0

Value

stage Life cycle stage, -

V Structure, cm^3

E Reserve density, J/cm^3

E_H Maturity, J

E_s Stomach energy content, J

E_R Reproduction buffer energy, J

E_B Reproduction batch energy, J

q Ageing acceleration, 1/time^2

hs Hazard rate

length Physical length, cm

wetmass Total wet mass, g

wetgonad Wet mass of gonad, g

wetgut Wet mass of food in gut, g

wetstorage Wet mass of reserve, g

p_surv Survival probability, -

fecundity Eggs produced at a given time point, #

clutches Clutches produced at a given time point

JMO2 Oxygen flux, mol/time

JMCO2 Carbon dioxide flux, mol/time

JMH2O metabolic water flux, mol/time

JMNWASTE nitrogenous waste flux, mol/time

O2ML Oxgen consumption rate, ml/hour

CO2ML Carbon dioxide production rate, ml/time

GH2OML Metabolic water production rate, ml/time

DEBQMET Metabolic heat generation, J/time

GDRYFOOD Dry food intake, g/time

GFAECES Faeces production, g/time

GNWASTE Nitrogenous waste production, g/time

p_A Assimilation power, J/time

p_C Catabolic power, J/time

p_M Somatic maintenance power, J/time

p_G Growth power, J/time

p_D Dissipation power, J/time

p_J Maturity power, J/time

p_R Reproduction power, J/time

p_B Reproduction batch power, J/time

L_b Structural length at birth, cm

L_j Structural length at end of metabolic acceleration (if occurring), cm

Examples

# simulate growth and reproduction at different constant body temperatures at
# constant food for a lizard (Tiliqua rugosa - default parameter values, starting
# as an egg)

n <- 3000 # time steps
step <- 1 # step size (days)

Tbs=seq(25, 35, 5) # sequence of body temperatures to use

for(j in 1:length(Tbs)){
  debout<-matrix(data = 0, nrow = n, ncol = 38)
  deb.names <- c("stage", "V", "E", "E_H", "E_s", "E_R", "E_B", "q", "hs", "length", "wetmass", "wetgonad", "wetgut", "wetstorage", "p_surv", "fecundity", "clutches", "JMO2", "JMCO2", "JMH2O", "JMNWASTE", "O2ML", "CO2ML", "GH2OMET", "DEBQMETW", "GDRYFOOD", "GFAECES", "GNWASTE", "p_A", "p_C", "p_M", "p_G", "p_D", "p_J", "p_R", "p_B", "L_b", "L_j")
  colnames(debout)<-deb.names

  # initialise
  debout[1,]<-DEB_euler(Tb = Tbs[j], step = step)

  for(i in 2:n){
    debout[i,] <- DEB_euler(Tb = Tbs[j], breeding = 1, step = step, E_pres = debout[(i - 1), 3], V_pres = debout[(i - 1), 2], E_H_pres = debout[(i - 1), 4], q = debout[(i - 1), 8], hs = debout[(i - 1), 9], p_surv = debout[(i - 1), 15], E_s_pres = debout[(i - 1), 5], E_R = debout[(i - 1), 6], E_B = debout[(i - 1), 7])
  }

  if(j == 1){
    plot((seq(1, n) / 365), debout[, 11], ylim = c(100, 1500), type = 'l', xlab = 'years', ylab = 'wet mass, g', col = j)
  }else{
    points((seq(1,n) / 365), debout[, 11], ylim = c(100, 1500), type = 'l', xlab = 'years', ylab = 'wet mass, g',col = j)
  }

} #end loop through body temperatures

mrke/NicheMapR documentation built on April 3, 2024, 10:05 a.m.