# Online algorithm for calculating EWMA at irregular intervals?

What is a fast online algorithm for calculating the EWMA (exponentially weighted moving average) of an input variable observed at irregular intervals?

I know the formula for when sampling at regular intervals:

Calculating alpha from halflife:

$$\alpha = 1 - e^{\frac{\ln{.5}}{H}}$$

Calculating the EWMA of x:

$$E = \alpha\cdot x + (1-\alpha)\cdot E_{-1}$$

What is an algorithm for doing the same where the sampling interval is irregular?

Edit:

I have found an algorithm online, which purports to achieve an irregular EWMA.

double operator()(double x)
{
if (isnan(prev_ewma_)) // we don't decay the first sample
{
prev_ewma_ = x;
prev_time_ = Time::now();
return x;
}

double time_decay = Time::now() - prev_time_;

double alpha = 1 - std::exp(-time_decay / halflife_);
double ewma  = alpha * x + (1 - alpha) * prev_ewma_;

prev_ewma_ = ewma;
prev_time_ = now;
return ewma;
}


Is this algorithm correct?

• Note that what is $\alpha$ in this code is $1-\alpha$ in your post (and vice versa). Other than that it looks OK to me. – noob2 Aug 31 '17 at 18:21
• @noob2 I think that's because alpha in the code is calculated as exp(...), whereas in the formula it is 1 - exp(...) – Steve Lorimer Aug 31 '17 at 18:24
• Besides this code you also need to decide what to do when you need EWMA value between observations. E.g. you observed a couple of values long ago and now need an up to date EWMA value. Do you decay them to zero or not. – LazyCat Aug 31 '17 at 18:44
• @LazyCat good point - I hadn't thought of that – Steve Lorimer Aug 31 '17 at 18:47
• @noob2 I've updated the code to reflect your comment – Steve Lorimer Aug 31 '17 at 18:49

You count the number of 'unit intervals' within that irregular duration between two events and repeat the update function by the count. Amortized time is the average number of unit intervals.

In practical use cases, the 'unit intervals' are larger (subsampling), so this is done in amortized constant time.

• Please could you elaborate - perhaps with a formula or pseudocode? – Steve Lorimer Aug 29 '17 at 19:53
• I presume you're asking this for a low latency use case and the irregular intervals are between, for example, order book events. There's an implicit time unit for your halflife parameter, say it's 1s. If the last time you had 'snapped' a value was 0.14s, you wait till the first event that is >=0.86s later to snap the next value for your update method. Since book events occur fairly frequently, most of the time you're not calling your update method. If a book event arrives 3 seconds later, then 3.14 seconds has elapsed since the last 'snap' and you call your update method 3 times. – madilyn Aug 30 '17 at 15:19
• If you're sensitive to the latency of method calls, there's 2-3 ways I know of that can further cut that down. 1 nontrivially depends on how you've architected the rest of your platform, the other trivially depends on the desired public interface of your EMA. – madilyn Aug 30 '17 at 15:27
• ^For clarity: 0.14s* ago. – madilyn Aug 30 '17 at 16:26
• If the order-book changes at 0.14s, I want the EWMA to reflect that latest "state-of-the-world" at that time, as it is off this event (the change in market) that my algorithm will react. As such, I don't want to use the "old" EWMA, as it won't have the most up-to-date information (the change that has just happened); what I want is to decay the new observation at "0.14s worth of the half-life" and fold that into the average. Not sure if I'm making sense?! – Steve Lorimer Aug 30 '17 at 16:31

The above code for an irregular EWMA doesn't quite give a half-life - the code is missing the $$e^{\ln(.5)}$$ term found in the preceding formula. To get a true half-life, the code should look like this:

double operator()(double x)
{
if (isnan(prev_ewma_)) // we don't decay the first sample
{
prev_ewma_ = x;
prev_time_ = Time::now();
return x;
}

double time_decay = Time::now() - prev_time_;

double alpha = 1 - std::exp(std::log(0.5) * time_decay / halflife_);
double ewma  = alpha * x + (1 - alpha) * prev_ewma_;

prev_ewma_ = ewma;
prev_time_ = now;
return ewma;
}


To show that this works, we look at how alpha is used in the EWMA formula.

$$E = \alpha \cdot x + (1-\alpha) \cdot E_{-1}$$

We expect that after the half-life has elapsed, exactly half of $$E_{-1}$$ will remain in our filtered value. For each filtering timestep, the remaining value of $$E_{-1}$$ will be multiplied by $$1-\alpha$$, meaning we want to solve for $$\alpha$$ such that

$$0.5 = (1-\alpha)^N$$

where $$N$$ is the the number of samples we filter on during our half-life. For a fixed timestep $$dt$$ (time_decay in the code), we calculate $$N$$ as

$$N = \frac{H}{dt}$$

where $$H$$ is the half-life. This gives us

$$0.5 = (1-\alpha)^{\frac{H}{dt}}$$

Plugging in our new formula for $$\alpha$$:

$$\alpha = 1 - e^{\ln(0.5) \cdot \frac {dt} {H}}$$

$$(e^{\ln(0.5)\cdot \frac {dt} {H}})^{\frac{H}{dt}}$$

$$e^{\ln(0.5)} = 0.5$$

Exactly half of the original value will remain.