import os, io, requests
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.formula.api as smf
from scipy import stats
from scipy.optimize import minimize
def import_csv(file, base_url="https://www.apradie.com/datos/", **kwargs):
"""Imports a csv file from base_url; downloads it only once and keeps a local copy."""
if not os.path.exists(file):
headers = {"User-Agent": "Mozilla/5.0 (Windows NT 10.0; Win64; x64) AppleWebKit/537.36 "
"(KHTML, like Gecko) Chrome/91.0.4472.124 Safari/537.36"}
response = requests.get(base_url + file, headers=headers)
response.raise_for_status()
with open(file, "wb") as f:
f.write(response.content)
return pd.read_csv(file, **kwargs)
prices = import_csv("fz2024_prices_monthly.csv", index_col="date", parse_dates=True)
rf = import_csv("fz2024_rf_monthly.csv", index_col="date", parse_dates=True)
R = prices.pct_change(fill_method=None).iloc[1:]
r = np.log(prices).diff().iloc[1:]
R_ex = R.sub(rf["rf"].reindex(R.index), axis=0)Workshop 4, Financial Modeling and Programming
Integration and exam warm-up: risk measures, the complete investment pipeline, and review
Abstract
This is an INDIVIDUAL workshop and your best preparation for the Evidence (practical exam) and the theory exam of week 5. We add three risk measures that complement volatility (maximum drawdown, historical Value at Risk, and the Information ratio), we put all the steps of Workshops 1-3 together in a pipeline of functions, and we run a complete dry run of the Evidence with a different ratio and different dates. We finish with interpretation drills, code reading and a bug hunt that cover the whole course.
1 General directions
Notebook W4-FinProg-YourFirstName-YourLastName, same rules as before. Grade: notebook 60% + Canvas Quiz W4 40%.
This workshop has two parts:
- Session 1: three new risk measures and the complete pipeline: a dry run of the Evidence with different data.
- Session 2: review for the theory exam: interpretation drills, code reading and a bug hunt with questions similar to the exam.
2 Risk beyond volatility
Volatility treats upside and downside surprises the same way, and it does not tell a client “how much can I lose?”. Three complementary measures:
2.1 Maximum drawdown
The drawdown at month t is the percentage loss from the highest value reached before t (the running peak). The maximum drawdown (MDD) is the worst drawdown of the period: the loss of an investor who bought at the worst possible peak and sold at the worst possible trough.
DD_t=\frac{V_t}{\max_{s\le t}V_s}-1 \qquad\qquad MDD=\min_t DD_t
def max_drawdown(port_R):
"""Maximum drawdown from monthly simple returns."""
value = (1 + port_R).cumprod()
peak = value.cummax() # running maximum
return (value / peak - 1).min()
spy = R.loc[:"2024-03-31", "SPY"] # all the months of our data up to March 2024
value = (1 + spy).cumprod()
dd = value / value.cummax() - 1
fig, ax = plt.subplots(2, 1, figsize=(9, 5), sharex=True)
ax[0].plot(value); ax[0].plot(value.cummax(), ls="--"); ax[0].set_title("SPY: value of $1 and running peak")
ax[1].fill_between(dd.index, dd, 0, color="red", alpha=0.4); ax[1].set_title("Drawdown")
for a in ax: a.grid(alpha=0.3)
plt.tight_layout(); plt.show()
print(f"Maximum drawdown of SPY, Jan 2015 - Mar 2024: {max_drawdown(spy):.2%}")
Maximum drawdown of SPY, Jan 2015 - Mar 2024: -23.93%cummax() is the accumulator pattern of the Algorithms course (keep the maximum seen so far), vectorized.
2.2 Historical Value at Risk (VaR)
The 95% monthly VaR is the loss that is exceeded only in 5% of the months. The historical method uses the 5th percentile of the past monthly returns, with no assumption about their distribution:
def hist_var(port_R, level=0.95):
"""Historical VaR (as a positive loss) of monthly simple returns."""
return -np.quantile(port_R.dropna(), 1 - level)
spy_all = R.loc[:"2024-03-31", "SPY"]
print(f"95% monthly historical VaR of SPY: {hist_var(spy_all):.2%}")95% monthly historical VaR of SPY: 7.41%Interpretation: “in 95% of the months, the loss was not larger than this percentage; in 1 of every 20 months it was larger.” VaR does not say how large the loss is in that bad month (the drawdown and the worst month help with that).
3 The complete pipeline
These are the steps of Workshops 1–3, written as functions. Read each one and identify the workshop where you programmed it.
def prepare_panel(uspanel, usfirms):
"""W1: merge, quarterly amounts from YTD, TTM figures and ratios (roa_ttm, roabit_ttm, oepspw, bmrw)."""
df = pd.merge(uspanel, usfirms[["empresa", "Nombre", "status", "partind", "naics1"]],
left_on="firm", right_on="empresa", how="left").drop(columns="empresa")
assert len(df) == len(uspanel)
df = df.sort_values(["firm", "q"]).reset_index(drop=True)
df["fq"] = np.where(df["fiscalmonth"].isin([3, 6, 9, 12]), df["fiscalmonth"] / 3, np.nan)
df["ebit_ytd"] = df["revenue"] - df["cogs"] - df["sgae"] - df["depreciation"]
df["netincome_ytd"] = (df["ebit_ytd"] + df["otherincome"] + df["extraordinaryitems"]
- df["finexp"] - df["incometax"])
g = df.groupby("firm")
for col, new in [("revenue", "revenueq"), ("ebit_ytd", "ebitq"), ("netincome_ytd", "netincomeq")]:
prev, prev_fq = g[col].shift(1), g["fq"].shift(1)
df[new] = np.select([df["fq"] == 1, prev_fq == df["fq"] - 1], [df[col], df[col] - prev], np.nan)
g = df.groupby("firm")
for q_col, ttm_col in [("revenueq", "revenue_ttm"), ("ebitq", "ebit_ttm"), ("netincomeq", "netincome_ttm")]:
df[ttm_col] = g[q_col].transform(lambda s: s.rolling(4, min_periods=4).sum())
df["totalassets_lag4"] = g["totalassets"].shift(4)
df["roa_ttm"] = df["netincome_ttm"] / df["totalassets_lag4"]
df["roabit_ttm"] = df["ebit_ttm"] / df["totalassets_lag4"]
df["marketvalue"] = df["originalprice"] * df["sharesoutstanding"]
df["oepsp"] = (df["ebit_ttm"] / df["sharesoutstanding"]) / df["originalprice"]
df["bmr"] = (df["totalassets"] - df["totalliabilities"]) / df["marketvalue"]
df["r"] = g["adjprice"].pct_change(fill_method=None)
df["r_next"] = df.groupby("firm")["r"].shift(-1)
return df.replace([np.inf, -np.inf], np.nan)
def winsorize(s, low=0.01, high=0.99):
return s.clip(lower=s.quantile(low), upper=s.quantile(high))
def screen(annual, ratio, n, last_q, first_q):
"""W1: top-n S&P 500 firms by ratio, last annual report with fiscal-year end in [first_q, last_q]."""
sel = annual[annual["partind"].notna() & annual["q"].between(first_q, last_q) & annual[ratio].notna()]
sel = sel.sort_values(["firm", "q"]).groupby("firm").tail(1)
return sel.nlargest(n, ratio)[["firm", "Nombre", "q", ratio]].reset_index(drop=True)
def capm(ticker, R_ex, start, end, market="SPY", min_obs=24):
"""W2: CAPM regression with simple excess returns."""
d = pd.DataFrame({"stock": R_ex.loc[start:end, ticker], "mkt": R_ex.loc[start:end, market]}).dropna()
if len(d) < min_obs:
return None
m = smf.ols("stock ~ mkt", data=d).fit()
t_b1 = (m.params["mkt"] - 1) / m.bse["mkt"]
return pd.Series({"alpha": m.params["Intercept"], "p_alpha": m.pvalues["Intercept"],
"beta": m.params["mkt"], "p_beta1": 2 * stats.t.sf(abs(t_b1), df=m.df_resid),
"n": int(m.nobs)})
def capm_many(tickers, R_ex, start, end):
"""W2: CAPM for many tickers; skips tickers without enough data."""
res = {tk: capm(tk, R_ex, start, end) for tk in tickers if tk in R_ex.columns}
return pd.DataFrame({k: v for k, v in res.items() if v is not None}).T
def opt_inputs(assets, start, end):
"""W3: annual expected simple returns (from cc means) and annual covariance of cc returns."""
rr = r.loc[start:end, assets].dropna()
return np.exp(12 * rr.mean()) - 1, 12 * rr.cov()
def gmv_long_only(cov):
"""W3: long-only Global Minimum Variance weights."""
n = cov.shape[0]
res = minimize(lambda w: w @ cov.values @ w, np.repeat(1 / n, n), method="SLSQP",
bounds=[(0, 1)] * n, constraints=[{"type": "eq", "fun": lambda w: w.sum() - 1}],
options={"ftol": 1e-12, "maxiter": 1000})
return pd.Series(res.x, index=cov.index)
def backtest(w, assets, start, end):
"""W3: buy-and-hold monthly simple returns."""
growth = (1 + R.loc[start:end, assets].fillna(0)).cumprod()
value = growth @ np.asarray(w)
return value.pct_change().fillna(value.iloc[0] - 1)
def report(port_R, bench_R, rf_annual):
"""W3 + W4: performance and risk summary."""
n = len(port_R)
hpr = (1 + port_R).prod() - 1
ann = (1 + hpr) ** (12 / n) - 1
vol = np.log(1 + port_R).std() * np.sqrt(12)
return pd.Series({"HPR": hpr, "annual return": ann, "annual volatility": vol,
"Sharpe": (ann - rf_annual) / vol, "max drawdown": max_drawdown(port_R),
"VaR 95% (monthly)": hist_var(port_R), "Information ratio": info_ratio(port_R, bench_R)})4 CHALLENGE 1: Dry run of the Evidence
You are a portfolio manager on March 31, 2022. Using the functions above:
- Import
usdataq.csvandfirmsus2024.csv, runprepare_panel, and createannual(rows withfiscalmonth == 12). - Screen: select the 20 S&P 500 firms with the highest operating ROA (
roabit_ttm= TTM EBIT / total assets one year before), using the last annual report with fiscal-year end between 2021q1 and 2021q4. - CAPM: estimate alpha and beta of the 20 firms with simple excess returns, window April 2019 – March 2022. Select the 10 with the highest alpha, and from them the 5 with the lowest beta.
- Optimize: calculate the long-only GMV portfolio of the 5 stocks (same window).
- Backtest (buy and hold) from April 2022 to March 2023 the GMV portfolio, the equally weighted portfolio of the 5 stocks, and SPY. Show the
reportof each one (use SPY as the benchmark for the Information ratio, and the annualized average T-Bill rate of the period asrf_annual). - Write a short paragraph for your investment committee. Mention one limitation of the process (look-ahead, survivorship, estimation risk or the length of the test).
# Write your code here5 CHALLENGE 2: A logistic regression in the style of the Evidence
Using df from Challenge 1:
- Create the event
beat_ind_next= 1 if the next-quarter returnr_nextis higher than the median next-quarter return of the firm’s industry (naics1) in the same quarter, 0 otherwise (missing ifr_nextis missing). Hint:groupby(["naics1", "q"])["r_next"].transform("median"). - Create
roaw=roa_ttmwinsorized by quarter at 1%–99%. Build the dataset with non-missingbeat_ind_nextandroawand report the number of observations and the percentage of ones. - Run the logit
beat_ind_next ~ roaw. Interpret sign and significance. - Calculate and interpret the odds ratio for an increase of +0.05 in
roaw. - Calculate the probability of beating the industry for a firm with
roaw = 0and for a firm withroaw = 0.10.
# Write your code here6 CHALLENGE 3: Interpretation drills
Answer in a text cell, without running code first. Then check your answers with Python where possible.
- A logit has b_1=2.0 for a ratio X. (a) What is the odds ratio for \Delta X=+0.10? (b) A firm has p=0.30; what is its new probability after \Delta X=+0.10?
- A logit gives b_0=-1.5. What is the probability of the event when all X=0?
- A CAPM regression gives \beta=1.35 with standard error 0.15. Is the stock significantly riskier than the market at 5%? (Use 2 as the critical value.)
- A monthly alpha is 0.006 with p-value 0.21. Write the interpretation in one sentence.
- A portfolio had an annual return of 14% and a volatility of 20%; the T-Bill rate was 4%. Its benchmark returned 12% with a tracking error of 5%. Calculate the Sharpe ratio and the Information ratio (assume the active return is 2%).
- Why does the GMV portfolio not need expected returns? Why is this an advantage?
- What is the difference between a backtest in sample and out of sample? If a screen uses the annual statements of fiscal 2023, why must the decision date be the end of March 2024 and the backtest start in April 2024?
Write your answers here
7 CHALLENGE 4: Code reading
Predict the output before running each piece of code. Then run it and explain any difference.
4.1
s = pd.Series([100, 120, 90, 108])
print(((s / s.cummax()) - 1).min())4.2
x = pd.DataFrame({"firm": ["A", "A", "B", "B"], "ta": [10, 12, 50, 55]})
print(x.groupby("firm")["ta"].shift(1).tolist())4.3
w = np.array([0.5, 0.5]); cov = np.array([[0.04, 0.0], [0.0, 0.09]])
print(round(float(np.sqrt(w @ cov @ w)), 4))4.4
R = pd.Series([0.10, -0.10])
print(round((1 + R).prod() - 1, 4), round(R.mean(), 4))Write your predictions and explanations here
8 CHALLENGE 5: Bug hunt (whole course)
5.1 ROA with the previous year’s assets, quarterly panel sorted by firm and q:
df["roa"] = df["netincome_ttm"] / df.groupby("firm")["totalassets"].shift(1)5.2 Odds ratio for +0.05 in ROA:
odds_ratio = np.exp(0.05) * model.params["roaw"]5.3 Excess returns:
R_ex = R.sub(rf["rf"], axis=0) # R starts in Jan 2015, rf starts in Dec 20145.4 The decision date is March 31, 2024 and the screen uses annual statements:
top = screen(annual, "roa_ttm", 20, "2024q1", "2023q2")5.5 GMV weights for the selected stocks:
mu, cov = opt_inputs(sorted(MY5), START, END)
w = gmv_long_only(cov)
port = R[MY5] @ w.valuesWrite your answers here
9 Preparing for the Evidence and the theory exam
Evidence (practical exam): you will run a complete pipeline like Challenge 1 plus a logistic regression like Challenge 2, with a different sample of firms. You can use your notebooks and an LLM, but you must explain every step with your own words; explanations written by an LLM are penalized. Keep only the code you need.
Theory exam: multiple-choice and short-answer questions about the concepts of the course: YTD vs quarterly vs TTM amounts, look-ahead and survivorship bias, probability/odds/log-odds and the interpretation of a logit, alpha and beta and their tests, simple vs cc returns, the GMV and tangency portfolios, the CML, estimation risk, backtesting, drawdown, VaR, and Sharpe vs Information ratio. Code-reading questions like Challenge 4 are included.
10 W4 submission
Submit your .ipynb file in Canvas and answer the Canvas Quiz W4.