JOPS (version 0.2.0)

pclm: Fit a composite link model

Description

Fit a smooth latent distribution using the penalized composite link model (PCLM).

Usage

pclm(y, C, B, lambda = 1, pord = 2, itmax = 50, show = FALSE)

Value

A list with the following items:

alpha

the estimated B-spline coefficients, length n.

gamma

the estimated latent distribution, length q.

mu

estimated values of y, length m.

dev

the deviance of the model.

ed

the effective model dimension.

aic

Akaike's Information Criterion.

Arguments

y

a vector of counts, length m.

C

a composition matrix, m by q.

B

a B-spline basis matrix, q by n.

lambda

the penalty parameter.

pord

the the order of the difference penalty (default = 2).

itmax

the maximum number of iterations (default = 50).

show

Set to TRUE or FALSE to display iteration history (default = FALSE).

Author

Paul Eilers and Jutta Gampe

Details

The composite link model assumes that \(E(y) = \mu = C\exp(B \alpha)\), where \(\exp(B\alpha)\) is a latent discrete distribution, usually on a finer grid than that for \(y\).

Note that sum(gamma) == sum(mu).

References

Eilers, P. H. C. (2007). III-posed problems with counts, the composite link model and penalized likelihood. Statistical Modelling, 7(3), 239–254.

Eilers, P.H.C. and Marx, B.D. (2021). Practical Smoothing, The Joys of P-splines. Cambridge University Press.

Examples

# Left and right boundaries, and counts, of wide intervals of the data
cb <- c( 0, 20, 30, 40, 50, 60)
ce <- c(20, 30, 40, 50, 60, 70)
y <- c(79, 54, 19, 1, 1, 0)

# Construct the composition matrix
m <- length(y)
n <- max(ce)
C <- matrix(0, m, n)
for (i in 1:m) C[i, cb[i]:ce[i]] <- 1

mids = (cb + ce) / 2 - 0.5
widths = ce - cb + 1
dens = y / widths / sum(y)
x = (1:n) - 0.5
B = bbase(x)
fit = pclm(y, C, B, lambda = 2, pord = 2, show = TRUE)
gamma = fit$gamma / sum(fit$gamma)
# Plot density estimate and data
plot(x, gamma, type = 'l', lwd = 2, xlab = "Lead Concentration", ylab = "Density")
rect(cb, 0, ce, dens, density = rep(10, 6), angle = rep(45, 6))