When creating regressors (or contrasts) for a parametric experimental design, it is often useful to model (or test) linear, quadratic and cubic effects separately. However, if you just use the usual x, x*x, x*x*x approach to create predictors, you wil create highly correlated variables. Specifically, the quadratic term will include some of the linear effect, and vice-versa. The solution is to use orthogonal polynomials.

The very simplest way to do this is with Matlabâ€™s `orth`

function. E.g. for 10-level parametric regressor

n=10; x=[1:n]'; One=ones(n,1); X0=[One x x.^2 x.^3]; X=orth(X0); plot(X)

But annoyingly, this standardizes the columns and produces a fairly uninterprtable basis. In particular, while a linear trend can be represented by a linear combination of those columns, it is not in any one single column.

Instead, we can carefully build up a basis that is mean zero, and uses successive othogonalization to ensure the interpretability of each column.

n=10; x=[1:n]'; One=ones(n,1); X0=[x x.^2 x.^3]; X=One; for i=1:size(X0,2) X=[X X0(:,i)-X*pinv(X)*X0(:,i)]; end plot(X)

*Now* we have three columns readily identifiable as intercept, linear, quadratic and cubic terms; after the first column, they are mean zero (so that they do not carry any of the grand mean), and all are orthogonal, as can be checked with `corrcoef(X)`

.

There are of course specialized functions to create orthogonal polynomials, but for this purpose, this code snippet is all you need.