diff --git a/docs/REF/Matrix_Profile_Derivation.ipynb b/docs/REF/Matrix_Profile_Derivation.ipynb index 2dd8deb33..2f456dfa1 100644 --- a/docs/REF/Matrix_Profile_Derivation.ipynb +++ b/docs/REF/Matrix_Profile_Derivation.ipynb @@ -6,7 +6,7 @@ "source": [ "# Computing the Z-Normalized Euclidean Distance from Dot Products\n", "\n", - "In the [Matrix Profile I](https://www.cs.ucr.edu/~eamonn/STOMP_GPU_final_submission_camera_ready.pdf) and [Matrix Profile II](https://www.cs.ucr.edu/~eamonn/PID4481997_extend_Matrix%20Profile_I.pdf) papers, the Z-normalized Euclidean distance between a query subsequence, $Q_{i,m}=(q_i, q_{i+1}, q_{i+2}\\ldots, q_{i+m-1})$, and the $i^{th}$ subsequence, $T_{i,m}=(t_i, t_{i+1}, t_{i+2}, \\ldots, t_{i+m-1})$, with window size, $m$, in the time series, $T$, can be computed following:" + "In the [Matrix Profile I](https://www.cs.ucr.edu/~eamonn/STOMP_GPU_final_submission_camera_ready.pdf) and [Matrix Profile II](https://www.cs.ucr.edu/~eamonn/PID4481997_extend_Matrix%20Profile_I.pdf) papers, the **Z-normalized Euclidean distance**, $D$, between a query subsequence, $Q_{i,m}=(q_i, q_{i+1}, q_{i+2}\\ldots, q_{i+m-1})$, with window length/size, $m$, and an arbitrary subsequence (also with size $m$) in the time series $T$, $T_{i,m}=(t_i, t_{i+1}, t_{i+2}, \\ldots, t_{i+m-1})$, can be computed from the following:" ] }, { @@ -14,7 +14,7 @@ "metadata": {}, "source": [ "\\begin{align}\n", - " D(Q_{i,m}, T_{i,m}) ={}&\n", + " D(Q_{i,m}, T_{i,m}) &=\n", " \\sqrt{\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", " \\left(\n", @@ -24,13 +24,13 @@ " \\right)^2\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", " \\left[\n", " \\left(\n", " \\frac{t_{i+j}-M_{T_{i,m}}}{\\Sigma_{T_{i,m}}}\n", - " \\right)\n", + " \\right)^2\n", " -\n", " 2\n", " \\left(\n", @@ -46,7 +46,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", " \\left(\n", @@ -59,7 +59,7 @@ " \\frac{t_{i+j}-M_{T_{i,m}}}{\\Sigma_{T_{i,m}}}\n", " \\right)\n", " \\left(\n", - " \\frac{q_j-\\mu_{Q_{i,m}}}{\\sigma_{Q_{i,m}}}\n", + " \\frac{q_{i+j}-\\mu_{Q_{i,m}}}{\\sigma_{Q_{i,m}}}\n", " \\right)\n", " +\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", @@ -68,7 +68,7 @@ " \\right)^2\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " m\n", " -\n", @@ -84,7 +84,7 @@ " m\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " -\n", @@ -94,11 +94,11 @@ " \\frac{t_{i+j}-M_{T_{i,m}}}{\\Sigma_{T_{i,m}}}\n", " \\right)\n", " \\left(\n", - " \\frac{q_{i+j}-\\mu_{Q_m}}{\\sigma_{Q_{i,m}}}\n", + " \\frac{q_{i+j}-\\mu_{Q_{i,m}}}{\\sigma_{Q_{i,m}}}\n", " \\right)\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -113,9 +113,10 @@ " \\frac{q_{i+j}-\\mu_{Q_{i,m}}}{\\sigma_{Q_{i,m}}}\n", " \\right)\n", " \\right]\n", - " }\n", + " } \n", + " \\tag{i}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -135,7 +136,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -143,17 +144,17 @@ " -\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", " \\frac{\n", - " t_{i+j}q_j\n", + " t_{i+j}q_{i+j}\n", " -t_{i+j}\\mu_{Q_{i,m}}\n", " -M_{T_{i,m}}q_{i+j}\n", - " +M_{T_{i,m}}\\mu_{Q_{i+m}}\n", + " +M_{T_{i,m}}\\mu_{Q_{i,m}}\n", " }{\n", - " m \\sigma_{Q_{i+m}} \\Sigma_{T_{i,m}}\n", + " m \\sigma_{Q_{i,m}} \\Sigma_{T_{i,m}}\n", " }\n", " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -177,7 +178,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -200,7 +201,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -220,7 +221,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -238,7 +239,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -247,7 +248,7 @@ " \\frac{\n", " Q_{i,m}\\cdot{T_{i,m}}\n", " -\n", - " \\mu_{Q{i,m}}m\n", + " \\mu_{Q_{i,m}}m\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", " \\frac{t_{i+j}}{m}\n", " }{\n", @@ -256,7 +257,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -275,6 +276,13 @@ "\\end{align}" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "where $\\mu_{Q_{i,m}}$, $\\sigma_{Q_{i,m}}$, $T_{i,m}$, and $\\Sigma_{T_{i,m}}$ are the corresponding population mean(s) and standard deviation(s) for subsequences $Q_{i,m}$ and $T_{i,m}$, respectively." + ] + }, { "cell_type": "markdown", "metadata": {}, @@ -295,7 +303,7 @@ " \\right]\n", " }{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", " \\\\\n", - " ={}& \n", + " &= \n", " \\frac{\n", " \\langle\n", " \\left(\n", @@ -308,7 +316,7 @@ " \\rangle\n", " }{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{1}{m}\n", " \\sum \\limits _{0 \\leq j \\lt m}\n", " \\frac{\n", @@ -320,7 +328,7 @@ " \\right)\n", " }{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{1}{m}\n", " \\sum \\limits _{0 \\leq j \\lt m}\n", " \\left(\n", @@ -341,7 +349,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Similar to above, the Z-normalized Euclidean distance can be computed from $\\rho$ following:" + "Substituting $\\rho(Q_{i,m},T_{i,m})$ into equation (i) above, the Z-normalized Euclidean distance becomes:" ] }, { @@ -349,7 +357,7 @@ "metadata": {}, "source": [ "\\begin{align}\n", - " D(Q_{i,m}, T_{i,m}) ={}&\n", + " D(Q_{i,m}, T_{i,m}) &=\n", " \\sqrt{\n", " \\sum \\limits _{0 \\leq {j} \\lt m}\n", " \\left(\n", @@ -361,7 +369,7 @@ " \\\\\n", " \\vdots\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -378,7 +386,7 @@ " \\right]\n", " }\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\sqrt{\n", " 2m\n", " \\left[\n", @@ -395,7 +403,9 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Thus, by employing the most efficient way to compute $\\rho(Q_{i,m},T_{i,m})$, then we'd also have an efficient way to directly compute $D(Q_{i,m},T_{i,m})$. Recall that:" + "Thus, by finding the most efficient way to compute $\\rho(Q_{i,m},T_{i,m})$, it follows that we would also have the most efficient way to directly compute $D(Q_{i,m},T_{i,m})$.\n", + "\n", + "Furthermore, recall the basic relationship between Pearson correlation and covariance:" ] }, { @@ -411,7 +421,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Thus, it follows that finding the most efficient way to compute the covariance matrix, $cov(Q_{i,m},T_{i,m})$ would result in the most efficient way to compute the distance. Also, remember that we would like to traverse our distance matrix along each diagonal rather than along each row/column." + "Therefore, it also follows that finding the most efficient way to compute the covariance, $cov(Q_{i,m},T_{i,m})$ would ultimately result in the most efficient way to compute the distance. Also, remember that we would like to traverse our distance matrix along each diagonal rather than along each row/column as this would allow us to leverage an iterative recurrence algorithm that can dramatically decrease the overall computational time." ] }, { @@ -433,20 +443,20 @@ "metadata": {}, "source": [ "\\begin{align}\n", - " cov(Q_{i,m},T_{i,m}) ={}& E\n", + " cov(Q_{i,m},T_{i,m}) &= E\n", " \\left[\n", " \\left(\n", - " Q-\\mu_{Q_{i,m}}\n", + " Q_{i,m}-\\mu_{Q_{i,m}}\n", " \\right)\n", " \\left(\n", " T_{i,m}-M_{T_{i,m}}\n", " \\right)\n", " \\right]\n", " \\\\\n", - " ={}& \n", + " &=\n", " \\langle\n", " \\left(\n", - " Q-\\mu_{Q_{i,m}}\n", + " Q_{i,m}-\\mu_{Q_{i,m}}\n", " \\right)\n", " ,\n", " \\left(\n", @@ -454,7 +464,7 @@ " \\right)\n", " \\rangle\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{1}{m}\n", " \\sum \\limits _{0 \\leq j \\lt m}\n", " \\left(\n", @@ -471,7 +481,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Note that we've explicitly called out the fact that the means, $\\mu_{Q_{{i,m}}}$ and $M_{T_{i,m}}$, are computed with the subsequences of length $m$. Additionally, according to Welford, we can express these means with respect to the means of the same subsequences that have their last elements removed (i.e., $\\mu_{Q_{i,m-1}}$ and $M_{T_{i,m-1}}$)." + "Note that we've explicitly called out the fact that the means, $\\mu_{Q_{{i,m}}}$ and $M_{T_{i,m}}$, are computed with the subsequences of length $m$. Additionally, we can express these means with respect to the means of the same subsequences that have their last elements removed (i.e., $\\mu_{Q_{i,m-1}}$ and $M_{T_{i,m-1}}$)." ] }, { @@ -480,7 +490,7 @@ "source": [ "\\begin{align}\n", " cov(Q_{i,m},T_{i,m}) \n", - " ={}&\n", + " &=\n", " \\frac{1}{m}\n", " \\sum \\limits _{0 \\leq j \\lt m}\n", " \\left(\n", @@ -490,7 +500,7 @@ " t_{i+j}-M_{T_{i,m}}\n", " \\right)\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{\n", " S(Q_{i,m-1}, T_{i,m-1})\n", " +\n", @@ -505,7 +515,7 @@ " \\right)\n", " }{m}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{\n", " \\frac{m-1}{m-1}S(Q_{i,m-1}, T_{i,m-1})\n", " +\n", @@ -520,7 +530,7 @@ " \\right)\n", " }{m}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{\n", " cov(Q_{i,m-1},T_{i,m-1}) (m-1)\n", " +\n", @@ -535,7 +545,7 @@ " \\right)\n", " }{m}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{m-1}{m} \n", " \\left[\n", " cov(Q_{i,m-1},T_{i,m-1})\n", @@ -566,7 +576,7 @@ "source": [ "\\begin{align}\n", " cov(Q_{i-1,m},T_{i-1,m}) \n", - " ={}&\n", + " &=\n", " \\frac{1}{m}\n", " \\sum \\limits _{0 \\leq j \\lt m}\n", " \\left(\n", @@ -576,7 +586,7 @@ " t_{i+j-1}-M_{T_{i-1,m}}\n", " \\right)\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{\n", " S(Q_{i,m-1},T_{i,m-1})\n", " +\n", @@ -593,7 +603,7 @@ " \\right)\n", " }{m}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{\n", " \\frac{m-1}{m-1}S(Q_{i,m-1},T_{i,m-1})\n", " +\n", @@ -610,7 +620,7 @@ " \\right)\n", " }{m}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{\n", " cov(Q_{i,m-1},T_{i,m-1}) (m-1)\n", " +\n", @@ -629,7 +639,7 @@ " \\right)\n", " }{m}\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{m-1}{m} \n", " \\left[\n", " cov(Q_{i,m-1},T_{i,m-1})\n", @@ -664,7 +674,7 @@ "source": [ "\\begin{align}\n", " cov(Q_{i-1,m},T_{i-1,m})\n", - " ={}&\n", + " &=\n", " \\frac{m-1}{m} \n", " \\left[\n", " cov(Q_{i,m-1},T_{i,m-1})\n", @@ -698,7 +708,7 @@ " M_{T_{i,m-1}} \n", " \\right)\n", " }{m}\n", - " ={}&\n", + " &=\n", " cov(Q_{i,m-1},T_{i,m-1})\n", " \\\\\n", "\\end{align}" @@ -717,7 +727,7 @@ "source": [ "\\begin{align}\n", " cov(Q_{i,m},T_{i,m})\n", - " ={}&\n", + " &=\n", " \\frac{m-1}{m} \n", " \\left[\n", " cov(Q_{i,m-1},T_{i,m-1})\n", @@ -732,7 +742,7 @@ " }{m}\n", " \\right]\n", " \\\\\n", - " ={}&\n", + " &=\n", " \\frac{m-1}{m} \n", " \\left[\n", " \\frac{m}{m-1}\n", @@ -761,7 +771,7 @@ " }{m}\n", " \\right]\n", " \\\\\n", - " ={}&\n", + " &=\n", " cov(Q_{i-1,m},T_{i-1,m})\n", " +\n", " \\frac{m-1}{m^2}\n", @@ -796,10 +806,10 @@ "\n", "\\begin{align}\n", " \\rho(Q_{i,m},T_{i,m}) \n", - " &{}= \n", + " &= \n", " \\frac{cov(Q_{i,m},T_{i,m})}{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", " \\\\\n", - " &{}=\n", + " &=\n", " \\frac{\n", " cov(Q_{i-1,m},T_{i-1,m})\n", " +\n", @@ -832,15 +842,674 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "# Z-Normalized Distance" + "# Catastrophic Cancellation\n", + "\n", + "However, if you look at the recurrence relationship, the covariance, $cov(Q_{i,m},T_{i,m})$, is computed based on the previously calculated covariance, $cov(Q_{i-1,m},T_{i-1,m})$, \"updated\" by a simple looking \"correction\" term $[(q_{i+m-1} - \\mu_{Q_{i,m-1}})(t_{i+m-1} - M_{T_{i,m-1}}) - (q_{i-1} - \\mu_{Q_{i,m-1}}) (t_{i-1} - M_{T_{i,m-1}})]$ and that is \"normalized\" by $\\frac{m-1}{m^2}$. Sadly, this correction term, which involves the subtraction of two nearby numbers, can lead to [catastrophic cancellation](https://en.wikipedia.org/wiki/Catastrophic_cancellation) and rounding errors introduced during the recurrence can accumulate over subsequent updates, potentially causing significant loss of accuracy. Luckily, Section 3.3 of the [Matrix Profile XIV: Scaling Time Series Motif Discovery with GPUs to Break a Quintillion Pairwise Comparisons a Day and Beyond](https://www.cs.ucr.edu/~eamonn/SCAMP-camera-ready-final1.pdf) paper recommends using a \"centered sum-of-products\" to improve the numerical stability (i.e., we can avoid the much larger catastrophic cancellation but we can't avoid the extremely small loss of precision in the recurrence)." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "> **Why are the two products \"nearby numbers?**\n", + "> \n", + "> 1. Same mean, similar data:\n", + ">\n", + ">> Both products use essentially the same window means, $\\mu_{Q_{i,m-1}}$ and $T_{T_{i,m-1}}$ and the window has only shifted by one time step. If the underlying time series doesn't change drastically point-to-point (which is typical for real-world sensor/financial/etc. data), then $q_{i+m-1}$ and $q_{i-1}$ are just two individual data points being compared against nearly the same mean and their deviations from that mean, are often similar in magnitude\n", + ">\n", + "> 2. It's essentially new-minus-old:\n", + ">\n", + ">> The term represents \"add the new point entering the window\" minus \"remove the old point leaving the window.\" When you slide a window by one sample, you're swapping out one point for another — the two points being subtracted are computed the same way, just at different time indices. Structurally this makes the whole bracket a difference of two similar quantities, not two independent large numbers.\n", + ">\n", + "> 2. Catastrophic cancellation risk.\n", + "> \n", + ">> Because the two products are numerically close, subtracting them in floating-point arithmetic can wipe out significant digits (i.e., you're left mostly with rounding error rather than meaningful signal). This is the classic \"subtracting two nearly-equal floats\" numerical stability problem (the same issue that motivates using Welford's method for streaming variance in the first place, rather than the naive $E[X^2] - E[X]^2$ formula).\n", + ">\n", + "> So \"nearby\" doesn't mean the two products are always numerically identical — it means that for realistic, not-wildly-discontinuous time series, consecutive/nearby data points tend to have deviations-from-mean of comparable size, making this update term prone to precision loss when computed naively in floating point.\n", + ">\n", + "> **Cumulative Loss of Precision**\n", + ">\n", + "> Since we are traversing the diagonal and using the previous value to update the next, the small loss of precision in the first iteration will affect the second iteration and so on and so forth up until the $n^{th}$ iteration. However, if the loss of precision is small (without catastrophic cancellation), the accumulated loss is actually imperceptible (i.e., the magnitude is much, much smaller than we care about and likely smaller than the precision of our raw data). When there is catastrophic cancellation, this really gets amplified quickly as the loss of precision is BIG. In other words, we can resolve the BIG catastrophic cancellation but we can (probably?) NEVER resolve the extremely small loss of precision being accumulated since we are dealing with floating point numbers, which must involve loss of precision." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Improving the Numerical Stability\n", + "\n", + "Starting with the original covariance equation above:\n", + "\n", + "\\begin{align}\n", + " cov(Q_{i,m},T_{i,m}) \n", + " &=\n", + " \\frac{1}{m}\n", + " \\sum \\limits _{0 \\leq j \\lt m}\n", + " \\left(\n", + " q_{i+j}-\\mu_{Q_{i,m}}\n", + " \\right)\n", + " \\left(\n", + " t_{i+j}-M_{T_{i,m}}\n", + " \\right)\n", + "\\end{align}\n", + "\n", + "the summation term on the right is called a \"centered sum-of-products\" because each term in the sum is a product of two deviations from their respective means, rather than a product of raw values. For simplicity, we can focus only on the summation and let:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}}\n", + " &=\n", + " \\sum \\limits_{0 \\leq j \\lt m}\n", + " \\left(\n", + " q_{i+j} - \\mu_{Q_{i,m}}\n", + " \\right)\n", + " \\left(\n", + " t_{i+j} - M_{T_{i,m}}\n", + " \\right)\n", + " \\qquad\\text{and} \\qquad \\therefore\n", + " cov(Q_{i,m},T_{i,m}) = \\frac{\\overline{Q_{i,m}T_{i,m}}}{m}\n", + " \\qquad\\text{and} \\qquad\n", + " \\rho(Q_{i,m},T_{i,m}) = \\frac{cov(Q_{i,m},T_{i,m})}{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", + "\\end{align}\n", + "\n", + "Thus, unlike the Welford approach earlier, if you can calculate $\\overline{Q_{i,m}T_{i,m}}$ in a numerically stable way, then $cov(Q_{i,m},T_{i,m})$ can be trivially computed by dividing $\\overline{Q_{i,m}T_{i,m}}$ by $m$. So, let's expand the product and see where we can get:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}}\n", + " &=\n", + " \\sum \\limits_{0 \\leq j \\lt m}\n", + " \\left(\n", + " q_{i+j} - \\mu_{Q_{i,m}}\n", + " \\right)\n", + " \\left(\n", + " t_{i+j} - M_{T_{i,m}}\n", + " \\right)\n", + " \\\\\n", + " &=\n", + " \\sum \\limits_{0 \\leq j \\lt m}\n", + " q_{i+j} t_{i+j}\n", + " - \\sum \\limits_{0 \\leq j \\lt m}\n", + " q_{i+j} M_{T_{i,m}}\n", + " - \\sum \\limits_{0 \\leq j \\lt m}\n", + " \\mu_{Q_{i,m}} t_{i+j}\n", + " + \\sum \\limits_{0 \\leq j \\lt m}\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\\\\n", + " &=\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - \\sum \\limits_{0 \\leq j \\lt m}\n", + " q_{i+j} M_{T_{i,m}}\n", + " - \\sum \\limits_{0 \\leq j \\lt m}\n", + " \\mu_{Q_{i,m}} t_{i+j}\n", + " + \\sum \\limits_{0 \\leq j \\lt m}\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\\\\n", + " &=\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - \\left(\n", + " M_{T_{i,m}} \n", + " \\sum \\limits_{0 \\leq j \\lt m} q_{i+j}\n", + " \\right)\n", + " - \\left(\n", + " \\mu_{Q_{i,m}}\n", + " \\sum \\limits_{0 \\leq j \\lt m} t_{i+j}\n", + " \\right)\n", + " + \\left(\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\sum \\limits_{0 \\leq j \\lt m} 1\n", + " \\right)\n", + " \\\\\n", + " &=\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - \\left(\n", + " M_{T_{i,m}}\n", + " \\frac{m}{m}\n", + " \\sum \\limits_{0 \\leq j \\lt m} q_{i+j}\n", + " \\right)\n", + " - \\left(\n", + " \\mu_{Q_{i,m}}\n", + " \\frac{m}{m}\n", + " \\sum \\limits_{0 \\leq j \\lt m} t_{i+j}\n", + " \\right)\n", + " + \\left(\n", + " \\mu_{Q_{i,m}}\n", + " M_{T_{i,m}}\n", + " m\n", + " \\right)\n", + " \\\\\n", + " &=\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - \\left(\n", + " m\n", + " M_{T_{i,m}}\n", + " \\frac{\\sum \\limits_{0 \\leq j \\lt m} q_{i+j}}{m}\n", + " \\right)\n", + " - \\left(\n", + " m\n", + " \\mu_{Q_{i,m}}\n", + " \\frac{\\sum \\limits_{0 \\leq j \\lt m} t_{i+j}}{m}\n", + " \\right)\n", + " + \\left(\n", + " m\n", + " \\mu_{Q_{i,m}}\n", + " M_{T_{i,m}}\n", + " \\right)\n", + " \\\\\n", + " &=\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - m M_{T_{i,m}} \\mu_{Q_{i,m}}\n", + " - m \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " + m \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\\\\n", + " \\\\\n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} ={}&\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - m \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\qquad\\text{and}\\qquad\n", + " cov(Q_{i,m},T_{i,m}) = \\frac{\\overline{Q_{i,m}T_{i,m}}}{m}\n", + " \\qquad\\text{and}\\qquad\n", + " \\rho(Q_{i,m},T_{i,m}) = \\frac{cov(Q_{i,m},T_{i,m})}{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", + " \\\\\n", + "\\end{align}\n", + "\n", + "So, if we can update $\\overline{Q_{i,m}T_{i,m}}$ directly along the diagonal, without ever propagating $cov(Q_{i,m},T_{i,m})$ itself, then we can sidestep the compounding cancellation described earlier! Unfortunately, \n", + "\n", + "1. $Q_{i,m} \\cdot T_{i,m}$ is still a sum of raw products $q_{i} t_{i}$, which can be large relative to the covariance itself, and\n", + "2. the final subtraction $cov(Q_{i,m},T_{i,m}) = \\frac{Q_{i,m} \\cdot T_{i,m}}{m} - \\mu_{Q_{i,m}} M_{T_{i,m}}$ is still a single cancellation-prone step per index (the classic $E[XY]-E[X]E[Y]$ instability)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}} &=\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - m \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - m \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " \\right)\n", + " - \\left(\n", + " Q_{i-1,m} \\cdot T_{i-1,m}\n", + " - m \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right) \n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " Q_{i,m} \\cdot T_{i,m}\n", + " - Q_{i-1,m} \\cdot T_{i-1,m}\n", + " \\right)\n", + " -m \\left(\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right) \n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " \\sum \\limits_{0 \\le j \\lt m} {q_{i+j}t_{i+j}}\n", + " - \\sum \\limits_{0 \\le j \\lt m} {q_{i-1+j}t_{i-1+j}}\n", + " \\right)\n", + " -m \\left(\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right) \n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i}t_{i} + \\ldots + q_{i+m-2}t_{i+m-2} + q_{i+m-1}t_{i+m-1}\n", + " - q_{i-1}t_{i-1} - q_{i}t_{i} - \\ldots - q_{i-1+m-1}t_{i-1+m-1}\n", + " \\right)\n", + " \\\\ \n", + " & \\qquad\n", + " -m \\left(\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right) \n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i}t_{i} - q_{i}t_{i} + \\ldots + q_{i+m-2}t_{i+m-2} - q_{i-1+m-1}t_{i-1+m-1}\n", + " + q_{i+m-1}t_{i+m-1} - q_{i-1}t_{i-1}\n", + " \\right)\n", + " \\\\\n", + " & \\qquad\n", + " -m \\left(\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right) \n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i+m-1}t_{i+m-1} - q_{i-1}t_{i-1}\n", + " \\right)\n", + " -m \\left(\n", + " \\mu_{Q_{i,m}} M_{T_{i,m}}\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right) \n", + " \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Recall, that the standard rolling mean recurrence can be computed with $\\mu_{Q_{i,m}} = \\mu_{Q_{i-1,m}} + \\frac{q_{i+m-1}-q_{i-1}}{m}$ and $M_{T_{i,m}} = M_{T_{i-1,m}} + \\frac{t_{i+m-1}-t_{i-1}}{m}$ and after substituing this into the first term of right-most parantheses in the final equation above, we obtain:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i+m-1}t_{i+m-1} - q_{i-1}t_{i-1}\n", + " \\right)\n", + " -m \\left[\n", + " \\left(\n", + " \\mu_{Q_{i-1,m}} + \\frac{q_{i+m-1}-q_{i-1}}{m}\n", + " \\right)\n", + " \\left(\n", + " M_{T_{i-1,m}} + \\frac{t_{i+m-1}-t_{i-1}}{m}\n", + " \\right)\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right]\n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i+m-1}t_{i+m-1} - q_{i-1}t_{i-1}\n", + " \\right)\n", + " \\\\\n", + " -m & \\left(\n", + " \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " + \\mu_{Q_{i-1,m}} \\frac{t_{i+m-1}-t_{i-1}}{m}\n", + " + \\frac{q_{i+m-1}-q_{i-1}}{m} M_{T_{i-1,m}}\n", + " + \\frac{q_{i+m-1}-q_{i-1}}{m} \\frac{t_{i+m-1}-t_{i-1}}{m}\n", + " - \\mu_{Q_{i-1,m}} M_{T_{i-1,m}}\n", + " \\right)\n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i+m-1}t_{i+m-1} - q_{i-1}t_{i-1}\n", + " \\right)\n", + " -m \\left(\n", + " \\frac{t_{i+m-1}-t_{i-1}}{m} \\mu_{Q_{i-1,m}}\n", + " + \\frac{q_{i+m-1}-q_{i-1}}{m} M_{T_{i-1,m}}\n", + " + \\frac{q_{i+m-1}-q_{i-1}}{m} \\frac{t_{i+m-1}-t_{i-1}}{m}\n", + " \\right)\n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left(\n", + " q_{i+m-1}t_{i+m-1} - q_{i-1}t_{i-1}\n", + " \\right)\n", + " - \\left( t_{i+m-1}-t_{i-1} \\right) \\mu_{Q_{i-1,m}}\n", + " - \\left( q_{i+m-1}-q_{i-1} \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left(q_{i+m-1}-q_{i-1}\\right)\\left(t_{i+m-1}-t_{i-1}\\right)}{m}\n", + " \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now, let's simplify the way that this looks and let:\n", + "\n", + "1. $a = q_{i+m-1}$\n", + "2. $b = q_{i-1}$\n", + "3. $c = t_{i+m-1}$\n", + "4. $d = t_{i-1}$\n", + "\n", + "Then, the above equation becomes a much more compact looking:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left( ac - bd \\right)\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now, we can manipulate the first term on the right, then eventually combine it with the last three terms, and then factor out like terms:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}} &=\n", + " \\left[\n", + " \\frac{2}{2} \\left( ac - bd \\right)\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{1}{2} \\left( 2ac - 2bd \\right)\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{1}{2} \\left( ac + ac - bd - bd + 0 + 0 \\right)\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{1}{2} \\left[ ac + ac - bd - bd + (ad - ad) + (bc - bc) \\right]\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{1}{2} \\left[ \\left( ac - ad + bc - bd \\right) + \\left( ac + ad - bd - bc \\right) \\right]\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{1}{2}\n", + " \\left( \n", + " \\left[\n", + " a\\left( c - d \\right) + b \\left( c - d \\right)\n", + " \\right] \n", + " +\n", + " \\left[\n", + " a\\left( c + d \\right) - b \\left( c + d \\right)\n", + " \\right] \n", + " \\right)\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{1}{2} \\left[ \\left( a + b \\right) \\left( c - d \\right) + \\left( a - b \\right) \\left( c + d \\right) \\right]\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{\\left( a + b \\right)}{2} \\left( c - d \\right)\n", + " + \\frac{\\left( c + d \\right)}{2} \\left( a - b \\right)\n", + " \\right]\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{m}\n", + " \\\\\n", + " &=\n", + " \\frac{\\left( c - d \\right)}{2} \\left( a + b \\right)\n", + " + \\frac{\\left( a - b \\right)}{2} \\left( c + d \\right)\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}} \n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{2m}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{2m}\n", + " \\\\\n", + " &=\n", + " \\left[\n", + " \\frac{\\left( c - d \\right)}{2} \\left( a + b \\right)\n", + " - \\left( c - d \\right) \\mu_{Q_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{2m}\n", + " \\right]\n", + " + \n", + " \\left[\n", + " \\frac{\\left( a - b \\right)}{2} \\left( c + d \\right)\n", + " - \\left( a - b \\right) M_{T_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) \\left( c - d \\right)}{2m}\n", + " \\right]\n", + " \\\\\n", + " &=\n", + " \\frac{\\left( c - d \\right)}{2} \n", + " \\left[\n", + " \\left( a + b \\right)\n", + " - 2 \\mu_{Q_{i-1,m}}\n", + " - \\frac{\\left( a - b \\right) }{m}\n", + " \\right]\n", + " + \\frac{\\left( a - b \\right)}{2}\n", + " \\left[\n", + " \\left( c + d \\right)\n", + " - 2 M_{T_{i-1,m}}\n", + " - \\frac{\\left( c - d \\right)}{m}\n", + " \\right]\n", + " \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Notice here that all of the major terms are written as sums/differences between values from the **same** subsequence (e.g., $(a+b)$, $(a-b)$, $(c+d)$, $(c-d)$), rather than products (e.g., $ac$, $bd$, or dot products) across different subsequences! This should help dramatically improve the numerical stability but how can we continue to make progress on this gnarly equation?" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Recall from the standard derivation of the rolling mean that:\n", + "\n", + "\\begin{align}\n", + " \\mu_{Q_{i,m}} &= \\mu_{Q_{i-1,m}} + \\frac{q_{i+m-1} - q_{i-1}}{m} \\\\\n", + " \\mu_{Q_{i,m}} - \\mu_{Q_{i-1,m}} &= \\frac{a-b}{m} \\\\\n", + "\\end{align}\n", + "\n", + "and, similarly,\n", + "\\begin{align}\n", + " M_{T_{i,m}} &= M_{T_{i-1,m}} + \\frac{t_{i+m-1} - t_{i-1}}{m} \\\\\n", + " M_{T_{i,m}} - M_{T_{i-1,m}} &= \\frac{c-d}{m} \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "So, substituting this into the equation above, we get:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}}\n", + " &=\n", + " \\frac{\\left( c - d \\right)}{2} \n", + " \\left[\n", + " \\left( a + b \\right)\n", + " - 2 \\mu_{Q_{i-1,m}}\n", + " - \\left( \\mu_{Q_{i,m}} - \\mu_{Q_{i-1,m}} \\right)\n", + " \\right]\n", + " + \\frac{\\left( a - b \\right)}{2}\n", + " \\left[\n", + " \\left( c + d \\right)\n", + " - 2 M_{T_{i-1,m}}\n", + " - \\left( M_{T_{i,m}} - M_{T_{i-1,m}} \\right)\n", + " \\right]\n", + " \\\\\n", + " &=\n", + " \\frac{\\left( c - d \\right)}{2} \n", + " \\left[\n", + " \\left( a + b \\right)\n", + " - \\mu_{Q_{i-1,m}}\n", + " - \\mu_{Q_{i,m}}\n", + " \\right]\n", + " + \\frac{\\left( a - b \\right)}{2}\n", + " \\left[\n", + " \\left( c + d \\right)\n", + " - M_{T_{i-1,m}}\n", + " - M_{T_{i,m}}\n", + " \\right]\n", + " \\\\\n", + " &=\n", + " \\frac{\\left( c - d \\right)}{2} \n", + " \\left[\n", + " \\left( a - \\mu_{Q_{i,m}} \\right)\n", + " +\n", + " \\left( b - \\mu_{Q_{i-1,m}} \\right)\n", + " \\right]\n", + " + \\frac{\\left( a - b \\right)}{2}\n", + " \\left[\n", + " \\left( c - M_{T_{i,m}} \\right)\n", + " +\n", + " \\left( d - M_{T_{i-1,m}} \\right)\n", + " \\right]\n", + " \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "And substituting back $a, b, c, d$ and moving $\\overline{Q_{i-1,m}T_{i-1,m}}$ over to the right hand side, we get the updated recurrence:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}} - \\overline{Q_{i-1,m}T_{i-1,m}}\n", + " &=\n", + " \\frac{\\left( t_{i+m-1} - t_{i-1} \\right)}{2} \n", + " \\left[\n", + " \\left( q_{i+m-1} - \\mu_{Q_{i,m}} \\right)\n", + " +\n", + " \\left( q_{i-1} - \\mu_{Q_{i-1,m}} \\right)\n", + " \\right]\n", + " + \\frac{\\left( q_{i+m-1} - q_{i-1} \\right)}{2}\n", + " \\left[\n", + " \\left( t_{i+m-1} - M_{T_{i,m}} \\right)\n", + " +\n", + " \\left( t_{i-1} - M_{T_{i-1,m}} \\right)\n", + " \\right]\n", + " \\\\\n", + " \\overline{Q_{i,m}T_{i,m}}\n", + " &=\n", + " \\overline{Q_{i-1,m}T_{i-1,m}}\n", + " \\\\\n", + " + & \\frac{ t_{i+m-1} - t_{i-1} }{2} \n", + " \\left[\n", + " \\left( q_{i+m-1} - \\mu_{Q_{i,m}} \\right)\n", + " +\n", + " \\left( q_{i-1} - \\mu_{Q_{i-1,m}} \\right)\n", + " \\right]\n", + " + \\frac{ q_{i+m-1} - q_{i-1} }{2}\n", + " \\left[\n", + " \\left( t_{i+m-1} - M_{T_{i,m}} \\right)\n", + " +\n", + " \\left( t_{i-1} - M_{T_{i-1,m}} \\right)\n", + " \\right]\n", + " \\\\\n", + "\\end{align}\n", + "\n", + "The recurrence avoids the potentially severe cancellation associated with forming large raw sum-of-products and subsequently subtracting $m\\mu_{Q_{i,m}}M_{T_{i,m}}$. It does not eliminate floating-point error entirely: the differences $q_{i+m-1}-q_{i-1}$ and $t_{i+m-1}-t_{i-1}$, as well as the recursive accumulation itself can still introduce some small amount of error. \n", + "\n", + "> An estimate of the accumulated error\n", + "> For a time series with $n = 100,000,000$ (with varying dynamic range), the old approach (with catastrophic cancellation) would be precise to $~10^{-4}-10^{-8}$ while the new approach would be (MORE) precise to $~10^{-8}-10^{-12}$. Note that double float precision is probably good to, at best, $10^{-16}$.\n", + " \n", + "\n", + "Additionally, all of these quantities can be precomputed! If we let:\n", + "\n", + "1. $df^{T}_{i,m} = \\frac{ t_{i+m-1} - t_{i-1} }{2}$\n", + "2. $df^{Q}_{i,m} = \\frac{ q_{i+m-1} - q_{i-1} }{2}$\n", + "3. $dg^{T}_{i,m} = \\left( t_{i+m-1} - M_{T_{i,m}} \\right) + \\left( t_{i-1} - M_{T_{i-1,m}} \\right)$\n", + "4. $dg^{Q}_{i,m} = \\left( q_{i+m-1} - \\mu_{Q_{i,m}} \\right) + \\left( q_{i-1} - \\mu_{Q_{i-1,m}} \\right)$\n", + "\n", + "Then our equation simply becomes:\n", + "\n", + "\\begin{align}\n", + " \\overline{Q_{i,m}T_{i,m}}\n", + " &=\n", + " \\overline{Q_{i-1,m}T_{i-1,m}}\n", + " + df^{T}_{i,m} dg^{Q}_{i,m}\n", + " + df^{Q}_{i,m} dg^{T}_{i,m}\n", + " \\\\\n", + "\\end{align}" ] }, { - "cell_type": "code", - "execution_count": null, + "cell_type": "markdown", "metadata": {}, - "outputs": [], - "source": [] + "source": [ + "Finally, we can compute the covariance via $cov(Q_{i,m},T_{i,m}) = \\frac{\\overline{Q_{i,m}T_{i,m}}}{m}$ and its corresponding Pearson correlation, $\\rho$ via:\n", + "\n", + "\\begin{align}\n", + " \\rho(Q_{i,m},T_{i,m}) \n", + " &= \n", + " \\frac{cov(Q_{i,m},T_{i,m})}{\\sigma_{Q_{i,m}}\\Sigma_{T_{i,m}}}\n", + " \\\\\n", + " \\\\\n", + " &=\n", + " \\frac{\\overline{Q_{i,m}T_{i,m}}}{m \\sigma_{Q_{i,m}} \\Sigma_{T_{i,m}}}\n", + " \\\\\n", + "\\end{align}\n", + "\n", + "However, recall that the (population) variances are simply:\n", + "\n", + "\\begin{align}\n", + " \\sigma^{2}_{Q_{i,m}} &= \\frac{1}{m} \\sum \\limits_{0 \\le j \\lt m} \\left( q_{i+j} - \\mu_{Q_{i,m}} \\right) ^{2}\n", + " \\\\\n", + " m \\sigma^{2}_{Q_{i,m}} &= \\sum \\limits_{0 \\le j \\lt m} \\left( q_{i+j} - \\mu_{Q_{i,m}} \\right) ^{2}\n", + " \\\\\n", + " \\sqrt{m \\sigma^{2}_{Q_{i,m}}} &= \\sqrt{\\sum \\limits_{0 \\le j \\lt m} \\left( q_{i+j} - \\mu_{Q_{i,m}} \\right) ^{2}}\n", + " \\\\\n", + " \\sqrt{m} \\sigma_{Q_{i,m}} &= \\left\\| Q_{i,m} - \\mu_{Q_{i,m}} \\right\\|\n", + " \\\\\n", + " \\\\ \\text{and}\n", + " \\\\\n", + " \\\\\n", + " \\Sigma^{2}_{T_{i,m}} &= \\frac{1}{m} \\sum \\limits_{0 \\le j \\lt m} \\left( t_{i+j} - M_{T_{i,m}} \\right) ^{2}\n", + " \\\\\n", + " m \\Sigma^{2}_{T_{i,m}} &= \\sum \\limits_{0 \\le j \\lt m} \\left( t_{i+j} - M_{T_{i,m}} \\right) ^{2}\n", + " \\\\\n", + " \\sqrt{m \\Sigma^{2}_{T_{i,m}}} &= \\sqrt{\\sum \\limits_{0 \\le j \\lt m} \\left( t_{i+j} - M_{T_{i,m}} \\right) ^{2}}\n", + " \\\\\n", + " \\sqrt{m} \\Sigma_{T_{i,m}} &= \\left\\| T_{i,m} - M_{T_{i,m}} \\right\\|\n", + "\\end{align}\n", + "\n", + "\n", + "Thus, this can be used to transform the denominator, $m \\sigma_{Q_{i,m}} \\Sigma_{T_{i,m}}$, in $\\rho(Q_{i,m},T_{i,m}) $:" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\\begin{align}\n", + " \\rho(Q_{i,m},T_{i,m}) \n", + " &=\n", + " \\frac{\\overline{Q_{i,m}T_{i,m}}}{m \\sigma_{Q_{i,m}} \\Sigma_{T_{i,m}}}\n", + " \\\\\n", + " \\\\\n", + " &=\n", + " \\frac{\\overline{Q_{i,m}T_{i,m}}}{\\sqrt{m} \\sigma_{Q_{i,m}} \\cdot \\sqrt{m} \\Sigma_{T_{i,m}}}\n", + " \\\\\n", + " \\\\\n", + " &=\n", + " \\frac{\\overline{Q_{i,m}T_{i,m}}}{\\left\\| Q_{i,m} - \\mu_{Q_{i,m}} \\right\\| \\left\\| T_{i,m} - M_{T_{i,m}} \\right\\|}\n", + " \\\\\n", + " \\\\\n", + " &=\n", + " \\overline{Q_{i,m}T_{i,m}} \\cdot\n", + " \\frac{1}{ \\left\\| Q_{i,m} - \\mu_{Q_{i,m}} \\right\\| } \\cdot\n", + " \\frac{1}{ \\left\\| T_{i,m} - M_{T_{i,m}} \\right\\|}\n", + " \\\\\n", + "\\end{align}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Notice that both inverse terms can be precomputed. Finally, the distance can be computed via $D(Q_{i,m},T_{i,m}) = \\sqrt{2m \\left[1 - \\rho(Q_{i,m}, T_{i,m}) \\right]}$." + ] } ], "metadata": { @@ -859,7 +1528,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.9.10" + "version": "3.14.6" } }, "nbformat": 4,