Summarizing model fits for stepmented regression

summary/print method for class stepmented.

Usage

# S3 method for class 'stepmented'
summary(object, short = FALSE, var.diff = FALSE, p.df="p", .vcov=NULL, ...)

# S3 method for class 'summary.stepmented'
print(x, short=x$short, var.diff=x$var.diff, 
    digits = max(3, getOption("digits") - 3),
    signif.stars = getOption("show.signif.stars"),...)

# S3 method for class 'stepmented'
print(x, digits = max(3, getOption("digits") - 3), 
    ...)

Arguments

object, x: Object of class "stepmented" or a summary.stepmented object produced by summary.stepmented().
short: logical indicating if the `short' summary should be printed.
var.diff: logical indicating if different error variances should be computed in each interval of the stepmented variable, see Details. If .vcov is provided, var.diff is set to FALSE.
p.df: A character as a function of 'p' (number of parameters) and 'K' (number of groups or segments) affecting computations of the group-specific variance (and the standard errors) if var.diff=TRUE, see Details.
.vcov: Optional. The full covariance matrix for the parameter estimates. If provided, standard errors are computed (and displayed) according to this matrix.

digits: controls number of digits printed in output.
signif.stars: logical, should stars be printed on summary tables of coefficients?
...: further arguments, notably type to be passed to vcov.stepmented to compute the standard errors. See vcov.stepmented.

Details

If short=TRUE only coefficients of the stepmented relationships are printed. If var.diff=TRUE and there is only one stepmented variable, different error variances are computed in the intervals defined by the estimated breakpoints of the stepmented variable. For the jth interval with $n_j$ observations, the error variance is estimated via $RSS_j/(n_j-p)$, where $RSS_j$ is the residual sum of squares in interval j, and $p$ is the number of model parameters. This number to be subtracted from $n_j$ can be changed via argument p.df. For instance p.df="0" uses $RSS_j/(n_j)$, and p.df="p/K" leads to $RSS_j/(n_j-p/K)$, where $K$ is the number of groups (segments), and $p/K$ can be interpreted as the average number of model parameter in that group.

Note var.diff=TRUE only affects the estimates covariance matrix. It does not affect the parameter estimates, neither the log likelihood and relevant measures, such as AIC or BIC. In other words, var.diff=TRUE just provides 'alternative' standard errors, probably appropriate when the error variances are different before/after the estimated breakpoints. Also $p-values$ are computed using the t-distribution with 'naive' degrees of freedom (as reported in object$df.residual).

If var.diff=TRUE the variance-covariance matrix of the estimates is computed via the sandwich formula, $$(X^TX)^{-1}X^TVX(X^TX)^{-1}$$ where V is the diagonal matrix including the different group-specific error variance estimates. Standard errors are the square root of the main diagonal of this matrix.

Value

A list (similar to one returned by stepmented.lm or stepmented.glm) with additional components:

psi: estimated break-points and relevant (approximate) standard errors
Ttable: estimates and standard errors of the model parameters. This is similar to the matrix coefficients returned by summary.lm or summary.glm, but without the rows corresponding to the breakpoints. Even the p-values relevant to the difference-in-slope parameters have been replaced by NA, since they are meaningless in this case, see davies.test.
cov.var.diff: if var.diff=TRUE, the covaraince matrix accounting for heteroscedastic errors.
sigma.new: if var.diff=TRUE, the square root of the estimated error variances in each interval.
df.new: if var.diff=TRUE, the residual degrees of freedom in each interval.

Author

Vito M.R. Muggeo

Warning

If type is not specified in ... (which means type="standard"), no standard error will be computed (and returned) for the jumpoint.

Examples

##continues example from stepmented()
# summary(stepmented.model,short=TRUE)

## an heteroscedastic example..
# set.seed(123)
# n<-100
# x<-1:n/n
# y<- -x+1.5*pmax(x-.5,0)+rnorm(n,0,1)*ifelse(x<=.5,.4,.1)
# o<-lm(y~x)
# oseg<-stepmented(o,seg.Z=~x,psi=.6)
# summary(oseg,var.diff=TRUE)$sigma.new