Seasonality: day of the week, month, turn of the month and intraday sessions.

Patterns in the calendar

A calendar effect is a systematic difference in return or volatility tied to a date: the day of the week, the month, the days around the turn of the month, the time of day, the days before a holiday. They are attractive because they are easy to code and easy to explain. They are also the simplest way to find a pattern in the data that does not exist: there are many possible calendars, each with little data, and it is enough to test a few of them for one to come out “significant”.

This lesson always follows the same flow:

flowchart TB
  H[Prior hypothesis] --> D[Data and calendar]
  D --> T[t-test or bootstrap]
  T --> M[Multiple-testing correction]
  M --> F[Out of sample]
  F --> C[Costs and decision]

The data is ISX and FIX, the two daily indices of level 1 (2019-2025), the 57 stocks of the universe from lesson 1.4.4 and the NOVA 5-minute bars for 2024 and 2025. A warning first, as in lesson 1.6.2: this is synthetic data and its generator does not include day-of-the-week, month or turn-of-the-month effects. Whatever comes out of those tests is chance or a side effect of the model, and that is why they make a good test bench: we know the right answer. The intraday profile of volatility and volume, on the other hand, is in the generator.

Day of the week

The t-statistic of a mean is the mean divided by its standard error:

With a few hundred data points, |t| > 1.96 corresponds to p < 0.05 (two-tailed).

import numpy as np
import pandas as pd
from scipy import stats

isx = pd.read_csv("etf_indice_diario.csv", parse_dates=["fecha"],
                  index_col="fecha")["indice"]
fix = pd.read_csv("fix_indice_diario.csv", parse_dates=["fecha"],
                  index_col="fecha")["indice"]
px = pd.DataFrame({"ISX": isx, "FIX": fix})
r = np.log(px).diff().dropna()
DIAS = ["Mon", "Tue", "Wed", "Thu", "Fri"]

def por_grupo(x, grupos):
    g = x.groupby(grupos)
    out = pd.DataFrame({"n": g.size(), "media_%": 100 * g.mean()})
    out["t"] = g.mean() / (g.std() / np.sqrt(g.size()))
    return out

tablas = {c: por_grupo(r[c], r.index.dayofweek) for c in r}
dia = pd.concat(tablas, axis=1)
dia.index = DIAS
dia.round(3)

     ISX                 FIX
       n media_%      t    n media_%      t
Mon  329   0.040  0.527  329   0.082  1.212
Tue  364   0.032  0.414  364  -0.017 -0.231
Wed  360   0.101  1.505  360   0.173  2.669
Thu  353   0.044  0.645  353   0.005  0.080
Fri  353   0.037  0.595  353   0.047  0.725

All the means are small and almost every t is below 2 in absolute value. The exception is FIX on Wednesday: a mean of 0.17% and t = 2.67. Before celebrating, two calculations. There are ten tests (two indices times five days): with ten independent tests at 5%, the probability that at least one comes out significant by chance is 1 − 0.95¹⁰ ≈ 40%. And the test is badly posed.

Against the rest, not against zero

The t in the table tests “the Wednesday mean is zero”. But an index that rises 14% a year has a positive daily mean on every day: the question for a calendar effect is whether Wednesday differs from the other days. With a Welch test (means of two groups with different variances):

def frente_al_resto(x, grupos):
    filas = {}
    for g in np.unique(grupos):
        dentro, fuera = x[grupos == g], x[grupos != g]
        t, p = stats.ttest_ind(dentro, fuera, equal_var=False)
        filas[g] = {"dif_%": 100 * (dentro.mean() - fuera.mean()),
                    "t": t, "p": p}
    return pd.DataFrame(filas).T

res = frente_al_resto(r["FIX"], r.index.dayofweek)
res.index = DIAS
res.round(3)

     dif_%      t      p
Mon  0.030  0.393  0.695
Tue -0.094 -1.173  0.241
Wed  0.145  1.981  0.048
Thu -0.066 -0.890  0.374
Fri -0.013 -0.181  0.856

