JOPS (version 0.1.19)

hist2dsm: Smooth a 2D histogram

Description

Fit a 2D smooth P-spline surface to a matrix of counts, assuming Poisson distributed observations.

Usage

hist2dsm(
  Y,
  nsegx = 10,
  nsegy = nsegx,
  bdeg = 3,
  lambdax = 10,
  lambday = lambdax,
  dx = 3,
  dy = dx,
  Mu = Y + 0.01,
  kappa = 1e-04,
  tol = 1e-05
)

Value

A list with elements:

ed

the effective dimension of the smooth 2D surface.

Mu

a matrix with the smooth estimates, with dimensions of dim(Y)

pen

the numerical value of the penalty.

Arguments

Y

a matrix of counts.

nsegx

the number of knots along x (default=10).

nsegy

the number of evenly spaced knots along y for Tensor product B-spline basis (default=10).

bdeg

the degree of the basis, default is 3.

lambdax

the positive number for the tuning parameter along x.

lambday

the positive number for the tuning parameter along y.

dx

the order of the difference penalty along x, default is 3.

dy

the order of the difference penalty along y, default is 3.

Mu

the initialization of the mean (default Y + 0.01).

kappa

a (small, positive) number for ridge tuning parameter to stabilize estimation (default 1e-4).

tol

the convergence criterion (default 1e-5).

Author

Paul Eilers

References

Eilers, P.H.C., Marx, B.D., and Durban, M. (2015). Twenty years of P-splines, SORT, 39(2): 149-186.

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

Examples

x = faithful$eruptions
y = faithful$waiting
h = hist2d(x, y, c(100, 100))
sm = hist2dsm(h$H, nsegx = 25, nsegy = 25, bdeg = 3, lambdax = 10, lambday = 10)
image(h$xgrid, h$ygrid, sm$Mu, xlab = 'Eruption length (min)',
      ylab = 'Waiting time (min)', main = 'Old Faithful')