1
vote
$\begingroup$

I want to create a model matrix for a quadratic regression in R. I have a matrix $X$, where each column is a covariate. Currently, I'm using the code:

Model = model.matrix(~ 1 + X[,1] + X[,2] + X[,3] + I(X[,1]^2/2) + 
I(X[,2]^2/2) + I(X[,3]^2/2)  + I(X[,1]*X[,2]) + I(X[,1]*X[,3])+I(X[,2]*X[,3]))

Then if I'm dealing with a problem with 4 or more covariates, I have to add the new terms by hand. I tried to create the model matrix myself using:

Model = cbind(1,X,X^2/2,mixed_terms)

I don't remember exactly how I was calculating the mixed terms (I think I used a permutation of the columns of $X$), but in any case my code was much slower in the model.matrix, while I need the model matrix to be created fast because it goes inside a loop in a simulation routine. Thanks!

$\endgroup$

2 Answers 2

Reset to default
1
vote
$\begingroup$

Just use the poly function in your formula, if you give it a matrix as an argument then it will give all the interactions and powers:

Model <- model.matrix( ~ poly(X, degree=2, raw=TRUE) )

This will work for any number of columns. The raw=TRUE prevents it from orthogonalizing the data so you can compare to your hand calculations, but in general your results will probably be better if you leave it as the default of FALSE.

$\endgroup$
2
  • $\begingroup$ That's great, and it's faster than my code. Maybe I'm asking too much: is there a way to have the polynomial terms ordered the way I want (i.e. first linear terms, then quadratic and finally mixed terms)? $\endgroup$
    – Matteo Fasiolo
    Nov 8, 2012 at 18:33
  • $\begingroup$ You can reorder the columns yourself after the fact (should be fairly quick). The linear terms are columns cumsum(1:3) and the quadratic terms are cumsum(1:3)+1 then everything else (assuming you only go to degree=2). $\endgroup$
    – Greg Snow
    Nov 8, 2012 at 19:26
1
vote
$\begingroup$

In the end I did it as follow:

comb <- combinations(n.par,2)
Model = cbind(1, X, (X^2)/2)
for(j in 1:nrow(comb)){Model = cbind(Model,X[,comb[j,1]]*X[,comb[j,2]])}

which turned out to be faster then using model matrix, and also I got the columns in the order that I wanted.

$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.