Against the rest of the days, FIX's Wednesday drops to t = 1.98 and p = 0.048: just under 0.05 and the only one of the five. With the hypothesis properly posed the result weakens but does not disappear, so the question becomes a different one: is a result like this to be expected among all the ones we have looked at? It is answered further down, with many more tests.

Month of the year: one data point per year

For the month effect, each observation is the return of a full month: with seven years there are seven Januaries, seven Februaries… A bootstrap interval (resample the seven values with replacement and recompute the mean thousands of times) shows how little they say:

mens = np.log(px.resample("ME").last()).diff().dropna()["ISX"]
rng = np.random.default_rng(42)
filas = {}
for mes, x in mens.groupby(mens.index.month):
    ic = stats.bootstrap((x.to_numpy(),), np.mean, n_resamples=4000,
                         method="percentile", rng=rng)
    filas[mes] = {"n": len(x), "media_%": 100 * x.mean(),
                  "ic_bajo_%": 100 * ic.confidence_interval.low,
                  "ic_alto_%": 100 * ic.confidence_interval.high,
                  "t": x.mean() / (x.std() / np.sqrt(len(x)))}
meses = pd.DataFrame(filas).T
meses.round(2)

      n  media_%  ic_bajo_%  ic_alto_%     t
1   6.0     2.84      -0.21       5.60  1.73
2   7.0     0.34      -0.87       1.71  0.46
3   7.0    -5.53     -18.92       6.68 -0.74
4   7.0     0.45      -4.68       5.86  0.15
5   7.0     5.19      -0.04      10.90  1.65
6   7.0     5.97       1.82      10.61  2.45
7   7.0     0.16      -4.76       4.91  0.06
8   7.0    -0.88      -3.68       2.03 -0.56
9   7.0    -0.44      -3.83       2.52 -0.26
10  7.0     2.90       0.69       5.40  2.21
11  7.0     1.01      -0.73       2.97  0.99
12  7.0     1.10      -0.84       3.06  1.02

The means range from −5.5% in March to +6.0% in June and the intervals are huge: March's goes from −18.9% to +6.7%. March is negative because of the 2020 crash, a single year. June and October come out with t above 2 and with intervals that exclude zero, but with seven data points each and twelve months tested. In this data, moreover, we know why: the generator places the big rallies and sell-offs on specific dates (March 2020, spring 2025), and that leaves a mark on the months in which they fell. It is an effect of history, not of the calendar. Real data behaves the same way: most “month effects” in an asset with ten or twenty years of history are the mark of a handful of events.

The turn of the month

The turn-of-the-month effect says that returns concentrate on the last trading day of a month and the first days of the next. A common definition: the last day and the first three.

mes = r.index.to_period("M")
orden = pd.Series(1, index=r.index).groupby(mes).cumsum()
total = pd.Series(1, index=r.index).groupby(mes).transform("sum")
cambio = (orden <= 3) | (orden == total)
print(f"Turn-of-the-month sessions: {cambio.sum()} of {len(r)}")
for c in r:
    a, b = r.loc[cambio, c], r.loc[~cambio, c]
    t, p = stats.ttest_ind(a, b, equal_var=False)
    print(f"{c}: {100 * a.mean():.3f}% versus "
          f"{100 * b.mean():.3f}%  (t = {t:.2f}, p = {p:.3f})")

Turn-of-the-month sessions: 336 of 1759
ISX: 0.096% versus 0.040%  (t = 0.78, p = 0.434)
FIX: 0.054% versus 0.059%  (t = -0.07, p = 0.947)

In ISX the turn-of-the-month sessions return more on average than the rest, but with t = 0.78; in FIX, practically the same. Nothing. Look at the definition: “last day and first three” is one of many possible ones (two and two, one and four…). Trying several and keeping the best is another way of multiplying the tests, even if it does not show in the final code.

Multiple testing: how many “significant” results chance gives

If there is no effect, the p-value of a well-designed test is uniformly distributed between 0 and 1: 5% of the tests come out below 0.05 by pure chance. With the 57 stocks of the universe and five days of the week there are 285 tests:

