1\name{predict.rqss}
2\alias{predict.rqss}
3\alias{predict.qss1}
4\alias{predict.qss2}
5\title{Predict from fitted nonparametric quantile regression smoothing spline models}
6\description{
7  Additive models for nonparametric quantile regression using total
8  variation penalty methods can be fit with the \code{\link{rqss}}
9  function.  Univarariate and bivariate components can be predicted
10  using these functions.
11}
12\usage{
13\method{predict}{rqss}(object, newdata, interval = "none", level = 0.95, ...)
14\method{predict}{qss1}(object, newdata, ...)
15\method{predict}{qss2}(object, newdata, ...)
16}
17\arguments{
18\item{object}{ is a fitted object produced by \code{\link{rqss}} }
19\item{newdata}{ a data frame describing the observations at which
20  prediction is to be made.  For qss components, newdata should
21  lie in strictly within the convex hull of the fitting data.  Newdata
22  corresponding to the partially linear component of the model
23  may require caution concerning the treatment of factor levels, if any.}
24\item{interval}{If set to \code{confidence} then a \code{level} confidence interval
25for the predictions is returned.}
26\item{level}{intended coverage probability for the confidence intervals}
27\item{\dots}{ optional arguments }
28}
29\details{
30For both univariate and bivariate prediction linear interpolation is
31done.  In the bivariate case, this involves computing barycentric
32coordinates of the new points relative to their enclosing triangles.
33It may be of interest to plot individual components of fitted rqss
34models:  this is usually best done by fixing the values of other
35covariates at reference values typical of the sample data and
36predicting the response at varying values of one qss term at a
37time.   Direct use of the \code{predict.qss1} and \code{predict.qss2} functions
38is discouraged since it usually corresponds to predicted values
39at absurd  reference values of the other covariates, i.e. zero.
40}
41\value{
42A vector of predictions, or in the case that \code{interval = "confidence")}
43a matrix whose first column is the vector of predictions and whose second and
44third columns are the lower and upper confidence limits for each prediction.
45}
46\author{ R. Koenker }
47\seealso{ \code{\link{rqss}}
48}
49\examples{
50n <- 200
51lam <- 2
52x <- sort(rchisq(n,4))
53z <- exp(rnorm(n)) + x
54y <- log(x)+ .1*(log(x))^2 + z/4 +  log(x)*rnorm(n)/4
55plot(x,y - z/4 + mean(z)/4)
56Ifit <- rqss(y ~ qss(x,constraint="I") + z)
57sfit <- rqss(y ~ qss(x,lambda = lam) + z)
58xz <- data.frame(z = mean(z),
59                 x = seq(min(x)+.01,max(x)-.01,by=.25))
60lines(xz[["x"]], predict(Ifit, xz), col=2)
61lines(xz[["x"]], predict(sfit, xz), col=3)
62legend(10,2,c("Increasing","Smooth"),lty = 1, col = c(2,3))
63title("Predicted Median Response at Mean Value of z")
64%%keep objects for inspection : do not rm(x,y,z,xz,fit)
65
66## Bivariate example -- loads pkg "tripack"
67require(tripack)
68require(akima)
69data(CobarOre)
70fit <- rqss(z ~ qss(cbind(x,y), lambda=.08),
71            data= CobarOre)
72plot(fit, col="grey",
73     main = "CobarOre data -- rqss(z ~ qss(cbind(x,y)))")
74T <- with(CobarOre, tri.mesh(x, y))
75set.seed(77)
76ndum <- 100
77xd <- with(CobarOre, runif(ndum, min(x), max(x)))
78yd <- with(CobarOre, runif(ndum, min(y), max(y)))
79table(s <- in.convex.hull(T, xd, yd))
80pred <- predict(fit, data.frame(x = xd[s], y = yd[s]))
81contour(interp(xd[s],yd[s], pred),
82        col="red", add = TRUE)
83}
84\keyword{regression}
85\keyword{smooth}
86\keyword{robust}
87