Time Series Anomaly Detection With Isolation Forest

Isolation Forest takes a table of rows and finds the rows that are easy to separate from the rest. It does not know that your rows are ordered, and every useful thing about applying it to a time series follows from arranging for the ordering to be in the columns instead.

What Isolation Forest can and cannot see

The algorithm, from Liu, Ting and Zhou’s 2008 ICDM paper, builds a forest of random trees. Each tree picks a random feature and a random split point inside that feature’s observed range, recursively, until points are isolated. Points in sparse regions get isolated after few splits; points in dense regions need many. The anomaly score is built from the average path length across the forest, and no distance or density is ever computed, which is why it scales to large tables.

Fed a single column of raw values, it will therefore find the globally extreme values, which is what a max and min already told you. It cannot find a contextual anomaly: a value of 50 that is unremarkable in December and impossible in June, a flat line at the daily mean where there should be a daily cycle, a spike that is within the global range but four standard deviations from the local level. Those are the anomalies worth detecting in a time series, and reaching them means giving the model columns in which they are extreme.

The score, and what 0.5 means

The score is s(x, n) = 2^(-E(h(x)) / c(n)), where E(h(x)) is the mean path length over the trees and c(n) normalises by the average path length of an unsuccessful search in a binary search tree of n points: c(n) = 2H(n−1) − 2(n−1)/n, with H the harmonic number, approximately ln(i) + 0.5772.

for the scikit-learn default sample size n = 256:

  H(255)  = ln(255) + 0.5772 = 5.5413 + 0.5772 = 6.1185
  c(256)  = 2(6.1185) - 2(255/256) = 12.2370 - 1.9922 = 10.2448

  a point isolated in 5 splits on average:
    s = 2^(-5 / 10.2448) = 2^(-0.4880) = 0.713   -> anomalous
  a point with average path length equal to c(n):
    s = 2^(-1) = 0.5                             -> the neutral value

So 0.5 is the pivot: shorter than average path means a score above 0.5 and an anomaly, longer means below. scikit-learn negates this, so its score_samples returns lower numbers for more anomalous points, and decision_function is score_samples minus an offset which is −0.5 when contamination is "auto", making the sign the decision: negative is an outlier, positive an inlier. If you set contamination to a number instead, the offset is moved so that exactly that fraction is flagged, which means you are asserting the anomaly rate rather than discovering it.

Turning the series into a table

Every feature below exists to make one kind of anomaly extreme in at least one column. Build them from a trailing window only, never from the whole series, or the features at time t will contain information from after t.

• Seasonal deviation. The value minus the typical value for its position in the cycle. This is what makes a June value of 50 anomalous while December’s is not. The remainder from an STL decomposition is the more careful version of the same column.
• Local z-score. The value minus a rolling mean, divided by a rolling standard deviation. Catches level shifts that are within the global range.
• First difference. Catches jumps regardless of level. A step change produces one extreme difference, then normality; a spike produces two, in opposite directions.
• Rolling standard deviation. Catches the flat-line failure, where a sensor freezes and the value is perfectly plausible but stops moving. This one is often the only column that fires on the most operationally serious fault.

One column is missing from that list on purpose: the raw value. Leaving it in gives the forest a feature in which every seasonal peak is extreme, and the resulting flags cluster on the top of the daily cycle every single day. If you want a magnitude check, express it as a deviation or a z-score so that the column’s extremes are the ones you actually consider abnormal. The general rule is that Isolation Forest will faithfully find whatever is rare in the columns you give it, so every column is a statement about what you think unusual means.

The script

1. Install numpy, pandas and scikit-learn.
2. Save the script below. It synthesises sixty days of hourly data with a daily and a weekly cycle, injects three spikes and one fifteen-hour flat stretch, and scores the result.
3. Run it. It prints the ten lowest-scoring timestamps and the number of rows flagged; compare that list against the indices it injected, which it also prints.
4. Remove the feature columns one at a time and re-run to see which injected fault each column was responsible for catching. The flat stretch is the interesting one.

import numpy as np
import pandas as pd
from sklearn.ensemble import IsolationForest

rng = np.random.default_rng(0)
n = 24 * 60                                   # 60 days, hourly
t = np.arange(n)
y = (100.0
     + 10 * np.sin(2 * np.pi * t / 24)        # daily cycle
     + 4 * np.sin(2 * np.pi * t / 168)        # weekly cycle
     + rng.normal(0, 1.5, n))

spikes = [300, 700, 1180]
y[spikes] += 18                               # inside the global range
flat = slice(900, 915)
y[flat] = 100.0                               # sensor freeze, plausible value

df = pd.DataFrame(
    {"y": y},
    index=pd.date_range("2026-01-01", periods=n, freq="h"),
)

# the ordering has to live in the columns, because the model ignores rows' order
df["hour"] = df.index.hour
df["dev"] = df["y"] - df.groupby("hour")["y"].transform("median")
df["diff1"] = df["y"].diff()
roll = df["y"].rolling(24, min_periods=24)
df["z"] = (df["y"] - roll.mean()) / roll.std()
df["rstd"] = roll.std()

feat = df[["dev", "diff1", "z", "rstd"]].dropna()

clf = IsolationForest(
    n_estimators=200,
    max_samples=256,          # the c(n) worked above assumes this
    contamination="auto",     # offset_ = -0.5; do not assert a rate you don't know
    random_state=0,
)
clf.fit(feat)

out = feat.assign(score=clf.score_samples(feat), flag=clf.predict(feat))
print("injected spikes at rows:", spikes, "flat stretch:", flat)
print(out.nsmallest(10, "score")[["dev", "z", "rstd", "score"]])
print("flagged:", int((out["flag"] == -1).sum()), "of", len(out))

Reading and tuning the output

The ranked list is the product, not the flags. contamination set to "auto" uses the fixed −0.5 offset, which flags whatever falls below it; that count is a property of the data and the features and is not a rate you chose. Ranking by score_samples and taking as many as a human can review is almost always the more useful operating mode, and it degrades gracefully when the true anomaly rate changes.

• max_samples is the main sensitivity knob. The default of 256 is small on purpose: the original paper’s argument is that subsampling reduces swamping and masking, where dense normal regions hide anomalies from each other. Raising it makes the trees deeper and often makes detection worse, which surprises people who expect more data to help.
• Scale does not matter but range does. Splits are drawn uniformly within each feature’s observed range, so a column with one enormous outlier has most of its range in empty space, and random splits there isolate quickly. That is the intended behaviour, but it means a badly scaled or heavy-tailed column can dominate the forest. Check which column is driving a flag before trusting it.
• Consecutive flags are one event. A fifteen-hour flat stretch produces up to fifteen flagged rows. Group adjacent flags into intervals before counting, or your alert volume is a function of your sampling rate.
• It will not tell you a level shift is permanent. Isolation Forest sees a regime change as a run of unusual rows and then, as the rolling features adapt, as nothing at all. If the question is where the process changed rather than which points are odd, that is change point detection instead.