Constrained Generalized Additive Model Fitting
Description
The partial linear generalized additive model is fitted using the method of maximum likelihood, where shape or order restrictions can be imposed on the nonparametrically modelled predictors with optional smoothing, and no restrictions are imposed on the optional parametrically modelled covariate.
Usage
1 2 
Arguments
formula 
A formula object which gives a symbolic description of the model to be fitted. It has the form "response ~ predictor". The response is a vector of length n. The specification of the model can be one of the three exponential families: gaussian, binomial and poisson. The systematic component η is E(y), the log odds of y = 1, and the logarithm of E(y) respectively. A predictor can be a nonparametrically modelled variable with or without a shape or order restriction, or a parametrically modelled unconstrained covariate. In terms of a nonparametrically modelled predictor, the user is supposed to indicate the relationship between the systematic component η and a predictor x in the following way: Assume that η is the systematic component and x is a predictor:

nsim 
The number of simulations used to get the cic parameter. The default is nsim = 1e+2. 
family 
A parameter indicating the error distribution and link function to be used in the model. It can be a character string naming a family function or the result of a call to a family function. This is borrowed from the glm routine in the stats package. There are three families used in cgam: gaussian, binomial and poisson. 
cpar 
A multiplier to estimate the model variance, which is defined as σ^2 = SSR / (n  d_0  cpar * edf). SSR is the sum of squared residuals for the full model, d_0 is the dimension of the linear space in the cone, and edf is the effective degrees of freedom. The default is cpar = 1.2. See Meyer, M. C. and M. Woodroofe (2000) for more details. 
data 
An optional data frame, list or environment containing the variables in the model. The default is data = NULL. 
weights 
An optional nonnegative vector of "replicate weights" which has the same length as the response vector. If weights are not given, all weights are taken to equal 1. The default is weights = NULL. 
sc_x 
A logical scalar showing if or not nonparametrically modelled predictors will be scaled. The default is sc_x = FALSE. 
sc_y 
A logical scalar showing if or not the response will be scaled. The default is sc_y = FALSE. 
Details
We consider generalized partial linear models with independent observations from an exponential family of the form p(y_i;θ,τ) = exp[\{y_iθ_i  b(θ_i)\}τ  c(y_i, τ)], i = 1,…,n, where the specifications of the functions b and c determine the subfamily of models. The mean vector μ = E(y) has values μ_i = b'(θ_i), and is related to a design matrix of predictor variables through a monotonically increasing link function g(μ_i) = η_i, i = 1,…,n, where η is the systematic component and describes the relationship with the predictors. The relationship between η and θ is determined by the link function b.
For the additive model, the systematic component is specified for each observation by η_i = f_1(x_{1i}) + … + f_L(x_{Li}) + z_i'β, where the functions f_l describe the relationships of the nonparametrically modelled predictors x_l, β is a parameter vector, and z_i contains the values of variables to be modelled parametrically. The nonparametric components are modelled with shape or order assumptions with optional smoothing, and the solution is obtained through an iteratively reweighted cone projection, with no backfitting of individual components.
Suppose that η is a n by 1 vector. The matrix form of the systematic component and the predictor is η = φ_1 + … + φ_L + Zβ, where φ_l is the individual component for the lth nonparametrically modelled predictor, l = 1, …, L, and Z is an n by p design matrix for the parametrically modelled covariate.
To model the component φ_l, smooth regression splines or nonsmooth ordinal basis functions can be used. The constraints for the component φ_l are in C_l. In the first case, C_l = \{φ_l \in R^n: φ_l = v_l+B_lβ_l, where β_l ≥ 0 and v_l\in V_l \}, where B_l has regression splines as columns and V_l is the linear space contained in C_l, and in the second case, C_l = \{φ \in R^n: A_lφ ≥ 0 and B_lφ = 0\}, for inequality constraint matrix A_l and equality constraint matrix B_l.
The set C_l is a convex cone and the set C = C_1 + … + C_p + Z is also a convex cone with a finite set of edges, where the edges are the generators of C, and Z is the column space of the design matrix Z for the parametrically modelled covariate.
An iteratively reweighted cone projection algorithm is used to fit the generalized regression model over the cone C.
See references cited in this section and the official manual (https://cran.rproject.org/package=coneproj) for the R package coneproj for more details.
Value
etahat 
The fitted systematic component η. 
muhat 
The fitted mean value, obtained by transforming the systematic component η by the inverse of the link function. 
vcoefs 
The estimated coefficients for the basis spanning the null space of the constraint set. 
xcoefs 
The estimated coefficients for the edges corresponding to the smooth predictors with no shape constraint and shaperestricted predictors. 
zcoefs 
The estimated coefficients for the parametrically modelled covariate, i.e., the estimation for the vector β. 
ucoefs 
The estimated coefficients for the edges corresponding to the predictors with an umbrellaordering constraint. 
tcoefs 
The estimated coefficients for the edges corresponding to the predictors with a treeordering constraint. 
coefs 
The estimated coefficients for the basis spanning the null space of the constraint set and edges corresponding to the shaperestricted and orderrestricted predictors. 
cic 
The cone information criterion proposed in Meyer(2013a). It uses the "null expected degrees of freedom" as a measure of the complexity of the model. See Meyer(2013a) for further details of cic. 
d0 
The dimension of the linear space contained in the cone generated by all constraint conditions. 
edf0 
The estimated "null expected degrees of freedom". It is a measure of the complexity of the model. See Meyer (2013a) and Meyer (2013b) for further details. 
etacomps 
The fitted systematic component value for nonparametrically modelled predictors. It is a matrix of which each row is the fitted systematic component value for a nonparametrically modelled predictor. If there are more than one such predictors, the order of the rows is the same as the order that the user defines such predictors in the formula argument of cgam. 
xmat 
A matrix whose columns represent the shaperestricted predictors and smooth predictors with no shape constraint. 
zmat 
A matrix whose columns represent the basis for the parametrically modelled covariate. The user can choose to include a constant vector in it or not. It must have full column rank. 
ztb 
A list keeping track of the order of the parametrically modelled covariate. 
tr 
A matrix whose columns represent the predictors with a treeordering constraint. 
umb 
A matrix whose columns represent the predictors with an umbrellaordering constraint. 
tree.delta 
A matrix whose rows are the edges corresponding to the predictors with a treeordering constraint. 
umbrella.delta 
A matrix whose rows are the edges corresponding to the predictors with an umbrellaordering constraint. 
bigmat 
A matrix whose rows are the basis spanning the null space of the constraint set and the edges corresponding to the shaperestricted and orderrestricted predictors. 
shapes 
A vector including the shape and partialordering constraints in a cgam fit. 
shapesx 
A vector including the shape constraints in a cgam fit. 
wt 
The weights in the final iteration of the iteratively reweighted cone projections. 
wt.iter 
A logical scalar indicating if or not iteratively reweighted cone projections may be used. If the response is gaussian, then wt.iter = FALSE; if the response is binomial or poisson, then wt.iter = TRUE. 
family 
The family parameter defined in a cgam formula. 
SSE0 
The sum of squared residuals for the linear part. 
SSE1 
The sum of squared residuals for the full model. 
pvals.beta 
The approximate pvalues for the estimation of the vector β. A tdistribution is used as the approximate distribution. 
se.beta 
The standard errors for the estimation of the vector β. 
null_df 
The degree of freedom for the null model of a cgam fit, i.e., the model only containing a constant vector. 
df 
The degree of freedom for the null space of a cgam fit. 
resid_df_obs 
The observed degree of freedom for the residuals of a cgam fit. 
null_deviance 
The deviance for the null model of a cgam fit, i.e., the model only containing a constant vector. 
deviance 
The residual deviance of a cgam fit. 
tms 
The terms objects extracted by the generic function terms from a cgam fit. See the official help page (http://stat.ethz.ch/Rmanual/Rpatched/library/stats/html/terms.html) of the terms function for more details. 
capm 
The number of edges corresponding to the shaperestricted predictors. 
capms 
The number of edges corresponding to the smooth predictors with no shape constraint. 
capk 
The number of nonconstant columns of zmat. 
capt 
The number of edges corresponding to the treeordering predictors. 
capu 
The number of edges corresponding to the umbrellaordering predictors. 
xid1 
A vector keeping track of the beginning position of the set of edges in bigmat for each shaperestricted predictor and smooth predictor with no shape constraint in xmat. 
xid2 
A vector keeping track of the end position of the set of edges in bigmat for each shaperestricted predictor and smooth predictor with no shape constraint in xmat. 
tid1 
A vector keeping track of the beginning position of the set of edges in bigmat for each treeordering factor in tr. 
tid2 
A vector keeping track of the end position of the set of edges in bigmat for each treeordering factor in tr. 
uid1 
A vector keeping track of the beginning position of the set of edges in bigmat for each umbrellaordering factor in umb. 
uid2 
A vector keeping track of the end position of the set of edges in bigmat for each umbrellaordering factor in umb. 
zid 
A vector keeping track of the positions of the parametrically modelled covariate. 
vals 
A vector storing the levels of each variable used as a factor. 
zid1 
A vector keeping track of the beginning position of the levels of each variable used as a factor. 
zid2 
A vector keeping track of the end position of the levels of each variable used as a factor. 
nsim 
The number of simulations used to get the cic parameter. 
xnms 
A vector storing the name of the shaperestricted predictor and the smooth predictor with no shape constraint in xmat. 
ynm 
The name of the response variable. 
znms 
A vector storing the name of the parametrically modelled covariate. 
is_param 
A logical scalar showing if or not a variable is a parametrically modelled covariate, which could be a linear term or a factor. 
is_fac 
A logical scalar showing if or not a variable is a factor. 
knots 
A list storing the knots used for each shaperestricted predictor and smooth predictor with no shape constraint. For a smooth, constrained and a smooth, unconstrainted predictor, knots is a vector of more than 1 elements, and for a shaperestricted predictor without smoothing, knots = 0. 
numknots 
A vector storing the number of knots for each shaperestricted predictor and smooth predictor with no shape constraint. For a smooth, constrained and a smooth, unconstrainted predictor, numknots > 1, and for a shaperestricted predictor without smoothing, numknots = 0. 
sps 
A character vector storing the space parameter to create knots for each shaperestricted predictor. 
ms 
The centering terms used to make edges for shaperestricted predictors. 
cpar 
The cpar argument in the cgam formula 
pl 
A vector of the placebo values for the nonparametrically modelled covariates with a treeordering constraint. 
call 
The matched call. 
Author(s)
Mary C. Meyer and Xiyue Liao
References
Meyer, M. C. (2013a) Semiparametric additive constrained regression. Journal of Nonparametric Statistics 25(3), 715
Meyer, M. C. (2013b) A simple new algorithm for quadratic programming with applications in statistics. Communications in Statistics 42(5), 1126–1139.
Meyer, M. C. and M. Woodroofe (2000) On the degrees of freedom in shaperestricted regression. Annals of Statistics 28, 1083–1104.
Meyer, M. C. (2008) Inference using shaperestricted regression splines. Annals of Applied Statistics 2(3), 1013–1033.
Mammen, E. and K. Yu (2007) Additive isotonic regression. IMS Lecture NotesMonograph Series Asymptotics: Particles, Process, and Inverse Problems 55, 179–195.
Huang, J. (2002) A note on estimating a partly linear model under monotonicity constraints. Journal of Statistical Planning and Inference 107, 343–351.
Cheng, G.(2009) Semiparametric additive isotonic regression. Journal of Statistical Planning and Inference 139, 1980–1991.
Bacchetti, P. (1989) Additive isotonic models. Journal of the American Statistical Association 84(405), 289–294.
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 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47  # Example 1.
data(cubic)
# extract x
x < cubic$x
# extract y
y < cubic$y
# regress y on x with no restriction with lm()
fit.lm < lm(y ~ x + I(x^2) + I(x^3))
# regress y on x under the restriction: "increasing and convex"
fit.cgam < cgam(y ~ incr.conv(x))
# make a plot to compare the two fits
par(mar = c(4, 4, 1, 1))
plot(x, y, cex = .7, xlab = "x", ylab = "y")
lines(x, fit.cgam$muhat, col = 2, lty = 2)
lines(x, fitted(fit.lm), col = 1, lty = 1)
legend("topleft", bty = "n", c("constrained cgam fit", "unconstrained lm fit"),
lty = c(2, 1), col = c(2, 1))
# Example 2.
## Not run:
library(gam)
data(kyphosis)
# regress Kyphosis on Age, Number, and Start under the restrictions:
# "concave", "increasing and concave", and "decreasing and concave"
fit < cgam(Kyphosis ~ conc(Age) + incr.conc(Number) + decr.conc(Start),
family = binomial(), data = kyphosis)
## End(Not run)
# Example 3.
library(MASS)
data(Rubber)
# regress loss on hard and tens under the restrictions:
# "decreasing" and "decreasing"
fit.cgam < cgam(loss ~ decr(hard) + decr(tens), data = Rubber)
# "smooth and decreasing" and "smooth and decreasing"
fit.cgam.s < cgam(loss ~ s.decr(hard) + s.decr(tens), data = Rubber)
# make a 3D plot based on fit.cgam and fit.cgam.s
plotpersp(fit.cgam, th = 120, main = "3D Plot of a Cgam Fit")
plotpersp(fit.cgam.s, tens, hard, data = Rubber, th = 120, main = "3D Plot of a Smooth Cgam Fit")
