From 80e170f96e96b240332d1458252ff6b78ea6475c Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 17:51:34 +0100 Subject: [PATCH 01/14] ARIMA: Cleanup inclu2 --- setup.py | 2 +- src/arima.cpp | 55 ++++++++++++++++++++++++++++++--------------------- 2 files changed, 34 insertions(+), 23 deletions(-) diff --git a/setup.py b/setup.py index 375563724..5091e63e9 100644 --- a/setup.py +++ b/setup.py @@ -58,7 +58,7 @@ name="statsforecast._lib", sources=glob.glob("src/*.cpp"), include_dirs=["include/statsforecast", "external_libs/eigen"], - cxx_std=17, + cxx_std=20, ) ] ParallelCompile("CMAKE_BUILD_PARALLEL_LEVEL", needs_recompile=naive_recompile).install() diff --git a/src/arima.cpp b/src/arima.cpp index 4d3879a49..daba827ef 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -1,6 +1,8 @@ #include #include +#include #include +#include #include #include @@ -326,33 +328,43 @@ arima_like(const py::array_t yv, const py::array_t phiv, return {ssq, sumlog, nu}; } -void inclu2(int np, const double *xnext, double *xrow, double ynext, double *d, - double *rbar, double *thetab) { - std::copy(xnext, xnext + np, xrow); - int ithisr = 0; - for (int i = 0; i < np; ++i) { +void inclu2(const std::vector &xnext, std::vector &xrow, + double ynext, std::span d, std::vector &rbar, + std::vector &thetab) { + const size_t np = xnext.size(); + assert(xrow.size() == np); + assert(thetab.size() == np); + + std::copy(xnext.begin(), xnext.end(), xrow.begin()); + size_t ithisr = 0; + + for (size_t i = 0; i < np; ++i) { if (xrow[i] != 0.0) { - double xi = xrow[i]; - double di = d[i]; - double dpi = di + xi * xi; + const double xi = xrow[i]; + const double di = d[i]; + const double dpi = di + xi * xi; d[i] = dpi; - double cbar, sbar; - if (dpi == 0) { - cbar = std::numeric_limits::infinity(); - sbar = std::numeric_limits::infinity(); - } else { - cbar = di / dpi; - sbar = xi / dpi; - } - for (int k = i + 1; k < np; ++k) { - double xk = xrow[k]; - double rbthis = rbar[ithisr]; + + const auto [cbar, sbar] = [dpi, di, xi]() -> std::pair { + if (dpi == 0) { + return {std::numeric_limits::infinity(), + std::numeric_limits::infinity()}; + } else { + return {di / dpi, xi / dpi}; + } + }(); + + for (size_t k = i + 1; k < np; ++k) { + const double xk = xrow[k]; + const double rbthis = rbar[ithisr]; xrow[k] = xk - xi * rbthis; rbar[ithisr++] = cbar * rbthis + sbar * xk; } - double xk = ynext; + + const double xk = ynext; ynext = xk - xi * thetab[i]; thetab[i] = cbar * thetab[i] + sbar * xk; + if (di == 0.0) { return; } @@ -430,8 +442,7 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, ind2 = 0; } xnext[ind2] += 1.0; - inclu2(np, xnext.data(), xrow.data(), ynext, res, rbar.data(), - thetab.data()); + inclu2(xnext, xrow, ynext, std::span(res, resv.size()), rbar, thetab); xnext[ind2] = 0.0; if (i != r - 1) { xnext[indi++] = 0.0; From c0843bb0f266812e7a5aabf66274fe0d5e7d4ad0 Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 18:24:17 +0100 Subject: [PATCH 02/14] ARIMA: Cleanup arima_gradtrans --- src/arima.cpp | 52 +++++++++++++++++++++++++++++++-------------------- 1 file changed, 32 insertions(+), 20 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index daba827ef..ebc1f0d5f 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -498,49 +498,61 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, py::array_t arima_gradtrans(const py::array_t xv, const py::array_t armav) { + assert(xv.ndim() == 1); + assert(armav.ndim() == 1); + constexpr double eps = 1e-3; - auto x = xv.data(); - auto arma = armav.data(); - int n = static_cast(xv.size()); - int mp = arma[0]; - int mq = arma[1]; - int msp = arma[2]; + const std::span arma(armav.data(), armav.size()); + const std::span x(xv.data(), xv.size()); + + const size_t n = x.size(); + const int mp = arma[0]; + const int mq = arma[1]; + const int msp = arma[2]; + + std::array w1; + std::array w2; + std::array w3; - auto w1 = std::array(); - auto w2 = std::array(); - auto w3 = std::array(); py::array_t outv({n, n}); auto out = outv.mutable_data(); - for (int i = 0; i < n; ++i) { - for (int j = 0; j < n; ++j) { + for (size_t i = 0; i < n; ++i) { + for (size_t j = 0; j < n; ++j) { out[i * n + j] = (i == j) ? 1.0 : 0.0; } } + if (mp > 0) { - std::copy(x, x + mp, w1.begin()); - partrans(mp, w1.data(), w2.data()); - for (int i = 0; i < mp; ++i) { + std::copy(x.begin(), x.begin() + mp, w1.begin()); + partrans(mp, w1, w2); + + for (size_t i = 0; i < mp; ++i) { w1[i] += eps; - partrans(mp, w1.data(), w3.data()); - for (int j = 0; j < mp; ++j) { + partrans(mp, w1, w3); + + for (size_t j = 0; j < mp; ++j) { out[n * i + j] = (w3[j] - w2[j]) / eps; } w1[i] -= eps; } } + if (msp > 0) { - int v = mp + mq; - std::copy(x + v, x + v + msp, w1.begin()); + const size_t v = mp + mq; + std::copy(x.begin() + v, x.begin() + v + msp, w1.begin()); partrans(msp, w1.data(), w2.data()); - for (int i = 0; i < msp; ++i) { + + for (size_t i = 0; i < msp; ++i) { w1[i] += eps; partrans(msp, w1.data(), w3.data()); - for (int j = 0; j < msp; ++j) { + + for (size_t j = 0; j < msp; ++j) { out[n * (i + v) + j + v] = (w3[j] - w2[j]) / eps; } w1[i] -= eps; } } + return outv; } From bd7923c78d16b66665507bdba5af75e4558819ce Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 19:24:49 +0100 Subject: [PATCH 03/14] ARIMA: Cleanup invpartrans --- src/arima.cpp | 54 ++++++++++++++++++++++++++++++++++++--------------- 1 file changed, 38 insertions(+), 16 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index ebc1f0d5f..d693d5140 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -10,6 +10,18 @@ namespace arima { namespace py = pybind11; +namespace { +template std::span make_span(py::array_t &array) { + // ssize_t to size_t is safe because the size of an array is non-negative + return {array.mutable_data(), static_cast(array.size())}; +} + +template +std::span make_cspan(const py::array_t &array) { + // ssize_t to size_t is safe because the size of an array is non-negative + return {array.data(), static_cast(array.size())}; +} +} // namespace void partrans(int p, const double *raw, double *newv) { std::transform(raw, raw + p, newv, [](double x) { return std::tanh(x); }); @@ -502,10 +514,10 @@ py::array_t arima_gradtrans(const py::array_t xv, assert(armav.ndim() == 1); constexpr double eps = 1e-3; - const std::span arma(armav.data(), armav.size()); - const std::span x(xv.data(), xv.size()); - + const auto arma = make_cspan(armav); + const auto x = make_cspan(xv); const size_t n = x.size(); + const int mp = arma[0]; const int mq = arma[1]; const int msp = arma[2]; @@ -513,9 +525,10 @@ py::array_t arima_gradtrans(const py::array_t xv, std::array w1; std::array w2; std::array w3; + assert(mp < 100); py::array_t outv({n, n}); - auto out = outv.mutable_data(); + const auto out = make_span(outv); for (size_t i = 0; i < n; ++i) { for (size_t j = 0; j < n; ++j) { out[i * n + j] = (i == j) ? 1.0 : 0.0; @@ -524,11 +537,11 @@ py::array_t arima_gradtrans(const py::array_t xv, if (mp > 0) { std::copy(x.begin(), x.begin() + mp, w1.begin()); - partrans(mp, w1, w2); + partrans(mp, w1.data(), w2.data()); for (size_t i = 0; i < mp; ++i) { w1[i] += eps; - partrans(mp, w1, w3); + partrans(mp, w1.data(), w3.data()); for (size_t j = 0; j < mp; ++j) { out[n * i + j] = (w3[j] - w2[j]) / eps; @@ -576,20 +589,29 @@ py::array_t arima_undopars(const py::array_t xv, return outv; } -void invpartrans(int p, const py::array_t phiv, +void invpartrans(const int p, const py::array_t phiv, py::array_t outv) { - auto phi = phiv.data(); - auto out = outv.mutable_data(); - std::copy(phi, phi + p, out); - std::vector work(phi, phi + p); - for (int j = p - 1; j > 0; --j) { - double a = out[j]; - for (int k = 0; k < j; ++k) { + assert(p >= 0); + assert(phiv.ndim() == 1); + assert(outv.ndim() == 1); + assert(p <= phiv.size()); + assert(p <= outv.size()); + + const auto phi = make_cspan(phiv); + const auto out = make_span(outv); + + std::copy(phi.begin(), phi.begin() + p, out.begin()); + + std::vector work(phi.begin(), phi.begin() + p); + for (size_t j = p - 1; j > 0; --j) { + const double a = out[j]; + for (size_t k = 0; k < j; ++k) { work[k] = (out[k] + a * out[j - k - 1]) / (1 - a * a); } - std::copy(work.begin(), work.begin() + j, out); + std::copy(work.begin(), work.begin() + j, out.begin()); } - for (int j = 0; j < p; ++j) { + + for (size_t j = 0; j < p; ++j) { out[j] = std::atanh(out[j]); } } From bd59a29cc9abf99596ab6224a425bad94b03e876 Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 19:43:44 +0100 Subject: [PATCH 04/14] ARIMA: Cleanup partrans --- src/arima.cpp | 123 ++++++++++++++++++++++++++++++++------------------ 1 file changed, 79 insertions(+), 44 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index d693d5140..fe324989c 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -23,64 +23,86 @@ std::span make_cspan(const py::array_t &array) { } } // namespace -void partrans(int p, const double *raw, double *newv) { - std::transform(raw, raw + p, newv, [](double x) { return std::tanh(x); }); - std::vector work(newv, newv + p); - for (int j = 1; j < p; ++j) { - for (int k = 0; k < j; ++k) { +void partrans(const size_t p, const std::span rawv, + const std::span newv) { + assert(p <= rawv.size()); + assert(p <= newv.size()); + + std::transform(rawv.begin(), rawv.begin() + p, newv.begin(), + [](double x) { return std::tanh(x); }); + + std::vector work(newv.begin(), newv.begin() + p); + for (size_t j = 1; j < p; ++j) { + for (size_t k = 0; k < j; ++k) { work[k] -= newv[j] * newv[j - k - 1]; } - std::copy(work.begin(), work.begin() + j, newv); + std::copy(work.begin(), work.begin() + j, newv.begin()); } } std::tuple, py::array_t> arima_transpar(const py::array_t params_inv, const py::array_t armav, bool trans) { - auto arma = armav.data(); - auto params_in = params_inv.data(); - int mp = arma[0]; - int mq = arma[1]; - int msp = arma[2]; - int msq = arma[3]; - int ns = arma[4]; - int p = mp + ns * msp; - int q = mq + ns * msq; - auto params = std::vector(params_in, params_in + params_inv.size()); + assert(params_inv.ndim() == 1); + assert(armav.ndim() == 1); + + const auto arma = make_cspan(armav); + const auto params_in = make_cspan(params_inv); + + const int mp = arma[0]; + const int mq = arma[1]; + const int msp = arma[2]; + const int msq = arma[3]; + const int ns = arma[4]; + assert(mp >= 0); + assert(mq >= 0); + assert(msp >= 0); + assert(msq >= 0); + assert(ns >= 0); + const int p = mp + ns * msp; + const int q = mq + ns * msq; + + std::vector params(params_in.begin(), params_in.end()); + py::array_t phiv(p); py::array_t thetav(q); - auto phi = phiv.mutable_data(); - auto theta = thetav.mutable_data(); + const auto phi = make_span(phiv); + const auto theta = make_span(thetav); + if (trans) { if (mp > 0) { - partrans(mp, params_in, params.data()); + partrans(mp, params_in, params); } int v = mp + mq; if (msp > 0) { - partrans(msp, params_in + v, params.data() + v); + partrans(msp, params_in.subspan(v), std::span(params).subspan(v)); } } + if (ns > 0) { - std::copy(params.begin(), params.begin() + mp, phi); - std::fill(phi + mp, phi + p, 0.0); - std::copy(params.begin() + mp, params.begin() + mp + mq, theta); - std::fill(theta + mq, theta + q, 0.0); - for (int j = 0; j < msp; ++j) { + std::copy(params.begin(), params.begin() + mp, phi.begin()); + std::fill(phi.begin() + mp, phi.begin() + p, 0.0); + std::copy(params.begin() + mp, params.begin() + mp + mq, theta.begin()); + std::fill(theta.begin() + mq, theta.begin() + q, 0.0); + + for (size_t j = 0; j < msp; ++j) { phi[(j + 1) * ns - 1] += params[j + mp + mq]; - for (int i = 0; i < mp; ++i) { + for (size_t i = 0; i < mp; ++i) { phi[(j + 1) * ns + i] -= params[i] * params[j + mp + mq]; } } - for (int j = 0; j < msq; ++j) { + + for (size_t j = 0; j < msq; ++j) { theta[(j + 1) * ns - 1] += params[j + mp + mq + msp]; - for (int i = 0; i < mq; ++i) { + for (size_t i = 0; i < mq; ++i) { theta[(j + 1) * ns + i] += params[i + mp] * params[j + mp + mq + msp]; } } } else { - std::copy(params.begin(), params.begin() + mp, phi); - std::copy(params.begin() + mp, params.begin() + mp + mq, theta); + std::copy(params.begin(), params.begin() + mp, phi.begin()); + std::copy(params.begin() + mp, params.begin() + mp + mq, theta.begin()); } + return {phiv, thetav}; } @@ -537,11 +559,11 @@ py::array_t arima_gradtrans(const py::array_t xv, if (mp > 0) { std::copy(x.begin(), x.begin() + mp, w1.begin()); - partrans(mp, w1.data(), w2.data()); + partrans(mp, w1, w2); for (size_t i = 0; i < mp; ++i) { w1[i] += eps; - partrans(mp, w1.data(), w3.data()); + partrans(mp, w1, w3); for (size_t j = 0; j < mp; ++j) { out[n * i + j] = (w3[j] - w2[j]) / eps; @@ -553,11 +575,11 @@ py::array_t arima_gradtrans(const py::array_t xv, if (msp > 0) { const size_t v = mp + mq; std::copy(x.begin() + v, x.begin() + v + msp, w1.begin()); - partrans(msp, w1.data(), w2.data()); + partrans(msp, w1, w2); for (size_t i = 0; i < msp; ++i) { w1[i] += eps; - partrans(msp, w1.data(), w3.data()); + partrans(msp, w1, w3); for (size_t j = 0; j < msp; ++j) { out[n * (i + v) + j + v] = (w3[j] - w2[j]) / eps; @@ -571,21 +593,34 @@ py::array_t arima_gradtrans(const py::array_t xv, py::array_t arima_undopars(const py::array_t xv, const py::array_t armav) { - auto x = xv.data(); - auto arma = armav.data(); - int mp = arma[0]; - int mq = arma[1]; - int msp = arma[2]; - py::array_t outv{xv.size()}; - auto out = outv.mutable_data(); - std::copy(xv.data(), xv.data() + xv.size(), out); + assert(xv.ndim() == 1); + assert(armav.ndim() == 1); + + const auto x = make_cspan(xv); + const auto arma = make_cspan(armav); + + const int mp = arma[0]; + const int mq = arma[1]; + const int msp = arma[2]; + assert(mp >= 0); + assert(mq >= 0); + assert(msp >= 0); + + py::array_t outv(xv.size()); + const auto out = make_span(outv); + + std::copy(x.begin(), x.end(), out.begin()); + if (mp > 0) { partrans(mp, x, out); } - int v = mp + mq; + + const size_t v = mp + mq; + if (msp > 0) { - partrans(msp, x + v, out + v); + partrans(msp, x.subspan(v), out.subspan(v)); } + return outv; } From d1aef03dade5130288d72a9ee9e4ce23bb19018f Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 19:59:19 +0100 Subject: [PATCH 05/14] ARIMA: convert ints to uints --- src/arima.cpp | 47 +++++++++++++++++++---------------------------- 1 file changed, 19 insertions(+), 28 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index fe324989c..4b410646d 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -23,7 +23,7 @@ std::span make_cspan(const py::array_t &array) { } } // namespace -void partrans(const size_t p, const std::span rawv, +void partrans(const uint32_t p, const std::span rawv, const std::span newv) { assert(p <= rawv.size()); assert(p <= newv.size()); @@ -42,25 +42,20 @@ void partrans(const size_t p, const std::span rawv, std::tuple, py::array_t> arima_transpar(const py::array_t params_inv, - const py::array_t armav, bool trans) { + const py::array_t armav, bool trans) { assert(params_inv.ndim() == 1); assert(armav.ndim() == 1); const auto arma = make_cspan(armav); const auto params_in = make_cspan(params_inv); - const int mp = arma[0]; - const int mq = arma[1]; - const int msp = arma[2]; - const int msq = arma[3]; - const int ns = arma[4]; - assert(mp >= 0); - assert(mq >= 0); - assert(msp >= 0); - assert(msq >= 0); - assert(ns >= 0); - const int p = mp + ns * msp; - const int q = mq + ns * msq; + const uint32_t mp = arma[0]; + const uint32_t mq = arma[1]; + const uint32_t msp = arma[2]; + const uint32_t msq = arma[3]; + const uint32_t ns = arma[4]; + const uint32_t p = mp + ns * msp; + const uint32_t q = mq + ns * msq; std::vector params(params_in.begin(), params_in.end()); @@ -73,7 +68,7 @@ arima_transpar(const py::array_t params_inv, if (mp > 0) { partrans(mp, params_in, params); } - int v = mp + mq; + const uint32_t v = mp + mq; if (msp > 0) { partrans(msp, params_in.subspan(v), std::span(params).subspan(v)); } @@ -531,7 +526,7 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, } py::array_t arima_gradtrans(const py::array_t xv, - const py::array_t armav) { + const py::array_t armav) { assert(xv.ndim() == 1); assert(armav.ndim() == 1); @@ -540,9 +535,9 @@ py::array_t arima_gradtrans(const py::array_t xv, const auto x = make_cspan(xv); const size_t n = x.size(); - const int mp = arma[0]; - const int mq = arma[1]; - const int msp = arma[2]; + const uint32_t mp = arma[0]; + const uint32_t mq = arma[1]; + const uint32_t msp = arma[2]; std::array w1; std::array w2; @@ -592,19 +587,16 @@ py::array_t arima_gradtrans(const py::array_t xv, } py::array_t arima_undopars(const py::array_t xv, - const py::array_t armav) { + const py::array_t armav) { assert(xv.ndim() == 1); assert(armav.ndim() == 1); const auto x = make_cspan(xv); const auto arma = make_cspan(armav); - const int mp = arma[0]; - const int mq = arma[1]; - const int msp = arma[2]; - assert(mp >= 0); - assert(mq >= 0); - assert(msp >= 0); + const uint32_t mp = arma[0]; + const uint32_t mq = arma[1]; + const uint32_t msp = arma[2]; py::array_t outv(xv.size()); const auto out = make_span(outv); @@ -624,9 +616,8 @@ py::array_t arima_undopars(const py::array_t xv, return outv; } -void invpartrans(const int p, const py::array_t phiv, +void invpartrans(const uint32_t p, const py::array_t phiv, py::array_t outv) { - assert(p >= 0); assert(phiv.ndim() == 1); assert(outv.ndim() == 1); assert(p <= phiv.size()); From 56c0e92324ec8787333f7f8e3efd84e862a581ff Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 20:27:31 +0100 Subject: [PATCH 06/14] ARIMA: Cleanup arima_css --- src/arima.cpp | 57 ++++++++++++++++++++++++++++++++------------------- 1 file changed, 36 insertions(+), 21 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index 4b410646d..d116aec46 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -102,50 +102,65 @@ arima_transpar(const py::array_t params_inv, } std::tuple> -arima_css(const py::array_t yv, const py::array_t armav, +arima_css(const py::array_t yv, const py::array_t armav, const py::array_t phiv, const py::array_t thetav) { - int n = static_cast(yv.size()); - int p = static_cast(phiv.size()); - int q = static_cast(thetav.size()); - auto y = yv.data(); - auto arma = armav.data(); - auto phi = phiv.data(); - auto theta = thetav.data(); - int ncond = arma[0] + arma[5] + arma[4] * (arma[2] + arma[6]); - int nu = 0; + assert(yv.ndim() == 1); + assert(armav.ndim() == 1); + assert(phiv.ndim() == 1); + assert(thetav.ndim() == 1); + + // ssize_t to size_t is safe because the size of an array is non-negative + const size_t n = static_cast(yv.size()); + const size_t p = static_cast(phiv.size()); + const size_t q = static_cast(thetav.size()); + + const auto y = make_cspan(yv); + const auto arma = make_cspan(armav); + const auto phi = make_cspan(phiv); + const auto theta = make_cspan(thetav); + + const uint32_t ncond = arma[0] + arma[5] + arma[4] * (arma[2] + arma[6]); + uint32_t nu = 0; double ssq = 0.0; - auto residv = py::array_t(n); - auto resid = residv.mutable_data(); - auto w = std::vector(y, y + yv.size()); - for (int _ = 0; _ < arma[5]; ++_) { - for (int l = n - 1; l > 0; --l) { + py::array_t residv(n); + const auto resid = make_span(residv); + std::vector w(y.begin(), y.end()); + + for (size_t _ = 0; _ < arma[5]; ++_) { + for (size_t l = n - 1; l > 0; --l) { w[l] -= w[l - 1]; } } - int ns = arma[4]; - for (int _ = 0; _ < arma[6]; ++_) { - for (int l = n - 1; l >= ns; --l) { + + const uint32_t ns = arma[4]; + for (size_t _ = 0; _ < arma[6]; ++_) { + for (size_t l = n - 1; l >= ns; --l) { w[l] -= w[l - ns]; } } - for (int l = ncond; l < n; ++l) { + + for (size_t l = ncond; l < n; ++l) { double tmp = w[l]; - for (int j = 0; j < p; ++j) { + for (size_t j = 0; j < p; ++j) { tmp -= phi[j] * w[l - j - 1]; } - for (int j = 0; j < std::min(l - ncond, q); ++j) { + + assert(l >= ncond); + for (size_t j = 0; j < std::min(static_cast(l - ncond), q); ++j) { if (l - j - 1 < 0) { continue; } tmp -= theta[j] * resid[l - j - 1]; } + resid[l] = tmp; if (!std::isnan(tmp)) { nu++; ssq += tmp * tmp; } } + return {ssq / nu, residv}; } From 7a035b8111ec448e85e7937988d5343b81bc4470 Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 20:44:35 +0100 Subject: [PATCH 07/14] ARIMA: Cleanup arima_like --- src/arima.cpp | 141 ++++++++++++++++++++++++++++++-------------------- 1 file changed, 84 insertions(+), 57 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index d116aec46..ce6dd4072 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -168,25 +168,34 @@ std::tuple arima_like(const py::array_t yv, const py::array_t phiv, const py::array_t thetav, const py::array_t deltav, py::array_t av, py::array_t Pv, - py::array_t Pnewv, int up, bool use_resid, + py::array_t Pnewv, uint32_t up, bool use_resid, py::array_t rsResidv) { - int n = static_cast(yv.size()); - int d = static_cast(deltav.size()); - int rd = static_cast(av.size()); - int p = static_cast(phiv.size()); - int q = static_cast(thetav.size()); - auto y = yv.data(); - auto phi = phiv.data(); - auto theta = thetav.data(); - auto delta = deltav.data(); - auto a = av.mutable_data(); - auto P = Pv.mutable_data(); - auto Pnew = Pnewv.mutable_data(); - auto rsResid = rsResidv.mutable_data(); + // ssize_t to size_t is safe because the size of an array is non-negative + const size_t n = static_cast(yv.size()); + const size_t d = static_cast(deltav.size()); + const size_t rd = static_cast(av.size()); + const size_t p = static_cast(phiv.size()); + const size_t q = static_cast(thetav.size()); + + const auto y = make_cspan(yv); + const auto phi = make_cspan(phiv); + const auto theta = make_cspan(thetav); + const auto delta = make_cspan(deltav); + const auto a = make_span(av); + const auto P = make_span(Pv); + const auto Pnew = make_span(Pnewv); + const auto rsResid = make_span(rsResidv); + + // asserts for copies at the end + assert(Pnew.size() >= rd * rd); + assert(P.size() >= Pnew.size()); + double ssq = 0.0; double sumlog = 0.0; - int nu = 0; - int r = rd - d; + + assert(rd >= d); + const size_t r = rd - d; + uint32_t nu = 0; std::vector anew(rd); std::vector M(rd); @@ -194,9 +203,10 @@ arima_like(const py::array_t yv, const py::array_t phiv, if (d > 0) { mm.resize(rd * rd); } + double tmp; - for (int l = 0; l < n; ++l) { - for (int i = 0; i < r; ++i) { + for (size_t l = 0; l < n; ++l) { + for (size_t i = 0; i < r; ++i) { if (i < r - 1) { tmp = a[i + 1]; } else { @@ -207,26 +217,32 @@ arima_like(const py::array_t yv, const py::array_t phiv, } anew[i] = tmp; } + if (d > 0) { - for (int i = r + 1; i < rd; ++i) { + for (size_t i = r + 1; i < rd; ++i) { anew[i] = a[i - 1]; } tmp = a[0]; - for (int i = 0; i < d; ++i) { + for (size_t i = 0; i < d; ++i) { tmp += delta[i] * a[r + i]; } anew[r] = tmp; } + if (l > up) { if (d == 0) { - for (int i = 0; i < r; ++i) { - double vi = 0.0; - if (i == 0) { - vi = 1.0; - } else if (i - 1 < q) { - vi = theta[i - 1]; - } - for (int j = 0; j < r; ++j) { + for (size_t i = 0; i < r; ++i) { + const double vi = [&]() { + if (i == 0) { + return 1.0; + } else if (i - 1 < q) { + return theta[i - 1]; + } else { + return 0.0; + } + }(); + + for (size_t j = 0; j < r; ++j) { tmp = 0.0; if (j == 0) { tmp = vi; @@ -248,9 +264,10 @@ arima_like(const py::array_t yv, const py::array_t phiv, Pnew[i + r * j] = tmp; } } + } else { - for (int i = 0; i < r; ++i) { - for (int j = 0; j < rd; ++j) { + for (size_t i = 0; i < r; ++i) { + for (size_t j = 0; j < rd; ++j) { tmp = 0.0; if (i < p) { tmp += phi[i] * P[rd * j]; @@ -261,20 +278,23 @@ arima_like(const py::array_t yv, const py::array_t phiv, mm[i + rd * j] = tmp; } } - for (int j = 0; j < rd; ++j) { + + for (size_t j = 0; j < rd; ++j) { tmp = P[rd * j]; - for (int k = 0; k < d; ++k) { + for (size_t k = 0; k < d; ++k) { tmp += delta[k] * P[r + k + rd * j]; } mm[r + rd * j] = tmp; } - for (int i = 1; i < d; ++i) { - for (int j = 0; j < rd; ++j) { + + for (size_t i = 1; i < d; ++i) { + for (size_t j = 0; j < rd; ++j) { mm[r + i + rd * j] = P[r + i - 1 + rd * j]; } } - for (int i = 0; i < r; ++i) { - for (int j = 0; j < rd; ++j) { + + for (size_t i = 0; i < r; ++i) { + for (size_t j = 0; j < rd; ++j) { tmp = 0.0; if (i < p) { tmp += phi[i] * mm[j]; @@ -285,26 +305,24 @@ arima_like(const py::array_t yv, const py::array_t phiv, Pnew[j + rd * i] = tmp; } } - for (int j = 0; j < rd; ++j) { + + for (size_t j = 0; j < rd; ++j) { tmp = mm[j]; - for (int k = 0; k < d; ++k) { + for (size_t k = 0; k < d; ++k) { tmp += delta[k] * mm[rd * (r + k) + j]; } Pnew[rd * r + j] = tmp; } - for (int i = 1; i < d; ++i) { - for (int j = 0; j < rd; ++j) { + + for (size_t i = 1; i < d; ++i) { + for (size_t j = 0; j < rd; ++j) { Pnew[rd * (r + i) + j] = mm[rd * (r + i - 1) + j]; } } - for (int i = 0; i < q + 1; ++i) { - double vi; - if (i == 0) { - vi = 1.0; - } else { - vi = theta[i - 1]; - } - for (int j = 0; j < q + 1; ++j) { + + for (size_t i = 0; i < q + 1; ++i) { + const double vi = i == 0 ? 1.0 : theta[i - 1]; + for (size_t j = 0; j < q + 1; ++j) { if (j == 0) { Pnew[i + rd * j] += vi; } else { @@ -314,22 +332,26 @@ arima_like(const py::array_t yv, const py::array_t phiv, } } } + if (!std::isnan(y[l])) { double resid = y[l] - anew[0]; - for (int i = 0; i < d; ++i) { + for (size_t i = 0; i < d; ++i) { resid -= delta[i] * anew[r + i]; } - for (int i = 0; i < rd; ++i) { + + for (size_t i = 0; i < rd; ++i) { tmp = Pnew[i]; - for (int j = 0; j < d; ++j) { + for (size_t j = 0; j < d; ++j) { tmp += Pnew[i + (r + j) * rd] * delta[j]; } M[i] = tmp; } + double gain = M[0]; for (int j = 0; j < d; ++j) { gain += delta[j] * M[r + j]; } + if (gain < 1e4) { nu++; if (gain == 0) { @@ -339,6 +361,7 @@ arima_like(const py::array_t yv, const py::array_t phiv, } sumlog += std::log(gain); } + if (use_resid) { if (gain == 0) { rsResid[l] = std::numeric_limits::infinity(); @@ -346,29 +369,33 @@ arima_like(const py::array_t yv, const py::array_t phiv, rsResid[l] = resid / std::sqrt(gain); } } + if (gain == 0) { - for (int i = 0; i < rd; ++i) { + for (size_t i = 0; i < rd; ++i) { a[i] = std::numeric_limits::infinity(); - for (int j = 0; j < rd; ++j) { + for (size_t j = 0; j < rd; ++j) { Pnew[i + j * rd] = std::numeric_limits::infinity(); } } + } else { - for (int i = 0; i < rd; ++i) { + for (size_t i = 0; i < rd; ++i) { a[i] = anew[i] + M[i] * resid / gain; - for (int j = 0; j < rd; ++j) { + for (size_t j = 0; j < rd; ++j) { P[i + j * rd] = Pnew[i + j * rd] - M[i] * M[j] / gain; } } } + } else { - std::copy(anew.begin(), anew.end(), a); - std::copy(Pnew, Pnew + rd * rd, P); + std::copy(anew.begin(), anew.end(), a.begin()); + std::copy(Pnew.begin(), Pnew.begin() + rd * rd, P.begin()); if (use_resid) { rsResid[l] = std::numeric_limits::quiet_NaN(); } } } + return {ssq, sumlog, nu}; } From ef1e02e30887b8dbdf1e7fe80ee22d3d8fb3fd5d Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Thu, 14 Nov 2024 23:14:01 +0100 Subject: [PATCH 08/14] ARIMA: Cleanup getQ0 --- src/arima.cpp | 145 +++++++++++++++++++++++++++++++------------------- 1 file changed, 89 insertions(+), 56 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index ce6dd4072..19e4f349c 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -447,34 +447,45 @@ void inclu2(const std::vector &xnext, std::vector &xrow, void getQ0(const py::array_t phiv, const py::array_t thetav, py::array_t resv) { - auto phi = phiv.data(); - auto theta = thetav.data(); - auto res = resv.mutable_data(); - int p = static_cast(phiv.size()); - int q = static_cast(thetav.size()); - int r = std::max(p, q + 1); - int np = r * (r + 1) / 2; - int nrbar = np * (np - 1) / 2; - int ind = 0; + assert(thetav.ndim() == 1); + assert(phiv.ndim() == 1); + assert(resv.ndim() == 1); + + const auto phi = make_cspan(phiv); + const auto theta = make_cspan(thetav); + const auto res = make_span(resv); + + const size_t p = phi.size(); + const size_t q = theta.size(); + const size_t r = std::max(p, q + 1); + const size_t np = r * (r + 1) / 2; + const size_t nrbar = np * (np - 1) / 2; std::vector V(np); - for (int j = 0; j < r; ++j) { - double vj = 0.0; - if (j == 0) { - vj = 1.0; - } else if (j - 1 < q) { - vj = theta[j - 1]; - } - for (int i = j; i < r; ++i) { - double vi = 0.0; - if (i == 0) { - vi = 1.0; - } else if (i - 1 < q) { - vi = theta[i - 1]; + + { + size_t ind = 0; + + for (size_t j = 0; j < r; ++j) { + double vj = 0.0; + if (j == 0) { + vj = 1.0; + } else if (j - 1 < q) { + vj = theta[j - 1]; + } + + for (size_t i = j; i < r; ++i) { + double vi = 0.0; + if (i == 0) { + vi = 1.0; + } else if (i - 1 < q) { + vi = theta[i - 1]; + } + V[ind++] = vi * vj; } - V[ind++] = vi * vj; } } + if (r == 1) { if (p == 0) { res[0] = 1.0; @@ -483,24 +494,31 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, } return; } + if (p > 0) { std::vector rbar(nrbar); std::vector thetab(np); std::vector xnext(np); std::vector xrow(np); - ind = 0; - int ind1 = -1; - int npr = np - r; - int npr1 = npr + 1; - int indj = npr; - int ind2 = npr - 1; - for (int j = 0; j < r; ++j) { - double phij = j < p ? phi[j] : 0.0; + + size_t ind = 0; + ssize_t ind1 = -1; + + const size_t npr = np - r; + const size_t npr1 = npr + 1; + size_t indj = npr; + size_t ind2 = npr - 1; + + for (size_t j = 0; j < r; ++j) { + const double phij = j < p ? phi[j] : 0.0; + size_t indi = npr1 + j; + xnext[indj++] = 0.0; - int indi = npr1 + j; - for (int i = j; i < r; ++i) { - double ynext = V[ind++]; - double phii = i < p ? phi[i] : 0.0; + + for (size_t i = j; i < r; ++i) { + const double ynext = V[ind++]; + const double phii = i < p ? phi[i] : 0.0; + if (j != r - 1) { xnext[indj] = -phii; if (i != r - 1) { @@ -508,44 +526,54 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, xnext[++ind1] = -1.0; } } + xnext[npr] = -phii * phij; + if (++ind2 >= np) { ind2 = 0; } + xnext[ind2] += 1.0; - inclu2(xnext, xrow, ynext, std::span(res, resv.size()), rbar, thetab); + inclu2(xnext, xrow, ynext, res, rbar, thetab); xnext[ind2] = 0.0; + if (i != r - 1) { xnext[indi++] = 0.0; xnext[ind1] = 0.0; } } } - int ithisr = nrbar - 1; - int im = np - 1; - for (int i = 0; i < np; ++i) { + + size_t ithisr = nrbar - 1; + size_t im = np - 1; + for (size_t i = 0; i < np; ++i) { + size_t jm = np - 1; double bi = thetab[im]; - int jm = np - 1; - for (int j = 0; j < i; ++j) { + + for (size_t j = 0; j < i; ++j) { bi -= rbar[ithisr--] * res[jm--]; } + res[im--] = bi; } + ind = npr; - for (int i = 0; i < r; ++i) { + + for (size_t i = 0; i < r; ++i) { xnext[i] = res[ind++]; } + ind = np - 1; ind1 = npr - 1; - for (int i = 0; i < npr; ++i) { + for (size_t i = 0; i < npr; ++i) { res[ind--] = res[ind1--]; } - std::copy(xnext.begin(), xnext.begin() + r, res); + std::copy(xnext.begin(), xnext.begin() + r, res.begin()); } else { - int indn = np; - ind = np; - for (int i = 0; i < r; ++i) { - for (int j = 0; j < i + 1; ++j) { + size_t ind = np; + size_t indn = np; + for (size_t i = 0; i < r; ++i) { + for (size_t j = 0; j < i + 1; ++j) { --ind; res[ind] = V[ind]; if (j != 0) { @@ -554,15 +582,20 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, } } } - ind = np; - for (int i = r - 1; i > 0; --i) { - for (int j = r - 1; j > i - 1; --j) { - res[r * i + j] = res[--ind]; + + { + size_t ind = np; + + for (size_t i = r - 1; i > 0; --i) { + for (size_t j = r - 1; j > i - 1; --j) { + res[r * i + j] = res[--ind]; + } } - } - for (int i = 0; i < r - 1; ++i) { - for (int j = i + 1; j < r; ++j) { - res[i + r * j] = res[j + r * i]; + + for (size_t i = 0; i < r - 1; ++i) { + for (size_t j = i + 1; j < r; ++j) { + res[i + r * j] = res[j + r * i]; + } } } } From 6c1f1acb26595d61ad04e14609a25fc073681636 Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Sat, 16 Nov 2024 22:35:39 +0100 Subject: [PATCH 09/14] ARIMA: Use py::ssize_t instead of posix size_t --- src/arima.cpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/arima.cpp b/src/arima.cpp index 19e4f349c..3f471f6f8 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -12,13 +12,13 @@ namespace arima { namespace py = pybind11; namespace { template std::span make_span(py::array_t &array) { - // ssize_t to size_t is safe because the size of an array is non-negative + // py::ssize_t to size_t is safe because the size of an array is non-negative return {array.mutable_data(), static_cast(array.size())}; } template std::span make_cspan(const py::array_t &array) { - // ssize_t to size_t is safe because the size of an array is non-negative + // py::ssize_t to size_t is safe because the size of an array is non-negative return {array.data(), static_cast(array.size())}; } } // namespace @@ -109,7 +109,7 @@ arima_css(const py::array_t yv, const py::array_t armav, assert(phiv.ndim() == 1); assert(thetav.ndim() == 1); - // ssize_t to size_t is safe because the size of an array is non-negative + // py::ssize_t to size_t is safe because the size of an array is non-negative const size_t n = static_cast(yv.size()); const size_t p = static_cast(phiv.size()); const size_t q = static_cast(thetav.size()); @@ -170,7 +170,7 @@ arima_like(const py::array_t yv, const py::array_t phiv, py::array_t av, py::array_t Pv, py::array_t Pnewv, uint32_t up, bool use_resid, py::array_t rsResidv) { - // ssize_t to size_t is safe because the size of an array is non-negative + // py::ssize_t to size_t is safe because the size of an array is non-negative const size_t n = static_cast(yv.size()); const size_t d = static_cast(deltav.size()); const size_t rd = static_cast(av.size()); @@ -502,7 +502,7 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, std::vector xrow(np); size_t ind = 0; - ssize_t ind1 = -1; + py::ssize_t ind1 = -1; const size_t npr = np - r; const size_t npr1 = npr + 1; From 5d61950d12eba0fba9c13742790cb38ab283ebcb Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Sat, 16 Nov 2024 22:36:51 +0100 Subject: [PATCH 10/14] ARIMA: Change arma array dtype from intc to uint32 --- python/statsforecast/arima.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/python/statsforecast/arima.py b/python/statsforecast/arima.py index 2e2b67e6a..ffcf35b26 100644 --- a/python/statsforecast/arima.py +++ b/python/statsforecast/arima.py @@ -313,7 +313,7 @@ def maInvert(ma): order[1], seasonal["order"][1], ], - dtype=np.intc, + dtype=np.uint32, ) narma = arma[:4].sum().item() From 63260204f19f2c9306e26ea2feb99bb8ed37e225 Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Sun, 17 Nov 2024 16:19:19 +0100 Subject: [PATCH 11/14] ARIMA: arima_like fix int iter --- src/arima.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/arima.cpp b/src/arima.cpp index 3f471f6f8..18e532676 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -348,7 +348,7 @@ arima_like(const py::array_t yv, const py::array_t phiv, } double gain = M[0]; - for (int j = 0; j < d; ++j) { + for (size_t j = 0; j < d; ++j) { gain += delta[j] * M[r + j]; } From b5bd7e3933a13ffd485425cac1d446a107d5a7e9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jos=C3=A9=20Morales?= Date: Tue, 19 Nov 2024 13:19:07 -0600 Subject: [PATCH 12/14] update nb --- nbs/src/arima.ipynb | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/nbs/src/arima.ipynb b/nbs/src/arima.ipynb index 658fc71fa..ae423ecc7 100644 --- a/nbs/src/arima.ipynb +++ b/nbs/src/arima.ipynb @@ -658,7 +658,7 @@ " order[1],\n", " seasonal['order'][1],\n", " ],\n", - " dtype=np.intc,\n", + " dtype=np.uint32,\n", " )\n", " narma = arma[:4].sum().item()\n", " \n", From bdef8f51c3065be5269fadafe8617df0ed33cc36 Mon Sep 17 00:00:00 2001 From: Filip Cacky Date: Wed, 20 Nov 2024 15:57:19 +0100 Subject: [PATCH 13/14] ARIMA: Change armav to int32 and add gt>0 asserts --- nbs/src/arima.ipynb | 2 +- python/statsforecast/adapters/prophet.py | 2 + python/statsforecast/arima.py | 2 +- python/statsforecast/core.py | 2 + .../statsforecast/distributed/multiprocess.py | 2 + python/statsforecast/feature_engineering.py | 2 + python/statsforecast/models.py | 2 + python/statsforecast/utils.py | 2 + src/arima.cpp | 46 +++++++++++++------ 9 files changed, 45 insertions(+), 17 deletions(-) diff --git a/nbs/src/arima.ipynb b/nbs/src/arima.ipynb index ae423ecc7..db37c8ea6 100644 --- a/nbs/src/arima.ipynb +++ b/nbs/src/arima.ipynb @@ -658,7 +658,7 @@ " order[1],\n", " seasonal['order'][1],\n", " ],\n", - " dtype=np.uint32,\n", + " dtype=np.int32,\n", " )\n", " narma = arma[:4].sum().item()\n", " \n", diff --git a/python/statsforecast/adapters/prophet.py b/python/statsforecast/adapters/prophet.py index ce3d3a667..01d307de5 100644 --- a/python/statsforecast/adapters/prophet.py +++ b/python/statsforecast/adapters/prophet.py @@ -1,3 +1,5 @@ +"""In 2017, Facebook open-sourced [Prophet](https://peerj.com/preprints/3190.pdf), with the promise of providing experts and non-experts the possibility of producing high-quality predictions. The forecasting community heavily adopted the solution, reaching millions of accumulated downloads. It became evident that its [quality is shadowed](https://www.reddit.com/r/MachineLearning/comments/wqrw8x/d_fool_me_once_shame_on_you_fool_me_twice_shame/) by simpler well-proven methods. This effort aims to provide an alternative to overcome the Prophet's memory.

