Description Usage Arguments Details Value Note Author(s) References Examples
There are two functions for fitting the expectile bundle model, one for estimating asymmetry parameters (fitasy
),
the other for estimating the amplitude function, fitampl
, this function. See the details below.
1 | fitampl(y, B, alpha, p, a, pord = 2, lambda)
|
y |
a response vector. |
B |
a proper B-spline basis matrix, see |
alpha |
a vector of B-spline coefficients. |
p |
a vector of asymmetries. |
a |
a vector of asymmetry parameters. |
pord |
the order of the difference penalty, default is 2. |
lambda |
the positive tuning parameter for the penalty. |
The expectile bundle model determines a set of expectile curves for a point cloud with data vectors x
and y
,
as ψ_j{x_i} = a_j g(x_i). Here a_j is the asymmetry parameter corresponding to a given asymmetry p_j.
A vector of asymmetries with all 0 <p_j < 1 is specified by the user.
The asymmetric least squares objective function is
∑_j ∑_i w_{ij}(y_i - ∑_j a_j g_j(x_i))^2.
The function g(\cdot) is called the amplitude. The weights depend on the residuals:
w_{ij} = p_j
if y_i > a_jg(x_i) and w_{ij} = 1- p_j otherwise.
The amplitude function is a sum of B-splines with coefficients alpha
. There is no direct solution, so alpha
and the asymmetry parameters a
must be updated alternatingly. See the example.
a vector of estimated B-spline coefficients.
This is a simplification of the model described in the reference. There is no explict term for the trend.
Paul Eilers
Schnabel, S.K. and Eilers, P.H.C. (2013) A location-scale model for non-crossing expectile curves. Stat 2: 171–183.
Eilers, P.H.C. and Marx, B.D. (2021). Practical Smoothing, The Joys of P-splines. Cambridge University Press.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 | # Get the data
data(bone_data)
x = bone_data$age
y = bone_data$spnbmd
m <- length(x)
# Set asymmetry levels
p = c(0.005, 0.01, 0.02, 0.05, 0.2, 0.5, 0.8, 0.9, 0.95, 0.98, 0.99, 0.995)
np <- length(p)
# Set P-spline parameters
x0 <- 5
x1 <- 30
ndx <- 20
bdeg <- 3
pord <- 2
# Compute bases
B <- bbase(x, x0, x1, ndx, bdeg)
xg <- seq(from = min(x), to = max(x), length = 100)
Bg <- clone_base(B, xg)
n <- ncol(B)
lambda = 1
alpha <- rep(1,n)
a = p
for (it in 1:20){
alpha <- fitampl(y, B, alpha, p, a, pord, lambda)
alpha <- alpha / sqrt(mean(alpha ^ 2))
anew <- fitasy(y, B, alpha, p, a)
da = max(abs(a - anew))
a = anew
cat(it, da, '\n')
if (da < 1e-6) break
}
# Compute bundle on grid
ampl <- Bg %*% alpha
Z <- ampl %*% a
# Plot data and bundle
plot(x, y, pch = 15, cex = 0.7, col = 'grey', xlab = 'Age', ylab = 'Density')
cols = colorspace::rainbow_hcl(np, start = 10, end = 350)
matlines(xg, Z, lty = 1, lwd = 2, col = cols)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.