from statsmodels.stats.multitest import multipletests

uni = pd.read_parquet("universo_precios.parquet").pivot(
    index="fecha", columns="ticker", values="cierre")
ru = np.log(uni).diff()
filas = []
for tk in ru:
    s = ru[tk].dropna()
    dias = s.index.dayofweek
    for d in range(5):
        t, p = stats.ttest_ind(s[dias == d], s[dias != d],
                               equal_var=False)
        filas.append({"ticker": tk, "dia": DIAS[d], "t": t, "p": p})
pr = pd.DataFrame(filas)
print(f"Tests: {len(pr)}; expected with p < 0.05: "
      f"{0.05 * len(pr):.1f}; observed: {(pr['p'] < 0.05).sum()}")
for metodo in ["bonferroni", "holm", "fdr_bh"]:
    rechaza = multipletests(pr["p"], alpha=0.05, method=metodo)[0]
    print(f"{metodo:>10}: {rechaza.sum()} significant")
pr.sort_values("p").head(3).round(4)

Tests: 285; expected with p < 0.05: 14.2; observed: 10
bonferroni: 0 significant
      holm: 0 significant
    fdr_bh: 0 significant
    ticker  dia       t       p
127   ZXLC  Wed  2.9605  0.0032
20    ZXBI  Mon  2.7957  0.0055
93    ZXHX  Thu -2.7190  0.0067

Ten tests below 0.05, fewer than the fourteen chance would give, and none survives a correction. The histogram of the 285 p-values is flat: the signature of nothing there. The three corrections answer different questions:

If FIX's Wednesdays were part of a search like this, their p = 0.048 would pass none of them. Lesson 2.1.8 comes back to this with the “factor zoo”: the literature on calendar effects and anomalies has the same problem at the scale of thousands of papers.

Intraday sessions

The finest calendar is the clock. With NOVA's 5-minute bars (502 sessions in 2024 and 2025), in New York time, you can compute for each time slot the mean return, its t, the volatility and the volume:

b = pd.read_parquet("nova_barras_5min.parquet")
ny = b["ts"].dt.tz_convert("America/New_York")
b["franja"] = ny.dt.strftime("%H:%M")
b["r"] = np.log(b["cierre"] / b["apertura"])
g = b.groupby("franja")
perfil = pd.DataFrame({
    "n": g.size(), "media_pb": 1e4 * g["r"].mean(),
    "vol_pb": 1e4 * g["r"].std(),
    "volumen_miles": g["volumen"].mean() / 1e3})
perfil["t"] = perfil["media_pb"] / (perfil["vol_pb"]
                                    / np.sqrt(perfil["n"]))
print(f"Slots: {len(perfil)}; |t| > 1.96: "
      f"{(perfil['t'].abs() > 1.96).sum()}")
perfil.iloc[[0, 1, 2, 30, 31, 76, 77]].round(1)

Slots: 78; |t| > 1.96: 2
          n  media_pb  vol_pb  volumen_miles    t
franja
09:30   502       2.1    47.8          351.6  1.0
09:35   502       2.1    42.2          282.4  1.1
09:40   502      -0.1    36.3          236.8 -0.0
12:00   502      -0.4    14.9          115.5 -0.6
12:05   502      -0.1    15.8          118.3 -0.2
15:50   496       0.5    29.2          205.2  0.4
15:55   496       2.2    33.2         1277.8  1.5

The figure at the top of the page draws the two columns that matter. The volatility of the first bar (09:30) is about three times that of midday and the last bar's rises again; volume traces the same U and shoots up in the 15:55 bar, because it includes the closing auction (lesson 1.1.2). The slots from 13:00 onward have 496 sessions instead of 502: on the six half-day sessions of the period the market closes at that time.

The mean returns, on the other hand, are noise: 2 slots of 78 with |t| > 1.96, when chance would give about 4.

perfil.loc[perfil["t"].abs() > 1.96, ["n", "media_pb", "t"]].round(2)

          n  media_pb     t
