PolynomialRegression.java

Below is the syntax highlighted version of PolynomialRegression.java from §9.7 Optimization.

```/******************************************************************************
*  Compilation:  javac -cp .:jama.jar PolynomialRegression.java
*  Execution:    java  -cp .:jama.jar PolynomialRegression
*
*  % java -cp .:jama.jar PolynomialRegression
*  0.01 N^3 + -1.64 N^2 + 168.92 N + -2113.73  (R^2 = 0.997)
*
******************************************************************************/

import Jama.Matrix;
import Jama.QRDecomposition;

public class PolynomialRegression {
private final int N;
private final int degree;
private final Matrix beta;
private double SSE;
private double SST;

public PolynomialRegression(double[] x, double[] y, int degree) {
this.degree = degree;
N = x.length;

// build Vandermonde matrix
double[][] vandermonde = new double[N][degree+1];
for (int i = 0; i < N; i++) {
for (int j = 0; j <= degree; j++) {
vandermonde[i][j] = Math.pow(x[i], j);
}
}
Matrix X = new Matrix(vandermonde);

// create matrix from vector
Matrix Y = new Matrix(y, N);

// find least squares solution
QRDecomposition qr = new QRDecomposition(X);
beta = qr.solve(Y);

// mean of y[] values
double sum = 0.0;
for (int i = 0; i < N; i++)
sum += y[i];
double mean = sum / N;

// total variation to be accounted for
for (int i = 0; i < N; i++) {
double dev = y[i] - mean;
SST += dev*dev;
}

// variation not accounted for
Matrix residuals = X.times(beta).minus(Y);
SSE = residuals.norm2() * residuals.norm2();

}

public double beta(int j) {
return beta.get(j, 0);
}

public int degree() {
return degree;
}

public double R2() {
return 1.0 - SSE/SST;
}

// predicted y value corresponding to x
public double predict(double x) {
// horner's method
double y = 0.0;
for (int j = degree; j >= 0; j--)
y = beta(j) + (x * y);
return y;
}

public String toString() {
String s = "";
int j = degree;

while (Math.abs(beta(j)) < 1E-5)
j--;

// create remaining terms
for (j = j; j >= 0; j--) {
if      (j == 0) s += String.format("%.2f ", beta(j));
else if (j == 1) s += String.format("%.2f N + ", beta(j));
else             s += String.format("%.2f N^%d + ", beta(j), j);
}
return s + "  (R^2 = " + String.format("%.3f", R2()) + ")";
}

public static void main(String[] args) {
double[] x = { 10, 20, 40, 80, 160, 200 };
double[] y = { 100, 350, 1500, 6700, 20160, 40000 };
PolynomialRegression regression = new PolynomialRegression(x, y, 3);
StdOut.println(regression);
}
}

```