Calculate the mean at time t, given X0, under a PCM model

Calculate the mean at time t, given X0, under a PCM model

PCMMeanAtTime(
  t,
  model,
  X0 = model$X0,
  regime = PCMRegimes(model)[1L],
  verbose = FALSE
)

Arguments

t	positive numeric denoting time
model	a PCM model object
X0	a numeric vector of length k, where k is the number of traits in the model (Defaults to model$X0).
regime	an integer or a character denoting the regime in model for which to do the calculation; Defaults to PCMRegimes(model)[1L], meaning the first regime in the model.
verbose	a logical indicating if (debug) messages should be written on the console (Defaults to FALSE).

Value

A numeric vector of length k

Examples

# a Brownian motion model with one regime
modelBM <- PCM(model = "BM", k = 2)
# print the model
modelBM
#> Brownian motion model
#> S3 class: BM, GaussianPCM, PCM; k=2; p=8; regimes: 1. Parameters/sub-models:
#> X0 (VectorParameter, _Global, numeric; trait values at the root):
#> [1] 0 0
#> Sigma_x (MatrixParameter, _UpperTriangularWithDiagonal, _WithNonNegativeDiagonal; factor of the unit-time variance rate):
#> , , 1
#> 
#>      [,1] [,2]
#> [1,]    0    0
#> [2,]    0    0
#> 
#> Sigmae_x (MatrixParameter, _UpperTriangularWithDiagonal, _WithNonNegativeDiagonal; factor of the non-heritable variance or the variance of the measurement error):
#> , , 1
#> 
#>      [,1] [,2]
#> [1,]    0    0
#> [2,]    0    0
#> 
#>  
# assign the model parameters at random: this will use uniform distribution
# with boundaries specified by PCMParamLowerLimit and PCMParamUpperLimit
# We do this in two steps:
# 1. First we generate a random vector. Note the length of the vector equals PCMParamCount(modelBM)
randomParams <- PCMParamRandomVecParams(modelBM, PCMNumTraits(modelBM), PCMNumRegimes(modelBM))
randomParams
#>          [,1]     [,2]     [,3]     [,4]     [,5]     [,6]      [,7]     [,8]
#> [1,] -3.29025 4.474519 3.376153 2.608282 8.406146 8.561317 -2.172814 3.804939
# 2. Then we load this random vector into the model.
PCMParamLoadOrStore(modelBM, randomParams, 0, PCMNumTraits(modelBM), PCMNumRegimes(modelBM), TRUE)
#> [1] 8

# PCMMeanAtTime(1, modelBM)

# note that the variance at time 0 is not the 0 matrix because the model has a non-zero
# environmental deviation
PCMMeanAtTime(0, modelBM)
#> [1] -3.290250  4.474519