[R] paraPen in gam [mgcv 1.4-1.1] and centering constraints

[Simon Wood]

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