Daniel Sabanés Bové
2009-Feb-07 21:00 UTC
[R] paraPen in gam [mgcv 1.4-1.1] and centering constraints
Dear Mr. Simon Wood, dear list members, I am trying to fit a similar model with gam from mgcv compared to what I did with BayesX, and have discovered the relatively new possibility of incorporating user-defined matrices for quadratic penalties on parametric terms using the "paraPen" argument. This was really a very good idea! However, I would like to constraint the coefficients belonging to one penalty matrix to sum to zero. So I would like to have the same centering constraint on user-defined penalized coefficient groups like it is implemented for the spline smoothing terms. The reason is that I have actually a factor coding different regions, and the penalty matrix results from the neighborhood structure in a Gaussian Markov Random Field (GMRF). So I can't choose one region as the reference category, because then the structure in the other regions would not contain the same information as before... Is there a way to constraint a group of coefficients to sum to zero? Thanks in advance, Daniel Sabanes
Simon Wood
2009-Feb-09 10:07 UTC
[R] paraPen in gam [mgcv 1.4-1.1] and centering constraints
Daniel, There's nothing built in to constraint the coefficients in this way, but it's not too difficult to reparameterize to impose any linear constraint. The key is to find an orthogonal basis for null space of the constraint. Then it's easy to re-parameterize to work in that space. Here's a simple linear model based example .... ## simulated problem: linear model subject to constraint ## i.e y = X b + e subject to C b=0 X <- matrix(runif(40),10,4) ## a model matrix y <- rnorm(10) ## a response C <- matrix(1,1,4) ## a constraint matrix so C b = 0 ## Get a basis for null space of constraint... qrc <- qr(t(C)) Z <- qr.Q(qrc,complete=TRUE)[,(nrow(C)+1):ncol(C)] ## absorb constraint into basis... XZ <- X%*%Z ## fit model.... b <- lm(y~XZ-1) ## back to original parameterization beta <- Z%*%coef(b) sum(beta) ## it works! ## If there had been a penalty matrix, S, as well, ## then it should have been transformed as follows.... S <- t(Z)%*%S%*%Z ... So provided that you have the model matrix for the term that you want to penalize in `gam' then it's just a matter of transforming that model matrix and corresponding penalty/ies, using a null space matrix like Z. Note that the explicit formation of Z is not optimally efficient here, but this won't have a noticeable impact on the total computational cost in this context anyway (given that mgcv is not able to make use of all that lovely sparcity in the MRF, :-( ). Hope this helps. best, Simon On Saturday 07 February 2009 21:00, Daniel Saban?s Bov? wrote:> Dear Mr. Simon Wood, dear list members, > > I am trying to fit a similar model with gam from mgcv compared to what I > did with BayesX, and have discovered the relatively new possibility of > incorporating user-defined matrices for quadratic penalties on > parametric terms using the "paraPen" argument. This was really a very > good idea! > > However, I would like to constraint the coefficients belonging to one > penalty matrix to sum to zero. So I would like to have the same > centering constraint on user-defined penalized coefficient groups like > it is implemented for the spline smoothing terms. The reason is that I > have actually a factor coding different regions, and the penalty matrix > results from the neighborhood structure in a Gaussian Markov Random > Field (GMRF). So I can't choose one region as the reference category, > because then the structure in the other regions would not contain the > same information as before... > > Is there a way to constraint a group of coefficients to sum to zero? > > Thanks in advance, > Daniel Sabanes--> Simon Wood, Mathematical Sciences, University of Bath, Bath, BA2 7AY UK > +44 1225 386603 www.maths.bath.ac.uk/~sw283