4
$\begingroup$

I'm trying to find the best polynomial regression for a dataset where the polynomial's power is between 2 and 10. So the regression can have an x10 term at most in it. The dataset itself is simply a set of x and y pairs as follows:

1,15.3
2,66.0
3,272.5
4,814.8

As I understand, the normal way to do a polynomial regression is to simply apply the power transformation to the x vector (i.e. put every element in the vector to the power of 6), add this vector to the dataset, and then treat this transformed vector as another independent variable.

However if I try this approach with high enough powers (usually 6 and above), my regression library (statsample for Ruby) tells me "Regressors are linearly dependent", and throws an error. I know that technically the x vectors are dependent on each other since they were derived from each other, but it's certainly not a linear dependency (where one is the same as the other multiplied by a scalar). What's going on? And what does this mean?

As an example, here's an example of my code in Ruby (I'm told that this library is a lot like R however, for all you R users):

# Read the dataset (like an R data frame??) containing just x and y (see the example above)
ds = CSV.read(file_name)

# For each possible polynomial power between 2 and 10, add a new vector which consists
# of the x vector to that power, and then run a regression
(2..10).each do |i|

  # Add a vector called x2, x3, x4 etc. and apply the power transformation
  ds.add_vector("x#{i}", ds["x"].map{ |x| x**i }.to_scale)

  #Run the regression
  reg = lr(ds,'y')
end
$\endgroup$
3
$\begingroup$

Recall from linear algebra that linearly dependent vectors are a set of vectors which can be expressed as a linear combination of each other. When performing regression, this creates problems because the matrix $X^TX$ is singular, so there is not a uniquely defined solution to estimating your regression coefficients. (The matrix $X^TX$ is important because of its role in estimating OLS regression coefficients $\hat{\beta}$: $\hat{\beta}=(X^TX)^{-1}X^Ty$.)

Linear dependence is a technical phenomenon distinct from the ordinary usage of "dependence" as you express it.

So with an understanding of linear dependence in hand, we can begin to examine likely sources of this problem.

  1. More features than rows It's unclear how much data you have, but it's possible that by adding polynomial terms, you're inadvertently creating a dimensionality problem: an ordinary regression with more columns than observations will fail! This is because the system of equations that you've defined has infinitely many solutions, i.e. there is not a pivot point in every column. This amounts to the same problem as having a singular $X^TX$.

  2. Duplicate data Even if you have more rows than features, it's important that each of those rows provides unique information. Duplicate rows are unhelpful in this context. By definition, the polynomial regression matrix $X$ (consisting of $[\mathbf{1} , x, x^2, ... x^m]$ as columns) is a Vandermonde matrix, so it will have a unique OLS solution provided that there are $m+1$ unique entries in your original $x$ vector. So despite having 100 observations, perhaps you have a smaller number of unique $x$ entries?

  3. Ill-conditioned matrix Even if neither (1) nor (2) is true, it's possible that the matrix is numerically singular, i.e. singular due to machine precision reasons. There are a few methods to dealing with this. One is to enforce 0 mean unit variance columns (just subtract the column mean from each column and divide by the column's standard deviation). This is an unusual problem with modern software, but it's possible that exponents of small/large numbers are the culprit here.

A general strategy to address an ill-conditioned matrix $X$ is called ridge regression, which works by finding an optimal amount of regularization to apply to your problem. Ridge regression is discussed all over this website. One place to start would be with this excellent answer.

AndyW points out that fitting very high-order polynomials is often ill-advised since it increases the risk of overfitting. In predictive settings, it's often advised that one use cross validation to assess the fitness of a given model. Depending on your application, you might care about different out-of-sample tests, but a typical metric for this type of problem is mean squared error.

$\endgroup$
  • $\begingroup$ Unfortunately I have very little stats/linear algebra background (university calculus is the extent of my knowledge), but I can understand that certain matrices are singular/have no inverse. What causes this to happen? Is it a fault in the data, or just some quirk of matrices? Also my actual data set has 100 cases (rows), so I suspect that isn't the problem. $\endgroup$ – Miguel Mar 26 '15 at 13:24
  • $\begingroup$ @Miguel Please see my expanded answer. Perhaps you are experiencing issue (2) or (3)? As for why some matrices have no inverse, it's essentially just one of those facts of linear algebra, related to how one computes determinants. But choices about how we represent our data have ramifications for matrix computations; some of these issues are addressed in my post. $\endgroup$ – Sycorax Mar 26 '15 at 14:23
  • $\begingroup$ +1, definitely without centering the design matrix can be ill-behaved for even the lower polynomials without even that high of values. That high of polynomials would be a problem for many practical applications I imagine, see this R code for one example set.seed(10);x <- 10 + rnorm(100);y <- 0.5 + 3*x + rnorm(100);summary(lm(y ~ poly(x,10,raw=TRUE))). Using orthogonal polynomials allows one to estimate the coefficients in that example, i.e. summary(lm(y ~ poly(x,10,raw=FALSE))), but still that many polynomials are unlikely to be reasonable in any real application. $\endgroup$ – Andy W Mar 26 '15 at 15:42
  • $\begingroup$ Thank you to everyone for your fantastic explanations. I'm still not sure that we've found the exact reason my code isn't working, but I've had a look at the stat library's source (here) and it seems that generates a matrix out of the x vertices (i.e. the entire dataset with y removed), and if the determinant is less than 1e-15, it throws this error. Does this make the error any clearer? $\endgroup$ – Miguel Mar 28 '15 at 6:22
  • 1
    $\begingroup$ Given your description of determinants, you are experiencing issue (3) from my answer. The reason the software checks if the determinant is small is because as the determinant approaches $0$, the matrix becomes more and more ill-conditioned until it's finally a singular matrix as $\text{det}(A)=0$. So it's still the same problem. In any case, a high-order polynomial regression is usually ill-advised, especially on such a small data. $\endgroup$ – Sycorax Mar 28 '15 at 12:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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