Online Covariance
by Joshua Burkholder
online_covariance.pdf
online_covariance.docx

Given the following set of two-dimensional inputs:

$\left\{\left(x_1,y_1\right),\left(x_2,y_2\right),\dots ,\left(x_{n-1},y_{n-1}\right),\left(x_n,y_n\right)\right\}$

Let $n$ be the number of two-dimensional inputs, $X$ represent the $x$ dimension, $Y$ represent the $y$ dimension, $Co{v}_n\left(X,Y\right)$ be the biased sample covariance of the $x$ and $y$ dimensions for the first $n$ two-dimensional inputs, $Co{v}_{n-1}\left(X,Y\right)$ be the biased sample covariance of the $x$ and $y$ dimensions for the first $n-1$ two-dimensional inputs, $x_n$ be the $x$ value of the $n$-th two-dimensional input, ${\overline{x}}_n$ be the sample mean of the $x$ values for the first $n$ two-dimensional inputs, $y_n$ be the $y$ value of the $n$-th two-dimensional input, and ${\overline{y}}_{n-1}$ be the sample mean of the $y$ values for the first $n-1$ two-dimensional inputs.  Then, the recurrence equation for the biased sample covariance (a.k.a. online covariance) is:

$Co{v}_n\left(X,Y\right)=Co{v}_{n-1}\left(X,Y\right)-\frac{Co{v}_{n-1}\left(X,Y\right)-\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}$

Note: The recurrence equation above also applies when computing the online covariance matrix:

${\Sigma }_n\left[j,k\right]={\Sigma }_{n-1}\left[j,k\right]-\frac{{\Sigma }_{n-1}\left[j,k\right]-\left(x_n\left[j\right]-{\overline{x\left[j\right]}}_n\right)\left(x_n\left[k\right]-{\overline{x\left[k\right]}}_{n-1}\right)}{n}$.

However, we will restrict ourselves to the online covariance computation of two-dimensional input in this post and explore the online covariance matrix computation of $m$-dimensional input in a later post.

Proof:

The definition of the biased sample covariance of the $x$ and $y$ dimensions for the first $n$ two-dimensional inputs is defined as:

$Co{v}_n\left(X,Y\right)=\frac{{\displaystyle \sum _{i=1}^n\left(x_i-{\overline{x}}_n\right)\left(y_i-{\overline{y}}_n\right)}}{n}$.

If we expand this definition, we have:

$Co{v}_n\left(X,Y\right)=\frac{{\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_n\right)\left(y_i-{\overline{y}}_n\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)}{n}$.

Since the recurrence equations for the sample mean of the $x$ and $y$ values are:

${\overline{x}}_n={\overline{x}}_{n-1}-\frac{{\overline{x}}_{n-1}-x_n}{n}$ and ${\overline{y}}_n={\overline{y}}_{n-1}-\frac{{\overline{y}}_{n-1}-y_n}{n}$,

then we have:

$\begin{array}{l}Co{v}_n\left(X,Y\right)=\frac{{\displaystyle \sum _{i=1}^{n-1}\left(x_i-\left({\overline{x}}_{n-1}-\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\right)\left(y_i-\left({\overline{y}}_{n-1}-\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{{\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_{n-1}+\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(y_i-{\overline{y}}_{n-1}+\frac{{\overline{y}}_{n-1}-y_n}{n}\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{{\displaystyle \sum _{i=1}^{n-1}\left(\begin{array}{l}{x}_i{y}_i-x_i{\overline{y}}_{n-1}+x_i\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)-{\overline{x}}_{n-1}{y}_i+{\overline{x}}_{n-1}{\overline{y}}_{n-1}-{\overline{x}}_{n-1}\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)y_i-\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right){\overline{y}}_{n-1}+\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}{\displaystyle \sum _{i=1}^{n-1}\left(x_i{y}_i-x_i{\overline{y}}_{n-1}-{\overline{x}}_{n-1}{y}_i+{\overline{x}}_{n-1}{\overline{y}}_{n-1}\right)}+{\displaystyle \sum _{i=1}^{n-1}\left(x_i\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)-{\overline{x}}_{n-1}\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\right)}\\ +{\displaystyle \sum _{i=1}^{n-1}\left(\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)y_i-\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right){\overline{y}}_{n-1}+\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}{\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_{n-1}\right)\left(y_i-{\overline{y}}_{n-1}\right)}+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right){\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_{n-1}\right)}\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right){\displaystyle \sum _{i=1}^{n-1}\left(y_i-{\overline{y}}_{n-1}+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\end{array}$

Since the biased sample covariance of the $x$ and $y$ dimensions for the first $n-1$ two-dimensional inputs is defined as:

$Co{v}_{n-1}\left(X,Y\right)=\frac{{\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_{n-1}\right)\left(y_i-{\overline{y}}_{n-1}\right)}}{n-1}$,

then we also have:

${\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_{n-1}\right)\left(y_i-{\overline{y}}_{n-1}\right)}=\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)$.

