The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

Getting started with the PCMBaseCpp R-package

2025-05-04

How to use the package?

There are two ways to use PCMBaseCpp:

Passing the function PCMInfoCpp as a metaI argument of PCMLik and/or PCMCreateLikelihood

library(PCMBase)
library(PCMBaseCpp)
## Loading required package: Rcpp
system.time(llR <- PCMLik(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab))
##    user  system elapsed 
##   0.071   0.000   0.072
system.time(llCpp <- PCMLik(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = PCMInfoCpp))
##    user  system elapsed 
##   0.005   0.000   0.004
print(llR)
## [1] -206.4146
## attr(,"X0")
## [1] 5 2 1
print(llCpp)
## [1] -206.4146
## attr(,"X0")
## [1] 5 2 1
logLikFunR <- PCMCreateLikelihood(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab)

logLikFunCpp <- PCMCreateLikelihood(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = PCMInfoCpp)

metaICpp <- PCMInfoCpp(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab)

logLikFunCpp2 <- PCMCreateLikelihood(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaICpp)

set.seed(1, kind = "Mersenne-Twister", normal.kind = "Inversion")
randParam <- PCMParamRandomVecParams(PCMBaseTestObjects$model_MixedGaussian_ab)

system.time(llR <- logLikFunR(randParam))
##    user  system elapsed 
##   0.059   0.000   0.060
system.time(llCpp <- logLikFunCpp(randParam))
##    user  system elapsed 
##   0.003   0.000   0.002
system.time(llCpp2 <- logLikFunCpp2(randParam))
##    user  system elapsed 
##   0.002   0.000   0.001
print(llR)
## [1] -598.092
## attr(,"X0")
## [1] -4.689827 -2.557522  1.457067
print(llCpp)
## [1] -598.092
## attr(,"X0")
## [1] -4.689827 -2.557522  1.457067
print(llCpp2)
## [1] -598.092
## attr(,"X0")
## [1] -4.689827 -2.557522  1.457067

Passing the meta-information object returned by PCMInfoCpp as a metaI argument of PCMLik and PCMCreateLikelihood

This is the recommended usage in the case of multiple likelihood evaluations, e.g. during model inference:

metaIR <- PCMInfo(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab)

metaICpp <- PCMInfoCpp(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab)

system.time(llR <- PCMLik(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaIR))
##    user  system elapsed 
##   0.059   0.000   0.058
system.time(llCpp <- PCMLik(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaICpp))
##    user  system elapsed 
##   0.001   0.000   0.001
print(llR)
## [1] -206.4146
## attr(,"X0")
## [1] 5 2 1
print(llCpp)
## [1] -206.4146
## attr(,"X0")
## [1] 5 2 1
logLikFunR <- PCMCreateLikelihood(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaIR)

logLikFunCpp <- PCMCreateLikelihood(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaICpp)

system.time(llR <- PCMLik(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaIR))
##    user  system elapsed 
##   0.059   0.000   0.059
system.time(llCpp <- PCMLik(
  X = PCMBaseTestObjects$traits.ab.123, 
  tree = PCMBaseTestObjects$tree.ab,
  model = PCMBaseTestObjects$model_MixedGaussian_ab, 
  metaI = metaICpp))
##    user  system elapsed 
##   0.001   0.000   0.001
print(llR)
## [1] -206.4146
## attr(,"X0")
## [1] 5 2 1
print(llCpp)
## [1] -206.4146
## attr(,"X0")
## [1] 5 2 1

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.