Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It only takes a minute to sign up.
Sign up to join this community
How can I compute the standard deviation in an incremental way (using the new value and the last computed mean and/or std deviation) ?
for the non incremental way, I just do something like: $$S_N=\sqrt{\frac1N\sum_{i=1}^N(x_i-\overline{x})^2}.$$
mean = Mean(list)
for i = 0 to list.size
stdev = stdev + (list[i] - mean)^2
stdev = sqrRoot( stdev / list.size )
I think the easiest way to do this is with an orthogonality trick. I'll show how to incrementally compute the variance instead of the standard deviation. Let $X_1, X_2, ...$ be an iid sequence of random variables with $\bar X = n^{-1} \sum_{j = 1} ^ n X_j$ and $s^2_n$ defined similarly as the $n$'th sample variance (I use a denominator of $n-1$ instead of $n$ in your picture to keep things unbiased for the variance, but you can use the same argument just adding $1$ to all the terms with $n$ in them). First write $$ s^2_n = \frac{\sum_{j = 1} ^ n (X_j - \bar X_n)^2}{n - 1} = \frac{\sum_{j = 1} ^ n (X_j - \bar X_{n - 1} + \bar X_{n - 1} - \bar X_n)^2}{n - 1}. $$ Expand this to get $$ s^2_n = \frac{(n - 2)s^2_{n - 1} + (n - 1) (\bar X_{n - 1} - \bar X_n)^2 + 2 \sum_{j = 1} ^ {n - 1} (X_j - \bar X_{n - 1})(\bar X_{n - 1} - \bar X_n) + (X_n - \bar X_{n})^2}{n - 1} $$ and it is easy to show that the summation term above is equal to $0$ which gives $$ s^2_n = \frac{(n - 2)s^2_{n - 1} + (n - 1)(\bar X_{n - 1} - \bar X_n)^2 + (X_n - \bar X_{n})^2}{n - 1}. $$
EDIT: I assumed you already have an incremental expression for the sample mean. It is much easier to get that: $\bar X_n = n^{-1}[X_n + (n-1)\bar X_{n-1}]$.
The standard deviation is a function of the totals $T_{\alpha}=\sum_{i=1}^{N}x_{i}^{\alpha}$ for $\alpha=0,1,2$, each of which can be calculated incrementally in an obvious way. In particular, $E[X]=T_{1}/T_{0}$ and $E[X^2]=T_{2}/T_{0}$, and the standard deviation is $$ \sigma = \sqrt{\text{Var}[X]} = \sqrt{E[X^2]-E[X]^2} = \frac{1}{T_0}\sqrt{T_{0}T_{2}-T_{1}^2}. $$ By maintaining totals of higher powers ($T_{\alpha}$ for $\alpha \ge 3$), you can derive similar "incremental" expressions for the skewness, kurtosis, and so on.
Another way of saying the above is that you need to keep a count, an incremental sum of values and an incremental sum of squares.
Let $N$ be the count of values seen so far, $S = \sum_{1,N} x_i$ and $Q = \sum_{1,N} x_i^2$ (where both $S$ and $Q$ are maintained incrementally).
Then any stage,
the mean is $\frac{S}{N}$
and the variance is $\frac{Q}{N} - \left( \frac{S}{N}\right)^2$
An advantage of this method is that it is less prone to rounding errors after a long series of calculations.
What you refer to as an incremental computation is very close to the computer scientist's notion of an online algorithm. There is in fact a well-known online algorithm for computing the variance and thus square root of a sequence of data, documented here.
Forgive my poor math background, what I need is detail!
I added my progress here for someone like me.
$$ s_n^2=\frac {\sum_{i=1}^{n}(x_i-\bar{x}_n)^2}{n-1} \\ = \frac {\sum_{i=1}^n(x_i - \bar{x}_{n-1} + \bar{x}_{n-1} - \bar{x}_n)^2}{n-1} \\ = \frac {\sum_{i=1}^{n}(x_i - \bar{x}_{n-1})^2 + 2\sum_{i=1}^n(x_i - \bar{x}_{n-1})(\bar{x}_{n-1} - \bar{x}_n) + \sum_{i=1}^n(\bar{x}_{n-1} - \bar{x}_n)^2} {n-1} \\ = \frac {(\sum_{i=1}^{n-1}(x_i - \bar{x}_{n-1})^2 + (x_n - \bar{x}_{n-1})^2) + (2\sum_{i=1}^{n-1}(x_i - \bar{x}_{n-1})(\bar{x}_{n-1} - \bar{x}_n) + 2(x_n - \bar{x}_{n-1})(\bar{x}_{n-1} - \bar{x}_n)) + \sum_{i=1}^n(\bar{x}_{n-1} - \bar{x}_n)^2} {n-1} \\ = \frac {(n-2)s_{n-1}^2 + (x_n - \bar{x}_{n-1})^2 + 0 + 2(x_n - \bar{x}_{n-1})(\bar{x}_{n-1} - \bar{x}_n)) + n(\bar{x}_{n-1} - \bar{x}_n)^2} {n-1} \\ = \frac {(n-2)s_{n-1}^2 + x_n^2 - 2x_n\bar{x}_{n-1} + \bar{x}_{n-1}^2 + 2x_n \bar{x}_{n-1} - 2x_n \bar{x}_n - 2\bar{x}_{n-1}^2 + 2\bar{x}_{n-1}\bar{x}_n + n\bar{x}_{n-1}^2 - 2n\bar{x}_{n-1}\bar{x}_n + n\bar{x}_n^2} {n-1} \\ = \frac {(n-2)s_{n-1}^2 + x_n^2 + \bar{x}_{n-1}^2 - 2x_n \bar{x}_n - 2\bar{x}_{n-1}^2 + 2\bar{x}_{n-1}\bar{x}_n + n\bar{x}_{n-1}^2 - 2n\bar{x}_{n-1}\bar{x}_n + n\bar{x}_n^2} {n-1} \\ = \frac {(n-2)s_{n-1}^2 + x_n^2 - 2x_n\bar{x}_{n-1} + \bar{x}_{n-1}^2 + 2x_n \bar{x}_{n-1} - 2x_n \bar{x}_n - 2\bar{x}_{n-1}^2 + 2\bar{x}_{n-1}\bar{x}_n + n\bar{x}_{n-1}^2 - 2n\bar{x}_{n-1}\bar{x}_n + n\bar{x}_n^2} {n-1} \\ = \frac {(n-2)s_{n-1}^2 + (x_n^2 - 2x_n\bar{x}_n + \bar{x}_n^2) + (n-1)(\bar{x}_{n-1}^2 - 2\bar{x}_{n-1}\bar{x}_n + \bar{x}_n^2)} {n-1} \\ = \frac {(n-2)s_{n-1}^2 + (n-1)(\bar{x}_{n-1} - \bar{x}_n)^2 + (x_n - \bar{x}_n)^2} {n-1} \\ = \frac {n-2}{n-1}s_{n-1}^2 + (\bar{x}_{n-1} - \bar{x}_n)^2 + \frac {(x_n - \bar{x}_n)^2}{n-1} $$
and
$$ (\bar{x}_{n-1} - \bar{x}_n)^2 + \frac {(x_n - \bar{x}_n)^2}{n-1} \\ = (\bar{x}_{n-1} - \frac {x_n + (n-1)\bar{x}_{n-1}}{n})^2 + \frac {(x_n - \frac {x_n + (n-1)\bar{x}_{n-1}}{n})^2}{n-1} \\ = \frac {1}{n} (x_n - \bar{x}_{n-1})^2 $$
so
$$ s_n^2 = \frac{n-2}{n-1}s_{n-1}^2 + \frac{1}{n}(x_n - \bar{x}_{n-1})^2 $$
In case someone has to "decrement" and not only "increment" the standard deviation $\sigma$ (for example, when a result $x_i$ in the set is incorrect and needs to be removed or recalculated), you can use this formula:
$ \sigma_{\text{without } x_i} = \sqrt{\frac{n}{n -1} \left[ \sigma_n^2 - \frac{(\bar{x}_n-x_i)^2}{n-1} \right]} $
Here is the derivation:
\begin{equation} \label{varianceDecrementale} \begin{split} \sigma^2_{\text{without } x_i} & = \sum_{j\in \{0, 1, \cdots, i-1, i+1, \cdots, n\}} \frac{x_j^2}{n-1} - \bar{x}_{n \text{ without } x_i}^2 \\ & = \frac{x_0^2 + x_1^2 + \cdots + x_{i - 1}^2 + x_{i + 1}^2 + \cdots + x_n^2}{n - 1} - \bar{x}_{n \text{ without } x_i}^2 \\ & = \frac{n}{n -1} \left[ \frac{x_0^2 + x_1^2 + \cdots + x_{i - 1}^2 + x_{i}^2 + x_{i + 1}^2 + \cdots + x_n^2}{n} - \frac{x_{i}^2}{n}\right] - \bar{x}_{n \text{ without } x_i}^2 \\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} \right] - \bar{x}_{n \text{ without } x_i}^2 \\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} \right] - \frac{n^2}{(n-1)^2} \left[ \bar{x}_n - \frac{x_i}{n} \right]^2 \\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} - \frac{n}{n-1} \left[ \bar{x}_n - \frac{x_i}{n} \right]^2 \right] \\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} - \frac{n}{n-1} \left[ \bar{x}_n^2 - 2 \bar{x}_n \frac{x_i}{n} + \frac{x_i^2}{n^2} \right] \right] \\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} - \frac{n}{n-1} \bar{x}_n^2 + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right] \end{split} \end{equation}
We can remark that
\begin{equation} \begin{split} -\frac{n}{n-1} \bar{x}_n^2 & = -\frac{n}{n-1} \bar{x}_n^2 + \bar{x}_n^2 - \bar{x}_n^2 \\ & = -\bar{x}_n^2 - \frac{n}{n-1} \bar{x}_n^2 + \frac{n - 1}{n - 1} \bar{x}_n^2 \\ & = -\bar{x}_n^2 + \frac{-n \bar{x}_n^2 + n \bar{x}_n^2 - \bar{x}_n^2}{n-1} \\ & = - \left (\bar{x}_n^2 + \frac{\bar{x}_n^2}{n-1}\right) \end{split} \end{equation}
thus
\begin{equation} \begin{split} \sigma^2_{\text{without } x_i} & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} - \frac{n}{n-1} \bar{x}_n^2 + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right] \\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \frac{x_{i}^2}{n} - \left (\bar{x}_n^2 + \frac{\bar{x}_n^2}{n-1}\right) + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right]\\ & = \frac{n}{n -1} \left[ \sum_i^n \frac{x_i^2}{n} - \bar{x}_n^2 - \frac{x_{i}^2}{n} - \frac{\bar{x}_n^2}{n-1} + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right]\\ & = \frac{n}{n -1} \left[ \sigma^2_n- \frac{x_{i}^2}{n} - \frac{\bar{x}_n^2}{n-1} + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right] \\ & = \frac{n}{n -1} \left[ \sigma^2_n - \frac{x_{i}^2}{n} - \frac{\bar{x}_n^2}{n-1} + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right] \\ & = \frac{n}{n -1} \left[ \sigma^2_n - \frac{x_{i}^2(n-1)}{n(n-1)} - \frac{\bar{x}_n^2}{n-1} + 2 \bar{x}_n \frac{x_i}{n - 1} - \frac{x_i^2}{n (n - 1)} \right] \\ & = \frac{n}{n -1} \left[ \sigma^2_n - \frac{x_{i}^2{n}}{{n}(n-1)} + {\frac{r^2}{n(n-1)}} - \frac{\bar{x}_n^2}{n-1} + 2 \bar{x}_n \frac{x_i}{n - 1} - {\frac{x_i^2}{n (n - 1)}} \right] \\ & = \frac{n}{n -1} \left[ \sigma^2_n - \frac{x_i^2}{n-1} - \frac{\bar{x}_n^2}{n-1} + 2\bar{x}_n\frac{x_i}{n-1} \right]\\ & = \frac{n}{n -1} \left[ \sigma^2_n - \frac{(\bar{x}_n-x_i)^2}{n-1} \right] \end{split} \end{equation}
If it is of any help I found what seems to be a much nicer way here wikipedia - online algorithm .
The following is how I have implemented it in java. It has work very well for me. Hope it is clear enough, Please feel free to ask me questions about it.
/**
* standardDeviation() - designed to calculate the standard deviation of a data set incrementally by taking the last entered value and the previous sum of differences to the mean recorded.
* (i.e; upon adding a value to the data set this function should immediately be called)
*
* NOTE: do not call this function if the data set size it less than 2 since standard deviation cannot be calculated on a single value
* NOTE: sum_avg, sum_sd and avg are all static variables
* NOTE: on attempting to use this on another set following previous use, the static values will have to be reset**
*
* @param vector - List<Double> - data with only one additional value from previous method call
* @return updated value for the Standard deviation
*/
public static double standardDeviation(List<Double> vector)
{
double N = (double) vector.size(); //size of the data set
double oldavg = avg; //save the old average
avg = updateAverage(vector); //update the new average
if(N==2.0) //if there are only two, we calculate the standard deviation using the standard formula
{
for(double d:vector) //cycle through the 2 elements of the data set - there is no need to use a loop here, the set is quite small to just do manually
{
sum_sd += Math.pow((Math.abs(d)-avg), 2); //sum the following according to the formula
}
}
else if(N>2) //once we have calculated the base sum_std
{
double newel = (vector.get(vector.size()-1)); //get the latest addition to the data set
sum_sd = sum_sd + (newel - oldavg)*(newel-avg); //https://en.wikipedia.org/wiki/Algorithms_for_calculating_variance#Online_algorithm
}
return Math.sqrt((sum_sd)/(N)); //N or N-1 depends on your choice of sample of population standard deviation
}
/**
* simplistic method for incrementally calculating the mean of a data set
*
* @param vector - List<Double> - data with only one additional value from previous method call
* @return updated value for the mean of the given data set
*/
public static double updateAverage(List<Double> vector)
{
if(vector.size()==2){
sum_avg = vector.get(vector.size()-1) + vector.get(vector.size()-2);
}
else{
sum_avg += vector.get(vector.size()-1);
}
return sum_avg/(double)vector.size();
}
I think it could be done in easier form
$$s_n^2 = \frac{\sum_{i=1}^n(x_i - \bar{x}_n)}{n-1}^2 = \frac{\sum_{i=1}^n(x_i^2-2x_i\bar{x}_n + \bar{x}_n^2)}{n-1} $$
expand numerator
$$\sum_{i=1}^n(x_i^2-2x_i\bar{x}_n + \bar{x}_n^2) = \sum_{i=1}^nx_i^2 - \sum_{i=1}^n(2x_i\frac{\sum_{j=1}^nx_j}{n}) + \sum_{i=1}^n(\frac{\sum_{j=1}^nx_j}{n})^2 = \sum_{i=1}^nx_i^2 - 2(\sum_{i=1}^nx_i)(\sum_{j=1}^n\frac{x_j}{n}) + n(\sum_{j=1}^n\frac{x_j}{n})^2 = \sum_{i=1}^nx_i^2 - \frac{2}{n}(\sum_{i=1}^nx_i)^2 + \frac{1}{n}(\sum_{j=1}^nx_j)^2$$
by this way
$$s_n^2 = \frac{\sum_{i=1}^n(x_i - \bar{x}_n)}{n-1} = \frac{\sum_{i=1}^nx_i^2 - \frac{1}{n}(\sum_{i=1}^nx_i)^2}{n-1}$$
Here example of how it can be implemented in Python
from math import sqrt
x = int(input())
n = 0
sum = 0
sumsq = 0
while x != 0: # or whenever you want
n += 1
sum += x
sumsq += x^^2
x = int(input())
print(sqrt((sumsq - sum^^2 / n)/(n-1)))