"It is important to note that false prophets sometimes prophesied accurately, ... "
(Deuteronomy 13:2,5)
""" + # AUTOGENERATED! DO NOT EDIT! File to edit: ../../../nbs/src/adapters.prophet.ipynb. # %% auto 0 diff --git a/python/statsforecast/arima.py b/python/statsforecast/arima.py index ffcf35b26..c149cea86 100644 --- a/python/statsforecast/arima.py +++ b/python/statsforecast/arima.py @@ -313,7 +313,7 @@ def maInvert(ma): order[1], seasonal["order"][1], ], - dtype=np.uint32, + dtype=np.int32, ) narma = arma[:4].sum().item() diff --git a/python/statsforecast/core.py b/python/statsforecast/core.py index c15af6e0b..ac9ef2b56 100644 --- a/python/statsforecast/core.py +++ b/python/statsforecast/core.py @@ -1,3 +1,5 @@ +"""Methods for Fit, Predict, Forecast (fast), Cross Validation and plotting""" + # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/core/core.ipynb. # %% auto 0 diff --git a/python/statsforecast/distributed/multiprocess.py b/python/statsforecast/distributed/multiprocess.py index 56f23f948..3f8607dab 100644 --- a/python/statsforecast/distributed/multiprocess.py +++ b/python/statsforecast/distributed/multiprocess.py @@ -1,3 +1,5 @@ +"""The computational efficiency of `StatsForecast` can be tracked to its two core components:
1. Its `models` written in NumBa that optimizes Python code to reach C speeds.
2. Its `core.StatsForecast` class that enables distributed computing.
This is a low-level class enabling other distribution methods.

