odin supports many functions that you’d expect to see
for constructing differential equation models; primarily mathematical
functions available through R’s “Rmath” library. These include all
mathematical operations, and many more obscure mathematical functions.
Special support is provided for working with arrays. A further set of
functions is available for working discrete time stochastic models.
+ – Plus: Both infix
(a + b) and prefix (+a) versions supported
(e.g., 1 + 2 → 3)
- – Minus: Both infix
(a - b) and prefix (-a) versions supported
(e.g., 10 - 1 → 9)
* – Multiply: Multiply two numbers
together (e.g., 2 * 6 → 12)
/ – Divide: Divide two numbers
(e.g., 12 / 6 → 2)
^ – Power: Raise the first number
to the power of the second. Either number may be a floating point number
(e.g., 2.3 ^ 1.2 → 2.716898)
( – Parenthesis: Group expressions
together (e.g., (1 + 5) * 2 → 12)
if – Conditional: Inline
conditional statement. This takes a form slightly different to typically
seen in R with the result of the statement directly assigned (e.g.,
if (9 > 10) 1 else 2 → 2)
Because general programming is not supported in odin and
because every line must contain an assignment, instead of writing
if (mycondition) {
a <- true_value
} else {
a <- false_value
}instead write
a <- if (mycondition) true_value else false_value(this works in normal R too!)
There are a group of functions for interacting with arrays
sum – Sum: Compute sum over an
array, or over some dimension of an array
length – Length: Total length of an
array
dim – Length of one of an array’s
dimensions: If an array x has 10 rows and 20
columns, then dim(x, 1) is 10 and
dim(x, 2) is 20. Note that this differs from
dim in R and dim(x, i) is closer to
dim(x)[[i]]
[ – Subset an array: See
below
interpolate – Interpolate an array over
time: See below
When working with arrays, use generally implies a “for loop” in the generated C code. For example, in the example in the main package vignette the derivatives are computed as
The indexes on the right hand side can be one of i,
j, k, l i5,
i6, i7 or i8 corresponding to the
index on the left hand side being iterated over
(odin supports arrays up to 8 dimensions). The
left-hand-side here contains no explicit entry (y[]) which
is equivalent to y[1:length(y)], which expands
(approximately) to the “for loop”
(except slightly different, and in C).
Similarly, the expression
ay[, ] <- a[i, j] * y[j]involves loops over two dimensions (ay[, ] becomes
ay[1:dim(ay, 1), 1:dim(ay, 2)] and so the loop becomes
Due to constraints with using C, few things can be used as an index; in particular the following will not work:
idx <- if (t > 5) 2 else 1
x <- vec[idx](or where idx is some general odin variable as the
result of a different assignment). You must use as.integer
to cast this to integer immediately before indexing:
idx <- if (t > 5) 2 else 1
x <- vec[as.integer(idx)]This will truncate the value (same behaviour as
truncate) so be warned if passing in things that may be
approximately integer - you may want to use
as.integer(round(x)) in that case.
The interpolation functions are described in more detail in the main package vignette
A number of logical-returning operators exist, primarily to support
the if statement; all the usual comparison operators exist
(though not vectorised | or &).
> – Greater than (e.g.,
1 > 2 → FALSE)
< – Less than (e.g.,
1 < 2 → TRUE)
>= – Greater than or equal to
(e.g., 1 >= 2 → FALSE)
<= – Less than or equal to
(e.g., 1 <= 2 → TRUE)
== – Is exactly equal to (e.g.,
1 == 1 → TRUE)
!= – Is not exactly equal to (e.g.,
1 != 2 → TRUE)
&& – Boolean AND (e.g.,
(1 == 1) && (2 > 1) →
TRUE)
|| – Boolean OR (e.g.,
(1 == 1) && (2 > 1) →
TRUE)
Be wary of strict equality with == or != as
numbers may be floating point numbers, which have some surprising
properties for the uninitiated, for example
sqrt(3)^2 == 3## [1] FALSE%% – Modulo: Finds the remainder
after division of one number by another (e.g., 123 %% 100 →
23)
%/% – Integer divide: Different to
floating point division, effectively the full number of times one number
divides into another (e.g., 20 %/% 7 →
2)
abs – Absolute value (e.g.,
abs(-1) → 1)
sign – Sign function: Returns the
sign of its argument as either -1, 0 or 1, which may be useful for
multiplying by another argument (e.g., sign(-100) →
-1)
round – Round a number (e.g.,
round(1.23) → 1; round(1.23, 1) →
1.2)
floor – Floor of a number: Largest
integer not greater than the provided number (e.g.,
floor(6.5) → 6)
ceiling – Ceiling of a number:
Smallest integer not less than the provided number (e.g.,
ceiling(6.5) → 7)
trunc – Truncate a number: Round a
number towards zero
max – Maximum: Returns maximum of
all arguments given (e.g., max(2, 6, 1) →
6)
min – Minimum (e.g.,
min(2, 6, 1) → 1)
exp – Exponential function (e.g.,
exp(1) → 2.718282)
expm1 – Computes exp(x) - 1 accurately for
small |x| (e.g., exp(1) →
1.718282)
log – Logarithmic function (e.g.,
log(1) → 0)
log2 – Logarithmic function in base
2 (e.g., log2(1024) → 10)
log10 – Logarithmic function in base
10 (e.g., log10(1000) → 3)
log1p – Computes log(x + 1) accurately for
small |x| (e.g., log1p(1) →
0.6931472)
sqrt – Square root function (e.g.,
sqrt(4) → 2)
beta – Beta function (e.g.,
beta(3, 5) → 0.00952381)
lbeta – Log beta function (e.g.,
lbeta(3, 5) → -4.65396)
choose – Binomial coefficients
(e.g., choose(60, 3) → 34220)
lchoose – Log binomial coefficients
(e.g., choose(60, 3) → 10.44057)
gamma – Gamma function (e.g.,
gamma(10) → 362880)
lgamma – Log gamma function (e.g.,
lgamma(10) → 12.80183)
The exact for %% and %/% for floating point
numbers and signed numbers are complicated - please see
?Arithmetic. The rules for operators in odin
are exactly those in R as the same underlying functions are used.
Similarly, for the differences between round,
floor, ceiling and truncate, see
the help page ?round. Note that R’s behaviour for rounding
away from 0.5 is exactly followed and that this slightly changed
behaviour at version 4.0.0
All the usual trig functions are also available:
cos – Cosine function
sin – Sine function
tan – Tangent function
acos – Arc-cosine function
asin – Arc-sin function
atan – Arc-tangent
function
atan2 – Two-arg arc-tangent
function
cosh – Hyperbolic cosine
function
sinh – Hyperbolic sine
function
tanh – Hyperbolic tangent
function
acosh – Hyperbolic arc-cosine
function
asinh – Hyperbolic arc-sine
function
atanh – Hyperbolic arc-tangent
function
For discrete time stochastic models, all of R’s normal stochastic distribution functions are available:
unif_rand – Standard uniform
distribution: Sample from the uniform distribution on [0, 1] -
more efficient than but equivalent to runif(0, 1)
norm_rand – Standard normal
distribution: Sample from the standard normal distribution -
more efficient than but equivalent to rnorm(0, 1)
exp_rand – Standard exponential
distribution: Sample from the exponential distribution with
rate 1 - more efficient than but equivalent to rexp(1)
rbeta – Beta distribution: With
parameters shape1 and shape2 (see ?rbeta for
details)
rbinom – Binomial distribution:
With parameters size (number of trials) and
prob (probability of success)
rcauchy – Cauchy distribution: With
parameters location and scale
rchisq – Chi-Squared distribution:
With parameter df
rexp – Exponential distribution:
With parameter rate
rf – F-distribution: With parameter
df1 and `df2
rgamma – Gamma distribution: With
parameters shape and rate
rgeom – Geometric distribution:
Distribution with parameters prob
rhyper – Hypergeometric
distribution: With parameters m (the number of
white balls in the urn), n (the number of black balls in
the urn) and k (the number of balls drawn from the
urn)
rlogis – Logistic distribution:
With parameters location and scale
rlnorm – Log-normal distribution:
With parameters meanlog and sdlog
rnbinom – Negative binomial
distribution: With parameters size,
prob and mu
rnorm – Normal distribution: With
parameters mean and sd
rpois – Poisson distribution: With
parameter lambda
rt – Student’s t distribution: With
parameter df
runif – uniform distribution: With
parameters min and max
rweibull – Weibull distribution:
With parameters shape and scale
rwilcox – Wilcoxon rank sum statistic
distribution: With parameters n and
m
rsignrank – Wilcoxon signed rank statistic
distribution: With parameter n
With random number functions we can write:
x <- runif(10, 20)which will generate a random number from the uniform distribution. If you write:
x[] <- runif(10, 20)then each element of x will be filled with a
different random number drawn from this distribution (which is
generally what you want). Random numbers are considered to be time
varying which means they will automatically generate each time
step, so if you write
then at each time step, each element of y will be
updated by the same random number from a normal distribution with a mean
of zero and a standard deviation of 10 - the number will change each
time step but be the same for each element of y in the
example above.
In addition, two functions that are vector returning and require some care to use:
rmultinom – multinomial
distribution: The first parameter is the number of samples and
the second is the per-class probability and must be a vector
rmhyper – Multivariate hypergeometric
distribution: The first parameter is the number of samples and
the second is the per-class count and must be a vector
Both these functions require a vector input (of probabilities for
rmultinom and of counts for rmhyper) and
return a vector the same length. So the expression
y[] <- rmultinom(10, p)will produce a vector y of samples from the multinomial
distribution with parameters size = 10 (so after wards
sum(y) is 10) and probabilities p. It is very
important that y and p have the same size.
At the moment it is not possible to use expressions like
y[1, ] <- rmultinom(10, p[i, ])but this is planned for implementation in the future. A full example
of using rmultinom is given in the discrete models vignette.