With this, we have:

$\begin{array}{c}Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right){\displaystyle \sum _{i=1}^{n-1}\left(x_i-{\overline{x}}_{n-1}\right)}\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right){\displaystyle \sum _{i=1}^{n-1}\left(y_i-{\overline{y}}_{n-1}+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\right)}+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left({\displaystyle \sum _{i=1}^{n-1}\left(x_i\right)}+{\displaystyle \sum _{i=1}^{n-1}\left(-{\overline{x}}_{n-1}\right)}\right)\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left({\displaystyle \sum _{i=1}^{n-1}\left(y_i\right)}+{\displaystyle \sum _{i=1}^{n-1}\left(-{\overline{y}}_{n-1}\right)}+{\displaystyle \sum _{i=1}^{n-1}\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)}\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left({\displaystyle \sum _{i=1}^{n-1}\left(x_i\right)}-{\overline{x}}_{n-1}{\displaystyle \sum _{i=1}^{n-1}\left(1\right)}\right)\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left({\displaystyle \sum _{i=1}^{n-1}\left(y_i\right)}-{\overline{y}}_{n-1}{\displaystyle \sum _{i=1}^{n-1}\left(1\right)}+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right){\displaystyle \sum _{i=1}^{n-1}\left(1\right)}\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left({\displaystyle \sum _{i=1}^{n-1}\left(x_i\right)}-{\overline{x}}_{n-1}\left(n-1\right)\right)\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left({\displaystyle \sum _{i=1}^{n-1}\left(y_i\right)}-{\overline{y}}_{n-1}\left(n-1\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left(n-1\right)\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\end{array}$

Since the sample mean for the first $n-1$ $x$ and $y$ values are defined as:

${\overline{x}}_{n-1}=\frac{{\displaystyle \sum _{i=1}^{n-1}{x}_i}}{n-1}$ and ${\overline{y}}_{n-1}=\frac{{\displaystyle \sum _{i=1}^{n-1}{y}_i}}{n-1}$,

then we also have:

${\displaystyle \sum _{i=1}^{n-1}{x}_i}={\overline{x}}_{n-1}\left(n-1\right)$ and ${\displaystyle \sum _{i=1}^{n-1}{y}_i}={\overline{y}}_{n-1}\left(n-1\right)$.

With that, we have:

$\begin{array}{l}Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left({\overline{x}}_{n-1}\left(n-1\right)-{\overline{x}}_{n-1}\left(n-1\right)\right)\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left({\overline{y}}_{n-1}\left(n-1\right)-{\overline{y}}_{n-1}\left(n-1\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left(n-1\right)\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left(\cancel{{\overline{x}}_{n-1}\left(n-1\right)}\cancel{-{\overline{x}}_{n-1}\left(n-1\right)}\right)\\ +\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\cancel{{\overline{y}}_{n-1}\left(n-1\right)}\cancel{-{\overline{y}}_{n-1}\left(n-1\right)}+\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left(n-1\right)\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\left(n-1\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_n\right)\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +x_n{y}_n-x_n{\overline{y}}_n-{\overline{x}}_n{y}_n+{\overline{x}}_n{\overline{y}}_n\end{array}\right)}{n}\end{array}$

Since the recurrence equation for the sample mean of the $y$ values is:

${\overline{y}}_n={\overline{y}}_{n-1}-\frac{{\overline{y}}_{n-1}-y_n}{n}$,

then we have:

$\begin{array}{c}Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +x_n{y}_n-x_n\left({\overline{y}}_{n-1}-\frac{{\overline{y}}_{n-1}-y_n}{n}\right)-{\overline{x}}_n{y}_n+{\overline{x}}_n\left({\overline{y}}_{n-1}-\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +x_n{y}_n-x_n{\overline{y}}_{n-1}+x_n\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)-{\overline{x}}_n{y}_n+{\overline{x}}_n{\overline{y}}_{n-1}-{\overline{x}}_n\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +x_n{y}_n-x_n{\overline{y}}_{n-1}-{\overline{x}}_n{y}_n+{\overline{x}}_n{\overline{y}}_{n-1}+\left(x_n-{\overline{x}}_n\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)+\left(x_n-{\overline{x}}_n\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\end{array}$

Since the recurrence equation for the sample mean of the $x$ values is:

$\begin{array}{l}{\overline{x}}_n={\overline{x}}_{n-1}-\frac{{\overline{x}}_{n-1}-x_n}{n}\\ {\overline{x}}_n=\frac{n{\overline{x}}_{n-1}}{n}+\frac{-{\overline{x}}_{n-1}+x_n}{n}\\ {\overline{x}}_n=\frac{\left(n-1\right){\overline{x}}_{n-1}+x_n}{n},\end{array}$

then we have:

$\begin{array}{l}Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)+\left(x_n-\left(\frac{\left(n-1\right){\overline{x}}_{n-1}+x_n}{n}\right)\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)+\left(\frac{n{x}_n}{n}+\frac{-\left(n-1\right){\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)+\left(\frac{-\left(n-1\right){\overline{x}}_{n-1}+\left(n-1\right)x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\\ +\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)-\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(\begin{array}{l}\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)\cancel{+\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)}\\ +\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)\cancel{-\left(n-1\right)\left(\frac{{\overline{x}}_{n-1}-x_n}{n}\right)\left(\frac{{\overline{y}}_{n-1}-y_n}{n}\right)}\end{array}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{\left(n-1\right)Co{v}_{n-1}\left(X,Y\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{nCo{v}_{n-1}\left(X,Y\right)-Co{v}_{n-1}\left(X,Y\right)+\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}\\ Co{v}_n\left(X,Y\right)=\frac{nCo{v}_{n-1}\left(X,Y\right)}{n}+\frac{-\left(Co{v}_{n-1}\left(X,Y\right)-\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)\right)}{n}\\ Co{v}_n\left(X,Y\right)=Co{v}_{n-1}\left(X,Y\right)-\frac{Co{v}_{n-1}\left(X,Y\right)-\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}\end{array}$

Therefore, the recurrence equation for the biased sample covariance (a.k.a. online covariance) is:

$Co{v}_n\left(X,Y\right)=Co{v}_{n-1}\left(X,Y\right)-\frac{Co{v}_{n-1}\left(X,Y\right)-\left(x_n-{\overline{x}}_n\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}$

Note: We can manipulate this recurrence equation such as that we also have:

$Co{v}_n\left(X,Y\right)=Co{v}_{n-1}\left(X,Y\right)-\frac{Co{v}_{n-1}\left(X,Y\right)-\left(x_n-{\overline{x}}_{n-1}\right)\left(y_n-{\overline{y}}_n\right)}{n}$,

$Co{v}_n\left(X,Y\right)=Co{v}_{n-1}\left(X,Y\right)-\frac{Co{v}_{n-1}\left(X,Y\right)-\left(\frac{n-1}{n}\right)\left(x_n-{\overline{x}}_{n-1}\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}$,

and

$Co{v}_n\left(X,Y\right)=\frac{\left(n-1\right)\left(Co{v}_{n-1}\left(X,Y\right)+\frac{\left(x_n-{\overline{x}}_{n-1}\right)\left(y_n-{\overline{y}}_{n-1}\right)}{n}\right)}{n}$.

Reference:
https://en.wikipedia.org/wiki/Algorithms_for_calculating_variance

Example of C++ code that computes the online covariance:
// Filename: main.cpp
#include <iostream>
#include <iomanip>

int main () {
    
    double x;
    double y;
    double n = 0;
    double mean_x = 0;  // mean of the x values
    double mean_y = 0;  // mean of the y values
    double cov = 0;     // covariance of the x and y values
    double prev_mean_x; // previous mean of the x values
    double prev_mean_y; // previous mean of the y values
    double prev_cov;    // previous covariance of the x and y values
    
    if ( std::cin >> x && std::cin >> y ) {
        ++n;
        mean_x = x;
        mean_y = y;
        cov = 0;
        while ( std::cin >> x && std::cin >> y ) {
            prev_mean_x = mean_x;
            prev_mean_y = mean_y;
            prev_cov = cov;
            ++n;
            mean_x = prev_mean_x - ( prev_mean_x - x ) / n;
            mean_y = prev_mean_y - ( prev_mean_y - y ) / n;
            cov = prev_cov - ( prev_cov - ( x - mean_x ) * ( y - prev_mean_y ) ) / n;
        }
    }
    
    std::cout << "n:      " << n << '\n';
    std::cout << "mean_x: " << std::setprecision( 17 ) << mean_x << '\n';
    std::cout << "mean_y: " << std::setprecision( 17 ) << mean_y << '\n';
    std::cout << "cov:    " << std::setprecision( 17 ) << cov    << '\n';
    
}

Example of data.txt:
-281.189       612.083
974.663        -24.0965
25.8526        401.539
.              .
.              .
.              .

Command Line:
g++ -o main.exe main.cpp -std=c++11 -march=native -O3 -Wall -Wextra -Werror -static
./main.exe < data.txt

Note: Mathematica's Covariance[] function computes the unbiased sample covariance matrix, not the biased sample covariance matrix; therefore, the biased sample covariance matrix is computed in Mathematica as:

( ( Length[ list ] - 1 ) / Length[ list ] ) * Covariance[ list ]