--- title: "3.2.6 Seasonality of futures spreads" description: "Multi-year seasonal analysis of futures spreads with fixed contracts, annual harmonic, Sharpe standard error and the traps: continuous series and chance." url: "https://pythonparatrading.com/en/python-trading-course/level-3/seasonality-of-futures-spreads/" lang: "en" image: "https://pythonparatrading.com/og/en/python-trading-course/level-3/seasonality-of-futures-spreads.png" --- # 3.2.6 · Seasonality of spreads: multi-year analysis, patterns and traps > Commodity spreads have a reputation for being seasonal: the gasoline crack before summer, gas before winter, soybeans after the harvest. This lesson measures seasonality in nine spreads from the data, five with a planted pattern and four controls with none, and shows why with 16 years of data a control can look better than a real effect. - Course: [Level 3: Derivatives: futures, spreads, forex and options](https://pythonparatrading.com/en/python-trading-course/level-3/) - Module: 3.2 · Futures spreads - Duration: 70 min - Prerequisites: [1.4.6 · Continuous futures series: roll rules, difference adjustment and ratio adjustment](https://pythonparatrading.com/en/python-trading-course/level-1/continuous-futures-series-and-roll-adjustment/), [1.6.3 · Seasonality: day of the week, month, turn of the month and intraday sessions](https://pythonparatrading.com/en/python-trading-course/level-1/seasonality-and-calendar-effects/), [2.1.8 · The multiple testing problem: Bonferroni, FDR and the “factor zoo”](https://pythonparatrading.com/en/python-trading-course/level-2/multiple-testing-bonferroni-fdr-factor-zoo/), [2.5.4 · Deflated Sharpe Ratio and Probability of Backtest Overfitting (PBO)](https://pythonparatrading.com/en/python-trading-course/level-2/deflated-sharpe-ratio-and-pbo/), [3.2.2 · Intercommodity spreads: crack, spark and crush](https://pythonparatrading.com/en/python-trading-course/level-3/intercommodity-spreads-crack-spark-crush/), [3.2.5 · Spread data and backtesting: listed, synthetic and legged](https://pythonparatrading.com/en/python-trading-course/level-3/spread-data-and-legged-backtesting/) ## By the end of this lesson - Build the daily change of a spread with fixed contracts, without roll jumps, for two- and three-leg baskets and for calendars. - Separate the seasonal shape of the curve (priced in) from a seasonal drift (the part that can be traded). - Measure seasonality with an annual harmonic and its F-test, and compare it with the known answer in planted spreads and in controls. - Compute the standard error of a Sharpe ratio with few years, apply multiple-testing corrections and validate the calendar out of sample. - Recognize the traps of seasonal analysis: the unadjusted continuous series, the chart of averages, the calendar chosen with the whole sample and the chance of 16 observations. ## Two different seasonalities That the gasoline crack “rises before summer” or that winter gas “is expensive” are statements with two meanings that are worth separating from the start. - **The seasonal shape of the curve.** Gasoline contracts for summer delivery are worth more than winter ones; January gas is worth more than October gas. It is a difference in **level** between delivery months, known to everyone and priced in. The U.S. Energy Information Administration explains the physical reasons: refining margins vary by season and gasoline demand usually rises in summer; gas storage covers the winter demand peaks (links in the resources). Buying the cheap contract and selling the expensive one does not earn that difference: each of the two contracts converges to its own delivery price. - **The seasonal drift.** What can be traded is something else: that the spread of **two fixed contracts** tends to rise in one part of the year and fall in another, more than the curve priced in. If it exists, it is the way to make money from seasonality. The data of the level separate the two on purpose. Its model has a seasonal shape of the curve (priced in) and, in five spreads, a planted seasonal drift whose cumulative sum is $$ A\,\cos\left(\frac{2\pi\,(d - d^{\ast})}{365.25}\right) $$ with $d$ the day of the year, $d^{\ast}$ the day of the peak and $A$ the amplitude in logarithms: the spread rises from trough to peak for half a year and falls the other half. The curve does not price it in, but it dampens it in the far contracts: a pair of contracts expiring in $\tau$ years receives it multiplied by $e^{-\kappa_Q \tau}$ ($\kappa_Q$ is the mean-reversion speed priced in by the curve). The other four spreads are **controls**: they have no seasonal drift at all. This is **synthetic** data: no figure in this lesson says that a real spread is seasonal or by how much. ```python import numpy as np import pandas as pd from scipy import stats ver = pd.read_csv("spreads3_verdad.csv", index_col="spread") NUEVE = ["crack_gasolina", "spark", "invierno_verano_gas", "crush", "trigo_maiz", "crack_gasoleo", "aceite_harina", "oro_plata", "calendario_crudo"] v = ver.loc[NUEVE, ["patas", "control", "amplitud_estacional_log", "fecha_pico", "kappa_q"]] v.columns = ["patas", "control", "A", "pico", "kappa_q"] v ``` ```output patas control A pico kappa_q spread crack_gasolina PEG/PEC False 0.030 05-31 4.0 spark ELB/GNA False 0.080 07-16 6.0 invierno_verano_gas GNA ene/GNA oct False 0.035 09-30 1.0 crush GRH/GRA/GRS False 0.025 12-15 2.5 trigo_maiz GRT/GRM False 0.018 01-14 1.0 crack_gasoleo PEH/PEC True NaN NaN 3.0 aceite_harina GRA/GRH True NaN NaN 0.8 oro_plata MTO/MTP True NaN NaN 0.0 calendario_crudo PEC 1.º/2.º True NaN NaN 1.8 ``` The gasoline crack has its planted peak on May 31; the spark spread (peak electricity versus gas), on July 16; the winter gas premium (January versus October), on September 30; the soybean crush, on December 15, and wheat versus corn, on January 14. The controls are the heating oil crack, the split of the crush between oil and meal, the gold/silver ratio (a random walk with a script of episodes) and the crude oil 1st-2nd calendar. All spreads are measured **in logarithms**: the two-leg one is $\ln F_1 - \ln F_2$, the equal-value bet on the two legs from lesson 3.2.4. The honest analysis has six steps, and each one avoids one of the traps shown below: ```mermaid flowchart TB A[Fixed contracts] --> B[Daily changes] B --> C[Harmonic and F] C --> D[Prior calendar] D --> E[Sharpe and its error] E --> F[Multiple testing] ``` ## Building the spread without jumps Seasonality is measured on the **daily change of a spread with the same contracts**, never on a series that changes contract. For an intercommodity spread (lesson 3.2.2), the two legs must have the same delivery month and change month at the same time: you take the first common month whose roll date has not arrived in any of the legs. The crush has three legs (meal and oil against soybeans, with the weights of the crushing margin), so each side is a basket. ```python cur = pd.read_parquet("futuros3_curvas.parquet") ctr = pd.read_csv("futuros3_contratos.csv", parse_dates=["vencimiento", "fecha_roll"]) px = cur.set_index(["fecha", "contrato"])["liquidacion"] nombre = ctr.set_index(["mercado", "mes_entrega"])["contrato"] fechas = pd.DatetimeIndex(np.sort(cur["fecha"].unique())) f0 = np.r_[fechas[:1], fechas[:-1]] # previous session LETRAS = "FGHJKMNQUVXZ" def mes_comun(mercados, letras=None): """Delivery month held in all the legs: the first one whose roll date has not arrived in any of the legs.""" t = ctr[ctr["mercado"].isin(mercados)] if letras: mes = t["mes_entrega"].str[5:].astype(int) t = t[mes.map(lambda k: LETRAS[k - 1]).isin(list(letras))] g = t.groupby("mes_entrega").agg(n=("mercado", "nunique"), roll=("fecha_roll", "min")) roll = g.loc[g["n"] == len(mercados), "roll"].sort_index() roll = roll[::-1].cummin()[::-1] i = np.searchsorted(roll.to_numpy(), fechas.to_numpy(), "right") meses = np.r_[roll.index.to_numpy(), [None]] return pd.Series(meses[i], index=fechas) def precio(f, contratos): return px.reindex(list(zip(f, contratos))).to_numpy() def cesta(patas, f, mes): """Value of a basket [(market, weight)] with month `mes`.""" v = np.zeros(len(f)) for mk, peso in patas: c = nombre.reindex(list(zip([mk] * len(mes), mes))).to_numpy() v = v + peso * precio(f, c) return v def simple(x): return x - 1 def cambio(largo, corto, letras=None, fn=np.log): """Daily change with the previous day's contracts (no roll jumps): of the log-spread ln(long) - ln(short) or, with fn=simple, simple return of the long side minus that of the short side.""" mes = mes_comun({m for m, _ in largo + corto}, letras).shift(1) d = (fn(cesta(largo, fechas, mes) / cesta(largo, f0, mes)) - fn(cesta(corto, fechas, mes) / cesta(corto, f0, mes))) return pd.Series(d, index=fechas).iloc[1:] PATAS = { "crack_gasolina": ([("PEG", 42)], [("PEC", 1)], None), "crack_gasoleo": ([("PEH", 42)], [("PEC", 1)], None), "spark": ([("ELB", 1)], [("GNA", 7)], None), "crush": ([("GRH", 0.022), ("GRA", 0.11)], [("GRS", 0.01)], "FHKNQU"), "aceite_harina": ([("GRA", 1)], [("GRH", 1)], "FHKNQUVZ"), "oro_plata": ([("MTO", 1)], [("MTP", 1)], None), "trigo_maiz": ([("GRT", 1)], [("GRM", 1)], "HKNUZ"), } D = {k: cambio(*v) for k, v in PATAS.items()} # log R = {k: cambio(*v, fn=simple) for k, v in PATAS.items()} # simple print({k: int(v.notna().sum()) for k, v in D.items()}) ``` ```output {'crack_gasolina': 4023, 'crack_gasoleo': 4023, 'spark': 4023, 'crush': 4023, 'aceite_harina': 4023, 'oro_plata': 4023, 'trigo_maiz': 4023} ``` The seven series have all 4,023 sessions with a change from 2010 to 2025. The common-month rule matters: with the **latest** roll date of the legs, instead of the earliest, in the cracks an already expired leg leaves about twenty sessions a year without a price (crude oil stops trading before gasoline of the same month). The table of measures in the data's README uses that other rule, and that is why its crack figures do not match those of this lesson: the construction of the spread already moves the result. The two calendars are built separately: the crude oil one with the listed pair from `spreads3_cotizados.parquet` (lesson 3.2.5), and the winter gas premium with the fixed pair of January versus the previous October, which is held until October expires: ```python cot = pd.read_parquet("spreads3_cotizados.parquet") pec = cot[cot["mercado"] == "PEC"].set_index("fecha") c1 = pec["pata_cerca"].shift(1).reindex(fechas) c2 = pec["pata_lejos"].shift(1).reindex(fechas) def par(a, b, fn): return (fn(precio(fechas, a) / precio(f0, a)) - fn(precio(fechas, b) / precio(f0, b))) # Gas: January of the following year versus October, until October expires ano = np.where(fechas.month <= 10, fechas.year, fechas.year + 1) ene = pd.Series([f"GNAF{(a + 1) % 100:02d}" for a in ano], fechas) octb = pd.Series([f"GNAV{a % 100:02d}" for a in ano], fechas) e1, o1 = ene.shift(1), octb.shift(1) vivo = fechas <= ctr.set_index("contrato")["vencimiento"].reindex(o1) for dic, fn in ((D, np.log), (R, simple)): dic["calendario_crudo"] = pd.Series(par(c1, c2, fn), index=fechas).iloc[1:] gas = np.where(vivo, par(e1, o1, fn), np.nan) dic["invierno_verano_gas"] = pd.Series(gas, index=fechas).iloc[1:] D = {k: v.replace([np.inf, -np.inf], np.nan) for k, v in D.items()} R = {k: v.replace([np.inf, -np.inf], np.nan) for k, v in R.items()} pd.Series({k: v.std() * np.sqrt(252) for k, v in D.items()}, name="vol_anual").round(3) ``` ```output crack_gasolina 0.163 crack_gasoleo 0.140 spark 0.287 crush 0.103 aceite_harina 0.272 oro_plata 0.180 trigo_maiz 0.166 calendario_crudo 0.040 invierno_verano_gas 0.156 Name: vol_anual, dtype: float64 ``` Volatilities range from 4% a year for the crude oil calendar to 29% for the spark. They matter for what comes next: a seasonal drift of a few points in a spread that moves 20% a year is hard to see with few years. ## The seasonal chart and its noise The classic tool of seasonal analysis is the **multi-year profile**: the cumulative change of the spread since January 1, year by year, and its average. It is what futures seasonality services publish. With the two cracks: ```python def perfil(serie): """Cumulative change of the log-spread since January 1 of each year, measured every week: mean and quartiles across years.""" s = serie.fillna(0) acum = s.groupby(s.index.year).cumsum() semana = np.minimum((s.index.dayofyear - 1) // 7 + 1, 52) tabla = acum.groupby([s.index.year, semana]).last().unstack(0) return pd.DataFrame({"media": tabla.mean(axis=1), "p25": tabla.quantile(0.25, axis=1), "p75": tabla.quantile(0.75, axis=1), "error": tabla.std(axis=1) / np.sqrt(16)}) PERF = {k: perfil(D[k]) for k in ["crack_gasolina", "crack_gasoleo"]} for k, p in PERF.items(): print(f"{k:15s} mean: minimum week {p['media'].idxmin():2d} " f"({p['media'].min():+.3f}), maximum week " f"{p['media'].idxmax():2d} ({p['media'].max():+.3f}); " f"standard error of the mean in week 22: " f"{p.at[22, 'error']:.3f}") ``` ```output crack_gasolina mean: minimum week 36 (-0.009), maximum week 16 (+0.068); standard error of the mean in week 22: 0.030 crack_gasoleo mean: minimum week 19 (-0.013), maximum week 26 (+0.013); standard error of the mean in week 22: 0.019 ``` *Fig. — Seasonal profile of two data cracks: cumulative log-spread change from 1 January, 16-year mean (solid line) and quartiles across years (dashed).* Two panels by week of year. In the gasoline crack, highlighted, the mean rises to 0.068 in week 16 and falls afterwards, but the quartiles run from −0.09 to 0.12. In the heating oil crack, with no planted seasonality, the mean stays between −0.013 and 0.013 and the quartiles run from −0.09 to 0.08: almost the same. The average of the gasoline crack rises 6.8% until week 16 (mid-April) and falls afterwards; the planted peak is in week 22 (May 31). That of heating oil, with no seasonality, stays between −1.3% and +1.3%. The profile points in the right direction, but two things force you not to trust it at first sight. The first is its error: the 16-year average has a standard error of 3 points in week 22, of the same order as the average rise itself in those weeks. The second is the quartiles in the figure: in both spreads, half of the years fall within a band of about ±10 points; individual years do not follow the profile, and the control has almost the same dispersion. A chart of averages shows a smooth curve even if the data is noise, because averaging smooths; that is why a test is needed. ## An annual harmonic with its test A natural test for a cosine-shaped drift is a regression of the daily change $\Delta s_t$ on the daily increments of an annual cosine and sine: $$ \Delta s_t = c + a\,\Delta\cos\left(\frac{2\pi d_t}{365.25}\right) + b\,\Delta\sin\left(\frac{2\pi d_t}{365.25}\right) + \varepsilon_t $$ The amplitude is $\sqrt{a^2 + b^2}$ and the peak, day $\operatorname{atan2}(b, a) \cdot 365.25 / 2\pi$. The F-test of $a = b = 0$ asks whether there is seasonality. It has two advantages over the chart: it uses every day (not 16 years) and it looks for a specific shape, the model's, instead of any curve. ```python W = 2 * np.pi / 365.25 def armonico(s): """Regression of the daily change on the daily increments of annual cos and sin: amplitude, day of the peak and F p-value.""" doy = s.index.dayofyear.to_numpy() dc = np.r_[np.nan, np.diff(np.cos(W * doy))] ds = np.r_[np.nan, np.diff(np.sin(W * doy))] ok = np.isfinite(s.to_numpy()) & np.isfinite(dc) X = np.c_[np.ones(ok.sum()), dc[ok], ds[ok]] y = s.to_numpy()[ok] b = np.linalg.lstsq(X, y, rcond=None)[0] rss1 = ((y - X @ b) ** 2).sum() rss0 = ((y - y.mean()) ** 2).sum() gl = len(y) - 3 p = stats.f.sf((rss0 - rss1) / 2 / (rss1 / gl), 2, gl) pico = (np.arctan2(b[2], b[1]) / W) % 365.25 return np.hypot(b[1], b[2]), pico, p def fecha(dia): return (pd.Timestamp("2025-01-01") + pd.Timedelta(days=round(dia) - 1)).strftime("%m-%d") filas = {} for k in NUEVE: a, pico, p = armonico(D[k]) filas[k] = {"control": ver.at[k, "control"], "A_plantada": ver.at[k, "amplitud_estacional_log"], "pico_plantado": ver.at[k, "fecha_pico"], "A_medida": a, "pico_medido": fecha(pico), "p_F": p} med = pd.DataFrame.from_dict(filas, orient="index") med.round(3) ``` ```output control A_plantada pico_plantado A_medida pico_medido p_F crack_gasolina False 0.030 05-31 0.032 05-21 0.000 spark False 0.080 07-16 0.052 08-17 0.001 invierno_verano_gas False 0.035 09-30 0.016 10-16 0.095 crush False 0.025 12-15 0.013 12-28 0.042 trigo_maiz False 0.018 01-14 0.021 01-12 0.033 crack_gasoleo True NaN NaN 0.007 07-23 0.592 aceite_harina True NaN NaN 0.010 07-01 0.767 oro_plata True NaN NaN 0.016 06-11 0.204 calendario_crudo True NaN NaN 0.001 09-02 0.915 ``` Against the known answer: - **The gasoline crack and the spark are recovered.** Measured peaks on May 21 and August 17 versus the planted May 31 and July 16, with p of 0.0002 and 0.0008. The measured amplitude of the spark (0.052) is smaller than the planted one (0.080): its $\kappa_Q = 6$ dampens it in the contracts held. That of the gasoline crack (0.032) is somewhat larger than the dampened planted one: noise adds up too. - **Wheat-corn and the crush come close.** Almost exact peaks (January 12 versus 14, December 28 versus 15) and p of 0.033 and 0.042: significant at 5% one by one, barely. - **The winter gas premium is not seen.** Peak on October 16 (planted: September 30), but p = 0.095. It is a real effect that the test does not detect with 16 years: a false negative. - **The controls give nothing.** p between 0.20 (gold-silver) and 0.91 (crude oil calendar), with amplitudes from 0.001 to 0.016. Gold-silver has a measured amplitude (0.016) larger than that of the crush (0.013), which does have an effect: its volatility keeps the test from taking it seriously. ## The continuous series trap The fastest way to build a seasonal profile is to take the unadjusted continuous series of the front contract of each leg (lesson 1.4.6) and subtract. It looks harmless. With the same data: ```python def frente(mk): """Price of the contract held, unadjusted (lesson 1.4.6).""" cm = ctr[ctr["mercado"] == mk].sort_values("fecha_roll") i = np.searchsorted(cm["fecha_roll"], fechas, side="right") return pd.Series(precio(fechas, cm["contrato"].to_numpy()[i]), index=fechas) ing = {} for k in ["crack_gasolina", "spark", "trigo_maiz", "crack_gasoleo", "aceite_harina"]: (a, _), (b, _) = PATAS[k][0][0], PATAS[k][1][0] nivel = np.log(frente(a)) - np.log(frente(b)) A1, p1, f1 = armonico(nivel.diff()) ing[k] = {"A_bien": med.at[k, "A_medida"], "pico_bien": med.at[k, "pico_medido"], "A_serie": A1, "pico_serie": fecha(p1), "p_serie": f1} pd.DataFrame.from_dict(ing, orient="index").round(3) ``` ```output A_bien pico_bien A_serie pico_serie p_serie crack_gasolina 0.032 05-21 0.090 04-30 0.000 spark 0.052 08-17 0.106 07-06 0.010 trigo_maiz 0.021 01-12 0.043 11-09 0.001 crack_gasoleo 0.007 07-23 0.004 07-18 0.911 aceite_harina 0.010 07-01 0.010 07-26 0.771 ``` *Fig. — Amplitude and peak date of the annual harmonic in three spreads with planted seasonality (planted peaks: 31 May, 16 July and 14 January), measured with the same contracts and with the unadjusted continuous series of the held contract in each leg.* Six horizontal bars in pairs. Gasoline crack: 0.032 peaking on 21 May with the same contracts and 0.090 peaking on 30 April with the continuous series, highlighted. Spark: 0.052 (17 August) against 0.106 (6 July). Wheat-corn: 0.021 (12 January) against 0.043 (9 November). The continuous series doubles or triples the amplitude and moves the peak by weeks or months. In the three spreads with seasonality, the continuous series **doubles or triples** the amplitude and moves the peak: the gasoline crack goes from 0.032 to 0.090 with the peak a month earlier (April 30), and wheat-corn, from January to November. The tests come out more significant than with the fixed contracts (p practically 0 in the crack). The reason is the seasonal shape of the curve: each roll changes from one delivery month to another with a different seasonal level, and since rolls fall on the same dates every year, their jumps repeat with the calendar. The harmonic takes them for a drift. In the two controls in the table the roll jumps do not have that annual pattern and the continuous series does not change the result (p of 0.91 and 0.77). It is the easiest mistake to make and the hardest to see: the continuous series gives a cleaner and more significant pattern than the true one, and on another date. The same confusion between level and drift appears with gas, where the shape of the curve is huge: ```python g = cur[cur["mercado"] == "GNA"].pivot(index="fecha", columns="contrato", values="liquidacion") nivel, subida = [], [] for a in range(2010, 2025): ene_a, oct_a = f"GNAF{(a + 1) % 100:02d}", f"GNAV{a % 100:02d}" par = np.log(g[ene_a] / g[oct_a]) nivel.append(par.dropna().mean()) tramo = D["invierno_verano_gas"].loc[f"{a}-04-10":f"{a}-09-20"] subida.append(tramo.sum()) print(f"Mean level of ln(January / October): {np.mean(nivel):+.3f}") print(f"Mean rise from 04-10 to 09-20: {np.mean(subida):+.3f} " f"(positive years: {np.sum(np.array(subida) > 0)} of 15)") ``` ```output Mean level of ln(January / October): +0.170 Mean rise from 04-10 to 09-20: +0.052 (positive years: 12 of 15) ``` The January contract is worth on average 17% more than the previous October (0.170 in log): it is the winter premium that any producer or storage operator knows and that the curve prices in. The average rise of the spread between April and September, the segment the strategy trades, is 0.052, with 12 positive years out of 15. Both figures are seasonal, but only the second is a possible gain, and it is a third of the first. ## A calendar strategy and the error of its Sharpe The simplest strategy on a seasonal drift is mechanical: long the spread from 10 days after the trough until 10 days before the peak, short from 10 days after the peak until 10 days before the trough, every year, with the previous day's position. The controls use the gasoline calendar (peak on May 31). Each day's result is that of a position with the same value in both legs: the simple return of the long leg minus that of the short one (the `R` series from the construction cell); the harmonic, on the other hand, is estimated with the changes of the log-spread (`D`). With 16 years of data there are, in practice, 16 observations of each calendar, and the standard error of an annual Sharpe, with approximately normal and independent returns (Lo, 2002; lesson 2.5.4), is $$ \mathrm{SE}(\widehat{SR}) \approx \sqrt{\frac{1 + \widehat{SR}^2/2}{T}} $$ with $T$ in years: about 0.25 with 16 years, whatever the spread. ```python def posicion(indice, pico): """+1 from trough + 10 days to peak - 10; -1 from peak + 10 to trough - 10; 0 otherwise.""" doy = indice.dayofyear.to_numpy() valle = (pico - 182.6) % 365.25 sube = ((doy - (valle + 10)) % 365.25) < 162.6 baja = ((doy - (pico + 10)) % 365.25) < 162.6 return pd.Series(np.where(sube, 1.0, np.where(baja, -1.0, 0.0)), index=indice) def pico_de(k): # the controls use the gasoline calendar p = ver.at[k, "dia_pico"] return p if p == p else ver.at["crack_gasolina", "dia_pico"] PNL = {} filas = {} for k in NUEVE: s = R[k] # simple returns of the legs pnl = (posicion(s.index, pico_de(k)).shift(1) * s).fillna(0) PNL[k] = pnl anos = len(pnl) / 252 sr = pnl.mean() / pnl.std() * np.sqrt(252) ee = np.sqrt((1 + sr**2 / 2) / anos) # Lo (2002), annual por_ano = pnl.groupby(pnl.index.year).sum() filas[k] = {"control": ver.at[k, "control"], "sharpe": sr, "error_tipico": ee, "t": sr / ee, "anos_positivos": f"{(por_ano > 0).sum()}/{len(por_ano)}"} estr = pd.DataFrame.from_dict(filas, orient="index") estr = estr.sort_values("sharpe", ascending=False) estr.round(2) ``` ```output control sharpe error_tipico t anos_positivos trigo_maiz False 0.70 0.28 2.52 13/16 invierno_verano_gas False 0.59 0.27 2.18 13/16 crack_gasolina False 0.54 0.27 2.03 11/16 spark False 0.51 0.27 1.91 13/16 oro_plata True 0.39 0.26 1.52 9/16 crush False 0.37 0.26 1.45 12/16 aceite_harina True 0.32 0.26 1.25 11/16 crack_gasoleo True 0.08 0.25 0.30 8/16 calendario_crudo True -0.08 0.25 -0.31 7/16 ``` The header figure draws the table with two-standard-error bars. The first four, wheat-corn (0.70), the winter gas premium (0.59), the gasoline crack (0.54) and the spark (0.51), have planted seasonality. The fifth is a control: gold-silver, with 0.39, **above the crush** (0.37), which does have it, and with only 9 positive years out of 16. The other control with a borrowed calendar, oil-meal, gives 0.32 with 11 positive years out of 16. Almost all the bars overlap: with 16 years, a Sharpe of 0.3 to 0.6 does not distinguish an effect from a coincidence. Notice also what is not seen in the F-test: the winter gas premium, which the test did not detect (p = 0.095), is here the second best; and gold-silver, which the test ruled out (p = 0.20), looks promising. The Sharpe and the test measure different things: the F, whether there is a specific annual shape; the Sharpe, whether a rule makes money, for whatever reason. With nine spreads tested at once, you also have to correct for multiple testing (lesson 2.1.8): ```python from statsmodels.stats.multitest import multipletests t = estr["t"] umbral = stats.norm.isf(0.05 / len(t)) # Bonferroni, one-sided p_holm = multipletests(med["p_F"], alpha=0.05, method="holm")[1] print(f"t required by Bonferroni with 9 tests: {umbral:.2f}, " f"a Sharpe of {umbral * estr['error_tipico'].mean():.2f}") pd.DataFrame({"control": estr["control"], "t": t.round(2), "pasa_bonferroni": t > umbral, "p_F": med["p_F"].reindex(estr.index).round(4), "p_F_holm": pd.Series(p_holm, med.index) .reindex(estr.index).round(3)}) ``` ```output t required by Bonferroni with 9 tests: 2.54, a Sharpe of 0.67 control t pasa_bonferroni p_F p_F_holm trigo_maiz False 2.52 False 0.0327 0.229 invierno_verano_gas False 2.18 False 0.0955 0.477 crack_gasolina False 2.03 False 0.0002 0.002 spark False 1.91 False 0.0008 0.006 oro_plata True 1.52 False 0.2039 0.816 crush False 1.45 False 0.0415 0.249 aceite_harina True 1.25 False 0.7670 1.000 crack_gasoleo True 0.30 False 0.5924 1.000 calendario_crudo True -0.31 False 0.9145 1.000 ``` Bonferroni with nine tests requires a t of 2.54 (one-sided): a Sharpe of about 0.67 with 16 years. None passes; wheat-corn stays at 2.52. With the Holm correction on the F-test p-values, which uses the daily data and the shape of the pattern, two do pass: the gasoline crack (corrected p 0.002) and the spark (0.006). The rest, including three real effects, cannot be distinguished from chance with this data. It is the honest conclusion from 16 years of data with effects of this size, and it is the one that seasonal spread analysis usually skips. ## Choosing the calendar without looking at the future In the previous table, the calendar was the true one, which nobody knows in real data. The usual thing is to estimate it: the day of the peak comes from the profile or from the harmonic. If it is estimated with the whole sample and the strategy is measured on the same sample, the day of the peak is a parameter fitted to the data (lesson 2.5.1). The honest way is to estimate it every year with the previous years, from 2014 to have at least four years: ```python def sharpe(x): return x.mean() / x.std() * np.sqrt(252) fuera = {} for k in NUEVE: s, lg = R[k], D[k] # simple P&L; calendar with the log trozos = [] for a in range(2014, 2026): _, pico, _ = armonico(lg[lg.index.year < a]) sa = s[s.index.year == a] trozos.append(posicion(sa.index, pico).shift(1).fillna(0) * sa) wf = pd.concat(trozos).fillna(0) s2 = s[s.index.year >= 2014] _, pico_todo, _ = armonico(lg) dentro = (posicion(s2.index, pico_todo).shift(1) * s2).fillna(0) verdad = (posicion(s2.index, pico_de(k)).shift(1) * s2).fillna(0) fuera[k] = {"control": ver.at[k, "control"], "calendario_verdadero": sharpe(verdad), "estimado_con_todo": sharpe(dentro), "estimado_antes": sharpe(wf)} pd.DataFrame.from_dict(fuera, orient="index").round(2) ``` ```output control calendario_verdadero estimado_con_todo estimado_antes crack_gasolina False 0.30 0.39 0.24 spark False 0.52 0.62 0.56 invierno_verano_gas False 0.51 0.47 0.68 crush False 0.24 0.65 0.38 trigo_maiz False 0.94 0.97 0.83 crack_gasoleo True 0.32 0.29 0.06 aceite_harina True 0.13 0.38 0.05 oro_plata True 0.24 0.26 0.26 calendario_crudo True -0.01 -0.20 -0.44 ``` The three columns cover the same years, 2014-2025. With the calendar estimated on the whole sample, the crush goes from 0.24 (true calendar) to 0.65, and the oil-meal control, from 0.13 to 0.38: the estimation looks for the dates that best fit the noise of those years. With the calendar estimated only with the previous years, the crush gives 0.38 and oil-meal, 0.05; the heating oil crack, 0.06. Out-of-sample validation removes most of the false edge, but not all of it: gold-silver gives between 0.24 and 0.26 in all three columns. It is not a flaw of the method. The gold/silver ratio in the data follows a script of episodes (long rises and falls on fixed dates) and, with 12 years, some fall in the same halves of the year. Out of sample, a control can win for a decade by calendar coincidence. ## What seasonality demands in live trading - **Fixed contracts or adjusted series.** The profile is computed with the pair held, without jumps. An unadjusted continuous series manufactures seasonality in markets with a seasonal curve shape, which are precisely the seasonal ones. - **Few observations.** Sixteen years are sixteen repetitions of each calendar. With a Sharpe of 0.5 and 16 years, the 95% interval goes from almost 0 to 1. More years help little if the market has changed (new specifications, liquefied natural gas exports, biofuels): the 1995 pattern need not be the 2025 one. - **Windows and buffers.** The 10-day buffers at each end are one more parameter; changing them changes the result (exercise 1). Choosing them by looking at the result is another form of fitting (lesson 2.5.6). - **Costs.** A seasonal strategy enters and exits twice a year, and it also rolls: in a calendar or intercommodity spread, each trade pays for both legs (lesson 3.2.5). - **Tail risk.** Energy and electricity spreads have jumps (hurricanes, cold snaps, refinery outages). The Sharpe of a seasonal strategy does not show its worst year: you have to look at it (exercise 2). ## Exercises 1. Repeat the gasoline crack strategy and the gold-silver control with buffers of 0, 10, 20 and 30 days at each end of the two halves of the year. Does the ranking between the two change? **Solution.** ```python def posicion_margen(indice, pico, margen): doy = indice.dayofyear.to_numpy() valle = (pico - 182.6) % 365.25 largo_tramo = 182.6 - 2 * margen sube = ((doy - (valle + margen)) % 365.25) < largo_tramo baja = ((doy - (pico + margen)) % 365.25) < largo_tramo return pd.Series(np.where(sube, 1.0, np.where(baja, -1.0, 0.0)), index=indice) tabla = {} for k in ["crack_gasolina", "oro_plata"]: s = R[k] tabla[k] = {m: sharpe((posicion_margen(s.index, pico_de(k), m) .shift(1) * s).fillna(0)) for m in [0, 10, 20, 30]} pd.DataFrame(tabla).T.round(2) ``` ```output 0 10 20 30 crack_gasolina 0.49 0.54 0.70 0.65 oro_plata 0.27 0.39 0.39 0.24 ``` The ranking does not change: the gasoline crack stays above gold-silver with any buffer. But its Sharpe goes from 0.49 to 0.70 just because of the buffer, a difference of almost one standard error; whoever publishes the 20-day one without saying they tried four is reporting the maximum of four attempts (lesson 2.5.4). The control has its best result with buffers of 10 and 20 days (0.39 in both) and drops to 0.24 with 30: it also has “an optimal buffer”, with no effect behind it. 2. For the nine spreads, compute the strategy's worst year (sum of its daily returns) and the Sharpe in two halves, 2010-2017 and 2018-2025. How many keep their sign? Is there a difference between planted spreads and controls? **Solution.** ```python mitades = {} for k in NUEVE: pnl = PNL[k] mitades[k] = {"control": ver.at[k, "control"], "peor_ano": pnl.groupby(pnl.index.year).sum().min(), "sharpe_2010_17": sharpe(pnl.loc[:"2017"]), "sharpe_2018_25": sharpe(pnl.loc["2018":])} pd.DataFrame.from_dict(mitades, orient="index").round(2) ``` ```output control peor_ano sharpe_2010_17 sharpe_2018_25 crack_gasolina False -0.19 0.80 0.30 spark False -0.53 0.67 0.34 invierno_verano_gas False -0.11 0.70 0.49 crush False -0.13 0.21 0.53 trigo_maiz False -0.22 0.38 1.06 crack_gasoleo True -0.25 -0.22 0.38 aceite_harina True -0.49 0.46 0.17 oro_plata True -0.15 0.65 0.11 calendario_crudo True -0.08 0.19 -0.38 ``` The five spreads with planted seasonality win in both halves. Of the four controls, two change sign (heating oil crack, from −0.22 to 0.38; crude oil calendar, from 0.19 to −0.38) and two stay positive (oil-meal, 0.46 and 0.17; gold-silver, 0.65 and 0.11): splitting the sample in two helps little with eight years per half. The worst year says what the Sharpe keeps quiet: the spark lost 0.53 in one year (sum of the position's daily returns) and the oil-meal control, 0.49. A seasonal strategy on a volatile spread needs a size that can withstand a year like that. 3. Design, without code, the study of the seasonality of the gasoline crack with real CME data (RBOB gasoline and WTI crude). What data would you request, how would you build the spread, what would you test and what would have to come out to trade it? **Solution.** Data: the daily settlements of all listed contracts of both legs, with their calendar (expiry and first notice day), not a continuous series; with Databento, the `statistics` schema for the settlements and `definition` for the expiries (lesson 3.2.5). Construction: pairs of the same delivery month with a fixed roll rule that comes before the first leg expires, in logarithms (the equal-value bet from lesson 3.2.4) or in barrels if what matters is the refining margin. Test: the annual harmonic (or a regression with several frequencies) on the daily changes with fixed contracts, with errors robust to autocorrelation (lesson 2.1.1), and its stability across subperiods; the calendar, estimated every year with the previous years. To trade it, what would have to come out is a pattern with a physical reason (the EIA describes summer demand and seasonal refining margins), a test that survives the correction for the number of spreads and windows tried, a positive out-of-sample result net of the costs of both legs and of the rolls, and a worst year acceptable for the position size. If all that comes out is a Sharpe of 0.4 with 15 years, the honest thing is to say that it cannot be distinguished from chance. ## Self-assessment 1. The gold-silver spread in the data, a control with no seasonality, gives a Sharpe of 0.39 with the seasonal strategy, above the crush, which does have it. With 16 years, the standard error is 0.26. What do you conclude? - A. That gold-silver has a seasonality the model does not capture - B. That with 16 years a Sharpe of 0.39 is 1.5 standard errors from zero: it is compatible with having nothing, and its F-test (p = 0.20) does not see seasonality either - C. That you should trade gold-silver and not the crush - D. That the standard error is computed wrong *Answer B: That with 16 years a Sharpe of 0.39 is 1.5 standard errors from zero: it is compatible with having nothing, and its F-test (p = 0.20) does not see seasonality either.* The standard error of an annual Sharpe is approximately √((1 + SR²/2) / years). With 16 years, 0.26: any Sharpe between −0.5 and +0.5 is compatible with zero at two standard errors. The ranking of Sharpes between 0.3 and 0.6 is noise. 2. With the unadjusted continuous series of the contract held in each leg, the gasoline crack has a seasonal amplitude of 0.090 and a peak on April 30; with fixed contracts, 0.032 and May 21. Where does the difference come from? - A. From the roll jumps: every year, on the same dates, the series moves from one delivery month to another with a different seasonal shape in the curve, and those jumps, repeated by the calendar, look like seasonality - B. From a rounding error - C. From the continuous series having more data - D. From gasoline being more volatile in April *Answer A: From the roll jumps: every year, on the same dates, the series moves from one delivery month to another with a different seasonal shape in the curve, and those jumps, repeated by the calendar, look like seasonality.* The gasoline curve has a seasonal shape by delivery month (summer contracts are worth more). When rolling from one month to another always on the same dates, the unadjusted series jumps the same way every year. The harmonic picks it up as if it were a drift: it triples the amplitude and moves the peak a month earlier. In wheat-corn it moves it from January to November. 3. In the data, the January gas contract trades on average 17% above the previous October contract. Is it an opportunity? - A. Yes: buying October and selling January earns 17% - B. Not by itself: that difference in level is the seasonal shape of the curve, which is already in the prices; what can be traded is how much the spread of two fixed contracts moves over the year, which is something else - C. Yes, but only in winter - D. No, because gas is not seasonal *Answer B: Not by itself: that difference in level is the seasonal shape of the curve, which is already in the prices; what can be traded is how much the spread of two fixed contracts moves over the year, which is something else.* That a winter contract is worth more than an autumn one is what you would expect with expensive storage and heating demand, and the market prices it in. The return of a position depends on how the spread changes, not on its level. In the data, the average rise of the spread between April and September is 0.052 in log, a third of the level. 4. You estimate the day of the seasonal peak with the 16 years and measure the strategy's Sharpe over those same years. What happens in the data? - A. Nothing: the day of the peak is not a parameter - B. The Sharpe always goes down - C. It only changes in the controls - D. The Sharpe gets inflated: in 2014-2025 the crush goes from 0.24 with the true calendar to 0.65 with the one estimated on the whole sample, and the oil-meal control, from 0.13 to 0.38; estimated only with the previous years, 0.38 and 0.05 *Answer D: The Sharpe gets inflated: in 2014-2025 the crush goes from 0.24 with the true calendar to 0.65 with the one estimated on the whole sample, and the oil-meal control, from 0.13 to 0.38; estimated only with the previous years, 0.38 and 0.05.* The day of the peak is chosen so that the strategy wins in the sample: it is a fitted parameter, and the in-sample result carries the selection bias. The honest way is to estimate it every year with the previous years (out of sample, lesson 2.5.2). 5. Nine spreads are tested with the seasonal strategy. With Bonferroni at 5% (one-sided), what t is needed and what Sharpe does it imply with 16 years? - A. t = 1.65, a Sharpe of 0.4 - B. t = 1.96, a Sharpe of 0.5 - C. t = 2.54, a Sharpe of about 0.67: none of the nine passes, not even wheat-corn (0.70, t = 2.52) - D. No correction is needed: each spread is independent *Answer C: t = 2.54, a Sharpe of about 0.67: none of the nine passes, not even wheat-corn (0.70, t = 2.52).* With 9 tests the threshold of each one drops to 0.05 / 9 and the required t rises to 2.54. With an average standard error of 0.26, that is a Sharpe of 0.67. The daily F-test, which uses the shape of the pattern, has more power: with the Holm correction, the gasoline crack and the spark pass. ## Key points - The seasonality of a spread is measured with fixed contracts: the daily change of a pair that does not change. With the unadjusted continuous series, rolls on fixed dates manufacture seasonality: the amplitude of the gasoline crack goes from 0.032 to 0.090. - The seasonal level of the curve (January 17% above October in gas) is priced in; what can be traded is the drift of a fixed pair, much smaller. - The annual harmonic with its F-test recovers the strongest planted effects (gasoline crack, spark: p < 0.001), comes close on the others and sees nothing in the four controls (p between 0.20 and 0.91). - With 16 years, the standard error of an annual Sharpe is about 0.26: the gold-silver control (0.39) ends up above the crush (0.37), which does have an effect, and no strategy passes Bonferroni with nine tests (a Sharpe of 0.67 is needed). - The calendar chosen with the whole sample inflates the result (crush from 0.24 to 0.65; the oil-meal control, from 0.13 to 0.38). It has to be estimated with the previous years. - Even out of sample, a control can win twelve years in a row by calendar coincidence (gold-silver, 0.26): seasonality needs a physical reason, not just a Sharpe. ## Resources ### Research - [Seasonality in Python II - Backtesting](https://pythonparatrading.com/en/seasonality-intraday-backtesting/) · Research-0025 · 8 min ### Documentation - [Lo (2002): The Statistics of Sharpe Ratios, Financial Analysts Journal 58(4)](https://doi.org/10.2469/faj.v58.n4.2453) · doi.org - [EIA: Factors affecting gasoline prices (refining margins and summer demand)](https://www.eia.gov/energyexplained/gasoline/factors-affecting-gasoline-prices.php) · eia.gov - [EIA: Factors affecting natural gas prices (storage and winter demand)](https://www.eia.gov/energyexplained/natural-gas/factors-affecting-natural-gas-prices.php) · eia.gov - [statsmodels: multipletests (Bonferroni, Holm, FDR)](https://www.statsmodels.org/stable/generated/statsmodels.stats.multitest.multipletests.html) · statsmodels.org ### Data - [Level 3 data: ground truth of the spreads (CSV)](https://pythonparatrading.com/curso/datos/nivel-3/spreads3_verdad.csv) · Download · spreads3_verdad.csv - [Level 3 data: full futures curves (Parquet)](https://pythonparatrading.com/curso/datos/nivel-3/futuros3_curvas.parquet) · Download · futuros3_curvas.parquet Previous lesson: [3.2.5 · Spread data and backtesting: listed, synthetic and legged](https://pythonparatrading.com/en/python-trading-course/level-3/spread-data-and-legged-backtesting/) Next lesson: [3.3.1 · FX microstructure: OTC, ECN, last look and sessions](https://pythonparatrading.com/en/python-trading-course/level-3/fx-microstructure-otc-ecn-last-look/)