func eval a b c x .
   return a + (b + c * x) * x
.
proc regression xa[] ya[] . .
   n = len xa[]
   for i = 1 to n
      xm = xm + xa[i]
      ym = ym + ya[i]
      x2m = x2m + xa[i] * xa[i]
      x3m = x3m + xa[i] * xa[i] * xa[i]
      x4m = x4m + xa[i] * xa[i] * xa[i] * xa[i]
      xym = xym + xa[i] * ya[i]
      x2ym = x2ym + xa[i] * xa[i] * ya[i]
   .
   xm = xm / n
   ym = ym / n
   x2m = x2m / n
   x3m = x3m / n
   x4m = x4m / n
   xym = xym / n
   x2ym = x2ym / n
   #
   sxx = x2m - xm * xm
   sxy = xym - xm * ym
   sxx2 = x3m - xm * x2m
   sx2x2 = x4m - x2m * x2m
   sx2y = x2ym - x2m * ym
   #
   b = (sxy * sx2x2 - sx2y * sxx2) / (sxx * sx2x2 - sxx2 * sxx2)
   c = (sx2y * sxx - sxy * sxx2) / (sxx * sx2x2 - sxx2 * sxx2)
   a = ym - b * xm - c * x2m
   print "y = " & a & " + " & b & "x + " & c & "x^2"
   numfmt 0 3
   for i = 1 to n
      print xa[i] & "  " & ya[i] & " " & eval a b c xa[i]
   .
.
xa[] = [ 0 1 2 3 4 5 6 7 8 9 10 ]
ya[] = [ 1 6 17 34 57 86 121 162 209 262 321 ]
regression xa[] ya[]