franja
13:40   496      1.76  2.43
15:35   496     -2.56 -2.07

A rule like “buy at 13:40 and sell at 13:45” comes from looking at 78 slots and keeping the best one. With a mean of 1.8 basis points in 5 minutes and a one-cent bid-ask spread on a stock of about 50 USD (2 basis points round trip on the spread alone), it would not even be profitable if it were real.

What is useful in the intraday profile is its shape. What share of a day's variance is concentrated in each segment?

def tramo(h):
    if h < "10:00":
        return "1 open 9:30-10:00"
    if h >= "15:30":
        return "3 close 15:30-16:00"
    return "2 rest 10:00-15:30"

completas = b.groupby(ny.dt.date)["franja"].transform("size") == 78
bc = b[completas]
var = bc.groupby(bc["franja"].map(tramo))["r"].apply(
    lambda x: (x ** 2).sum())
reparto = pd.DataFrame({"barras": bc.groupby(
    bc["franja"].map(tramo)).size() / completas.sum() * 78,
    "varianza_%": 100 * var / var.sum()})
reparto.round(1)

                     barras  varianza_%
franja
1 open 9:30-10:00       6.0        23.9
2 rest 10:00-15:30     66.0        61.1
3 close 15:30-16:00     6.0        15.0

The six bars of the first half hour and the six of the last (15% of the session's time) account for 38.9% of the day's variance: 23.9% the open and 15.0% the close. For a system this has three direct uses: splitting a large order along the volume profile (as a VWAP algorithm does, lesson 1.7.1), placing wider stops at the open than at midday, and normalizing any intraday signal by the typical volatility of its slot before comparing it with other times of day.

How not to fool yourself with the calendar

• Hypothesis before data. Write down which effect you are going to test, with which definition and over which period, before looking. Lesson 2.2.4 proposes a research journal for that.
• The right test. One day against the rest, one month against the rest; not against zero.
• Count the tests. The ones you ran and the ones you discarded in your head (other definitions, other assets, other periods). Correct for all of them.
• Out of sample and with costs. An effect that survives in a period you did not use to find it and that pays the bid-ask spread and the commissions is a different matter. The site's research on intraday seasonality in futures follows that path all the way to the backtest and its statistics.

Exercises

1. Repeat the 285 tests on the universe comparing each day's mean with zero (stats.ttest_1samp) instead of with the rest of the days. How many come out with p < 0.05 and how many survive Bonferroni? Which stocks are they?

filas = []
for tk in ru:
    s = ru[tk].dropna()
    for d in range(5):
        x = s[s.index.dayofweek == d]
        t, p = stats.ttest_1samp(x, 0.0)
        filas.append({"ticker": tk, "dia": DIAS[d], "t": t, "p": p})
p0 = pd.DataFrame(filas)
bonf = multipletests(p0["p"], alpha=0.05, method="bonferroni")[0]
print(f"p < 0.05: {(p0['p'] < 0.05).sum()};"
      f" after Bonferroni: {bonf.sum()}")
p0[bonf].round(5)

p < 0.05: 23; after Bonferroni: 1
    ticker  dia        t        p
100   ZXIR  Mon -4.27665  0.00006

2. Compute NOVA's volatility by slot using only 2024 and only 2025, and the correlation between the two profiles. Is the shape stable?

ano = ny.dt.year
v24 = b[ano == 2024].groupby("franja")["r"].std()
v25 = b[ano == 2025].groupby("franja")["r"].std()
print(f"Correlation of the profiles: {v24.corr(v25):.3f}")
print(f"Open/midday 2024: {v24['09:30'] / v24['12:00']:.2f}")
print(f"Open/midday 2025: {v25['09:30'] / v25['12:00']:.2f}")

Correlation of the profiles: 0.974
Open/midday 2024: 3.12
Open/midday 2025: 3.26

3. An article claims that “the market goes up in December” and backs it with 30 years of data: 21 positive Decembers out of 30. What else do you need to take it seriously?