| joint.marginal {INLA} | R Documentation |
Sample and evalue from from a joint marginal approximation
as returned using argument selection in inla.
inla.rjmarginal(n, jmarginal, constr) inla.rjmarginal.eval(fun, samples, ...)
n |
The number of samples |
jmarginal |
A marginal object given either by a |
constr |
Optional linear constraints;
see |
fun |
A function which is evaluated for each sample, similar to
|
samples |
The samples, as in the form of the output from |
inla.rjmarginal returns a list with the samples in samples
(matrix) and the corresponding log-densities
in log.density (vector). Each column in samples contains one sample.
inla.rjmarginal.eval returns a matrix, where each row is the (vector) function
evaluated at each sample.
Havard Rue hrue@r-inla.org
n = 10
x = 1+rnorm(n)
xx = 3 + rnorm(n)
y = 1 + x + xx + rnorm(n)
selection = list(xx=1, Predictor = 3:4, x=1)
r = inla(y ~ 1 + x + xx,
data = data.frame(y, x, xx),
selection = selection)
ns = 100
xx = inla.rjmarginal(ns, r)
print(cbind(mean = r$selection$mean, sample.mean = rowMeans(xx$samples)))
print("cov matrix")
print(round(r$selection$cov.matrix, dig=3))
print("sample cov matrix")
print(round(cov(t(xx$samples)), dig=3))
skew = function(z) mean((z-mean(z))^3)/var(z)^1.5
print(round(cbind(skew = r$selection$skewness,
sample.skew = apply(xx$sample, 1, skew)), dig=3))
## illustrating the eval function
n = 10
x = rnorm(n)
eta = 1 + x
y = eta + rnorm(n, sd=0.1)
selection = list(x = 1, Predictor = c(1, 2, 4, 5), '(Intercept)' = 1)
r = inla(y ~ 1 + x,
data = data.frame(y, x),
selection = selection)
xx = inla.rjmarginal(100, r)
xx.eval = inla.rjmarginal.eval(function() c(x, Predictor, Intercept), xx)
print(cbind(xx$samples[, 1]))
print(cbind(xx.eval[, 1]))