code: func PolynomialRegression(xs, ys, weights []float64, degree int) PolynomialRegressionResult {
  	terms := make([]func(xs, termOut []float64), degree+1)
  	terms[0] = func(xs, termsOut []float64) {
  		for i := range termsOut {
  			termsOut[i] = 1
  		}
  	}
  	if degree >= 1 {
  		terms[1] = func(xs, termOut []float64) {
  			copy(termOut, xs)
  		}
  	}
  	if degree >= 2 {
  		terms[2] = func(xs, termOut []float64) {
  			for i, x := range xs {
  				termOut[i] = x * x
  			}
  		}
  	}
  	for d := 3; d < len(terms); d++ {
  		d := d
  		terms[d] = func(xs, termOut []float64) {
  			for i, x := range xs {
  				termOut[i] = math.Pow(x, float64(d+1))
  			}
  		}
  	}
  
  	coeffs := LinearLeastSquares(xs, ys, weights, terms...)
  	f := func(x float64) float64 {
  		y := coeffs[0]
  		xp := x
  		for _, c := range coeffs[1:] {
  			y += xp * c
  			xp *= x
  		}
  		return y
  	}
  	return PolynomialRegressionResult{coeffs, f}
  }
docs: // PolynomialRegression performs a least squares regression with a
  // polynomial of the given degree. If weights is non-nil, it is used
  // to weight the residuals.