Skip to contents

Given renewal rates for different isotopic probes over time, a model estimates a normal distribution of the isotopic values for each time interval. Out of the given rates of renewal, first the influence of each time interval on the final isotopic value is calculated. This proportion is used as a regressor in a fully Bayesian linear model. The variance is estimated as rbf-Kernel-Matrix, where the hyperparameters are chosen such that the correlation between time intervals close to each other is about 0.5. The estimation is implemented using link[rstan]{rstan}.

Usage

estimateIntervals(
  renewalRates,
  timeVars,
  boneVars = NULL,
  indVar = "",
  isoMean,
  isoSigma,
  renewalRatesSD = NULL,
  mu_df = 3,
  mu_mean = 0,
  mu_sd = 20,
  rho_mean = 1,
  rho_sd = 0.25,
  alpha_mean = 2,
  alpha_sd = 0.5,
  mc = TRUE,
  adapt_delta = 0.98,
  max_treedepth = 15,
  chains = 4,
  iter = 2000,
  burnin = 500,
  ...
)

Arguments

renewalRates

A dataframe specifying the renewal rates of different probes for each time interval. The renewalRates should be between 0 and 100 (percentages). The dataframe should include a column specifying a time-index (e.g. 1, 2, 3, ...) as well as columns for the different bones. The renewal rates should contain the times of origin, containing 100.

timeVars

A character string specifying the name of the column indicating the time.

boneVars

A vector of character strings indicating the relevant variables containing the renewal rates of bones and teeth. If not given, the variables with bones are all variables of the dataframe apart from the time variables.

indVar

A character string specifying the name of the column indicating the individual.

isoMean

A numeric number indicating the mean of the isotopic values measured.

isoSigma

A numeric, positive number indicating the standard deviation of the isotopic values measured.

renewalRatesSD

An optional dataframe specifying the renewal rates uncertainties in standard deviations of different probes for each time interval.

mu_df

Hyperparameter for the mean of the interval estimates: degrees of freedom of a non-standardized t-Student distribution. Defaults to 1.

mu_mean

Hyperparameter for the mean of the interval estimates: mean of a non-standardized t-Student distribution. Defaults to 0.

mu_sd

Hyperparameter for the mean of the interval estimates: standard deviation of a non-standardized t-Student distribution. Defaults to 10.

rho_mean

Hyperparameter for the length scale of the rbf-kernel: mean of a normal. Defaults to 1.

rho_sd

Hyperparameter for the length scale of the rbf-kernel: standard deviation of a normal. Defaults to 0.25.

alpha_mean

Hyperparameter for the signal variance of the rbf-kernel: mean of a normal. Defaults to 2.

alpha_sd

Hyperparameter for the signal variance of the rbf-kernel: standard deviation of a normal. Defaults to 0.5.

mc

A boolean indicating if multiple cores should be used. If TRUE, which is the default, 4 cores are used.

adapt_delta

A numeric value between 0 and 1 controlling the samplers behavior. Defaults to 0.9999. See stan for more details.

max_treedepth

A numeric, positive value controling the NUTS sampler. Defaults to 25. See stan for more details.

chains

Number of chains for mcmc. Defaults to 4

iter

Number of iterations per chain for mcmc. Defaults to 2000

burnin

Number of burnin iterations per chain for mcmc. Defaults to 500

...

Arguments passed to rstand sampling

Value

A fitted object of class TemporalIso.

See also

Examples


if (FALSE) {
#dataset with renewal rates
data <- data.frame(
 intStart = 0:5,
 intEnd = 1:6,
 bone1 = c(100, 50, 20, 10, 5, 2),
 bone2 = c(100, 10, 5, 1, 1, 1),
 tooth1 = c(0, 100, 0, 0, 0, 0),
 tooth2 = c(0, 0, 100, 0, 0, 0)
)
# isotope-values and sd
y_mean <- c(-10, -7, -12, -9)
y_sigma <- c(2, 1.5, 2.5, 2.5)


###############################
#estimate values
###############################
fit <- estimateIntervals(renewalRates = data,
                        timeVars = "intStart",
                        boneVars = c("bone1", "bone2", "tooth1", "tooth2"),
                        isoMean = y_mean,
                      isoSigma = y_sigma)
print(fit)
# implemented Methods
plotTime(fit)
#get estimates for specific time points
estimateTimePoint(fit, time = seq(0,5, by = 0.5))

###############################
#shift point estimation
###############################
plotTime(fit, plotShifts = TRUE, threshold = 0.5)
plotTime(fit, plotShifts = TRUE, threshold = 0.5,
         type = TRUE, probability = 0.5)

getShiftTime(fit, threshold = 0.5)


###############################
#Staying time estimation
###############################

estimatedStayTimes <- getSiteStayTimes(object = fit,
                                      siteMeans = c(-8, -10),
                                      siteSigma = c(1, 1.5))


###############################
#compute resulting isotopic values out of historic ones
###############################

testDat <- data.frame(
 t = 1:3,
 bone1 = c(100, 50, 10, 10, 10, 10),
 bone2 = c(0, 20, 20, 20, 20, 20),
 teeth1 = c(0, 0, 0, 10, 50, 60),
 mean = c(1, 3, 3,0, -1, - 5),
 sd = c(1, 3, 1, 1.5, 2, 3)
)

computeResult(
 data = testDat,
 timeVars = "t",
 boneVars = c("bone1","bone2","teeth1"),
 meanVar = "mean",
 sdVar = "sd"
)

shiny::runApp(paste0(system.file(package = "OsteoBioR"),"/app"))
}