d3-regression.esm.js 20 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869
  1. // https://github.com/HarryStevens/d3-regression#readme Version 1.3.10. Copyright 2022 Harry Stevens.
  2. function _slicedToArray(arr, i) {
  3. return _arrayWithHoles(arr) || _iterableToArrayLimit(arr, i) || _nonIterableRest();
  4. }
  5. function _arrayWithHoles(arr) {
  6. if (Array.isArray(arr)) return arr;
  7. }
  8. function _iterableToArrayLimit(arr, i) {
  9. var _arr = [];
  10. var _n = true;
  11. var _d = false;
  12. var _e = undefined;
  13. try {
  14. for (var _i = arr[Symbol.iterator](), _s; !(_n = (_s = _i.next()).done); _n = true) {
  15. _arr.push(_s.value);
  16. if (i && _arr.length === i) break;
  17. }
  18. } catch (err) {
  19. _d = true;
  20. _e = err;
  21. } finally {
  22. try {
  23. if (!_n && _i["return"] != null) _i["return"]();
  24. } finally {
  25. if (_d) throw _e;
  26. }
  27. }
  28. return _arr;
  29. }
  30. function _nonIterableRest() {
  31. throw new TypeError("Invalid attempt to destructure non-iterable instance");
  32. }
  33. // Adapted from vega-statistics by Jeffrey Heer
  34. // License: https://github.com/vega/vega/blob/f058b099decad9db78301405dd0d2e9d8ba3d51a/LICENSE
  35. // Source: https://github.com/vega/vega/blob/f058b099decad9db78301405dd0d2e9d8ba3d51a/packages/vega-statistics/src/regression/points.js
  36. function points(data, x, y, sort) {
  37. data = data.filter(function (d, i) {
  38. var u = x(d, i),
  39. v = y(d, i);
  40. return u != null && isFinite(u) && v != null && isFinite(v);
  41. });
  42. if (sort) {
  43. data.sort(function (a, b) {
  44. return x(a) - x(b);
  45. });
  46. }
  47. var n = data.length,
  48. X = new Float64Array(n),
  49. Y = new Float64Array(n); // extract values, calculate means
  50. var ux = 0,
  51. uy = 0,
  52. xv,
  53. yv,
  54. d;
  55. for (var i = 0; i < n;) {
  56. d = data[i];
  57. X[i] = xv = +x(d, i, data);
  58. Y[i] = yv = +y(d, i, data);
  59. ++i;
  60. ux += (xv - ux) / i;
  61. uy += (yv - uy) / i;
  62. } // mean center the data
  63. for (var _i = 0; _i < n; ++_i) {
  64. X[_i] -= ux;
  65. Y[_i] -= uy;
  66. }
  67. return [X, Y, ux, uy];
  68. }
  69. function visitPoints(data, x, y, cb) {
  70. var iterations = 0;
  71. for (var i = 0, n = data.length; i < n; i++) {
  72. var d = data[i],
  73. dx = +x(d, i, data),
  74. dy = +y(d, i, data);
  75. if (dx != null && isFinite(dx) && dy != null && isFinite(dy)) {
  76. cb(dx, dy, iterations++);
  77. }
  78. }
  79. }
  80. // return the coefficient of determination, or R squared.
  81. function determination(data, x, y, uY, predict) {
  82. var SSE = 0,
  83. SST = 0;
  84. visitPoints(data, x, y, function (dx, dy) {
  85. var sse = dy - predict(dx),
  86. sst = dy - uY;
  87. SSE += sse * sse;
  88. SST += sst * sst;
  89. });
  90. return 1 - SSE / SST;
  91. }
  92. // Returns the angle of a line in degrees.
  93. function angle(line) {
  94. return Math.atan2(line[1][1] - line[0][1], line[1][0] - line[0][0]) * 180 / Math.PI;
  95. } // Returns the midpoint of a line.
  96. function midpoint(line) {
  97. return [(line[0][0] + line[1][0]) / 2, (line[0][1] + line[1][1]) / 2];
  98. }
  99. // returns a smooth line.
  100. function interpose(xmin, xmax, predict) {
  101. var l = Math.log(xmax - xmin) * Math.LOG10E + 1 | 0;
  102. var precision = 1 * Math.pow(10, -l / 2 - 1),
  103. maxIter = 1e4;
  104. var points = [px(xmin), px(xmax)],
  105. iter = 0;
  106. while (find(points) && iter < maxIter) {
  107. }
  108. return points;
  109. function px(x) {
  110. return [x, predict(x)];
  111. }
  112. function find(points) {
  113. iter++;
  114. var n = points.length;
  115. var found = false;
  116. for (var i = 0; i < n - 1; i++) {
  117. var p0 = points[i],
  118. p1 = points[i + 1],
  119. m = midpoint([p0, p1]),
  120. mp = px(m[0]),
  121. a0 = angle([p0, m]),
  122. a1 = angle([p0, mp]),
  123. a = Math.abs(a0 - a1);
  124. if (a > precision) {
  125. points.splice(i + 1, 0, mp);
  126. found = true;
  127. }
  128. }
  129. return found;
  130. }
  131. }
  132. // Ordinary Least Squares from vega-statistics by Jeffrey Heer
  133. // License: https://github.com/vega/vega/blob/f058b099decad9db78301405dd0d2e9d8ba3d51a/LICENSE
  134. // Source: https://github.com/vega/vega/blob/f058b099decad9db78301405dd0d2e9d8ba3d51a/packages/vega-statistics/src/regression/ols.js
  135. function ols(uX, uY, uXY, uX2) {
  136. var delta = uX2 - uX * uX,
  137. slope = Math.abs(delta) < 1e-24 ? 0 : (uXY - uX * uY) / delta,
  138. intercept = uY - slope * uX;
  139. return [intercept, slope];
  140. }
  141. function exponential () {
  142. var x = function x(d) {
  143. return d[0];
  144. },
  145. y = function y(d) {
  146. return d[1];
  147. },
  148. domain;
  149. function exponential(data) {
  150. var n = 0,
  151. Y = 0,
  152. YL = 0,
  153. XY = 0,
  154. XYL = 0,
  155. X2Y = 0,
  156. xmin = domain ? +domain[0] : Infinity,
  157. xmax = domain ? +domain[1] : -Infinity;
  158. visitPoints(data, x, y, function (dx, dy) {
  159. var ly = Math.log(dy),
  160. xy = dx * dy;
  161. ++n;
  162. Y += (dy - Y) / n;
  163. XY += (xy - XY) / n;
  164. X2Y += (dx * xy - X2Y) / n;
  165. YL += (dy * ly - YL) / n;
  166. XYL += (xy * ly - XYL) / n;
  167. if (!domain) {
  168. if (dx < xmin) xmin = dx;
  169. if (dx > xmax) xmax = dx;
  170. }
  171. });
  172. var _ols = ols(XY / Y, YL / Y, XYL / Y, X2Y / Y),
  173. _ols2 = _slicedToArray(_ols, 2),
  174. a = _ols2[0],
  175. b = _ols2[1];
  176. a = Math.exp(a);
  177. var fn = function fn(x) {
  178. return a * Math.exp(b * x);
  179. },
  180. out = interpose(xmin, xmax, fn);
  181. out.a = a;
  182. out.b = b;
  183. out.predict = fn;
  184. out.rSquared = determination(data, x, y, Y, fn);
  185. return out;
  186. }
  187. exponential.domain = function (arr) {
  188. return arguments.length ? (domain = arr, exponential) : domain;
  189. };
  190. exponential.x = function (fn) {
  191. return arguments.length ? (x = fn, exponential) : x;
  192. };
  193. exponential.y = function (fn) {
  194. return arguments.length ? (y = fn, exponential) : y;
  195. };
  196. return exponential;
  197. }
  198. function linear () {
  199. var x = function x(d) {
  200. return d[0];
  201. },
  202. y = function y(d) {
  203. return d[1];
  204. },
  205. domain;
  206. function linear(data) {
  207. var n = 0,
  208. X = 0,
  209. // sum of x
  210. Y = 0,
  211. // sum of y
  212. XY = 0,
  213. // sum of x * y
  214. X2 = 0,
  215. // sum of x * x
  216. xmin = domain ? +domain[0] : Infinity,
  217. xmax = domain ? +domain[1] : -Infinity;
  218. visitPoints(data, x, y, function (dx, dy) {
  219. ++n;
  220. X += (dx - X) / n;
  221. Y += (dy - Y) / n;
  222. XY += (dx * dy - XY) / n;
  223. X2 += (dx * dx - X2) / n;
  224. if (!domain) {
  225. if (dx < xmin) xmin = dx;
  226. if (dx > xmax) xmax = dx;
  227. }
  228. });
  229. var _ols = ols(X, Y, XY, X2),
  230. _ols2 = _slicedToArray(_ols, 2),
  231. intercept = _ols2[0],
  232. slope = _ols2[1],
  233. fn = function fn(x) {
  234. return slope * x + intercept;
  235. },
  236. out = [[xmin, fn(xmin)], [xmax, fn(xmax)]];
  237. out.a = slope;
  238. out.b = intercept;
  239. out.predict = fn;
  240. out.rSquared = determination(data, x, y, Y, fn);
  241. return out;
  242. }
  243. linear.domain = function (arr) {
  244. return arguments.length ? (domain = arr, linear) : domain;
  245. };
  246. linear.x = function (fn) {
  247. return arguments.length ? (x = fn, linear) : x;
  248. };
  249. linear.y = function (fn) {
  250. return arguments.length ? (y = fn, linear) : y;
  251. };
  252. return linear;
  253. }
  254. // Returns the medium value of an array of numbers.
  255. function median(arr) {
  256. arr.sort(function (a, b) {
  257. return a - b;
  258. });
  259. var i = arr.length / 2;
  260. return i % 1 === 0 ? (arr[i - 1] + arr[i]) / 2 : arr[Math.floor(i)];
  261. }
  262. var maxiters = 2,
  263. epsilon = 1e-12;
  264. function loess () {
  265. var x = function x(d) {
  266. return d[0];
  267. },
  268. y = function y(d) {
  269. return d[1];
  270. },
  271. bandwidth = .3;
  272. function loess(data) {
  273. var _points = points(data, x, y, true),
  274. _points2 = _slicedToArray(_points, 4),
  275. xv = _points2[0],
  276. yv = _points2[1],
  277. ux = _points2[2],
  278. uy = _points2[3],
  279. n = xv.length,
  280. bw = Math.max(2, ~~(bandwidth * n)),
  281. yhat = new Float64Array(n),
  282. residuals = new Float64Array(n),
  283. robustWeights = new Float64Array(n).fill(1);
  284. for (var iter = -1; ++iter <= maxiters;) {
  285. var interval = [0, bw - 1];
  286. for (var i = 0; i < n; ++i) {
  287. var dx = xv[i],
  288. i0 = interval[0],
  289. i1 = interval[1],
  290. edge = dx - xv[i0] > xv[i1] - dx ? i0 : i1;
  291. var W = 0,
  292. X = 0,
  293. Y = 0,
  294. XY = 0,
  295. X2 = 0,
  296. denom = 1 / Math.abs(xv[edge] - dx || 1); // Avoid singularity
  297. for (var k = i0; k <= i1; ++k) {
  298. var xk = xv[k],
  299. yk = yv[k],
  300. w = tricube(Math.abs(dx - xk) * denom) * robustWeights[k],
  301. xkw = xk * w;
  302. W += w;
  303. X += xkw;
  304. Y += yk * w;
  305. XY += yk * xkw;
  306. X2 += xk * xkw;
  307. } // Linear regression fit
  308. var _ols = ols(X / W, Y / W, XY / W, X2 / W),
  309. _ols2 = _slicedToArray(_ols, 2),
  310. a = _ols2[0],
  311. b = _ols2[1];
  312. yhat[i] = a + b * dx;
  313. residuals[i] = Math.abs(yv[i] - yhat[i]);
  314. updateInterval(xv, i + 1, interval);
  315. }
  316. if (iter === maxiters) {
  317. break;
  318. }
  319. var medianResidual = median(residuals);
  320. if (Math.abs(medianResidual) < epsilon) break;
  321. for (var _i = 0, arg, _w; _i < n; ++_i) {
  322. arg = residuals[_i] / (6 * medianResidual); // Default to epsilon (rather than zero) for large deviations
  323. // Keeping weights tiny but non-zero prevents singularites
  324. robustWeights[_i] = arg >= 1 ? epsilon : (_w = 1 - arg * arg) * _w;
  325. }
  326. }
  327. return output(xv, yhat, ux, uy);
  328. }
  329. loess.bandwidth = function (bw) {
  330. return arguments.length ? (bandwidth = bw, loess) : bandwidth;
  331. };
  332. loess.x = function (fn) {
  333. return arguments.length ? (x = fn, loess) : x;
  334. };
  335. loess.y = function (fn) {
  336. return arguments.length ? (y = fn, loess) : y;
  337. };
  338. return loess;
  339. } // Weighting kernel for local regression
  340. function tricube(x) {
  341. return (x = 1 - x * x * x) * x * x;
  342. } // Advance sliding window interval of nearest neighbors
  343. function updateInterval(xv, i, interval) {
  344. var val = xv[i],
  345. left = interval[0],
  346. right = interval[1] + 1;
  347. if (right >= xv.length) return; // Step right if distance to new right edge is <= distance to old left edge
  348. // Step when distance is equal to ensure movement over duplicate x values
  349. while (i > left && xv[right] - val <= val - xv[left]) {
  350. interval[0] = ++left;
  351. interval[1] = right;
  352. ++right;
  353. }
  354. } // Generate smoothed output points
  355. // Average points with repeated x values
  356. function output(xv, yhat, ux, uy) {
  357. var n = xv.length,
  358. out = [];
  359. var i = 0,
  360. cnt = 0,
  361. prev = [],
  362. v;
  363. for (; i < n; ++i) {
  364. v = xv[i] + ux;
  365. if (prev[0] === v) {
  366. // Average output values via online update
  367. prev[1] += (yhat[i] - prev[1]) / ++cnt;
  368. } else {
  369. // Add new output point
  370. cnt = 0;
  371. prev[1] += uy;
  372. prev = [v, yhat[i]];
  373. out.push(prev);
  374. }
  375. }
  376. prev[1] += uy;
  377. return out;
  378. }
  379. function logarithmic () {
  380. var x = function x(d) {
  381. return d[0];
  382. },
  383. y = function y(d) {
  384. return d[1];
  385. },
  386. base = Math.E,
  387. domain;
  388. function logarithmic(data) {
  389. var n = 0,
  390. X = 0,
  391. Y = 0,
  392. XY = 0,
  393. X2 = 0,
  394. xmin = domain ? +domain[0] : Infinity,
  395. xmax = domain ? +domain[1] : -Infinity,
  396. lb = Math.log(base);
  397. visitPoints(data, x, y, function (dx, dy) {
  398. var lx = Math.log(dx) / lb;
  399. ++n;
  400. X += (lx - X) / n;
  401. Y += (dy - Y) / n;
  402. XY += (lx * dy - XY) / n;
  403. X2 += (lx * lx - X2) / n;
  404. if (!domain) {
  405. if (dx < xmin) xmin = dx;
  406. if (dx > xmax) xmax = dx;
  407. }
  408. });
  409. var _ols = ols(X, Y, XY, X2),
  410. _ols2 = _slicedToArray(_ols, 2),
  411. intercept = _ols2[0],
  412. slope = _ols2[1],
  413. fn = function fn(x) {
  414. return slope * Math.log(x) / lb + intercept;
  415. },
  416. out = interpose(xmin, xmax, fn);
  417. out.a = slope;
  418. out.b = intercept;
  419. out.predict = fn;
  420. out.rSquared = determination(data, x, y, Y, fn);
  421. return out;
  422. }
  423. logarithmic.domain = function (arr) {
  424. return arguments.length ? (domain = arr, logarithmic) : domain;
  425. };
  426. logarithmic.x = function (fn) {
  427. return arguments.length ? (x = fn, logarithmic) : x;
  428. };
  429. logarithmic.y = function (fn) {
  430. return arguments.length ? (y = fn, logarithmic) : y;
  431. };
  432. logarithmic.base = function (n) {
  433. return arguments.length ? (base = n, logarithmic) : base;
  434. };
  435. return logarithmic;
  436. }
  437. function quad () {
  438. var x = function x(d) {
  439. return d[0];
  440. },
  441. y = function y(d) {
  442. return d[1];
  443. },
  444. domain;
  445. function quadratic(data) {
  446. var _points = points(data, x, y),
  447. _points2 = _slicedToArray(_points, 4),
  448. xv = _points2[0],
  449. yv = _points2[1],
  450. ux = _points2[2],
  451. uy = _points2[3],
  452. n = xv.length;
  453. var X2 = 0,
  454. X3 = 0,
  455. X4 = 0,
  456. XY = 0,
  457. X2Y = 0,
  458. i,
  459. dx,
  460. dy,
  461. x2;
  462. for (i = 0; i < n;) {
  463. dx = xv[i];
  464. dy = yv[i++];
  465. x2 = dx * dx;
  466. X2 += (x2 - X2) / i;
  467. X3 += (x2 * dx - X3) / i;
  468. X4 += (x2 * x2 - X4) / i;
  469. XY += (dx * dy - XY) / i;
  470. X2Y += (x2 * dy - X2Y) / i;
  471. }
  472. var Y = 0,
  473. n0 = 0,
  474. xmin = domain ? +domain[0] : Infinity,
  475. xmax = domain ? +domain[1] : -Infinity;
  476. visitPoints(data, x, y, function (dx, dy) {
  477. n0++;
  478. Y += (dy - Y) / n0;
  479. if (!domain) {
  480. if (dx < xmin) xmin = dx;
  481. if (dx > xmax) xmax = dx;
  482. }
  483. });
  484. var X2X2 = X4 - X2 * X2,
  485. d = X2 * X2X2 - X3 * X3,
  486. a = (X2Y * X2 - XY * X3) / d,
  487. b = (XY * X2X2 - X2Y * X3) / d,
  488. c = -a * X2,
  489. fn = function fn(x) {
  490. x = x - ux;
  491. return a * x * x + b * x + c + uy;
  492. };
  493. var out = interpose(xmin, xmax, fn);
  494. out.a = a;
  495. out.b = b - 2 * a * ux;
  496. out.c = c - b * ux + a * ux * ux + uy;
  497. out.predict = fn;
  498. out.rSquared = determination(data, x, y, Y, fn);
  499. return out;
  500. }
  501. quadratic.domain = function (arr) {
  502. return arguments.length ? (domain = arr, quadratic) : domain;
  503. };
  504. quadratic.x = function (fn) {
  505. return arguments.length ? (x = fn, quadratic) : x;
  506. };
  507. quadratic.y = function (fn) {
  508. return arguments.length ? (y = fn, quadratic) : y;
  509. };
  510. return quadratic;
  511. }
  512. // Source: https://github.com/Tom-Alexander/regression-js/blob/master/src/regression.js#L246
  513. // License: https://github.com/Tom-Alexander/regression-js/blob/master/LICENSE
  514. // ...with ideas from vega-statistics by Jeffrey Heer
  515. // Source: https://github.com/vega/vega/blob/f21cb8792b4e0cbe2b1a3fd44b0f5db370dbaadb/packages/vega-statistics/src/regression/poly.js
  516. // License: https://github.com/vega/vega/blob/f058b099decad9db78301405dd0d2e9d8ba3d51a/LICENSE
  517. function polynomial () {
  518. var x = function x(d) {
  519. return d[0];
  520. },
  521. y = function y(d) {
  522. return d[1];
  523. },
  524. order = 3,
  525. domain;
  526. function polynomial(data) {
  527. // Use more efficient methods for lower orders
  528. if (order === 1) {
  529. var o = linear().x(x).y(y).domain(domain)(data);
  530. o.coefficients = [o.b, o.a];
  531. delete o.a;
  532. delete o.b;
  533. return o;
  534. }
  535. if (order === 2) {
  536. var _o = quad().x(x).y(y).domain(domain)(data);
  537. _o.coefficients = [_o.c, _o.b, _o.a];
  538. delete _o.a;
  539. delete _o.b;
  540. delete _o.c;
  541. return _o;
  542. }
  543. var _points = points(data, x, y),
  544. _points2 = _slicedToArray(_points, 4),
  545. xv = _points2[0],
  546. yv = _points2[1],
  547. ux = _points2[2],
  548. uy = _points2[3],
  549. n = xv.length,
  550. lhs = [],
  551. rhs = [],
  552. k = order + 1;
  553. var Y = 0,
  554. n0 = 0,
  555. xmin = domain ? +domain[0] : Infinity,
  556. xmax = domain ? +domain[1] : -Infinity;
  557. visitPoints(data, x, y, function (dx, dy) {
  558. ++n0;
  559. Y += (dy - Y) / n0;
  560. if (!domain) {
  561. if (dx < xmin) xmin = dx;
  562. if (dx > xmax) xmax = dx;
  563. }
  564. });
  565. var i, j, l, v, c;
  566. for (i = 0; i < k; ++i) {
  567. for (l = 0, v = 0; l < n; ++l) {
  568. v += Math.pow(xv[l], i) * yv[l];
  569. }
  570. lhs.push(v);
  571. c = new Float64Array(k);
  572. for (j = 0; j < k; ++j) {
  573. for (l = 0, v = 0; l < n; ++l) {
  574. v += Math.pow(xv[l], i + j);
  575. }
  576. c[j] = v;
  577. }
  578. rhs.push(c);
  579. }
  580. rhs.push(lhs);
  581. var coef = gaussianElimination(rhs),
  582. fn = function fn(x) {
  583. x -= ux;
  584. var y = uy + coef[0] + coef[1] * x + coef[2] * x * x;
  585. for (i = 3; i < k; ++i) {
  586. y += coef[i] * Math.pow(x, i);
  587. }
  588. return y;
  589. },
  590. out = interpose(xmin, xmax, fn);
  591. out.coefficients = uncenter(k, coef, -ux, uy);
  592. out.predict = fn;
  593. out.rSquared = determination(data, x, y, Y, fn);
  594. return out;
  595. }
  596. polynomial.domain = function (arr) {
  597. return arguments.length ? (domain = arr, polynomial) : domain;
  598. };
  599. polynomial.x = function (fn) {
  600. return arguments.length ? (x = fn, polynomial) : x;
  601. };
  602. polynomial.y = function (fn) {
  603. return arguments.length ? (y = fn, polynomial) : y;
  604. };
  605. polynomial.order = function (n) {
  606. return arguments.length ? (order = n, polynomial) : order;
  607. };
  608. return polynomial;
  609. }
  610. function uncenter(k, a, x, y) {
  611. var z = Array(k);
  612. var i, j, v, c; // initialize to zero
  613. for (i = 0; i < k; ++i) {
  614. z[i] = 0;
  615. } // polynomial expansion
  616. for (i = k - 1; i >= 0; --i) {
  617. v = a[i];
  618. c = 1;
  619. z[i] += v;
  620. for (j = 1; j <= i; ++j) {
  621. c *= (i + 1 - j) / j; // binomial coefficent
  622. z[i - j] += v * Math.pow(x, j) * c;
  623. }
  624. } // bias term
  625. z[0] += y;
  626. return z;
  627. } // Given an array for a two-dimensional matrix and the polynomial order,
  628. // solve A * x = b using Gaussian elimination.
  629. function gaussianElimination(matrix) {
  630. var n = matrix.length - 1,
  631. coef = [];
  632. var i, j, k, r, t;
  633. for (i = 0; i < n; ++i) {
  634. r = i; // max row
  635. for (j = i + 1; j < n; ++j) {
  636. if (Math.abs(matrix[i][j]) > Math.abs(matrix[i][r])) {
  637. r = j;
  638. }
  639. }
  640. for (k = i; k < n + 1; ++k) {
  641. t = matrix[k][i];
  642. matrix[k][i] = matrix[k][r];
  643. matrix[k][r] = t;
  644. }
  645. for (j = i + 1; j < n; ++j) {
  646. for (k = n; k >= i; k--) {
  647. matrix[k][j] -= matrix[k][i] * matrix[i][j] / matrix[i][i];
  648. }
  649. }
  650. }
  651. for (j = n - 1; j >= 0; --j) {
  652. t = 0;
  653. for (k = j + 1; k < n; ++k) {
  654. t += matrix[k][j] * coef[k];
  655. }
  656. coef[j] = (matrix[n][j] - t) / matrix[j][j];
  657. }
  658. return coef;
  659. }
  660. function power () {
  661. var x = function x(d) {
  662. return d[0];
  663. },
  664. y = function y(d) {
  665. return d[1];
  666. },
  667. domain;
  668. function power(data) {
  669. var n = 0,
  670. X = 0,
  671. Y = 0,
  672. XY = 0,
  673. X2 = 0,
  674. YS = 0,
  675. xmin = domain ? +domain[0] : Infinity,
  676. xmax = domain ? +domain[1] : -Infinity;
  677. visitPoints(data, x, y, function (dx, dy) {
  678. var lx = Math.log(dx),
  679. ly = Math.log(dy);
  680. ++n;
  681. X += (lx - X) / n;
  682. Y += (ly - Y) / n;
  683. XY += (lx * ly - XY) / n;
  684. X2 += (lx * lx - X2) / n;
  685. YS += (dy - YS) / n;
  686. if (!domain) {
  687. if (dx < xmin) xmin = dx;
  688. if (dx > xmax) xmax = dx;
  689. }
  690. });
  691. var _ols = ols(X, Y, XY, X2),
  692. _ols2 = _slicedToArray(_ols, 2),
  693. a = _ols2[0],
  694. b = _ols2[1];
  695. a = Math.exp(a);
  696. var fn = function fn(x) {
  697. return a * Math.pow(x, b);
  698. },
  699. out = interpose(xmin, xmax, fn);
  700. out.a = a;
  701. out.b = b;
  702. out.predict = fn;
  703. out.rSquared = determination(data, x, y, YS, fn);
  704. return out;
  705. }
  706. power.domain = function (arr) {
  707. return arguments.length ? (domain = arr, power) : domain;
  708. };
  709. power.x = function (fn) {
  710. return arguments.length ? (x = fn, power) : x;
  711. };
  712. power.y = function (fn) {
  713. return arguments.length ? (y = fn, power) : y;
  714. };
  715. return power;
  716. }
  717. export { exponential as regressionExp, linear as regressionLinear, loess as regressionLoess, logarithmic as regressionLog, polynomial as regressionPoly, power as regressionPow, quad as regressionQuad };