""" + # AUTOGENERATED! DO NOT EDIT! File to edit: ../../../nbs/src/distributed.multiprocess.ipynb. # %% auto 0 diff --git a/python/statsforecast/feature_engineering.py b/python/statsforecast/feature_engineering.py index a3f04c14e..5d7dbf5ff 100644 --- a/python/statsforecast/feature_engineering.py +++ b/python/statsforecast/feature_engineering.py @@ -1,3 +1,5 @@ +"""Generate features for downstream models""" + # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/feature_engineering.ipynb. # %% auto 0 diff --git a/python/statsforecast/models.py b/python/statsforecast/models.py index dad5285ce..c4dc0bd5f 100644 --- a/python/statsforecast/models.py +++ b/python/statsforecast/models.py @@ -1,3 +1,5 @@ +"""Models currently supported by StatsForecast""" + # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/core/models.ipynb. # %% auto 0 diff --git a/python/statsforecast/utils.py b/python/statsforecast/utils.py index af12326f0..c528d1f39 100644 --- a/python/statsforecast/utils.py +++ b/python/statsforecast/utils.py @@ -1,3 +1,5 @@ +"""The `core.StatsForecast` class allows you to efficiently fit multiple `StatsForecast` models for large sets of time series. It operates with pandas DataFrame `df` that identifies individual series and datestamps with the `unique_id` and `ds` columns, and the `y` column denotes the target time series variable. To assist development, we declare useful datasets that we use throughout all `StatsForecast`'s unit tests.""" + # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/utils.ipynb. # %% auto 0 diff --git a/src/arima.cpp b/src/arima.cpp index 18e532676..698a8fa17 100644 --- a/src/arima.cpp +++ b/src/arima.cpp @@ -42,18 +42,22 @@ void partrans(const uint32_t p, const std::span rawv, std::tuple, py::array_t> arima_transpar(const py::array_t params_inv, - const py::array_t armav, bool trans) { + const py::array_t armav, bool trans) { assert(params_inv.ndim() == 1); assert(armav.ndim() == 1); const auto arma = make_cspan(armav); const auto params_in = make_cspan(params_inv); - const uint32_t mp = arma[0]; - const uint32_t mq = arma[1]; - const uint32_t msp = arma[2]; - const uint32_t msq = arma[3]; - const uint32_t ns = arma[4]; + assert(arma.size() == 7); + assert(arma[0] >= 0 && arma[1] >= 0 && arma[2] >= 0 && arma[3] >= 0 && + arma[4] >= 0 && arma[5] >= 0 && arma[6] >= 0); + + const int32_t mp = arma[0]; + const int32_t mq = arma[1]; + const int32_t msp = arma[2]; + const int32_t msq = arma[3]; + const int32_t ns = arma[4]; const uint32_t p = mp + ns * msp; const uint32_t q = mq + ns * msq; @@ -102,7 +106,7 @@ arima_transpar(const py::array_t params_inv, } std::tuple> -arima_css(const py::array_t yv, const py::array_t armav, +arima_css(const py::array_t yv, const py::array_t armav, const py::array_t phiv, const py::array_t thetav) { assert(yv.ndim() == 1); assert(armav.ndim() == 1); @@ -119,6 +123,10 @@ arima_css(const py::array_t yv, const py::array_t armav, const auto phi = make_cspan(phiv); const auto theta = make_cspan(thetav); + assert(arma.size() == 7); + assert(arma[0] >= 0 && arma[1] >= 0 && arma[2] >= 0 && arma[3] >= 0 && + arma[4] >= 0 && arma[5] >= 0 && arma[6] >= 0); + const uint32_t ncond = arma[0] + arma[5] + arma[4] * (arma[2] + arma[6]); uint32_t nu = 0; double ssq = 0.0; @@ -601,7 +609,7 @@ void getQ0(const py::array_t phiv, const py::array_t thetav, } py::array_t arima_gradtrans(const py::array_t xv, - const py::array_t armav) { + const py::array_t armav) { assert(xv.ndim() == 1); assert(armav.ndim() == 1); @@ -610,9 +618,13 @@ py::array_t arima_gradtrans(const py::array_t xv, const auto x = make_cspan(xv); const size_t n = x.size(); - const uint32_t mp = arma[0]; - const uint32_t mq = arma[1]; - const uint32_t msp = arma[2]; + assert(arma.size() == 7); + assert(arma[0] >= 0 && arma[1] >= 0 && arma[2] >= 0 && arma[3] >= 0 && + arma[4] >= 0 && arma[5] >= 0 && arma[6] >= 0); + + const int32_t mp = arma[0]; + const int32_t mq = arma[1]; + const int32_t msp = arma[2]; std::array w1; std::array w2; @@ -662,16 +674,20 @@ py::array_t arima_gradtrans(const py::array_t xv, } py::array_t arima_undopars(const py::array_t xv, - const py::array_t armav) { + const py::array_t armav) { assert(xv.ndim() == 1); assert(armav.ndim() == 1); const auto x = make_cspan(xv); const auto arma = make_cspan(armav); - const uint32_t mp = arma[0]; - const uint32_t mq = arma[1]; - const uint32_t msp = arma[2]; + assert(arma.size() == 7); + assert(arma[0] >= 0 && arma[1] >= 0 && arma[2] >= 0 && arma[3] >= 0 && + arma[4] >= 0 && arma[5] >= 0 && arma[6] >= 0); + + const int32_t mp = arma[0]; + const int32_t mq = arma[1]; + const int32_t msp = arma[2]; py::array_t outv(xv.size()); const auto out = make_span(outv); From 6d2af42266e1be31f1cc2cbc25086bd0495805ba Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jos=C3=A9=20Morales?= Date: Wed, 20 Nov 2024 10:55:50 -0600 Subject: [PATCH 14/14] remove docstrings --- python/statsforecast/adapters/prophet.py | 2 -- python/statsforecast/core.py | 2 -- python/statsforecast/distributed/multiprocess.py | 2 -- python/statsforecast/feature_engineering.py | 2 -- python/statsforecast/models.py | 2 -- python/statsforecast/utils.py | 2 -- 6 files changed, 12 deletions(-) diff --git a/python/statsforecast/adapters/prophet.py b/python/statsforecast/adapters/prophet.py index 01d307de5..ce3d3a667 100644 --- a/python/statsforecast/adapters/prophet.py +++ b/python/statsforecast/adapters/prophet.py @@ -1,5 +1,3 @@ -"""In 2017, Facebook open-sourced [Prophet](https://peerj.com/preprints/3190.pdf), with the promise of providing experts and non-experts the possibility of producing high-quality predictions. The forecasting community heavily adopted the solution, reaching millions of accumulated downloads. It became evident that its [quality is shadowed](https://www.reddit.com/r/MachineLearning/comments/wqrw8x/d_fool_me_once_shame_on_you_fool_me_twice_shame/) by simpler well-proven methods. This effort aims to provide an alternative to overcome the Prophet's memory.

