Including Interactions and Splines in Piecewise Constant Hazard Regression

Description

This function can be used within a call to pchreg to automatically include spline functions in the linear predictor of the model.

Usage

1
  splinex(method = c("ns", "bs"), df = 2, degree = 2, v = 0.98, ...)

Arguments

method

a character string indicating whether natural splines (ns) or B-splines (bs) should be used. Default is “ns”.

df

the degrees of freedom of the spline basis.

degree

the degree of the polynomial (only for bs).

v

a value between 0 and 1 determining how many principal components of the design matrix must be used in model fitting (see “Details”).

...

for future arguments.

Details

The piecewise constant hazard model implemented by pchreg can be used as a nonparametric maximum likelihood estimator, in which the number of parameters is allowed to increase with the sample size in order to achieve any desired flexibility. Modeling the effect of covariates is as important as setting a sufficiently large number of breaks.

By letting splinex = splinex(...), each column of the original design matrix is automatically replaced by the corresponding spline basis, defined by method, df, and degree.

This modeling approach has the drawback of generating a potentially large design matrix: to reduce its dimension, you can choose v < 1. In this case, only the principal components explaining at least a proportion v of the variance are used to fit the model (see “Examples”).

Value

The function returns its arguments, to be passed to an internal function build.splinex that actually computes the design matrix.

Note

A multidimensional spline can be created by including a term like, for example ns(x1,df)*ns(x2,df). This is not supported by splinex, as it may generate a very large design matrix.

Author(s)

Paolo Frumento <paolo.frumento@ki.se>

See Also

pchreg

Examples

 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
  require(splines)
  n <- 1000
  x1 <- runif(n,-2,2)
  x2 <- runif(n,-2,2)
  t <- rexp(n, exp(-abs(x1 - x2)))

  # a simple model
  model1 <- pchreg(Surv(t) ~ x1 + x2)

  # using splinex: the same as ~ ns(x1, df = 2) + ns(x2, df = 2)
  model2 <- pchreg(Surv(t) ~ x1 + x2, splinex = splinex("ns", v = 1))

  # include interaction: ~ ns(x1, df = 2) + ns(x2, df = 2) + ns(x1*x2, df = 2)
  model3 <- pchreg(Surv(t) ~ x1 * x2, splinex = splinex("ns", v = 1))

  # the same as model 3, only keep the PCs explaining at least 95 percent of the variance
  model4 <- pchreg(Surv(t) ~ x1 * x2, splinex = splinex("ns", v = 0.95))

  # true CDF vs fitted

  trueF <- pexp(t, exp(-abs(x1 - x2)))
  par(mfrow = c(2,2))
  plot(trueF, 1 - predict(model1)$Surv); abline(0,1, col = "red") # does not fit
  plot(trueF, 1 - predict(model2)$Surv); abline(0,1, col = "red") # neither
  plot(trueF, 1 - predict(model3)$Surv); abline(0,1, col = "red") # great!
  plot(trueF, 1 - predict(model4)$Surv); abline(0,1, col = "red") # almost as good