Skip to contents

Function that runs the eDITH model for a given parameter set

Usage

run_eDITH_single(param, river, covariates,  Z.normalize = TRUE, 
        no.det = FALSE, ll.type = NULL, 
        data = NULL, source.area = "AG",
                tau.prior = list(spec="lnorm",a=0,b=Inf, 
          meanlog=log(5), sd=sqrt(log(5)-log(4))),
                log_p0.prior = list(spec="unif",min=-20, max=0),
                    beta.prior = list(spec="norm",sd=1),
                sigma.prior = list(spec="unif",min=0, 
          max=max(data$values, na.rm = TRUE)),
                omega.prior = list(spec="unif",min=1, 
          max=10*max(data$values, na.rm = TRUE)),
                 Cstar.prior = list(spec="unif",min=0, 
          max=max(data$values, na.rm = TRUE)))

Arguments

param

Parameter set. It has to be a named vector, with names:

tau

Decay time (expressed in h).

log_p0

Natural logarithm of the baseline production rate.

beta_X

Effect size of covariate X. There must be as many beta_X as columns in covariates. X must be the name of the corresponding column in covariates.

omega, sigma, Cstar

Parameters for estimation of the log-likelihood and detection probability. Only required if ll.type is provided.

river

A river object generated via aggregate_river.

covariates

Data frame containing covariate values for all river reaches.

Z.normalize

Logical. Should covariates be Z-normalized?

no.det

Logical. Should a probability of non-detection be included in the model?

ll.type

Character. String defining the error distribution used in the log-likelihood formulation. Allowed values are norm (for normal distribution), lnorm (for lognormal distribution), nbinom (for negative binomial distribution) and geom (for geometric distribution). The two latter choices are suited when eDNA data are expressed as read numbers, while norm and lnorm are better suited to eDNA concentrations.

data

eDNA data. Data frame containing columns ID (index of the AG node/reach where the eDNA sample was taken) and values (value of the eDNA measurement, expressed as concentration or number of reads).

source.area

Defines the extent of the source area of a node. Possible values are "AG" (if the source area is the reach surface, i.e. length*width), "SC" (if the source area is the subcatchment area), or, alternatively, a vector with length river$AG$nodes.

tau.prior, log_p0.prior,beta.prior,sigma.prior,omega.prior,Cstar.prior

Prior distribution for the relevant parameters of the eDITH model. Only used if both ll.type and data are provided.

Value

A list with objects:

p

Vector of eDNA production rates corresponding to the parameter set param. It has length equal to river$AG$nNodes.

C

Vector of eDNA values (concentrations or read numbers) corresponding to the parameter set param. It has length equal to river$AG$nNodes.

probDet

Vector of detection probabilities corresponding to the parameter set param. It is only computed if ll.type is provided. It has length equal to river$AG$nNodes.

logprior

Value of the log-prior distribution (computed only if ll.type and data are provided).

loglik

Value of the log-likelihood distribution (computed only if ll.type and data are provided).

logpost

Value of the log-posterior distribution (computed only if ll.type and data are provided).

See also

See run_eDITH_BT, run_eDITH_optim for details on parameters names and log-likelihood specification.

Examples

library(rivnet)
data(wigger)

# calculate AEMs and use the first 10 as covariates
ae <- river_to_AEM(wigger)
covariates <- data.frame(ae$vectors[,1:10])
names(covariates) <- paste0("AEM",1:10) 
# covariates names must correspond to param names
set.seed(1); param <- c(3,-15, runif(10,-1,1))
names(param) <- c("tau", "log_p0", paste0("beta_AEM",1:10))
# param names must correspond to covariates names

out <- run_eDITH_single(param, wigger, covariates)

# add parameter sigma and compute detection probability
param <- c(param, 5e-12) 
names(param)[length(param)] <- "sigma"
# note that the value of sigma has to be within the range indicated by sigma.prior
out2 <- run_eDITH_single(param, wigger, covariates, ll.type="norm")

# include data and compute logprior, loglikelihood, logposterior
data(dataC)
out3 <- run_eDITH_single(param, wigger, covariates, 
    ll.type="norm", data=dataC)