"It is important to note that false prophets sometimes prophesied accurately, ... "
(Deuteronomy 13:2,5)
""" - # AUTOGENERATED! DO NOT EDIT! File to edit: ../../../nbs/src/adapters.prophet.ipynb. # %% auto 0 diff --git a/python/statsforecast/core.py b/python/statsforecast/core.py index ac9ef2b56..c15af6e0b 100644 --- a/python/statsforecast/core.py +++ b/python/statsforecast/core.py @@ -1,5 +1,3 @@ -"""Methods for Fit, Predict, Forecast (fast), Cross Validation and plotting""" - # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/core/core.ipynb. # %% auto 0 diff --git a/python/statsforecast/distributed/multiprocess.py b/python/statsforecast/distributed/multiprocess.py index 3f8607dab..56f23f948 100644 --- a/python/statsforecast/distributed/multiprocess.py +++ b/python/statsforecast/distributed/multiprocess.py @@ -1,5 +1,3 @@ -"""The computational efficiency of `StatsForecast` can be tracked to its two core components:
1. Its `models` written in NumBa that optimizes Python code to reach C speeds.
2. Its `core.StatsForecast` class that enables distributed computing.
This is a low-level class enabling other distribution methods.

""" - # AUTOGENERATED! DO NOT EDIT! File to edit: ../../../nbs/src/distributed.multiprocess.ipynb. # %% auto 0 diff --git a/python/statsforecast/feature_engineering.py b/python/statsforecast/feature_engineering.py index 5d7dbf5ff..a3f04c14e 100644 --- a/python/statsforecast/feature_engineering.py +++ b/python/statsforecast/feature_engineering.py @@ -1,5 +1,3 @@ -"""Generate features for downstream models""" - # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/feature_engineering.ipynb. # %% auto 0 diff --git a/python/statsforecast/models.py b/python/statsforecast/models.py index c4dc0bd5f..dad5285ce 100644 --- a/python/statsforecast/models.py +++ b/python/statsforecast/models.py @@ -1,5 +1,3 @@ -"""Models currently supported by StatsForecast""" - # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/core/models.ipynb. # %% auto 0 diff --git a/python/statsforecast/utils.py b/python/statsforecast/utils.py index c528d1f39..af12326f0 100644 --- a/python/statsforecast/utils.py +++ b/python/statsforecast/utils.py @@ -1,5 +1,3 @@ -"""The `core.StatsForecast` class allows you to efficiently fit multiple `StatsForecast` models for large sets of time series. It operates with pandas DataFrame `df` that identifies individual series and datestamps with the `unique_id` and `ds` columns, and the `y` column denotes the target time series variable. To assist development, we declare useful datasets that we use throughout all `StatsForecast`'s unit tests.""" - # AUTOGENERATED! DO NOT EDIT! File to edit: ../../nbs/src/utils.ipynb. # %% auto 0