summaryrefslogtreecommitdiff
path: root/graphics/asymptote/examples/linearregression.asy
diff options
context:
space:
mode:
Diffstat (limited to 'graphics/asymptote/examples/linearregression.asy')
-rw-r--r--graphics/asymptote/examples/linearregression.asy64
1 files changed, 64 insertions, 0 deletions
diff --git a/graphics/asymptote/examples/linearregression.asy b/graphics/asymptote/examples/linearregression.asy
new file mode 100644
index 0000000000..a1cc6ff6bc
--- /dev/null
+++ b/graphics/asymptote/examples/linearregression.asy
@@ -0,0 +1,64 @@
+import graph3;
+import math; // for the leastsquares routine
+
+Billboard.targetsize = true; // Perspective should not affect the labels.
+currentprojection = perspective(60 * (5, 2, 3));
+
+file duncan = input("linearregression.dat");
+
+string headers = duncan;
+
+real[][] independentvars;
+real[] dependentvars;
+
+while (!eof(duncan)) {
+ string line = duncan;
+ string[] entries = split(line);
+ if (entries.length < 5) continue;
+ string type = entries[1];
+ real income = (real)(entries[2]);
+ real education = (real)(entries[3]);
+ real prestige = (real)(entries[4]);
+
+ // include 1.0 for the residue
+ independentvars.push(new real[] {income, education, 1.0});
+ dependentvars.push(prestige);
+}
+
+real[] coeffs = leastsquares(independentvars, dependentvars, warn=false);
+if (coeffs.length == 0) {
+ abort("Unable to find regression: independent variables are "
+ + "linearly dependent.");
+}
+
+real f(pair xy) {
+ return coeffs[0] * xy.x // income
+ + coeffs[1] * xy.y // education
+ + coeffs[2]; // residue
+}
+
+real xmin = infinity, xmax = -infinity, ymin = infinity, ymax = -infinity;
+for (real[] row : independentvars) {
+ if (row[0] < xmin) xmin = row[0];
+ if (row[0] > xmax) xmax = row[0];
+ if (row[1] < ymin) ymin = row[1];
+ if (row[1] > ymax) ymax = row[1];
+}
+
+// Draw the plane
+draw(surface(f, (xmin, ymin), (xmax, ymax)),
+ surfacepen=emissive(blue + opacity(0.6)),
+ meshpen = blue);
+
+for (int ii = 0; ii < independentvars.length; ++ii) {
+ triple pt = (independentvars[ii][0], independentvars[ii][1],
+ dependentvars[ii]);
+ draw(shift(pt) * unitsphere, material(yellow, emissivepen=0.2*yellow));
+ real z = f((pt.x, pt.y));
+ if (pt.z > z) draw (pt -- (pt.x, pt.y, z), green);
+ else draw(pt -- (pt.x, pt.y, z), red);
+}
+
+xaxis3("income", Bounds(Min, Min), InTicks);
+yaxis3("education", Bounds(Min, Min), InTicks);
+zaxis3("prestige", Bounds(Min, Min), InTicks);