From 4bf19126ecaad032294a0a99a9a4dbc6679c06b0 Mon Sep 17 00:00:00 2001 From: Yuning598 <1147936698@qq.com> Date: Sat, 8 Aug 2026 16:43:04 +0800 Subject: [PATCH] Remove monthly pipeline from main --- chars_ciz_monthly/abr.py | 230 --- chars_ciz_monthly/accounting.py | 2511 ----------------------- chars_ciz_monthly/download_data.py | 284 --- chars_ciz_monthly/functions.py | 478 ----- chars_ciz_monthly/iclink_ciz.sas | 281 --- chars_ciz_monthly/impute_rank_output.py | 271 --- chars_ciz_monthly/merge_chars.py | 209 -- chars_ciz_monthly/myre.py | 133 -- chars_ciz_monthly/rolling_chars.py | 489 ----- chars_ciz_monthly/sue.py | 158 -- 10 files changed, 5044 deletions(-) delete mode 100644 chars_ciz_monthly/abr.py delete mode 100644 chars_ciz_monthly/accounting.py delete mode 100644 chars_ciz_monthly/download_data.py delete mode 100644 chars_ciz_monthly/functions.py delete mode 100644 chars_ciz_monthly/iclink_ciz.sas delete mode 100644 chars_ciz_monthly/impute_rank_output.py delete mode 100644 chars_ciz_monthly/merge_chars.py delete mode 100644 chars_ciz_monthly/myre.py delete mode 100644 chars_ciz_monthly/rolling_chars.py delete mode 100644 chars_ciz_monthly/sue.py diff --git a/chars_ciz_monthly/abr.py b/chars_ciz_monthly/abr.py deleted file mode 100644 index fd29afb..0000000 --- a/chars_ciz_monthly/abr.py +++ /dev/null @@ -1,230 +0,0 @@ -# Calculate HSZ Replicating Anomalies -# ABR: Cumulative abnormal stock returns around earnings announcements -# Optimized version with Polars + DuckDB - -import polars as pl -import duckdb -import os -from functions import INPUT_PATH, OUTPUT_PATH - -num_threads = max(1, os.cpu_count() // 2) - -################### -# Compustat Block # -################### -comp = ( - pl.scan_parquet(INPUT_PATH + "comp_fundq.parquet") - .select(["gvkey", "datadate", "rdq", "fyearq", "fqtr"]) - .with_columns( - pl.col("datadate").cast(pl.Date), - pl.col("rdq").cast(pl.Date), - ) - .collect() -) - -################### -# CCM Block # -################### -ccm = ( - pl.scan_parquet(INPUT_PATH + "ccm.parquet") - .select([ - "gvkey", - pl.col("permno").cast(pl.Int64), - pl.col("linkdt").cast(pl.Date), - pl.col("linkenddt").cast(pl.Date), - ]) - .with_columns( - pl.col("linkenddt").fill_null(pl.date(2099, 12, 31)) - ) - .collect() -) - -# Join comp with CCM and filter by link dates -ccm_merged = ( - comp.join(ccm, on="gvkey", how="left") - .filter( - (pl.col("datadate") >= pl.col("linkdt")) - & (pl.col("datadate") <= pl.col("linkenddt")) - ) - .select(["gvkey", "datadate", "rdq", "fyearq", "fqtr", "permno"]) -) - -################### -# CRSP Block # -################### -# Map RDQ to first trading day on or after -# Get trading days from CRSP daily data -trading_days = ( - pl.scan_parquet(INPUT_PATH + "crsp_dsf.parquet") - .select(pl.col("dlycaldt").cast(pl.Date).alias("date")) - .unique() - .sort("date") - .collect() -) - -# Find the closest trading day (rdq_trad) on or after rdq -ccm1 = ( - ccm_merged - .sort("rdq") - .join_asof( - trading_days, - left_on="rdq", - right_on="date", - strategy="forward", - ) - .rename({"date": "rdq_trad"}) - .select(["gvkey", "permno", "datadate", "fyearq", "fqtr", "rdq", "rdq_trad"]) -) - -############################# -# CRSP abnormal return # -############################# -# Load CRSP daily data with Fama-French market return (mktrf + rf) -crsp_d = ( - pl.scan_parquet(INPUT_PATH + "crsp_dsf.parquet") - .select([ - pl.col("permno").cast(pl.Int64), - pl.col("dlycaldt").cast(pl.Date).alias("date"), - pl.col("dlyret").alias("ret"), - (pl.col("mktrf") + pl.col("rf")).alias("mkt"), - ]) - .collect() -) - -# Abnormal return = stock total return - CRSP VW market total return -crsp_d = ( - crsp_d.with_columns( - (pl.col("ret") - pl.col("mkt")).alias("abrd") - ) - .select(["date", "permno", "ret", "mkt", "abrd"]) -) - -################################ -# Event window range join # -################################ -# Add window bounds to nearest trading day of rdq -ccm1 = ccm1.filter(pl.col("rdq_trad").is_not_null()).with_columns([ - (pl.col("rdq_trad") - pl.duration(days=10)).alias("minus10d"), - (pl.col("rdq_trad") + pl.duration(days=5)).alias("plus5d"), -]) - -# Make sure the trading day version of rdq is within the window bounds -con = duckdb.connect(":memory:") -try: - con.execute(f"SET threads TO {num_threads};") - - con.register("ccm_data", ccm1.to_arrow()) - con.register("crsp_data", crsp_d.to_arrow()) - - df = con.execute(""" - SELECT a.gvkey, a.permno, a.datadate, a.fyearq, a.fqtr, - a.rdq, a.rdq_trad, b.date, b.abrd - FROM ccm_data a - LEFT JOIN crsp_data b - ON a.permno = b.permno - AND a.minus10d <= b.date - AND b.date <= a.plus5d - ORDER BY a.permno, a.rdq_trad, b.date - """).pl() - - # Filter out missing returns - df = df.filter(pl.col("abrd").is_not_null()) - - ############################### - # Count trading days # - ############################### - df = df.sort(["permno", "rdq_trad", "date"]) - - # Assign direction indicator: 0 = rdq, positive = after, negative = before - # This is used to count trading days before and after rdq_trad - df = df.with_columns( - pl.when(pl.col("date") == pl.col("rdq_trad")).then(pl.lit(0)) - .when(pl.col("date") > pl.col("rdq_trad")).then(pl.lit(1)) - .when(pl.col("date") < pl.col("rdq_trad")).then(pl.lit(-1)) - .alias("c_1") - ) - - # Trading days before rdq_trad (count descending: -1, -2, -3, ...) - df_before = ( - df.filter(pl.col("c_1") == -1) - .sort(["permno", "rdq_trad", "date"], descending=[False, False, True]) - .with_columns( - (-(pl.col("date").cum_count().over(["permno", "rdq_trad"]).cast(pl.Int64))).alias("count") - ) - .sort(["permno", "rdq_trad", "date"]) - ) - - # Trading days on or after rdq_trad (count: 0, 1, 2, ...) - df_after = ( - df.filter(pl.col("c_1") >= 0) - .with_columns( - (pl.col("date").cum_count().over(["permno", "rdq_trad"]).cast(pl.Int64) - 1).alias("count") - ) - ) - - df = pl.concat([df_before, df_after]) - - ############################### - # Calculate ABR # - ############################### - # Filter to event window [-2, +1] - df = df.filter((pl.col("count") >= -2) & (pl.col("count") <= 1)) - - # Sum abnormal returns by group - df_abr = ( - df.group_by(["permno", "rdq_trad"]) - .agg(pl.col("abrd").sum().alias("abr")) - ) - - # Join ABR back and keep only count == 1 rows (rdq + 1 day) - df = ( - df.join(df_abr, on=["permno", "rdq_trad"], how="left") - .filter(pl.col("count") == 1) - .rename({"date": "rdq_plus_1d"}) - .select(["gvkey", "permno", "datadate", "rdq", "rdq_plus_1d", "abr"]) - ) - - ############################### - # Populate to monthly # - ############################### - # Make sure the abr is used between rdq_plus_1d and plus12m - # Get monthly dates from CRSP monthly file - crsp_msf = ( - pl.scan_parquet(INPUT_PATH + "crsp_msf.parquet") - .select(pl.col("mthcaldt").cast(pl.Date).alias("date")) - .unique() - .collect() - ) - - # Add 12-month forward bound - df = df.with_columns( - pl.col("datadate").dt.offset_by("12mo").dt.month_end().alias("plus12m") - ) - - # Use DuckDB for the range join (reuse connection) - con.register("df_data", df.to_arrow()) - con.register("msf_data", crsp_msf.to_arrow()) - - df = con.execute(""" - SELECT a.gvkey, a.permno, a.datadate, a.rdq, a.rdq_plus_1d, a.abr, b.date - FROM df_data a - LEFT JOIN msf_data b - ON a.rdq_plus_1d < b.date - AND a.plus12m >= b.date - ORDER BY a.permno, b.date, a.datadate DESC - """).pl() -finally: - con.close() - -# Drop duplicates keeping first (most recent datadate per permno-month) -df = ( - df.unique(["permno", "date"], keep="first", maintain_order=True) - .filter(pl.col("date").is_not_null()) - .select(["gvkey", "permno", "datadate", "rdq", "rdq_plus_1d", "abr", "date"]) -) - -############################### -# Write output # -############################### -df.write_parquet(OUTPUT_PATH + "abr.parquet") -print(f"ABR data written to abr.parquet") diff --git a/chars_ciz_monthly/accounting.py b/chars_ciz_monthly/accounting.py deleted file mode 100644 index 034e837..0000000 --- a/chars_ciz_monthly/accounting.py +++ /dev/null @@ -1,2511 +0,0 @@ -import polars as pl -from functions import * -import datetime -import os - -# Configuration -INPUT_PATH = "../data/raw/" -OUTPUT_PATH = "../data/processed/" - -# Create output directory if it doesn't exist -os.makedirs(OUTPUT_PATH, exist_ok=True) - - -####################################################################################################################### -# Compustat Block # -####################################################################################################################### -comp = pl.read_parquet(INPUT_PATH + 'comp_funda.parquet') -# comp = pl.scan_parquet(INPUT_PATH + 'comp_funda.parquet') - -# cast all Decimal columns to Float64 (Decimal causes division by zero errors) -comp = comp.with_columns([ - pl.col(c).cast(pl.Float64) - for c in comp.columns - if str(comp[c].dtype).startswith('Decimal') -]) - -# convert datadate to date fmt and sort/clean up -comp = (comp - .with_columns([ - pl.col('datadate').cast(pl.Date) - ]) - .sort(['gvkey', 'datadate']) - .unique() -) - -# (fixed)2026-03-06: split into two with_columns to ensure mve_f uses cleaned csho -# clean up csho and calculate market equity -comp = comp.with_columns([ - # Replace 0 with null in csho - pl.when(pl.col('csho') == 0) - .then(None) - .otherwise(pl.col('csho')) - .alias('csho'), -]) - -comp = comp.with_columns([ - # Calculate Compustat market equity (now uses cleaned csho) - (pl.col('csho') * pl.col('prcc_f')).alias('mve_f') -]) - -# do some clean up for dr -comp = comp.with_columns([ - pl.when(pl.col('drc').is_not_null() & pl.col('drlt').is_not_null()) - .then(pl.col('drc') + pl.col('drlt')) - .when(pl.col('drc').is_not_null() & pl.col('drlt').is_null()) - .then(pl.col('drc')) - .when(pl.col('drlt').is_not_null() & pl.col('drc').is_null()) - .then(pl.col('drlt')) - .otherwise(None) - .alias('dr') -]) - -# do some clean up for dc -comp = comp.with_columns([ - pl.when(pl.col('dcvt').is_null() & - pl.col('dcpstk').is_not_null() & - pl.col('pstk').is_not_null() & - (pl.col('dcpstk') > pl.col('pstk'))) - .then(pl.col('dcpstk') - pl.col('pstk')) - .when(pl.col('dcvt').is_null() & - pl.col('dcpstk').is_not_null() & - pl.col('pstk').is_null()) - .then(pl.col('dcpstk')) - .otherwise(None) - .cast(pl.Float64) - .alias('dc') -]) - -# Fill dc with dcvt if dc is null -comp = comp.with_columns([ - pl.when(pl.col('dc').is_null()) - .then(pl.col('dcvt')) - .otherwise(pl.col('dc')) - .alias('dc') -]) - -# (removed) xint0/xsga0: moved to _ANNUAL_FILL_ZERO unified block (xint/xsga filled, aliased later) - -# Replace 0 with null in ceq and at, then filter out null at -comp = (comp - .with_columns([ - pl.when(pl.col('ceq') == 0).then(None).otherwise(pl.col('ceq')).alias('ceq'), - pl.when(pl.col('at') == 0).then(None).otherwise(pl.col('at')).alias('at') - ]) - .filter(pl.col('at').is_not_null()) -) - -comp = comp.rename({'cusip': 'cusip_comp'}) -####################################################################################################################### -# CRSP Block # -####################################################################################################################### -# Create a CRSP Subsample with Monthly Stock and Event Variables -# Restrictions will be applied later -# Select variables from the CRSP monthly stock and event datasets -crsp = pl.read_parquet(INPUT_PATH + 'crsp_msf.parquet') - -# rename cusip as cusip_crsp -crsp = crsp.rename({'cusip': 'cusip_crsp'}) - -# filter exchcd & shrcd -# equivalent to legacy code exchcd = 1, 2 or 3 -crsp = crsp.filter( - (pl.col('primaryexch').is_in(['N', 'A', 'Q'])) & - (pl.col('conditionaltype') == 'RW') & - (pl.col('tradingstatusflg') == 'A') -) -# crsp['exchcd'] = crsp['primaryexch'].map({'N': 1, 'A': 2, 'Q': 3}) - -# (fixed): control shrcd -crsp = crsp.filter( - (pl.col('sharetype') == 'NS') & - (pl.col('securitytype') == 'EQTY') & - (pl.col('securitysubtype') == 'COM') & - (pl.col('usincflg') == 'Y') & - (pl.col('issuertype').is_in(['ACOR', 'CORP'])) -) - - -# (fixed): usincflg='Y' (above) already restricts to US-incorporated issuers, -# which excludes China-incorporated ADRs (Baidu, JD, etc.) at source. -# If CIZ StkMthSecurityData ever exposes a 'primaryissue' or 'foreigncommonflg' field, -# add: crsp = crsp.filter(pl.col('foreigncommonflg') != 'Y') — but verify field availability first. - -# Mapping CIZ variables to SIZ varialbles -crsp = crsp.rename({ - 'mthprc': 'prc', - 'mthret': 'ret', - 'mthretx': 'retx', - 'mthvol': 'vol', - 'mthcumfacpr': 'cfacpr', - 'mthcumfacshr': 'cfacshr', - 'mthcaldt': 'date', - 'issuernm': 'comnam' -}) - -# change variable format to int -crsp = crsp.with_columns([ - pl.col('permco').cast(pl.Int64), - pl.col('permno').cast(pl.Int64) -]) - -# Line up date to be end of month -# set all the date to the standard end date of month -crsp = crsp.with_columns([ - pl.col('date').dt.month_end().alias('monthend') -]) - -# Drop nulls and calculate market equity -crsp = (crsp - .filter(pl.col('prc').is_not_null()) - .with_columns([ - (pl.col('prc').abs() * pl.col('shrout')).alias('me') - ]) -) - -# Unified fill_null(0) for CRSP data -# ret/retx: missing return treated as 0 (no trading / delisting handled separately) -_CRSP_FILL_ZERO = [ - 'ret', # monthly return - 'retx', # ex-dividend return -] -crsp = crsp.with_columns([ - pl.col(c).fill_null(0) for c in _CRSP_FILL_ZERO - if c in crsp.columns -]) - -# impute me - sort and deduplicate -crsp = crsp.sort(['permno', 'date']).unique() - -# Forward fill me within each permno group (only forward fill when same permno) -crsp = crsp.with_columns([ - pl.when(pl.col('permno') == pl.col('permno').shift(1)) - .then(pl.col('me').forward_fill().over('permno')) - .otherwise(pl.col('me')) - .alias('me') -]) - -# Aggregate Market Cap -''' -There are cases when the same firm (permco) has two or more securities (permno) at same date. -For the purpose of ME for the firm, we aggregated all ME for a given permco, date. -This aggregated ME will be assigned to the permno with the largest ME. -''' - -# sum of me across different permno belonging to same permco a given date -crsp_summe = (crsp - .group_by(['monthend', 'permco']) - .agg(pl.col('me').sum()) -) - -# largest mktcap within a permco/date -crsp_maxme = (crsp - .group_by(['monthend', 'permco']) - .agg(pl.col('me').max()) -) - -# join by monthend/maxme to find the permno -crsp1 = crsp.join(crsp_maxme, on=['monthend', 'permco', 'me'], how='inner') - -# join with sum of me to get the correct market cap info -# (no need to drop 'me' column first since we're joining and the sum will replace it) -crsp2 = (crsp1 - .drop('me') - .join(crsp_summe, on=['monthend', 'permco'], how='inner') - .sort(['permno', 'monthend']) - .unique() -) - -# Save full CRSP data for later use (momentum + ME-dependent characteristics) -# This avoids needing to reload CRSP later -crsp_full = crsp2.clone() - -# Create permno-only subset for initial CCM merge -# Only need permno and monthend (as jdate) for filtering the Compustat sample -crsp_permno_only = (crsp2 - .select(['permno', 'monthend']) - .rename({'monthend': 'jdate'}) - .unique() -) - -####################################################################################################################### -# CCM Block # - -####################################################################################################################### -# merge CRSP and Compustat -# reference: https://wrds-www.wharton.upenn.edu/pages/support/applications/linking-databases/linking-crsp-and-compustat/ -ccm = pl.read_parquet(INPUT_PATH + 'ccm.parquet') - -# convert the permno to int64 -ccm = ccm.with_columns([ - pl.col('permno').cast(pl.Int64) -]) - -# if linkenddt is missing then set to today date -ccm = ccm.with_columns([ - pl.col('linkenddt').fill_null(pl.lit(datetime.date.today())) -]) - -# merge ccm and comp -ccm1 = comp.join(ccm, on='gvkey', how='left') - -# we can only get the accounting data after the firm public their report -# for annual data, we use 4, 5 or 6 months lagged data, now we follow Hou, Xue and Zhang (2015) use 4 months lag -ccm1 = ccm1.with_columns([ - # Year end: December 31st of the same year - pl.date(pl.col('datadate').dt.year(), 12, 31).alias('yearend'), - # jdate: 4 months after datadate, then month end - pl.col('datadate').dt.offset_by('4mo').dt.month_end().alias('jdate') -]) - -# set link date bounds -ccm2 = ccm1.filter( - (pl.col('jdate') >= pl.col('linkdt')) & - (pl.col('jdate') <= pl.col('linkenddt')) -) - -# link comp and crsp (using permno only for initial sample restriction) -# Full CRSP data (me, ret, etc.) will be merged later for ME-dependent characteristics -data_rawa = ccm2.join(crsp_permno_only, on=['permno', 'jdate'], how='inner') - -# filter exchcd & shrcd and at least more than 1 year data -# Already filtered earlier in crsp - -# count single stock years -data_rawa = data_rawa.with_columns([ - (pl.col('gvkey').cum_count().over('gvkey')).alias('count') -]) - -# (fixed-20260316) check数据比例，对比去重逻辑执行前后的差别 -# 结果: 267524 行 → 267124 行，删去 400 行 (0.15%)。 -# 结论: 去重影响极小，但属于必要操作（去除同一 permno/datadate 的多条 link 重复行），保留不动。 - -# deal with the duplicates (align with data_rawq dedup logic) -# Keep first occurrence for each group of ['datadate', 'permno', 'linkprim'] -data_rawa = data_rawa.with_row_index('_temp_idx') -temp_first = (data_rawa - .group_by(['datadate', 'permno', 'linkprim'], maintain_order=True) - .agg(pl.col('_temp_idx').first()) -) -data_rawa = data_rawa.join(temp_first, on=['datadate', 'permno', 'linkprim', '_temp_idx'], how='semi').drop('_temp_idx') - -# Keep last occurrence for each group of ['permno', 'yearend', 'datadate'] -data_rawa = data_rawa.with_row_index('_temp_idx') -temp_last = (data_rawa - .group_by(['permno', 'yearend', 'datadate'], maintain_order=True) - .agg(pl.col('_temp_idx').last()) -) -data_rawa = data_rawa.join(temp_last, on=['permno', 'yearend', 'datadate', '_temp_idx'], how='semi').drop('_temp_idx') - -# Sort -data_rawa = data_rawa.sort(['permno', 'jdate']) - -# data_rawa.filter(data_rawa.is_duplicated)) - -# # Keep first occurrence within each group -# data_rawa = (data_rawa -# .with_columns([ -# pl.lit(1).alias('temp').over(['datadate', 'permno', 'linkprim']).cum_count() -# ]) -# .filter(pl.col('temp') == 1) -# .drop('temp') -# ) - -# # Keep last occurrence within each group -# data_rawa = (data_rawa -# .sort(['permno', 'yearend', 'datadate']) -# .with_columns([ -# pl.lit(1).alias('temp').over(['permno', 'yearend', 'datadate']).cum_count() -# ]) -# .with_columns([ -# pl.col('temp').max().over(['permno', 'yearend', 'datadate']).alias('max_temp') -# ]) -# .filter(pl.col('temp') == pl.col('max_temp')) -# .drop(['temp', 'max_temp']) -# ) - -# Sort -data_rawa = data_rawa.sort(['permno', 'jdate']) - -# Unified fill_null(0) for annual data -# These columns are used in formulas where null should be treated as 0 -_ANNUAL_FILL_ZERO = [ - 'ps', 'txditc', 'cogs', 'xint', 'xsga', - 'ivao', 'dlc', 'dltt', 'mib', 'pstk', - 'gdwl', 'intan', 'che', 'act', 'at', - 'dp', 'txp', 'aco', 'ao', 'ap', 'lco', - 'lo', 'rect','invt', 'ppent' -] -data_rawa = data_rawa.with_columns([ - pl.col(c).fill_null(0) for c in _ANNUAL_FILL_ZERO - if c in data_rawa.columns -]) - -# fama-french 49 industry -data_rawa = data_rawa.with_columns([ - pl.col('sic').cast(pl.Int64) -]) - -# Apply ffi49 function (assuming it returns a Series/column) -data_rawa = data_rawa.with_columns([ - ffi49().alias('ffi49') -]) - -data_rawa = data_rawa.with_columns([ - pl.col('ffi49').fill_null(49).cast(pl.Int64).alias('ffi49') -]) - -####################################################################################################################### -# Annual Variables # -####################################################################################################################### -# preferrerd stock -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('pstkrv').is_null()) - .then(pl.col('pstkl')) - .otherwise(pl.col('pstkrv')) - .alias('ps') -]) - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('ps').is_null()) - .then(pl.col('pstk')) - .otherwise(pl.col('ps')) - .alias('ps') -]) - -# (HXZ): "Stockholders' equity is the value reported by Compustat (item SEQ), if it is available. If not, we measure stockholders' equity as the book value of common equity (item CEQ) plus the par value of preferred stock (item PSTK), or the book value of assets (item AT) minus total liabilities (item LT)." - -# book equity -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('seq').is_not_null()).then(pl.col('seq')) - .when(pl.col('ceq').is_not_null() & pl.col('pstk').is_not_null()) - .then(pl.col('ceq') + pl.col('pstk')) - .otherwise(pl.col('at') - pl.col('lt')) - .alias('seq') -]) -data_rawa = data_rawa.with_columns([ - (pl.col('seq') + pl.col('txditc') - pl.col('ps')).alias('be') -]) - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('be') > 0) - .then(pl.col('be')) - .otherwise(None) - .alias('be') -]) - -# acc - lagged variables -data_rawa = data_rawa.with_columns([ - pl.col('act').shift(1).over('permno').alias('act_l1'), - pl.col('lct').shift(1).over('permno').alias('lct_l1'), - pl.col('at').shift(1).over('permno').alias('at_l1') -]) - -# #################### Add np lag (also fixed row 272 below) on 2025.02.23 #################### -# data_rawa['np_l1'] = data_rawa.groupby(['permno'])['np'].shift(1) - -# condlist = [data_rawa['np'].isnull(), -# data_rawa['act'].isnull() | data_rawa['lct'].isnull()] -# choicelist = [((data_rawa['act'] - data_rawa['lct']) - (data_rawa['act_l1'] - data_rawa['lct_l1']) / (data_rawa['be'])), -# (data_rawa['ib'] - data_rawa['oancf']) / (data_rawa['be'])] ##### Delete "10*" on 2025.02.26 ##### -# data_rawa['acc'] = np.select(condlist, -# choicelist, -# default=((data_rawa['act'] - data_rawa['lct'] + data_rawa['np']) - -# (data_rawa['act_l1'] - data_rawa['lct_l1'] + data_rawa['np_l1'])) / (data_rawa['be'])) - -#################### Add Sloan(1996) or HXZ and GHZ operating accruals on 2025.02.28 #################### -# More lagged variables -data_rawa = data_rawa.with_columns([ - pl.col('che').shift(1).over('permno').alias('che_l1'), - pl.col('dlc').shift(1).over('permno').alias('dlc_l1'), - pl.col('txp').shift(1).over('permno').alias('txp_l1') -]) -# txp is 0-filled; fill its lag to 0 as well (handles first obs per firm) -data_rawa = data_rawa.with_columns([ - pl.col('txp_l1').fill_null(0) -]) - -# acc calculation -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('oancf').is_null()) - .then( - (((pl.col('act') - pl.col('act_l1')) - (pl.col('che') - pl.col('che_l1')) - - (pl.col('lct') - pl.col('lct_l1')) + (pl.col('dlc') - pl.col('dlc_l1')) + - (pl.col('txp') - pl.col('txp_l1')) - pl.col('dp')) / - ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - ) - .otherwise( - (pl.col('ni') - pl.col('oancf')) / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None) - ) - .alias('acc') -]) - -# absacc -data_rawa = data_rawa.with_columns([ - pl.col('acc').abs().alias('absacc') -]) - -# agr -data_rawa = data_rawa.with_columns([ - ((pl.col('at') - pl.col('at_l1')) / pl.col('at_l1').replace(0, None)).alias('agr') -]) - -# bm -# data_rawa['bm'] = data_rawa['be'] / data_rawa['me'] - -# cfp -# condlist = [data_rawa['dp'].isnull(), -# data_rawa['ib'].isnull()] -# choicelist = [data_rawa['ib']/data_rawa['me'], -# np.nan] -# data_rawa['cfp'] = np.select(condlist, choicelist, default=(data_rawa['ib']+data_rawa['dp'])/data_rawa['me']) - -# ep -# data_rawa['ep'] = data_rawa['ib']/data_rawa['me'] - -# ni -data_rawa = data_rawa.with_columns([ - pl.col('csho').shift(1).over('permno').alias('csho_l1'), - pl.col('ajex').shift(1).over('permno').alias('ajex_l1') -]) - -# log() result: fill_nan(0) handles log(0)→−inf/nan, fill_null(0) handles null input -# order: fill_nan first then fill_null -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('gvkey') != pl.col('gvkey').shift(1).over('permno')) # 2026-02-12: add over - .then(None) - .otherwise( - (pl.col('csho') * pl.col('ajex')).log() - .fill_nan(0) - .fill_null(0) - - (pl.col('csho_l1') * pl.col('ajex_l1')).log() - .fill_nan(0) - .fill_null(0) - ) - .alias('ni') -]) - -# op -# cogs / xint / xsga are already 0-filled via _ANNUAL_FILL_ZERO; no alias needed -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('revt').is_null()) - .then(None) - .when(pl.col('be').is_null()) - .then(None) - .otherwise( - (pl.col('revt') - pl.col('cogs') - pl.col('xsga') - pl.col('xint')) / pl.col('be').replace(0, None) - ) - .alias('op') -]) - -# rsup -data_rawa = data_rawa.with_columns([ - pl.col('sale').shift(1).over('permno').alias('sale_l1') -]) -# data_rawa['rsup'] = (data_rawa['sale']-data_rawa['sale_l1'])/data_rawa['me'] - -# cash -data_rawa = data_rawa.with_columns([ - (pl.col('che') / pl.col('at').replace(0, None)).alias('cash') -]) - -# lev -# data_rawa['lev'] = data_rawa['lt']/data_rawa['me'] - -# sp -# data_rawa['sp'] = data_rawa['sale']/data_rawa['me'] - -# rd_sale -data_rawa = data_rawa.with_columns([ - (pl.col('xrd') / pl.col('sale').replace(0, None)).alias('rd_sale') -]) - -# rdm -# data_rawa['rdm'] = data_rawa['xrd']/data_rawa['me'] - -# adm hxz adm -# data_rawa['adm'] = data_rawa['xad']/data_rawa['me'] - -# gma -data_rawa = data_rawa.with_columns([ - ((pl.col('revt') - pl.col('cogs')) / pl.col('at_l1').replace(0, None)).alias('gma') -]) - -# chcsho -data_rawa = data_rawa.with_columns([ - ((pl.col('csho') / pl.col('csho_l1').replace(0, None)) - 1).alias('chcsho') -]) - -# lgr -data_rawa = data_rawa.with_columns([ - pl.col('lt').shift(1).over('permno').alias('lt_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('lt') / pl.col('lt_l1').replace(0, None)) - 1).alias('lgr') -]) - -#################### Follow Hafzalla, Lundholm, and Van Winkle (2011) and GHZ on 2025.02.28 #################### -# pctacc (follow HXZ-A3.22) -data_rawa = data_rawa.with_columns([ - (pl.col('acc') / pl.col('ni').abs().replace(0, None)).alias('pctacc') -]) - -# sgr -data_rawa = data_rawa.with_columns([ - ((pl.col('sale') / pl.col('sale_l1').replace(0, None)) - 1).alias('sgr') -]) - -# chato -# (fixed-20260316) `OnlineAppendixOPCSAP.pdf` (pp. 40-41) defines AssetTurnover as sales divided by two-year average net operating assets, and ChAssetTurnover as the annual change in that ratio. HXZ decomposes RNA into PM × ATO using NOA-based turnover, so using NOA keeps chato, noa, rna, and ato internally consistent. -# (fixed-20260330) align noa_raw with HXZ-style operating assets definition by excluding ivao/ivaoq from operating assets. This keeps noa/chato/rna/ato on the same base and matches the legacy benchmark more closely. -data_rawa = data_rawa.with_columns([ - ((pl.col('at') - pl.col('che') - pl.col('ivao')) - - (pl.col('at') - pl.col('dlc') - pl.col('dltt') - - pl.col('mib') - pl.col('pstk') - pl.col('ceq'))).alias('noa_raw') -]) -data_rawa = data_rawa.with_columns([ - pl.col('noa_raw').shift(1).over('permno').alias('noa_raw_l1'), - pl.col('noa_raw').shift(2).over('permno').alias('noa_raw_l2') -]) -# ATO_t = sale_t / avg(noa_raw_t, noa_raw_{t-1}); -# ATO_{t-1} = sale_{t-1} / avg(noa_raw_{t-1}, noa_raw_{t-2}) -data_rawa = data_rawa.with_columns([ - ((pl.col('sale') / ((pl.col('noa_raw') + pl.col('noa_raw_l1')) / 2).replace(0, None)) - - (pl.col('sale_l1') / ((pl.col('noa_raw_l1') + pl.col('noa_raw_l2')) / 2).replace(0, None))).alias('chato') -]) - -# chtx -data_rawa = data_rawa.with_columns([ - pl.col('txt').shift(1).over('permno').alias('txt_l1') -]) -data_rawa = data_rawa.with_columns([ - ((pl.col('txt') - pl.col('txt_l1')) / pl.col('at_l1').replace(0, None)).alias('chtx') -]) - -# noa -# delete fill_null(0) -# noa_raw is computed above so that chato, rna, and ato all use the same -# NOA definition from the reference documents. -data_rawa = data_rawa.with_columns([ - (pl.col('noa_raw') / pl.col('at_l1').replace(0, None)).alias('noa') -]) - -# rna -# (fix)2026-02-27: use noa_raw (unscaled) as denominator instead of noa (scaled by at_l1) -data_rawa = data_rawa.with_columns([ - ((pl.col('oiadp') / pl.col('noa_raw_l1').replace(0, None))).alias('rna') -]) - -# pm -data_rawa = data_rawa.with_columns([ - (pl.col('oiadp') / pl.col('sale').replace(0, None)).alias('pm') -]) - -# ato -# (fix)2026-02-27: use noa_raw (unscaled) as denominator instead of noa (scaled by at_l1) -data_rawa = data_rawa.with_columns([ - ((pl.col('sale') / pl.col('noa_raw_l1').replace(0, None))).alias('ato') -]) - -# depr -data_rawa = data_rawa.with_columns([ - (pl.col('dp') / pl.col('ppent').replace(0, None)).alias('depr') -]) - -# invest -data_rawa = data_rawa.with_columns([ - pl.col('ppent').shift(1).over('permno').alias('ppent_l1'), - pl.col('invt').shift(1).over('permno').alias('invt_l1'), - pl.col('ppegt').shift(1).over('permno').alias('ppegt_l1') # add ppegt_l1 -]) -# (fixed)2026-03-06: replace ppent_l1 with ppegt_l1 for consistency. -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('ppegt').is_null()) - .then( - ((pl.col('ppent') - pl.col('ppent_l1')) + - (pl.col('invt') - pl.col('invt_l1'))) / pl.col('at_l1').replace(0, None) - ) - .otherwise( - ((pl.col('ppegt') - pl.col('ppegt_l1')) + - (pl.col('invt') - pl.col('invt_l1'))) / pl.col('at_l1').replace(0, None) - ) - .alias('invest') -]) - -# egr -data_rawa = data_rawa.with_columns([ - pl.col('ceq').shift(1).over('permno').alias('ceq_l1') -]) -data_rawa = data_rawa.with_columns([ - ((pl.col('ceq') - pl.col('ceq_l1')) / pl.col('ceq_l1').replace(0, None)).alias('egr') -]) - -# cashdebt -data_rawa = data_rawa.with_columns([ - ((pl.col('ib') + pl.col('dp')) / - ((pl.col('lt') + pl.col('lt_l1')) / 2).replace(0, None)).alias('cashdebt') -]) - -# rd -data_rawa = data_rawa.with_columns([ - (pl.col('xrd') / pl.col('at_l1').replace(0, None)).alias('xrd/at_l1') -]) -data_rawa = data_rawa.with_columns([ - pl.col('xrd/at_l1').shift(1).over('permno').alias('xrd/at_l1_l1') -]) -data_rawa = data_rawa.with_columns([ - pl.when( - (((pl.col('xrd') / pl.col('at').replace(0, None)) - pl.col('xrd/at_l1_l1')) / - pl.col('xrd/at_l1_l1').replace(0, None)) > 0.05 - ) - .then(1) - .otherwise(0) - .alias('rd') -]) - -# roa -data_rawa = data_rawa.with_columns([ - (pl.col('ib') / pl.col('at_l1').replace(0, None)).alias('roa') -]) - -# roe -data_rawa = data_rawa.with_columns([ - (pl.col('ib') / pl.col('ceq_l1').replace(0, None)).alias('roe') -]) - -# dy -# data_rawa['dy'] = data_rawa['dvt']/data_rawa['me'] - -################## Added on 2020.07.28 ################## - -# roic -data_rawa = data_rawa.with_columns([ - ((pl.col('ebit') - pl.col('nopi')) / (pl.col('ceq') + pl.col('lt') - pl.col('che')).replace(0, None)).alias('roic') -]) - -# chinv -data_rawa = data_rawa.with_columns([ - ((pl.col('invt') - pl.col('invt_l1')) / - ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None) # HXZ(A.3.15) - ).alias('chinv') -]) - -# pchsale_pchinvt -data_rawa = data_rawa.with_columns([ - (((pl.col('sale') - pl.col('sale_l1')) / pl.col('sale_l1').replace(0, None)) - - ((pl.col('invt') - pl.col('invt_l1')) / pl.col('invt_l1').replace(0, None)) - ).alias('pchsale_pchinvt') -]) - -# pchsale_pchrect -data_rawa = data_rawa.with_columns([ - pl.col('rect').shift(1).over('permno').alias('rect_l1') -]) - -data_rawa = data_rawa.with_columns([ - (((pl.col('sale') - pl.col('sale_l1')) / pl.col('sale_l1').replace(0, None)) - - ((pl.col('rect') - pl.col('rect_l1')) / pl.col('rect_l1').replace(0, None)) - ).alias('pchsale_pchrect') -]) - -# pchgm_pchsale -data_rawa = data_rawa.with_columns([ - pl.col('cogs').shift(1).over('permno').alias('cogs_l1') -]) -# (fixed)2026-03-06: replace sale with sale_l1 for consistency. -data_rawa = data_rawa.with_columns([ - ((((pl.col('sale') - pl.col('cogs')) - (pl.col('sale_l1') - pl.col('cogs_l1'))) / - (pl.col('sale_l1') - pl.col('cogs_l1')).replace(0, None)) - - ((pl.col('sale') - pl.col('sale_l1')) / pl.col('sale_l1').replace(0, None)) - ).alias('pchgm_pchsale') -]) - -# pchsale_pchxsga -data_rawa = data_rawa.with_columns([ - pl.col('xsga').shift(1).over('permno').alias('xsga_l1') -]) - -data_rawa = data_rawa.with_columns([ - (((pl.col('sale') - pl.col('sale_l1')) / pl.col('sale_l1').replace(0, None)) - - ((pl.col('xsga') - pl.col('xsga_l1')) / pl.col('xsga_l1').replace(0, None)) - ).alias('pchsale_pchxsga') -]) - -# pchdepr -data_rawa = data_rawa.with_columns([ - pl.col('dp').shift(1).over('permno').alias('dp_l1') -]) -data_rawa = data_rawa.with_columns([ - (((pl.col('dp') / pl.col('ppent').replace(0, None)) - - (pl.col('dp_l1') / pl.col('ppent_l1').replace(0, None))) / - (pl.col('dp_l1') / pl.col('ppent_l1').replace(0, None)).replace(0, None)) - .alias('pchdepr') -]) - -# chadv: (Lou, 2014) https://academic.oup.com/rfs/article/27/6/1797/1596985#114323634 -data_rawa = data_rawa.with_columns([ - pl.col('xad').shift(1).over('permno').alias('xad_l1') -]) -# (fixed-20260322) first when, then check both current and lagged xad for null and >= 0.1, then compute log difference; otherwise set to null. -data_rawa = data_rawa.with_columns([ - pl.when( - pl.col('xad').is_not_null() & pl.col('xad_l1').is_not_null() & - (pl.col('xad') > 0.1) & (pl.col('xad_l1') > 0.1) - ) - .then(pl.col('xad').log() - pl.col('xad_l1').log()) - .otherwise(None) - .alias('chadv') -]) - -# pchcapx -data_rawa = data_rawa.with_columns([ - pl.col('capx').shift(1).over('permno').alias('capx_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('capx') - pl.col('capx_l1')) / pl.col('capx_l1').replace(0, None)) - .alias('pchcapx') -]) - -# grcapx (GHZ method) -data_rawa = data_rawa.with_columns([ - pl.col('capx').shift(2).over('permno').alias('capx_l2') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('capx') - pl.col('capx_l2')) / pl.col('capx_l2').replace(0, None)) - .alias('grcapx') -]) - -# grGW -data_rawa = data_rawa.with_columns([ - pl.col('gdwl').shift(1).over('permno').alias('gdwl_l1') -]) -data_rawa = data_rawa.with_columns([ - ((pl.col('gdwl') - pl.col('gdwl_l1')) / pl.col('gdwl_l1').replace(0, None)) - .alias('grGW') -]) - -data_rawa = data_rawa.with_columns([ - pl.when((pl.col('gdwl') == 0) | pl.col('gdwl').is_null()) - .then(0) - .when(pl.col('gdwl').is_not_null() & (pl.col('gdwl') != 0) & pl.col('grGW').is_null()) - .then(1) - .otherwise(pl.col('grGW')) - .alias('grGW') -]) - -# currat -data_rawa = data_rawa.with_columns([ - (pl.col('act') / pl.col('lct').replace(0, None)).alias('currat') -]) - -# pchcurrat -data_rawa = data_rawa.with_columns([ - (((pl.col('act') / pl.col('lct').replace(0, None)) - - (pl.col('act_l1') / pl.col('lct_l1').replace(0, None))) / - (pl.col('act_l1') / pl.col('lct_l1').replace(0, None)).replace(0, None)) - .alias('pchcurrat') -]) - -# quick -data_rawa = data_rawa.with_columns([ - ((pl.col('act') - pl.col('invt')) / pl.col('lct').replace(0, None)).alias('quick') -]) - -# pchquick -data_rawa = data_rawa.with_columns([ - ((((pl.col('act') - pl.col('invt')) / pl.col('lct').replace(0, None)) - - ((pl.col('act_l1') - pl.col('invt_l1')) / pl.col('lct_l1').replace(0, None))) / - ((pl.col('act_l1') - pl.col('invt_l1')) / pl.col('lct_l1').replace(0, None)).replace(0, None)) - .alias('pchquick') -]) - -# salecash -data_rawa = data_rawa.with_columns([ - (pl.col('sale') / pl.col('che').replace(0, None)).alias('salecash') -]) - -# salerec -data_rawa = data_rawa.with_columns([ - (pl.col('sale') / pl.col('rect').replace(0, None)).alias('salerec') -]) - -# saleinv -data_rawa = data_rawa.with_columns([ - (pl.col('sale') / pl.col('invt').replace(0, None)).alias('saleinv') -]) - -# pchsaleinv -data_rawa = data_rawa.with_columns([ - (((pl.col('sale') / pl.col('invt').replace(0, None)) - (pl.col('sale_l1') / pl.col('invt_l1').replace(0, None))) / - (pl.col('sale_l1') / pl.col('invt_l1').replace(0, None)).replace(0, None)).alias('pchsaleinv') -]) - -# realestate -data_rawa = data_rawa.with_columns([ - ((pl.col('fatb') + pl.col('fatl')) / pl.col('ppegt').replace(0, None)).alias('realestate') -]) - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('ppegt').is_null()) - .then((pl.col('fatb') + pl.col('fatl')) / pl.col('ppent').replace(0, None)) - .otherwise(pl.col('realestate')) - .alias('realestate') -]) - -# obklg -data_rawa = data_rawa.with_columns([ - (pl.col('ob') / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)).alias('obklg') -]) - -# chobklg -data_rawa = data_rawa.with_columns([ - pl.col('ob').shift(1).over('permno').alias('ob_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('ob') - pl.col('ob_l1')) / - ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)).alias('chobklg') -]) - -# grltnoa -data_rawa = data_rawa.with_columns([ - pl.col('aco').shift(1).over('permno').alias('aco_l1'), - pl.col('intan').shift(1).over('permno').alias('intan_l1'), - pl.col('ao').shift(1).over('permno').alias('ao_l1'), - pl.col('ap').shift(1).over('permno').alias('ap_l1'), - pl.col('lco').shift(1).over('permno').alias('lco_l1'), - pl.col('lo').shift(1).over('permno').alias('lo_l1'), - pl.col('rect').shift(1).over('permno').alias('rect_l1'), - pl.col('invt').shift(1).over('permno').alias('invt_l1'), - pl.col('ppent').shift(1).over('permno').alias('ppent_l1') -]) -# LTNOA_t and LTNOA_{t-1} -data_rawa = data_rawa.with_columns([ - ( - pl.col('rect') + pl.col('invt') + pl.col('ppent') + - pl.col('aco') + pl.col('intan') + pl.col('ao') - - pl.col('ap') - pl.col('lco') - pl.col('lo') - ).alias('ltnoa_t'), - ( - pl.col('rect_l1') + pl.col('invt_l1') + pl.col('ppent_l1') + - pl.col('aco_l1') + pl.col('intan_l1') + pl.col('ao_l1') - - pl.col('ap_l1') - pl.col('lco_l1') - pl.col('lo_l1') - ).alias('ltnoa_l1') -]) -# Working-capital operating accrual component -data_rawa = data_rawa.with_columns([ - ( - (pl.col('rect') - pl.col('rect_l1')) + - (pl.col('invt') - pl.col('invt_l1')) + - (pl.col('aco') - pl.col('aco_l1')) - - ( - (pl.col('ap') - pl.col('ap_l1')) + - (pl.col('lco') - pl.col('lco_l1')) - ) - ).alias('wcnoa_change') -]) -# Fairfield et al. style grltnoa: -# LTNOA_t / at_t - LTNOA_{t-1} / at_{t-1} - WCNOA_change / avg(at_t, at_{t-1}) -data_rawa = data_rawa.with_columns([ - ( - (pl.col('ltnoa_t') / pl.col('at').replace(0, None)) - - (pl.col('ltnoa_l1') / pl.col('at_l1').replace(0, None)) - - (pl.col('wcnoa_change') / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - ).alias('grltnoa') -]) - -# conv -data_rawa = data_rawa.with_columns([ - (pl.col('dc') / pl.col('dltt').replace(0, None)).alias('conv') -]) - -# convind -data_rawa = data_rawa.with_columns([ - pl.when( - ((pl.col('dc').is_not_null()) & (pl.col('dc') != 0)) | - ((pl.col('cshrc').is_not_null()) & (pl.col('cshrc') != 0)) - ) - .then(1) - .otherwise(0) - .alias('convind') -]) - -# chdrc -data_rawa = data_rawa.with_columns([ - pl.col('dr').shift(1).over('permno').alias('dr_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('dr') - pl.col('dr_l1')) / - ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - .alias('chdrc') -]) - -# rdbias -data_rawa = data_rawa.with_columns([ - pl.col('xrd').shift(1).over('permno').alias('xrd_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('xrd') / pl.col('xrd_l1').replace(0, None)) - 1 - - (pl.col('ib') / pl.col('ceq_l1').replace(0, None))) - .alias('rdbias') -]) - -# operprof -# cogs / xint / xsga are already 0-filled via _ANNUAL_FILL_ZERO -data_rawa = data_rawa.with_columns([ - ((pl.col('revt') - pl.col('cogs') - pl.col('xsga') - pl.col('xint')) / - pl.col('ceq_l1').replace(0, None)) - .alias('operprof') -]) - -# cfroa -data_rawa = data_rawa.with_columns([ - (pl.col('oancf') / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - .alias('cfroa') -]) - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('oancf').is_null()) - .then( - (pl.col('ib') + pl.col('dp')) / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None) - ) - .otherwise(pl.col('cfroa')) - .alias('cfroa') -]) - -# xrdint (HXZ-A.5.4, but denominator is market equity) -data_rawa = data_rawa.with_columns([ - (pl.col('xrd') / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - .alias('xrdint') -]) - -# capxint -data_rawa = data_rawa.with_columns([ - (pl.col('capx') / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - .alias('capxint') -]) - -# xadint (HXZ-A.5.2) -data_rawa = data_rawa.with_columns([ - (pl.col('xad') / ((pl.col('at') + pl.col('at_l1')) / 2).replace(0, None)) - .alias('xadint') -]) - -# chpm -data_rawa = data_rawa.with_columns([ - pl.col('ib').shift(1).over('permno').alias('ib_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('ib') / pl.col('sale').replace(0, None)) - - (pl.col('ib_l1') / pl.col('sale_l1').replace(0, None))).alias('chpm') -]) - -# ala -data_rawa = data_rawa.with_columns([ - (pl.col('che') + 0.75 * (pl.col('act') - pl.col('che')) - - 0.5 * (pl.col('at') - pl.col('act') - pl.col('gdwl') - pl.col('intan'))).alias('ala') -]) - -# alm -data_rawa = data_rawa.with_columns([ - (pl.col('ala') / - (pl.col('at') + pl.col('prcc_f') * pl.col('csho') - pl.col('ceq')).replace(0, None)) - .alias('alm') -]) - -# hire -data_rawa = data_rawa.with_columns([ - pl.col('emp').shift(1).over('permno').alias('emp_l1') -]) - -data_rawa = data_rawa.with_columns([ - ((pl.col('emp') - pl.col('emp_l1')) / pl.col('emp_l1').replace(0, None)).alias('hire') -]) - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('emp').is_null() | pl.col('emp_l1').is_null()) - .then(0) - .otherwise(pl.col('hire')) - .alias('hire') -]) - -# herf -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('sale').sum().alias('indsale')) -) - -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - ((pl.col('sale') / pl.col('indsale').replace(0, None)) * - (pl.col('sale') / pl.col('indsale').replace(0, None))) - .alias('herf') -]) - -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('herf').sum()) -) - -data_rawa = data_rawa.drop('herf') -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') -################## Added on 2022.09.06 ################## -# age -data_rawa = data_rawa.with_columns([ - pl.col('count').alias('age') -]) - -# cashpr -# data_rawa['cashpr'] = ((data_rawa['me'] + data_rawa['dltt'] - data_rawa['at']) / data_rawa['che']) - -# chempia -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('hire').mean().alias('hire_ind')) -) - -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - (pl.col('hire') - pl.col('hire_ind')).alias('chempia') -]) - -# chpmia -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('chpm').mean().alias('chpm_ind')) -) - -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - (pl.col('chpm') - pl.col('chpm_ind')).alias('chpmia') -]) - -# chatoia -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('chato').mean().alias('chato_ind')) -) - -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - (pl.col('chato') - pl.col('chato_ind')).alias('chatoia') -]) - -# divi -data_rawa = data_rawa.with_columns([ - pl.col('dvt').shift(1).over('permno').alias('dvt_l1') -]) - -data_rawa = data_rawa.with_columns([ - pl.when( - (pl.col('dvt').is_not_null()) & (pl.col('dvt') > 0) & - ((pl.col('dvt_l1') == 0) | pl.col('dvt_l1').is_null()) - ) - .then(1) - .otherwise(0) - .alias('divi') -]) - -# divo -# (dvt=0 or null) dvt_l1>0--> divo=1 -# (fix)2026-02-27: if dvt_l1=0, dvt>0, divo should be 0. The previous version was wrong since it treated dvt_l1=0 as dvt_l1 is null, which caused divo to be 1 when dvt_l1=0 and dvt>0, which is not correct since divo should be 0 in this case. -data_rawa = data_rawa.with_columns([ - pl.when( - (pl.col('dvt').is_null() | (pl.col('dvt') == 0)) & - ((pl.col('dvt_l1') > 0) & pl.col('dvt_l1').is_not_null()) - ) - .then(1) - .otherwise(0) - .alias('divo') -]) - -# Mohanram (2005) score (Annual Related) -df_temp = (data_rawa - .group_by(['fyear', 'ffi49']) - .agg(pl.col('roa').median().alias('md_roa')) -) -data_rawa = data_rawa.join(df_temp, on=['fyear', 'ffi49'], how='left') - -df_temp = (data_rawa - .group_by(['fyear', 'ffi49']) - .agg(pl.col('cfroa').median().alias('md_cfroa')) -) -data_rawa = data_rawa.join(df_temp, on=['fyear', 'ffi49'], how='left') - -df_temp = (data_rawa - .group_by(['fyear', 'ffi49']) - .agg(pl.col('oancf').median().alias('md_oancf')) -) -data_rawa = data_rawa.join(df_temp, on=['fyear', 'ffi49'], how='left') - -df_temp = (data_rawa - .group_by(['fyear', 'ffi49']) - .agg(pl.col('xrdint').median().alias('md_xrdint')) -) -data_rawa = data_rawa.join(df_temp, on=['fyear', 'ffi49'], how='left') - -df_temp = (data_rawa - .group_by(['fyear', 'ffi49']) - .agg(pl.col('capxint').median().alias('md_capxint')) -) -data_rawa = data_rawa.join(df_temp, on=['fyear', 'ffi49'], how='left') - -df_temp = (data_rawa - .group_by(['fyear', 'ffi49']) - .agg(pl.col('xadint').median().alias('md_xadint')) -) -data_rawa = data_rawa.join(df_temp, on=['fyear', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('roa') > pl.col('md_roa')).then(1).otherwise(0).alias('m1'), - pl.when(pl.col('cfroa') > pl.col('md_cfroa')).then(1).otherwise(0).alias('m2'), - pl.when(pl.col('oancf') > pl.col('md_oancf')).then(1).otherwise(0).alias('m3'), - pl.when(pl.col('xrdint') > pl.col('md_xrdint')).then(1).otherwise(0).alias('m4'), - pl.when(pl.col('capxint') > pl.col('md_capxint')).then(1).otherwise(0).alias('m5'), - pl.when(pl.col('xadint') > pl.col('md_xadint')).then(1).otherwise(0).alias('m6') -]) - -# pchcapx_ia -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('pchcapx').mean().alias('pchcapx_ind')) -) - -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - (pl.col('pchcapx') - pl.col('pchcapx_ind')).alias('pchcapx_ia') -]) - -# secured -data_rawa = data_rawa.with_columns([ - (pl.col('dm') / pl.col('dltt').replace(0, None)).alias('secured') -]) - -# securedind -data_rawa = data_rawa.with_columns([ - pl.when((pl.col('dm').is_not_null()) & (pl.col('dm') != 0)) - .then(1) - .otherwise(0) - .alias('securedind') -]) - -# sin -data_rawa = data_rawa.with_columns([ - pl.when( - ((pl.col('sic') >= 2100) & (pl.col('sic') <= 2199)) | - ((pl.col('sic') >= 2080) & (pl.col('sic') <= 2085)) | - (pl.col('naics') == '7132') | - (pl.col('naics') == '71312') | - (pl.col('naics') == '713210') | - (pl.col('naics') == '71329') | - (pl.col('naics') == '713290') | - (pl.col('naics') == '72112') | - (pl.col('naics') == '721120') - ) - .then(1) - .otherwise(0) - .alias('sin') -]) - -# tang -data_rawa = data_rawa.with_columns([ - ((pl.col('che') + pl.col('rect') * 0.715 + - pl.col('invt') * 0.547 + pl.col('ppent') * 0.535) / pl.col('at').replace(0, None)) - .alias('tang') -]) - -# tb, Lev and Nissim (2004) -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('fyear') <= 1978) - .then(0.48) - .when((pl.col('fyear') >= 1979) & (pl.col('fyear') <= 1986)) - .then(0.46) - .when(pl.col('fyear') == 1987) - .then(0.4) - .when((pl.col('fyear') >= 1988) & (pl.col('fyear') <= 1992)) - .then(0.34) - .when(pl.col('fyear') >= 1993) - .then(0.35) - .otherwise(None) - .alias('tr') -]) - -data_rawa = data_rawa.with_columns([ - (((pl.col('txfo') + pl.col('txfed').replace(0, None)) / - pl.col('tr').replace(0, None)) / - pl.col('ib').replace(0, None)) - .alias('tb_1') -]) - -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('txfo').is_null() | pl.col('txfed').is_null()) - .then( - ((pl.col('txt') - pl.col('txdi')) / pl.col('tr').replace(0, None)) / - pl.col('ib').replace(0, None) - ) - .otherwise(pl.col('tb_1')) - .alias('tb_1') -]) - -data_rawa = data_rawa.with_columns([ - pl.when( - (((pl.col('txfo') + pl.col('txfed') > 0) | (pl.col('txt') > pl.col('txdi'))) & - (pl.col('ib') <= 0)) - ) - .then(1) - .otherwise(pl.col('tb_1')) - .alias('tb_1') -]) - -df_temp = (data_rawa - .group_by(['datadate', 'ffi49']) - .agg(pl.col('tb_1').mean().alias('tb_1_ind')) -) - -data_rawa = data_rawa.join(df_temp, on=['datadate', 'ffi49'], how='left') - -data_rawa = data_rawa.with_columns([ - (pl.col('tb_1') - pl.col('tb_1_ind')).alias('tb') -]) - -print("Finish Annual Variables Calculation! \n") - -####################################################################################################################### -# Compustat Quarterly Raw Info # -####################################################################################################################### -comp = pl.read_parquet(INPUT_PATH + 'comp_fundq.parquet') - -# cast all Decimal columns to Float64 -comp = comp.with_columns([ - pl.col(c).cast(pl.Float64) - for c in comp.columns - if str(comp[c].dtype).startswith('Decimal') -]) - -# rename cusip as cusip_comp -comp = comp.rename({'cusip': 'cusip_comp'}) - -# comp['cusip6'] = comp['cusip'].str.strip().str[0:6] -comp = comp.filter(pl.col('ibq').is_not_null()) - -# sort and clean up -comp = comp.sort(['gvkey', 'datadate']).unique() -comp = comp.with_columns([ - pl.when(pl.col('cshoq') == 0).then(None).otherwise(pl.col('cshoq')).alias('cshoq'), - pl.when(pl.col('ceqq') == 0).then(None).otherwise(pl.col('ceqq')).alias('ceqq'), - pl.when(pl.col('atq') == 0).then(None).otherwise(pl.col('atq')).alias('atq') -]) -comp = comp.filter(pl.col('atq').is_not_null()) - -# convert datadate to date fmt -comp = comp.with_columns([ - pl.col('datadate').cast(pl.Date).alias('datadate') -]) - -# merge ccm and comp -# Lag rule: Following Hou, Xue and Zhang (2015), We use earnings immediately after the announcement day -# For those data with missing announcement date record, we straightly let the data available after 4 month -ccm1 = comp.join(ccm, on='gvkey', how='left') -ccm1 = ccm1.with_columns([ - # Year end: December 31st of the same year as datadate - pl.date(pl.col('datadate').dt.year(), 12, 31).alias('yearend'), - pl.col('datadate').dt.offset_by('4mo').dt.month_end().alias('jdate') -]) - -# deal with ibq to make it as up-to-date as possible -ccm1 = ccm1.with_columns([ - pl.col('rdq').cast(pl.Date).dt.month_end().alias('rdq') -]) -ccm1 = ccm1.with_columns([ - pl.when(pl.col('rdq').is_null()).then(pl.col('jdate')).otherwise(pl.col('rdq')).alias('rdq') -]) -# IMPORTANT: enforce chronological order within permno before lead/lag logic -ccm1 = ccm1.sort(['permno', 'datadate', 'rdq', 'jdate']) -# compare next quarter's announcement date with jdate -ccm1 = ccm1.with_columns([ - pl.col('rdq').shift(-1).over('permno').alias('rdq_temp'), - pl.col('ibq').shift(-1).over('permno').alias('ibq_new') -]) -ccm1 = ccm1.with_columns([ - pl.when(pl.col('rdq_temp').is_null()).then(pl.col('jdate')).otherwise(pl.col('rdq_temp')).alias('rdq_temp') -]) -ccm1 = ccm1.with_columns([ - (pl.col('jdate') - pl.col('rdq_temp')).dt.total_days().alias('ibq_diff') -]) -ccm1 = ccm1.rename({'ibq': 'ibq_old'}) # original ibq -''' -if the announcement date is same or in front of jdate, we can use the up-to-date ibq. -otherwise, we consider the up-to-date ibq is not available and still use the lag-4-months ibq -''' -ccm1 = ccm1.with_columns([ - pl.when(pl.col('ibq_diff') >= 0).then(pl.col('ibq_new')).otherwise(pl.col('ibq_old')).alias('ibq') -]) -# for most recent record we can only use the lag-4-months ibq -ccm1 = ccm1.with_columns([ - pl.when(pl.col('ibq').is_null()).then(pl.col('ibq_old')).otherwise(pl.col('ibq')).alias('ibq') -]) - -# set link date bounds -ccm2 = ccm1.filter( - (pl.col('jdate') >= pl.col('linkdt')) & (pl.col('jdate') <= pl.col('linkenddt')) -) - -# merge ccm2 and crsp (using permno only for initial sample restriction) -# Full CRSP data (me, ret, etc.) will be merged later for ME-dependent characteristics -data_rawq = ccm2.join(crsp_permno_only, on=['permno', 'jdate'], how='inner') - -# # filter exchcd & shrcd and at least one year data after the IPO -# data_rawq = data_rawq[((data_rawq['exchcd'] == 1) | (data_rawq['exchcd'] == 2) | (data_rawq['exchcd'] == 3)) & -# ((data_rawq['shrcd'] == 10) | (data_rawq['shrcd'] == 11))].reset_index(drop=True) - -# deal with the duplicates -# Keep first occurrence for each group of ['datadate', 'permno', 'linkprim'] -data_rawq = data_rawq.with_row_index('_temp_idx') -temp_first = (data_rawq - .group_by(['datadate', 'permno', 'linkprim'], maintain_order=True) - .agg(pl.col('_temp_idx').first()) -) -data_rawq = data_rawq.join(temp_first, on=['datadate', 'permno', 'linkprim', '_temp_idx'], how='semi').drop('_temp_idx') - -# Keep last occurrence for each group of ['permno', 'yearend', 'datadate'] -data_rawq = data_rawq.with_row_index('_temp_idx') -temp_last = (data_rawq - .group_by(['permno', 'yearend', 'datadate'], maintain_order=True) - .agg(pl.col('_temp_idx').last()) -) -data_rawq = data_rawq.join(temp_last, on=['permno', 'yearend', 'datadate', '_temp_idx'], how='semi').drop('_temp_idx') - -data_rawq = data_rawq.sort(['permno', 'jdate']) - -# Unified fill_null(0) for quarterly data -_QUARTERLY_FILL_ZERO = [ - 'ivaoq', 'dlcq', 'dlttq', 'mibq', 'pstkq', - 'gdwlq', 'intanq', 'xintq', 'xsgaq', 'cheq', - 'actq', 'lctq', 'dpq', 'txditcq', 'acoq', - 'aoq', 'apq', 'lcoq', 'loq', 'txpq', 'cogsq', - 'rectq', 'invtq', 'ppentq' -] -data_rawq = data_rawq.with_columns([ - pl.col(c).fill_null(0) for c in _QUARTERLY_FILL_ZERO - if c in data_rawq.columns -]) - -# add industry code for quarterly data -data_rawq = data_rawq.filter(pl.col('sic').is_not_null()) # gvkey 039750 does not have sic -data_rawq = data_rawq.with_columns([ - pl.col('sic').cast(pl.Int64).alias('sic') -]) - -data_rawq = data_rawq.with_columns([ - ffi49().alias('ffi49') -]) -data_rawq = data_rawq.with_columns([ - pl.col('ffi49').fill_nan(None).fill_null(49).cast(pl.Int64).alias('ffi49') -]) -####################################################################################################################### -# Quarterly Variables # -####################################################################################################################### -# prepare be -# (fixed): be(ps) beq(pstkq)具体的差异 -# data_rawq = data_rawq.with_columns([ -# pl.when(pl.col('seqq') > 0) -# .then(pl.col('seqq') + pl.col('txditcq') - pl.col('pstkq')) -# .otherwise(None) -# .alias('beq') -# ]) -# data_rawq = data_rawq.with_columns([ -# pl.when(pl.col('beq') <= 0).then(None).otherwise(pl.col('beq')).alias('beq') -# ]) - - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('seqq').is_not_null()).then(pl.col('seqq')) - .when(pl.col('ceqq').is_not_null() & pl.col('pstkq').is_not_null()) - .then(pl.col('ceqq') + pl.col('pstkq')) - .otherwise(pl.col('atq') - pl.col('ltq')) - .alias('seqq') -]) -# (fixed-20260316) compute beq = seqq + txditcq - pstkq from the fallback seqq above, -# consistent with annual be = seq + txditc - ps (HXZ A.1; GHZ Online Appendix p.7). -# Note: quarterly ps hierarchy is limited to pstkq only — pstkrv and pstkl are annual -# Compustat fields not available in comp_fundq, so we cannot replicate the full annual -# coalesce(pstkrv, pstkl, pstk) hierarchy. pstkq and txditcq are already 0-filled via -# _QUARTERLY_FILL_ZERO. Previously this block only filtered the Compustat beq by > 0 -# without applying the seqq fallback, silently leaving beq null for firms where seqq was -# missing but ceqq + pstkq or atq - ltq was available. -data_rawq = data_rawq.with_columns([ - (pl.col('seqq') + pl.col('txditcq') - pl.col('pstkq')).alias('beq') -]) -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('beq') > 0).then(pl.col('beq')).otherwise(None).alias('beq') -]) - -# (fixed): dy_a dy_q的计算方式不同，follow quarterly -# dy -# data_rawq['me_l1'] = data_rawq.groupby(['permno'])['me'].shift(1) -# data_rawq['retdy'] = data_rawq['ret'] - data_rawq['retx'] -# data_rawq['mdivpay'] = data_rawq['retdy']*data_rawq['me_l1'] -# -# data_rawq['dy'] = ttm12(series='mdivpay', df=data_rawq)/data_rawq['me'] - -# chtx -data_rawq = data_rawq.with_columns([ - pl.col('txtq').shift(4).over('permno').alias('txtq_l4'), - pl.col('atq').shift(4).over('permno').alias('atq_l4') -]) -data_rawq = data_rawq.with_columns([ - ((pl.col('txtq') - pl.col('txtq_l4')) / pl.col('atq_l4').replace(0, None)).alias('chtx') -]) - -# roa -data_rawq = data_rawq.with_columns([ - pl.col('atq').shift(1).over('permno').alias('atq_l1') -]) -data_rawq = data_rawq.with_columns([ - (pl.col('ibq') / pl.col('atq_l1').replace(0, None)).alias('roa') -]) - -# cash -data_rawq = data_rawq.with_columns([ - (pl.col('cheq') / pl.col('atq').replace(0, None)).alias('cash') -]) - -# acc -data_rawq = data_rawq.with_columns([ - pl.col('actq').shift(4).over('permno').alias('actq_l4'), - pl.col('lctq').shift(4).over('permno').alias('lctq_l4'), - pl.col('cheq').shift(4).over('permno').alias('cheq_l4'), - pl.col('dlcq').shift(4).over('permno').alias('dlcq_l4'), - pl.col('txpq').shift(4).over('permno').alias('txpq_l4') -]) - -data_rawq = data_rawq.with_columns([ - pl.col('oancfy').shift(1).over('permno').alias('oancfy_l1'), - pl.col('fyearq').shift(1).over('permno').alias('fyearq_l1'), - pl.col('fqtr').shift(1).over('permno').alias('fqtr_l1'), -]) - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('fqtr') == 1) - .then(pl.col('oancfy')) - .when( - (pl.col('fyearq') == pl.col('fyearq_l1')) & - (pl.col('fqtr') == pl.col('fqtr_l1') + 1) - ) - .then(pl.col('oancfy') - pl.col('oancfy_l1')) - .otherwise(None) - .alias('oancfq') -]) - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('oancfq').is_null()) - .then( - ((pl.col('actq') - pl.col('actq_l4')) - (pl.col('cheq') - pl.col('cheq_l4')) - - (pl.col('lctq') - pl.col('lctq_l4')) + (pl.col('dlcq') - pl.col('dlcq_l4')) + - (pl.col('txpq') - pl.col('txpq_l4')) - pl.col('dpq')) / - ((pl.col('atq') + pl.col('atq_l4')) / 2).replace(0, None) - ) - .otherwise( - (pl.col('niq') - pl.col('oancfq')) / ((pl.col('atq') + pl.col('atq_l4')) / 2).replace(0, None) - ) - .alias('acc') -]) - -# absacc -data_rawq = data_rawq.with_columns([ - pl.col('acc').abs().alias('absacc') -]) - -# bm -# data_rawq['bm'] = data_rawq['beq']/data_rawq['me'] - -# cfp -data_rawq = data_rawq.with_columns([ - ttm4('ibq', data_rawq).alias('ibq4'), - ttm4('dpq', data_rawq).alias('dpq4') -]) -# data_rawq['cfp'] = np.where(data_rawq['dpq'].isnull(), -# data_rawq['ibq4']/data_rawq['me'], -# (data_rawq['ibq4']+data_rawq['dpq4'])/data_rawq['me']) - -# ep -# data_rawq['ep'] = data_rawq['ibq4']/data_rawq['me'] - -# agr -data_rawq = data_rawq.with_columns([ - ((pl.col('atq') - pl.col('atq_l4')) / pl.col('atq_l4').replace(0, None)).alias('agr') -]) - -# ni -data_rawq = data_rawq.with_columns([ - pl.col('cshoq').shift(4).over('permno').alias('cshoq_l4'), - pl.col('ajexq').shift(4).over('permno').alias('ajexq_l4') -]) -# log() result: fill_nan(0) handles log(0)→−inf/nan, fill_null(0) handles null input -# order: fill_nan first then fill_null -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('cshoq').is_null()) - .then(None) - .otherwise( - (pl.col('cshoq') * pl.col('ajexq')).log().fill_nan(0).fill_null(0) - - (pl.col('cshoq_l4') * pl.col('ajexq_l4')).log().fill_nan(0).fill_null(0) - ) - .alias('ni') -]) - -# op: HXZ(A.4.12) scaled by book equity (current, not lagged) -# data_rawq = data_rawq.with_columns([ -# pl.col('beq').shift(4).over('permno').alias('beq_l4') -# ]) -data_rawq = data_rawq.with_columns([ - ((ttm4('revtq', data_rawq) - ttm4('cogsq', data_rawq) - ttm4('xsgaq', data_rawq) - ttm4('xintq', data_rawq)) / pl.col('beq').replace(0, None)).alias('op') -]) - -# chcsho -data_rawq = data_rawq.with_columns([ - ((pl.col('cshoq') / pl.col('cshoq_l4').replace(0, None)) - 1).alias('chcsho') -]) - -# cashdebt -data_rawq = data_rawq.with_columns([ - pl.col('ltq').shift(4).over('permno').alias('ltq_l4') -]) -data_rawq = data_rawq.with_columns([ - ((ttm4('ibq', data_rawq) + ttm4('dpq', data_rawq)) / ((pl.col('ltq') + pl.col('ltq_l4')) / 2).replace(0, None)).alias('cashdebt') -]) - -# rd -data_rawq = data_rawq.with_columns([ttm4('xrdq', data_rawq).alias('xrdq4')]) -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('xrdq4').is_null()).then(pl.col('xrdy')).otherwise(pl.col('xrdq4')).alias('xrdq4') -]) - -data_rawq = data_rawq.with_columns([ - (pl.col('xrdq4') / pl.col('atq_l4').replace(0, None)).alias('xrdq4/atq_l4') -]) -data_rawq = data_rawq.with_columns([ - pl.col('xrdq4/atq_l4').shift(4).over('permno').alias('xrdq4/atq_l4_l4') -]) -data_rawq = data_rawq.with_columns([ - pl.when( - ((pl.col('xrdq4') / pl.col('atq').replace(0, None)) - pl.col('xrdq4/atq_l4_l4')) / - pl.col('xrdq4/atq_l4_l4').replace(0, None) > 0.05 - ) - .then(1) - .otherwise(0) - .alias('rd') -]) - -#################### Follow Hafzalla, Lundholm, and Van Winkle (2011) and GHZ on 2025.02.28 #################### - -# pctacc (follow HXZ-A3.22) -data_rawq = data_rawq.with_columns([ - (pl.col('acc') / pl.col('niq').abs().replace(0, None)).alias('pctacc') -]) - - -# gma -data_rawq = data_rawq.with_columns([ - ttm4('revtq', data_rawq).alias('revtq4'), - ttm4('cogsq', data_rawq).alias('cogsq4') -]) -data_rawq = data_rawq.with_columns([ - ((pl.col('revtq4') - pl.col('cogsq4')) / pl.col('atq_l4').replace(0, None)).alias('gma') -]) - -# lev -# data_rawq['lev'] = data_rawq['ltq']/data_rawq['me'] - -# rdm -# data_rawq['rdm'] = data_rawq['xrdq4']/data_rawq['me'] - -# sgr -data_rawq = data_rawq.with_columns([ttm4('saleq', data_rawq).alias('saleq4')]) -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('saleq4').is_null()).then(pl.col('saley')).otherwise(pl.col('saleq4')).alias('saleq4') -]) - -data_rawq = data_rawq.with_columns([ - pl.col('saleq4').shift(4).over('permno').alias('saleq4_l4') -]) -data_rawq = data_rawq.with_columns([ - ((pl.col('saleq4') / pl.col('saleq4_l4').replace(0, None)) - 1).alias('sgr') -]) - -# sp -# data_rawq['sp'] = data_rawq['saleq4']/data_rawq['me'] - -# invest -data_rawq = data_rawq.with_columns([ - pl.col('ppentq').shift(4).over('permno').alias('ppentq_l4'), - pl.col('invtq').shift(4).over('permno').alias('invtq_l4'), - pl.col('ppegtq').shift(4).over('permno').alias('ppegtq_l4') -]) - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('ppegtq').is_null()) - .then( - ((pl.col('ppentq') - pl.col('ppentq_l4')) + - (pl.col('invtq') - pl.col('invtq_l4'))) / pl.col('atq_l4').replace(0, None) - ) - .otherwise( - ((pl.col('ppegtq') - pl.col('ppegtq_l4')) + - (pl.col('invtq') - pl.col('invtq_l4'))) / pl.col('atq_l4').replace(0, None) - ) - .alias('invest') -]) - -# rd_sale -data_rawq = data_rawq.with_columns([ - (pl.col('xrdq4') / pl.col('saleq4').replace(0, None)).alias('rd_sale') -]) - -# lgr -data_rawq = data_rawq.with_columns([ - ((pl.col('ltq') / pl.col('ltq_l4').replace(0, None)) - 1).alias('lgr') -]) - -# depr -data_rawq = data_rawq.with_columns([ - (ttm4('dpq', data_rawq) / pl.col('ppentq').replace(0, None)).alias('depr') -]) - -# egr -data_rawq = data_rawq.with_columns([ - pl.col('ceqq').shift(4).over('permno').alias('ceqq_l4') -]) -data_rawq = data_rawq.with_columns([ - ((pl.col('ceqq') - pl.col('ceqq_l4')) / pl.col('ceqq_l4').replace(0, None)).alias('egr') -]) - -# chpm -data_rawq = data_rawq.with_columns([ - pl.col('ibq4').shift(1).over('permno').alias('ibq4_l1'), - pl.col('saleq4').shift(1).over('permno').alias('saleq4_l1') -]) - -data_rawq = data_rawq.with_columns([ - ((pl.col('ibq4') / pl.col('saleq4').replace(0, None)) - - (pl.col('ibq4_l1') / pl.col('saleq4_l1').replace(0, None))).alias('chpm') -]) - -# noa -# (fixed-20260316) document quarterly noa_raw explicitly as operating -# assets minus operating liabilities before scaling. -# (fixed-20260330) align quarterly noa_raw with HXZ-style operating assets -# definition by excluding ivaoq from operating assets. -data_rawq = data_rawq.with_columns([ - ((pl.col('atq') - pl.col('cheq') - pl.col('ivaoq')) - - (pl.col('atq') - pl.col('dlcq') - pl.col('dlttq') - pl.col('mibq') - - pl.col('pstkq') - pl.col('ceqq'))).alias('noa_raw') -]) - -data_rawq = data_rawq.with_columns([ - pl.col('noa_raw').shift(4).over('permno').alias('noa_raw_l4'), - pl.col('noa_raw').shift(8).over('permno').alias('noa_raw_l8') -]) - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('atq_l4') != 0) - .then(pl.col('noa_raw') / pl.col('atq_l4')) - .otherwise(None) - .alias('noa') -]) - -# chato -# (fixed-20260316) document quarterly noa_raw explicitly as operating -# assets minus operating liabilities before scaling. Reason: -# `OnlineAppendixOPCSAP.pdf` (p. 24) defines NOA as operating assets minus -# operating liabilities, scaled by lagged total assets; `Global Factor Data -# Documentation.pdf` (p. 14) likewise defines NOA as OA - OL. Here, -# operating assets = atq - cheq, and operating liabilities = -# atq - dlcq - dlttq - mibq - pstkq - ceqq, so noa_raw is the unscaled -# dollar NOA used by noa, rna, and ato before scaling by atq_l4. -# (fixed-20260325) use NOA-based turnover change for quarterly chato so the -# denominator is consistent with annual chato and with quarterly rna / ato. -data_rawq = data_rawq.with_columns([ - ((pl.col('saleq4') / ((pl.col('noa_raw') + pl.col('noa_raw_l4')) / 2).replace(0, None)) - - (pl.col('saleq4_l4') / ((pl.col('noa_raw_l4') + pl.col('noa_raw_l8')) / 2).replace(0, None))).alias('chato') -]) - -# chatoia -df_temp = (data_rawq - .group_by(['datadate', 'ffi49']) - .agg(pl.col('chato').mean().alias('chato_ind')) -) -data_rawq = data_rawq.join(df_temp, on=['datadate', 'ffi49'], how='left') -data_rawq = data_rawq.with_columns([ - (pl.col('chato') - pl.col('chato_ind')).alias('chatoia') -]) - -# rna -# (fix)2026-02-27: use noa_raw (unscaled) as denominator instead of noa (scaled by atq_l4) -data_rawq = data_rawq.with_columns([ - pl.col('noa_raw').shift(4).over('permno').alias('noa_raw_l4') -]) -data_rawq = data_rawq.with_columns([ - ((pl.col('oiadpq') / pl.col('noa_raw_l4').replace(0, None))).alias('rna') -]) - -# pm -data_rawq = data_rawq.with_columns([ - (pl.col('oiadpq') / pl.col('saleq').replace(0, None)).alias('pm') -]) - -# ato -# (fix)2026-02-27: use noa_raw (unscaled) as denominator instead of noa (scaled by atq_l4) -data_rawq = data_rawq.with_columns([ - ((pl.col('saleq') / pl.col('noa_raw_l4').replace(0, None))).alias('ato') -]) - -# roe -data_rawq = data_rawq.with_columns([ - pl.col('ceqq').shift(1).over('permno').alias('ceqq_l1') -]) -data_rawq = data_rawq.with_columns([ - (pl.col('ibq') / pl.col('ceqq_l1').replace(0, None)).alias('roe') -]) - -################################## New Added ################################## - -# grltnoa -data_rawq = data_rawq.with_columns([ - pl.col('rectq').shift(4).over('permno').alias('rectq_l4'), - pl.col('acoq').shift(4).over('permno').alias('acoq_l4'), - pl.col('apq').shift(4).over('permno').alias('apq_l4'), - pl.col('lcoq').shift(4).over('permno').alias('lcoq_l4'), - pl.col('loq').shift(4).over('permno').alias('loq_l4'), - pl.col('aoq').shift(4).over('permno').alias('aoq_l4'), - pl.col('invtq').shift(4).over('permno').alias('invtq_l4'), - pl.col('ppentq').shift(4).over('permno').alias('ppentq_l4'), - pl.col('intanq').shift(4).over('permno').alias('intanq_l4') -]) - -# LTNOA_t and LTNOA_{t-1} using quarterly fields -data_rawq = data_rawq.with_columns([ - ( - pl.col('rectq') + pl.col('invtq') + pl.col('ppentq') + - pl.col('acoq') + pl.col('intanq') + pl.col('aoq') - - pl.col('apq') - pl.col('lcoq') - pl.col('loq') - ).alias('ltnoaq_t'), - ( - pl.col('rectq_l4') + pl.col('invtq_l4') + pl.col('ppentq_l4') + - pl.col('acoq_l4') + pl.col('intanq_l4') + pl.col('aoq_l4') - - pl.col('apq_l4') - pl.col('lcoq_l4') - pl.col('loq_l4') - ).alias('ltnoaq_l4') -]) -# Working-capital operating accrual component -data_rawq = data_rawq.with_columns([ - ( - (pl.col('rectq') - pl.col('rectq_l4')) + - (pl.col('invtq') - pl.col('invtq_l4')) + - (pl.col('acoq') - pl.col('acoq_l4')) - - ( - (pl.col('apq') - pl.col('apq_l4')) + - (pl.col('lcoq') - pl.col('lcoq_l4')) - ) - ).alias('wcnoaq_change') -]) -# Quarterly version: same structure, using t and t-4 assets -data_rawq = data_rawq.with_columns([ - ( - (pl.col('ltnoaq_t') / pl.col('atq').replace(0, None)) - - (pl.col('ltnoaq_l4') / pl.col('atq_l4').replace(0, None)) - - (pl.col('wcnoaq_change') / ((pl.col('atq') + pl.col('atq_l4')) / 2).replace(0, None)) - ).alias('grltnoa') -]) - -# scal -# condlist = [data_rawq['seqq'].isnull(), -# data_rawq['seqq'].isnull() & (data_rawq['ceqq'].isnull() | data_rawq['pstk'].isnull())] -# choicelist = [data_rawq['ceqq']+data_rawq['pstk'], -# data_rawq['atq']-data_rawq['ltq']] -# data_rawq['scal'] = np.select(condlist, choicelist, default=data_rawq['seqq']) - -# ala -# data_rawq = data_rawq.with_columns([ -# pl.when(pl.col('gdwlq').is_null()).then(0).otherwise(pl.col('gdwlq')).alias('gdwlq'), -# pl.when(pl.col('intanq').is_null()).then(0).otherwise(pl.col('intanq')).alias('intanq') -# ]) - -# (fix)2026-02-25: error in +0.5*..., should be -0.5*... -data_rawq = data_rawq.with_columns([ - (pl.col('cheq') + 0.75 * (pl.col('actq') - pl.col('cheq')) - - 0.5 * (pl.col('atq') - pl.col('actq') - pl.col('gdwlq') - pl.col('intanq'))).alias('ala') -]) - -# alm -# data_rawq['alm'] = data_rawq['ala']/(data_rawq['atq']+data_rawq['me']-data_rawq['ceqq']) - -# rsup -data_rawq = data_rawq.with_columns([ - pl.col('saleq').shift(4).over('permno').alias('saleq_l4') -]) -# data_rawq['rsup'] = (data_rawq['saleq'] - data_rawq['saleq_l4'])/data_rawq['me'] - -# (fixed-20260316) sacc, stdacc, scf, stdcf dropped — not needed in output. -# sacc = scaled 1-quarter accrual (Bandyopadhyay et al. 2010); stdacc = std of sacc over 16 qtrs; -# scf = cash-flow-to-sales proxy (ibq/saleq - sacc); stdcf = std of scf over 16 qtrs. -# roavol and sgrvol (used for Mohanram m7/m8) are retained. -# data_rawq = data_rawq.with_columns([ -# pl.col('actq').shift(1).over('permno').alias('actq_l1'), -# pl.col('cheq').shift(1).over('permno').alias('cheq_l1'), -# pl.col('lctq').shift(1).over('permno').alias('lctq_l1'), -# pl.col('dlcq').shift(1).over('permno').alias('dlcq_l1') -# ]) -# data_rawq = data_rawq.with_columns([ -# (((pl.col('actq') - pl.col('actq_l1')) - (pl.col('cheq') - pl.col('cheq_l1'))) - -# ((pl.col('lctq') - pl.col('lctq_l1')) - (pl.col('dlcq') - pl.col('dlcq_l1')))).alias('sacc_temp') -# ]) -# data_rawq = data_rawq.with_columns([ -# pl.when(pl.col('saleq') <= 0) -# .then(pl.col('sacc_temp') / 0.01) -# .otherwise(pl.col('sacc_temp') / pl.col('saleq').replace(0, None)) -# .alias('sacc') -# ]).drop('sacc_temp') - - -def chars_std(start, end, df, chars): - """ - Calculate rolling standard deviation across multiple lags using polars - - :param start: Order of starting lag - :param end: Order of ending lag - :param df: Polars DataFrame - :param chars: column name for which to calculate std - :return: polars Series with std of factor - """ - # Create list of lagged columns - lag_exprs = [pl.col(chars).shift(i).over('permno').alias(f'chars_l{i}') for i in range(start, end)] - - # Add all lag columns temporarily - df_temp = df.select(lag_exprs) - - # Calculate std across all lag columns (row-wise) - result = df_temp.select( - pl.concat_list([f'chars_l{i}' for i in range(start, end)]).list.std().alias('std_result') - )['std_result'] - - return result - -# stdacc — removed, see sacc block above -# data_rawq = data_rawq.with_columns([ -# pl.Series('stdacc', chars_std(0, 16, data_rawq, 'sacc')) -# ]) - -# roavol -data_rawq = data_rawq.with_columns([ - pl.Series('roavol', chars_std(0, 16, data_rawq, 'roa')) -]) - -# scf / stdcf — removed, see sacc block above -# data_rawq = data_rawq.with_columns([ -# ((pl.col('ibq') / pl.col('saleq').replace(0, None)) - pl.col('sacc')).alias('scf') -# ]) -# data_rawq = data_rawq.with_columns([ -# pl.when(pl.col('saleq') <= 0) -# .then((pl.col('ibq') / 0.01) - pl.col('sacc')) -# .otherwise(pl.col('scf')) -# .alias('scf') -# ]) -# data_rawq = data_rawq.with_columns([ -# pl.Series('stdcf', chars_std(0, 16, data_rawq, 'scf')) -# ]) - -# cinvest -data_rawq = data_rawq.with_columns([ - pl.col('ppentq').shift(1).over('permno').alias('ppentq_l1'), - pl.col('ppentq').shift(2).over('permno').alias('ppentq_l2'), - pl.col('ppentq').shift(3).over('permno').alias('ppentq_l3'), - pl.col('ppentq').shift(4).over('permno').alias('ppentq_l4'), - pl.col('saleq').shift(1).over('permno').alias('saleq_l1'), - pl.col('saleq').shift(2).over('permno').alias('saleq_l2'), - pl.col('saleq').shift(3).over('permno').alias('saleq_l3') -]) - -# Calculate temp columns for normal case (saleq > 0) -data_rawq = data_rawq.with_columns([ - ((pl.col('ppentq_l1') - pl.col('ppentq_l2')) / pl.col('saleq_l1').replace(0, None)).alias('c_temp1'), - ((pl.col('ppentq_l2') - pl.col('ppentq_l3')) / pl.col('saleq_l2').replace(0, None)).alias('c_temp2'), - ((pl.col('ppentq_l3') - pl.col('ppentq_l4')) / pl.col('saleq_l3').replace(0, None)).alias('c_temp3') -]) - -# Calculate cinvest for normal case -data_rawq = data_rawq.with_columns([ - (((pl.col('ppentq') - pl.col('ppentq_l1')) / pl.col('saleq').replace(0, None)) - - pl.concat_list(['c_temp1', 'c_temp2', 'c_temp3']).list.mean()).alias('cinvest') -]) - -# Recalculate temp columns for saleq <= 0 case -data_rawq = data_rawq.with_columns([ - ((pl.col('ppentq_l1') - pl.col('ppentq_l2')) / 0.01).alias('c_temp1_alt'), - ((pl.col('ppentq_l2') - pl.col('ppentq_l3')) / 0.01).alias('c_temp2_alt'), - ((pl.col('ppentq_l3') - pl.col('ppentq_l4')) / 0.01).alias('c_temp3_alt') -]) - -# Update cinvest for saleq <= 0 case -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('saleq') <= 0) - .then( - ((pl.col('ppentq') - pl.col('ppentq_l1')) / 0.01) - - pl.concat_list(['c_temp1_alt', 'c_temp2_alt', 'c_temp3_alt']).list.mean() - ) - .otherwise(pl.col('cinvest')) - .alias('cinvest') -]) - -data_rawq = data_rawq.drop(['c_temp1', 'c_temp2', 'c_temp3', 'c_temp1_alt', 'c_temp2_alt', 'c_temp3_alt']) - -# nincr -# (fixed-20260316) confirmed nincr should use YoY comparisons against the -# same quarter in the prior year, not QoQ comparisons. Reason: `Green 等 - -# 2017 - The Characteristics that Provide Independent Infor.pdf` (p. 41) -# defines nincr as the number of consecutive quarters (up to eight) with an -# increase in IBQ over the same quarter in the prior year; `Hou 等 - 2020 - -# Replicating Anomalies.pdf` (p. 56) uses the same definition for Nei. -data_rawq = data_rawq.with_columns([ - pl.col('ibq').shift(1).over('permno').alias('ibq_l1'), - pl.col('ibq').shift(2).over('permno').alias('ibq_l2'), - pl.col('ibq').shift(3).over('permno').alias('ibq_l3'), - pl.col('ibq').shift(4).over('permno').alias('ibq_l4'), - pl.col('ibq').shift(5).over('permno').alias('ibq_l5'), - pl.col('ibq').shift(6).over('permno').alias('ibq_l6'), - pl.col('ibq').shift(7).over('permno').alias('ibq_l7'), - pl.col('ibq').shift(8).over('permno').alias('ibq_l8'), - pl.col('ibq').shift(9).over('permno').alias('ibq_l9'), - pl.col('ibq').shift(10).over('permno').alias('ibq_l10'), - pl.col('ibq').shift(11).over('permno').alias('ibq_l11') -]) - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('ibq') > pl.col('ibq_l4')).then(1).otherwise(0).alias('nincr_temp1'), - pl.when(pl.col('ibq_l1') > pl.col('ibq_l5')).then(1).otherwise(0).alias('nincr_temp2'), - pl.when(pl.col('ibq_l2') > pl.col('ibq_l6')).then(1).otherwise(0).alias('nincr_temp3'), - pl.when(pl.col('ibq_l3') > pl.col('ibq_l7')).then(1).otherwise(0).alias('nincr_temp4'), - pl.when(pl.col('ibq_l4') > pl.col('ibq_l8')).then(1).otherwise(0).alias('nincr_temp5'), - pl.when(pl.col('ibq_l5') > pl.col('ibq_l9')).then(1).otherwise(0).alias('nincr_temp6'), - pl.when(pl.col('ibq_l6') > pl.col('ibq_l10')).then(1).otherwise(0).alias('nincr_temp7'), - pl.when(pl.col('ibq_l7') > pl.col('ibq_l11')).then(1).otherwise(0).alias('nincr_temp8') -]) - -data_rawq = data_rawq.with_columns([ - (pl.col('nincr_temp1') + - (pl.col('nincr_temp1') * pl.col('nincr_temp2')) + - (pl.col('nincr_temp1') * pl.col('nincr_temp2') * pl.col('nincr_temp3')) + - (pl.col('nincr_temp1') * pl.col('nincr_temp2') * pl.col('nincr_temp3') * pl.col('nincr_temp4')) + - (pl.col('nincr_temp1') * pl.col('nincr_temp2') * pl.col('nincr_temp3') * pl.col('nincr_temp4') * pl.col('nincr_temp5')) + - (pl.col('nincr_temp1') * pl.col('nincr_temp2') * pl.col('nincr_temp3') * pl.col('nincr_temp4') * pl.col('nincr_temp5') * pl.col('nincr_temp6')) + - (pl.col('nincr_temp1') * pl.col('nincr_temp2') * pl.col('nincr_temp3') * pl.col('nincr_temp4') * pl.col('nincr_temp5') * pl.col('nincr_temp6') * pl.col('nincr_temp7')) + - (pl.col('nincr_temp1') * pl.col('nincr_temp2') * pl.col('nincr_temp3') * pl.col('nincr_temp4') * pl.col('nincr_temp5') * pl.col('nincr_temp6') * pl.col('nincr_temp7') * pl.col('nincr_temp8')) - ).alias('nincr') -]) - -data_rawq = data_rawq.drop(['ibq_l1', 'ibq_l2', 'ibq_l3', 'ibq_l4', 'ibq_l5', 'ibq_l6', 'ibq_l7', 'ibq_l8', - 'nincr_temp1', 'nincr_temp2', 'nincr_temp3', 'nincr_temp4', 'nincr_temp5', 'nincr_temp6', 'nincr_temp7', 'nincr_temp8']) - - -# performance score -data_rawq = data_rawq.with_columns([ttm4('niq', data_rawq).alias('niq4')]) -data_rawq = data_rawq.with_columns([ - pl.col('niq4').shift(4).over('permno').alias('niq4_l4'), - pl.col('dlttq').shift(4).over('permno').alias('dlttq_l4'), - pl.col('cogsq4').shift(4).over('permno').alias('cogsq4_l4'), -]) -data_rawq = data_rawq.with_columns([ttm4('oancfq', data_rawq).alias('oancfq4')]) - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('niq4') > 0).then(1).otherwise(0).alias('p_temp1'), - pl.when(pl.col('oancfq4') > 0).then(1).otherwise(0).alias('p_temp2'), - pl.when( - (pl.col('niq4') / pl.col('atq').replace(0, None)) > - (pl.col('niq4_l4') / pl.col('atq_l4').replace(0, None)) - ).then(1).otherwise(0).alias('p_temp3'), - pl.when(pl.col('oancfq4') > pl.col('niq4')).then(1).otherwise(0).alias('p_temp4'), - pl.when( - (pl.col('dlttq') / pl.col('atq').replace(0, None)) < - (pl.col('dlttq_l4') / pl.col('atq_l4').replace(0, None)) - ).then(1).otherwise(0).alias('p_temp5'), - pl.when( - (pl.col('actq') / pl.col('lctq').replace(0, None)) > - (pl.col('actq_l4') / pl.col('lctq_l4').replace(0, None)) - ).then(1).otherwise(0).alias('p_temp6'), - pl.when( - ((pl.col('saleq4') - pl.col('cogsq4')) / pl.col('saleq4').replace(0, None)) > - ((pl.col('saleq4_l4') - pl.col('cogsq4_l4')) / pl.col('saleq4_l4').replace(0, None)) - ).then(1).otherwise(0).alias('p_temp7'), - pl.when( - (pl.col('saleq4') / pl.col('atq').replace(0, None)) > - (pl.col('saleq4_l4') / pl.col('atq_l4').replace(0, None)) - ).then(1).otherwise(0).alias('p_temp8'), - pl.when(pl.col('scstkcy') == 0).then(1).otherwise(0).alias('p_temp9') -]) - -data_rawq = data_rawq.with_columns([ - (pl.col('p_temp1') + pl.col('p_temp2') + pl.col('p_temp3') + pl.col('p_temp4') + - pl.col('p_temp5') + pl.col('p_temp6') + pl.col('p_temp7') + pl.col('p_temp8') + - pl.col('p_temp9')).alias('pscore') -]) - -data_rawq = data_rawq.drop(['p_temp1', 'p_temp2', 'p_temp3', 'p_temp4', 'p_temp5', 'p_temp6', 'p_temp7', 'p_temp8', 'p_temp9']) - -################## Added on 2022.09.06 ################## -# cashpr -# data_rawq['cashpr'] = ((data_rawq['me'] + data_rawq['dlttq'] - data_rawq['atq']) / data_rawq['cheq']) - -print("Finish Quarterly Variables Calculation! \n") - -####################################################################################################################### -# Momentum # -####################################################################################################################### -# Use crsp_full that was prepared at the beginning (no need to reload CRSP) -crsp_mom = crsp_full.clone() - -# Rename monthend to jdate for consistency -crsp_mom = crsp_mom.rename({'monthend': 'jdate'}) - -# Convert ME to millions -crsp_mom = crsp_mom.with_columns([ - (pl.col('me') / 1000).alias('me') # CRSP ME in million unit -]) - -crsp_mom = crsp_mom.sort(['permno', 'date']) - -def mom(start, end, df): - """ - - :param start: Order of starting lag - :param end: Order of ending lag - :param df: Dataframe - :return: Momentum factor - """ - # Calculate cumulative product: (1 + ret_lag_start) * (1 + ret_lag_start+1) * ... - 1 - # Build the expression without adding columns to the dataframe - result_expr = pl.lit(1) - for i in range(start, end): - result_expr = result_expr * (1 + pl.col('ret').shift(i).over('permno')) - - return result_expr - 1 - - -def chmom(start, end, df): - """ - - :param start: Order of starting lag - :param end: Order of ending lag - :param df: Dataframe - :return: Momentum factor - """ - # Calculate cumulative product for first half (without adding columns) - result_first_half = pl.lit(1) - for i in range(start, end): - result_first_half = result_first_half * (1 + pl.col('ret').shift(i).over('permno')) - result_first_half = result_first_half - 1 - - # Calculate cumulative product for second half (6 months later) - result_second_half = pl.lit(1) - for i in range(start + 6, end + 6): - result_second_half = result_second_half * (1 + pl.col('ret').shift(i).over('permno')) - result_second_half = result_second_half - 1 - - return result_first_half - result_second_half - -# (checked)2026-03-13: mom(1,12) = lags 1..11, skipping lag 0 (short-term reversal) -crsp_mom = crsp_mom.with_columns([ - chmom(1, 7, crsp_mom).alias('chmom'), - pl.col('ret').alias('mom1m'), - mom(1, 7, crsp_mom).alias('mom6m'), - mom(1, 13, crsp_mom).alias('mom12m'), - mom(13, 37, crsp_mom).alias('mom36m'), - mom(13, 61, crsp_mom).alias('mom60m'), - pl.col('ret').shift(11).over('permno').alias('seas1a'), - pl.col('vol').shift(1).over('permno').alias('vol_l1'), - pl.col('vol').shift(2).over('permno').alias('vol_l2'), - pl.col('vol').shift(3).over('permno').alias('vol_l3'), - pl.col('prc').shift(1).over('permno').alias('prc_l1'), - pl.col('prc').shift(2).over('permno').alias('prc_l2'), - pl.col('cfacshr').shift(1).over('permno').alias('cfacshr_l1') -]) -# crsp_mom['dolvol'] = np.log((crsp_mom['vol_l2']*100)*crsp_mom['prc_l2']).replace([np.inf, -np.inf], np.nan) ##### Added "*100" on 2025.02.23 (change "vol" unit from hundreds to one unit) ##### -# crsp_mom['turn'] = ((crsp_mom['vol_l1']+crsp_mom['vol_l2']+crsp_mom['vol_l3'])/3/10)/crsp_mom['shrout'] ##### Added "/10" on 2025.02.23 (change "vol" unit from hundreds to thousand unit, same as shrout) ##### - -# 2025-07-06 updates: In SIZ version, vol(daily-1 units, monthly-100 units), shrout-1000 units. -# In CIZ version, vol(monthly-1 units) -crsp_mom = crsp_mom.with_columns([ - (pl.col('vol_l2') * pl.col('prc_l2')).log() - .replace([float('inf'), float('-inf')], None).alias('dolvol'), - # (fixed-20260325) cast turn to Float64 at creation time so monthly share - # turnover keeps fractional precision and is not written as an integer-scale Decimal. - (((pl.col('vol_l1') + pl.col('vol_l2') + pl.col('vol_l3')).cast(pl.Float64) / 3.0 / 1000.0) / - pl.col('shrout').cast(pl.Float64).replace(0, None)).cast(pl.Float64).alias('turn'), - pl.col('me').shift(1).over('permno').alias('me_l1'), - (pl.col('ret') - pl.col('retx')).alias('retdy') -]) - - -# (fixed-20260316) JKP: div1m_me = (ret - retx) * prc_l1 * (cfacshr_t / cfacshr_l1) * shares_t -# Uses current shrout adjusted by cfacshr ratio so that share issuances and stock -# splits/dividends between t-1 and t are correctly accounted for (vs naive shrout_l1). -# Units: prc($/sh) * shrout(thousands) / 1000 = $M, matching me units. -crsp_mom = crsp_mom.with_columns([ - (pl.col('retdy') * pl.col('prc_l1') * - (pl.col('cfacshr') / pl.col('cfacshr_l1').replace(0, None)) * - pl.col('shrout') / 1000).alias('mdivpay') -]) - -crsp_mom = crsp_mom.with_columns([ - (ttm12('mdivpay', crsp_mom) / pl.col('me').replace(0, None)).alias('dy') -]) - -# 2026-02-11 updates: Add size group classification -# NYSE monthly size cutoffs and size group classification -nyse_cutoffs = (crsp_mom - .filter( - (pl.col('primaryexch') == 'N') & - pl.col('me').is_not_null() - ) - .group_by('jdate') - .agg([ - pl.len().alias('n'), - pl.col('me').quantile(0.01, interpolation='higher').alias('nyse_p1'), - pl.col('me').quantile(0.20, interpolation='higher').alias('nyse_p20'), - pl.col('me').quantile(0.50, interpolation='higher').alias('nyse_p50'), - pl.col('me').quantile(0.80, interpolation='higher').alias('nyse_p80') - ]) -) - -crsp_mom = (crsp_mom - .join(nyse_cutoffs, on='jdate', how='left') - .with_columns([ - pl.when(pl.col('me').is_null()) - .then(None) - .when(pl.col('nyse_p80').is_null()) - .then(pl.lit('mega')) - .when(pl.col('me') >= pl.col('nyse_p80')) - .then(pl.lit('mega')) - .when(pl.col('me') >= pl.col('nyse_p50')) - .then(pl.lit('large')) - .when(pl.col('me') >= pl.col('nyse_p20')) - .then(pl.lit('small')) - .when(pl.col('me') >= pl.col('nyse_p1')) - .then(pl.lit('micro')) - .otherwise(pl.lit('nano')) - .alias('size_grp') - ]) - .drop(['n', 'nyse_p1', 'nyse_p20', 'nyse_p50', 'nyse_p80']) -) - -# populate the chars to monthly - -# data_rawa -data_rawa = data_rawa.drop(['date', 'ret', 'retx', 'me', 'vol', 'permco', 'prc', 'shrout'], strict=False) -data_rawa = crsp_mom.join(data_rawa, on=['permno', 'jdate'], how='left') -data_rawa = data_rawa.sort(['permno', 'jdate']) - -# (fixed-20260316) forward_fill datadate over permno is correct. -# After left-join with monthly crsp_mom, most monthly rows have null datadate (no matching -# fiscal-year report). Forward_fill carries the most recent datadate forward within each -# permno until the next report arrives — standard stale-accounting practice. -data_rawa = data_rawa.with_columns([ - pl.col('datadate').forward_fill().over('permno') -]) -# (fixed): check-处理pandas才加入的datadate1和permno1，polars不需要，可以直接用datadate和permno -# data_rawa = data_rawa.with_columns([ -# pl.col('permno').alias('permno1'), -# pl.col('datadate').alias('datadate1') -# ]) -data_rawa = data_rawa.with_columns([ - pl.all().forward_fill().over(['permno', 'datadate']) -]) -# (fixed): check是否重复筛选 -# data_rawa = data_rawa.filter( -# (pl.col('primaryexch').is_in(['N', 'A', 'Q'])) & -# (pl.col('conditionaltype') == 'RW') & -# (pl.col('tradingstatusflg') == 'A') -# ) - -# data_rawq -data_rawq = data_rawq.drop(['date', 'ret', 'retx', 'me', 'vol', 'permco', 'prc', 'shrout'], strict=False) -data_rawq = crsp_mom.join(data_rawq, on=['permno', 'jdate'], how='left') -data_rawq = data_rawq.sort(['permno', 'jdate']) -data_rawq = data_rawq.with_columns([ - pl.col('datadate').forward_fill().over('permno') -]) - - -# (fixed): check-处理pandas才加入的datadate1和permno1，polars不需要，可以直接用datadate和permno -# data_rawq = data_rawq.with_columns([ -# pl.col('permno').alias('permno1'), -# pl.col('datadate').alias('datadate1') -# ]) -data_rawq = data_rawq.with_columns([ - pl.all().forward_fill().over(['permno', 'datadate']) -]) - -####################################################################################################################### -# Monthly ME # -####################################################################################################################### - -######################################## -# Annual # -######################################## - -# bm -data_rawa = data_rawa.with_columns([ - (pl.col('be') / pl.col('me').replace(0, None)).alias('bm') -]) - -# bm_ia -# (fixed): 用date还是datadate -df_temp = data_rawa.group_by(['jdate', 'ffi49']).agg(pl.col('bm').mean().alias('bm_ind')) -data_rawa = data_rawa.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawa = data_rawa.with_columns([ - (pl.col('bm') - pl.col('bm_ind')).alias('bm_ia') -]) - -# me_ia -df_temp = data_rawa.group_by(['jdate', 'ffi49']).agg(pl.col('me').mean().alias('me_ind')) -data_rawa = data_rawa.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawa = data_rawa.with_columns([ - (pl.col('me') - pl.col('me_ind')).alias('me_ia') -]) - -# cfp -data_rawa = data_rawa.with_columns([ - pl.when(pl.col('dp').is_null()) - .then(pl.col('ib') / pl.col('me').replace(0, None)) - .when(pl.col('ib').is_null()) - .then(None) - .otherwise((pl.col('ib') + pl.col('dp')) / pl.col('me').replace(0, None)) - .alias('cfp') -]) - -# cfp_ia -# (fix)2026-03-13: use jdate (not datadate) for industry adjustment to avoid look-ahead bias -df_temp = data_rawa.group_by(['jdate', 'ffi49']).agg(pl.col('cfp').mean().alias('cfp_ind')) -data_rawa = data_rawa.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawa = data_rawa.with_columns([ - (pl.col('cfp') - pl.col('cfp_ind')).alias('cfp_ia') -]) - -# ep -data_rawa = data_rawa.with_columns([ - (pl.col('ib') / pl.col('me').replace(0, None)).alias('ep') -]) - -# rsup -data_rawa = data_rawa.with_columns([ - ((pl.col('sale') - pl.col('sale_l1')) / pl.col('me').replace(0, None)).alias('rsup') -]) - -# lev -data_rawa = data_rawa.with_columns([ - (pl.col('lt') / pl.col('me').replace(0, None)).alias('lev') -]) - -# sp -data_rawa = data_rawa.with_columns([ - (pl.col('sale') / pl.col('me').replace(0, None)).alias('sp') -]) - -# rdm -data_rawa = data_rawa.with_columns([ - (pl.col('xrd') / pl.col('me').replace(0, None)).alias('rdm') -]) - -# adm hxz adm -data_rawa = data_rawa.with_columns([ - (pl.col('xad') / pl.col('me').replace(0, None)).alias('adm') -]) - -# dy — REMOVED: monthly dy from crsp_mom (TTM ret-retx method) is preferred -# data_rawa = data_rawa.with_columns([ -# (pl.col('dvt') / pl.col('me').replace(0, None)).alias('dy') -# ]) - -# cashpr -data_rawa = data_rawa.with_columns([ - ((pl.col('me') + pl.col('dltt') - pl.col('at')) / pl.col('che').replace(0, None)).alias('cashpr') -]) - -# indmom -# (fixed-20260316) keep indmom grouped on jdate and compute it from prior -# 6-month value-weighted industry returns. Reason: `OnlineAppendixOPCSAP.pdf` -# (p. 22) defines Industry Momentum as the weighted average of firm-level -# 6-month buy-and-hold return within each industry using market equity -# weights; `Hou 等 - 2020 - Replicating Anomalies.pdf` (p. 55) further uses -# Fama-French 49 industries and prior 6-month value-weighted returns from -# t-6 to t-1. Grouping on jdate avoids look-ahead and matches the predictor -# timestamp used elsewhere in this file. -data_rawa = data_rawa.with_columns([ - (pl.col('mom6m') * pl.col('me')).alias('_vw_mom6m') -]) -df_temp = (data_rawa - .filter(pl.col('me').is_not_null() & pl.col('mom6m').is_not_null()) - .group_by(['jdate', 'ffi49']) - .agg([ - pl.col('_vw_mom6m').sum().alias('_sum_vw'), - pl.col('me').sum().alias('_sum_me') - ]) - .with_columns([ - (pl.col('_sum_vw') / pl.col('_sum_me')).alias('indmom') - ]) - .select(['jdate', 'ffi49', 'indmom']) -) -data_rawa = data_rawa.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawa = data_rawa.drop('_vw_mom6m') - -# Annual Accounting Variables -# replace 'exchcd','shrcd' with 'primaryexch', 'conditionaltype', 'tradingstatusflg', 'sharetype', 'securitytype', 'securitysubtype', 'usincflg', 'issuertype' -chars_a = data_rawa.select(['cusip_comp', 'cusip_crsp', 'hdrcusip', 'gvkey', 'permno', 'primaryexch', 'conditionaltype', 'tradingstatusflg', - 'sharetype', 'securitytype', 'securitysubtype', 'usincflg', 'issuertype', - 'datadate', 'jdate', 'ticker', 'conm', 'comnam', 'prc', 'shrout', - 'sic', 'ret', 'retx', 'acc', 'agr', 'bm', 'cfp', 'ep', 'ni', 'op', 'rsup', 'cash', 'chcsho', - 'rd', 'cashdebt', 'pctacc', 'gma', 'lev', 'rdm', 'adm', 'sgr', 'sp', 'invest', 'roe', - 'rd_sale', 'lgr', 'roa', 'depr', 'egr', 'chato', 'chtx', 'noa', 'rna', 'pm', 'ato', - 'roic', 'chinv', 'pchsale_pchinvt', 'pchsale_pchrect', 'pchgm_pchsale', 'pchsale_pchxsga', - 'pchdepr', 'chadv', 'pchcapx', 'grcapx', 'grGW', 'currat', 'pchcurrat', 'quick', 'pchquick', - 'salecash', 'salerec', 'saleinv', 'pchsaleinv', 'realestate', 'obklg', 'chobklg', 'grltnoa', - 'conv', 'chdrc', 'rdbias', 'operprof', 'capxint', 'xadint', 'chpm', 'ala', 'alm', - 'mom1m', 'mom6m', 'mom12m', 'mom60m', 'mom36m', 'seas1a', 'me', 'size_grp', 'hire', 'herf', 'bm_ia', - 'me_ia', 'turn', 'dolvol', 'dy', 'absacc', 'age', 'cashpr', 'chatoia', 'chempia', 'chmom', 'chpmia', - 'convind', 'divi', 'divo', 'secured', 'securedind', 'sin', 'cfp_ia', 'indmom', 'pchcapx_ia', - 'tang', 'tb', 'm1', 'm2', 'm3', 'm4', 'm5', 'm6']) - -######################################## -# Quarterly # -######################################## -# bm -data_rawq = data_rawq.with_columns([ - (pl.col('beq') / pl.col('me').replace(0, None)).alias('bm') -]) - -# bm_ia -df_temp = data_rawq.group_by(['jdate', 'ffi49']).agg(pl.col('bm').mean().alias('bm_ind')) -data_rawq = data_rawq.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawq = data_rawq.with_columns([ - (pl.col('bm') - pl.col('bm_ind')).alias('bm_ia') -]) - -# me_ia -df_temp = data_rawq.group_by(['jdate', 'ffi49']).agg(pl.col('me').mean().alias('me_ind')) -data_rawq = data_rawq.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawq = data_rawq.with_columns([ - (pl.col('me') - pl.col('me_ind')).alias('me_ia') -]) - -# cfp -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('dpq').is_null()) - .then(pl.col('ibq4') / pl.col('me').replace(0, None)) - .otherwise((pl.col('ibq4') + pl.col('dpq4')) / pl.col('me').replace(0, None)) - .alias('cfp') -]) - -# cfp_ia -# (fix)2026-03-13: use jdate (not datadate) for industry adjustment to avoid look-ahead bias -df_temp = data_rawq.group_by(['jdate', 'ffi49']).agg(pl.col('cfp').mean().alias('cfp_ind')) -data_rawq = data_rawq.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawq = data_rawq.with_columns([ - (pl.col('cfp') - pl.col('cfp_ind')).alias('cfp_ia') -]) - -# ep -data_rawq = data_rawq.with_columns([ - (pl.col('ibq4') / pl.col('me').replace(0, None)).alias('ep') -]) - -# lev -data_rawq = data_rawq.with_columns([ - (pl.col('ltq') / pl.col('me').replace(0, None)).alias('lev') -]) - -# rdm -data_rawq = data_rawq.with_columns([ - (pl.col('xrdq4') / pl.col('me').replace(0, None)).alias('rdm') -]) - -# sp -data_rawq = data_rawq.with_columns([ - (pl.col('saleq4') / pl.col('me').replace(0, None)).alias('sp') -]) - -# alm -data_rawq = data_rawq.with_columns([ - (pl.col('ala') / (pl.col('atq') + pl.col('me') - pl.col('ceqq')).replace(0, None)).alias('alm') -]) - -# rsup -data_rawq = data_rawq.with_columns([ - ((pl.col('saleq') - pl.col('saleq_l4')) / pl.col('me').replace(0, None)).alias('rsup') -]) - -# (checked): check 0-15/0-16 -# sgrvol: 为什么用rsup计算sgrvol(reference) -data_rawq = data_rawq.with_columns([ - chars_std(0, 16, data_rawq, 'rsup').alias('sgrvol') -]) - -# cashpr -data_rawq = data_rawq.with_columns([ - ((pl.col('me') + pl.col('dlttq') - pl.col('atq')) / pl.col('cheq').replace(0, None)).alias('cashpr') -]) - -# indmom -# (fix)2026-03-13: per Grinblatt & Moskowitz (1999) / GHZ: value-weighted 6-month return -data_rawq = data_rawq.with_columns([ - (pl.col('mom6m') * pl.col('me')).alias('_vw_mom6m') -]) -df_temp = (data_rawq - .filter(pl.col('me').is_not_null() & pl.col('mom6m').is_not_null()) - .group_by(['jdate', 'ffi49']) - .agg([ - pl.col('_vw_mom6m').sum().alias('_sum_vw'), - pl.col('me').sum().alias('_sum_me') - ]) - .with_columns([ - (pl.col('_sum_vw') / pl.col('_sum_me')).alias('indmom') - ]) - .select(['jdate', 'ffi49', 'indmom']) -) -data_rawq = data_rawq.join(df_temp, on=['jdate', 'ffi49'], how='left') -data_rawq = data_rawq.drop('_vw_mom6m') - -# Mohanram (2005) score (Quarterly Related) -df_temp = data_rawq.group_by(['fyearq', 'fqtr', 'ffi49']).agg(pl.col('roavol').median().alias('md_roavol')) -data_rawq = data_rawq.join(df_temp, on=['fyearq', 'fqtr', 'ffi49'], how='left') - -df_temp = data_rawq.group_by(['fyearq', 'fqtr', 'ffi49']).agg(pl.col('sgrvol').median().alias('md_sgrvol')) -data_rawq = data_rawq.join(df_temp, on=['fyearq', 'fqtr', 'ffi49'], how='left') - -data_rawq = data_rawq.with_columns([ - pl.when(pl.col('roavol') < pl.col('md_roavol')).then(1).otherwise(0).alias('m7'), - pl.when(pl.col('sgrvol') < pl.col('md_sgrvol')).then(1).otherwise(0).alias('m8') -]) - -# Quarterly Accounting Variables -# replace 'exchcd','shrcd' with 'primaryexch', 'conditionaltype', 'tradingstatusflg', 'sharetype', 'securitytype', 'securitysubtype', 'usincflg', 'issuertype' -chars_q = data_rawq.select(['gvkey', 'permno', 'datadate', 'jdate', 'sic', 'primaryexch', 'conditionaltype', 'tradingstatusflg', - 'sharetype', 'securitytype', 'securitysubtype', 'usincflg', 'issuertype', 'ticker', 'conm', 'comnam', 'prc', 'shrout', - 'ret', 'retx', - 'acc', 'bm', 'cfp', 'ep', 'agr', 'ni', 'op', 'cash', 'chcsho', 'rd', - 'cashdebt', 'pctacc', 'gma', 'lev', 'rdm', 'sgr', 'sp', 'invest', 'rd_sale', 'lgr', - 'roa', 'depr', 'egr', 'roe', 'chato', 'chpm', 'chtx', 'noa', 'rna', 'pm', 'ato', - 'grltnoa', 'ala', 'alm', 'rsup', 'sgrvol', 'roavol', 'cinvest', - 'mom1m', 'mom6m', 'mom12m', 'mom60m', 'mom36m', 'seas1a', 'me', 'size_grp', 'pscore', 'nincr', - 'cfp_ia', 'bm_ia', 'me_ia', 'chatoia', 'chmom', - 'turn', 'dolvol', 'cashpr', 'dy', 'indmom', 'm7', 'm8']) - -chars_a.write_parquet(OUTPUT_PATH + 'chars_a_accounting.parquet') - -chars_q.write_parquet(OUTPUT_PATH + 'chars_q_accounting.parquet') diff --git a/chars_ciz_monthly/download_data.py b/chars_ciz_monthly/download_data.py deleted file mode 100644 index 73a42d6..0000000 --- a/chars_ciz_monthly/download_data.py +++ /dev/null @@ -1,284 +0,0 @@ -import os -import duckdb -import polars as pl -import time - - -# Configuration -OUTPUT_PATH = "../data/raw/" - -os.makedirs(OUTPUT_PATH, exist_ok=True) - - -def measure_time(func): - """ - Decorator to time a function and print start/end timestamps and elapsed minutes:seconds. - """ - def wrapper(*args, **kwargs): - start_time = time.time() - print(f"Function : {func.__name__.upper()}", flush=True) - print( - f"Start : {time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(start_time))}", - flush=True, - ) - result = func(*args, **kwargs) - end_time = time.time() - print( - f"End : {time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(end_time))}", - flush=True, - ) - total_seconds = end_time - start_time - minutes = int(total_seconds // 60) - seconds = total_seconds % 60 - print( - f"Execution time : {minutes} minutes and {seconds:.2f} seconds", flush=True - ) - print() - return result - return wrapper - - -def gen_wrds_connection_info(user, password): - """Generate WRDS PostgreSQL connection string for DuckDB.""" - return ( - f"host=wrds-pgdata.wharton.upenn.edu " - f"port=9737 dbname=wrds " - f"user={user} password={password} sslmode=require" - ) - - -def _execute_download(con, table_name, query, output_file): - """ - Helper function to execute a single table download. - - Args: - con: DuckDB connection (already attached to WRDS) - table_name: Name of table for logging - query: SQL query to execute on WRDS - output_file: Path to save parquet file - """ - print(f"Downloading {table_name}...", flush=True) - - con.execute(f""" - COPY ( - SELECT * FROM postgres_query('wrds', '{query}') - ) - TO '{output_file}' (FORMAT PARQUET); - """) - - file_size_mb = os.path.getsize(output_file) / (1024 * 1024) - print(f"{table_name} saved to {output_file} ({file_size_mb:.1f} MB)", flush=True) - - -# ==================================================================================================== -# TABLE DEFINITIONS - Edit these queries to modify what data to download -# ==================================================================================================== - -def get_tables_config(start_date='2020-01-01'): - """ - Define all tables to download with their queries. - Edit the queries in this function to modify what data to download. - - Args: - start_date: Start date for filtering (format: 'YYYY-MM-DD') - - Returns: - Dictionary of table configurations with 'output' and 'query' keys - """ - return { - 'comp_funda': { - 'output': os.path.join(OUTPUT_PATH, 'comp_funda.parquet'), - 'query': f""" - SELECT - f.gvkey, f.cusip, f.datadate, f.fyear, c.cik, substr(c.sic,1,2) as sic2, c.sic, c.naics, - - /* income statement */ - f.sale, f.revt, f.cogs, f.xsga, f.dp, f.xrd, f.xad, f.ib, f.ebitda, - f.ebit, f.nopi, f.spi, f.pi, f.txp, f.ni, f.txfed, f.txfo, f.txt, f.xint, - - /* CF statement and others */ - f.capx, f.oancf, f.dvt, f.ob, f.gdwlia, f.gdwlip, f.gwo, f.mib, f.oiadp, f.ivao, - - /* assets */ - f.rect, f.act, f.che, f.ppegt, f.invt, f.at, f.aco, f.intan, f.ao, f.ppent, f.gdwl, f.fatb, f.fatl, - - /* liabilities */ - f.lct, f.dlc, f.dltt, f.lt, f.dm, f.dcvt, f.cshrc, - f.dcpstk, f.pstk, f.ap, f.lco, f.lo, f.drc, f.drlt, f.txdi, - - /* equity and other */ - f.ceq, f.scstkc, f.emp, f.csho, f.seq, f.txditc, f.pstkrv, f.pstkl, f.np, f.txdc, f.dpc, f.ajex, f.conm, - - /* market */ - ABS(f.prcc_f) AS prcc_f - FROM comp.funda AS f - LEFT JOIN comp.company AS c - ON f.gvkey = c.gvkey - WHERE f.indfmt = ''INDL'' - AND f.datafmt = ''STD'' - AND f.popsrc = ''D'' - AND f.consol = ''C'' - AND f.datadate >= ''{start_date}'' - """ - }, -# (fixed): add cfacshr(crsp_msf: mthcumfacshr), wcaptq (is wcapq) - 'comp_fundq': { - 'output': os.path.join(OUTPUT_PATH, 'comp_fundq.parquet'), - 'query': f""" - SELECT - /*header info*/ - c.gvkey, f.cusip, f.datadate, f.fyearq, substr(c.sic,1,2) as sic2, c.sic, f.fqtr, f.rdq, - - /*income statement*/ - f.ibq, f.saleq, f.txtq, f.revtq, f.cogsq, f.xsgaq, f.revty, f.cogsy, f.saley, - - /*balance sheet items*/ - f.atq, f.actq, f.cheq, f.lctq, f.dlcq, f.ppentq, f.ppegtq, f.txpq, - - /*others*/ - abs(f.prccq) as prccq, abs(f.prccq)*f.cshoq as mveq_f, f.ceqq, f.seqq, f.pstkq, f.ltq, - f.pstkrq, f.gdwlq, f.intanq, f.mibq, f.oiadpq, f.ivaoq, f.conm, - - /* v3 my formula add*/ - f.ajexq, f.cshoq, f.txditcq, f.npq, f.xrdy, f.xrdq, f.dpq, f.xintq, f.invtq, f.scstkcy, f.niq, - f.oancfy, f.wcapq, f.dlttq, f.rectq, f.acoq, f.apq, f.lcoq, f.loq, f.aoq, - - /* SUE calculation */ - f.epspxq - - FROM comp.fundq as f - LEFT JOIN comp.company as c - ON f.gvkey = c.gvkey - - /*get consolidated, standardized, industrial format statements*/ - WHERE f.indfmt = ''INDL'' - AND f.datafmt = ''STD'' - AND f.popsrc = ''D'' - AND f.consol = ''C'' - AND f.datadate >= ''{start_date}'' - """ - }, -# cfacshr - is mthcumfacshr - 'crsp_msf': { - 'output': os.path.join(OUTPUT_PATH, 'crsp_msf.parquet'), - 'query': f""" - SELECT - mthprc, mthret, mthretx, mthvol, - shrout, mthcumfacpr, mthcumfacshr, - permno, permco, mthcaldt, ticker, cusip, hdrcusip, - issuernm, issuertype, securitytype, securitysubtype, sharetype, usincflg, - primaryexch, conditionaltype, TradingStatusFlg - FROM crspm.msf_v2 - WHERE mthcaldt >= ''{start_date}'' - """ - }, - - # 2026-02-10 updates: confirm and use ccmxpf_lnkhist - 'ccm': { - 'output': os.path.join(OUTPUT_PATH, 'ccm.parquet'), - 'query': """ - SELECT - gvkey, lpermno as permno, linktype, linkprim, - linkdt, linkenddt - FROM crsp.ccmxpf_lnkhist - WHERE linktype IN (''LC'', ''LU'', ''LS'') - """ - }, - - 'crsp_dsf': { - 'output': os.path.join(OUTPUT_PATH, 'crsp_dsf.parquet'), - 'query': f""" - SELECT - a.permno, a.permco, a.dlycaldt, a.dlyret, a.dlyvol, a.dlyprc, a.dlyhigh, a.dlylow, - a.shrout, a.dlydelflg, a.dlycumfacpr, a.dlycumfacshr, - a.primaryexch, a.conditionaltype, a.tradingstatusflg, - a.cusip, a.hdrcusip, a.siccd, - b.rf, b.mktrf, b.smb, b.hml, b.umd, b.rmw, b.cma - FROM crspm.dsf_v2 as a - LEFT JOIN ff_all.fivefactors_daily as b - ON a.dlycaldt = b.date - WHERE a.dlycaldt >= ''{start_date}'' - """ - }, - - 'crsp_ind': { - 'output': os.path.join(OUTPUT_PATH, 'crsp_ind.parquet'), - 'query': f""" - SELECT - dlycaldt, dlyprcret - FROM crspm.inddlyseriesdata - WHERE indno = 1000502 /* S&P 500 Composite */ - AND dlycaldt >= ''{start_date}'' - """ - }, - - 'ibes': { - 'output': os.path.join(OUTPUT_PATH, 'ibes.parquet'), - 'query': """ - SELECT - ticker, statpers, meanest, fpedats, anndats_act, curr_act, fpi, medest - FROM ibes.statsum_epsus - WHERE - /* filtering IBES */ - statpers < ANNDATS_ACT /*only keep summarized forecasts prior to earnings annoucement*/ - AND measure=''EPS'' - AND (fpedats-statpers)>=0 - AND CURCODE=''USD'' - AND fpi in (''1'',''2'') - """ - } - } - -# ==================================================================================================== - - -@measure_time -def download_all_tables(username, password, start_date='2023-01-01'): - """ - Download all required tables from WRDS with fresh connection per table to avoid timeouts. - - Creates a new DuckDB/WRDS connection for each table download to prevent - connection timeout issues during long multi-table downloads. - - Args: - username: WRDS username - password: WRDS password - start_date: Filter for datadate >= start_date (default '2020-01-01') - - Output: - Parquet files for comp_funda, crsp_msf, and ccm. - """ - os.makedirs(OUTPUT_PATH, exist_ok=True) - - wrds_conninfo = gen_wrds_connection_info(username, password) - - # Get table configurations - tables = get_tables_config(start_date) - - # Download each table with a fresh connection to avoid WRDS timeout - for table_name, config in tables.items(): - # Create fresh connection for each table - con = duckdb.connect(":memory:") - con.execute("INSTALL postgres; LOAD postgres;") - con.execute(f"ATTACH '{wrds_conninfo}' AS wrds (TYPE postgres, READ_ONLY)") - - try: - _execute_download( - con=con, - table_name=table_name, - query=config['query'], - output_file=config['output'] - ) - finally: - con.close() - - print("All tables downloaded successfully!", flush=True) - - -if __name__ == "__main__": - # Prompt for WRDS credentials - username = input("Enter WRDS username: ") - password = input("Enter WRDS password: ") - - # Download all tables in one session (recommended to avoid connection timeouts) - download_all_tables(username, password, start_date='1940-01-01') diff --git a/chars_ciz_monthly/functions.py b/chars_ciz_monthly/functions.py deleted file mode 100644 index 9c33e83..0000000 --- a/chars_ciz_monthly/functions.py +++ /dev/null @@ -1,478 +0,0 @@ -import polars as pl -import numpy as np -from tqdm import tqdm -import re - -INPUT_PATH = "../data/raw/" -OUTPUT_PATH = "../data/processed/" - -def ffi49(): - """ - Returns a Polars expression that classifies SIC codes into 49 Fama-French industries. - Usage: df.with_columns(ffi49_polars().alias('ffi49')) - """ - sic = pl.col('sic') - - return ( - pl.when(sic.is_between(100, 199) | sic.is_between(200, 299) | sic.is_between(700, 799) | - sic.is_between(910, 919) | (sic == 2048)).then(1) - .when(sic.is_between(2000, 2009) | sic.is_between(2010, 2019) | sic.is_between(2020, 2029) | - sic.is_between(2030, 2039) | sic.is_between(2040, 2046) | sic.is_between(2050, 2059) | - sic.is_between(2060, 2063) | sic.is_between(2070, 2079) | sic.is_between(2090, 2092) | - (sic == 2095) | sic.is_between(2098, 2099)).then(2) - .when(sic.is_between(2064, 2068) | (sic == 2086) | (sic == 2087) | (sic == 2096) | (sic == 2097)).then(3) - .when(sic.is_between(2080, 2080) | (sic == 2082) | (sic == 2083) | (sic == 2084) | (sic == 2085)).then(4) - .when(sic.is_between(2100, 2199)).then(5) - .when(sic.is_between(920, 999) | sic.is_between(3650, 3652) | (sic == 3732) | - sic.is_between(3930, 3931) | sic.is_between(3940, 3949)).then(6) - .when(sic.is_between(7800, 7829) | sic.is_between(7830, 7833) | sic.is_between(7840, 7841) | - (sic == 7900) | sic.is_between(7910, 7911) | sic.is_between(7920, 7929) | - sic.is_between(7930, 7933) | sic.is_between(7940, 7949) | (sic == 7980) | sic.is_between(7990, 7999)).then(7) - .when(sic.is_between(2700, 2709) | sic.is_between(2710, 2719) | sic.is_between(2720, 2729) | - sic.is_between(2730, 2739) | sic.is_between(2740, 2749) | sic.is_between(2770, 2771) | - sic.is_between(2780, 2789) | sic.is_between(2790, 2799)).then(8) - .when((sic == 2047) | sic.is_between(2391, 2392) | sic.is_between(2510, 2519) | sic.is_between(2590, 2599) | - sic.is_between(2840, 2844) | sic.is_between(3160, 3161) | sic.is_between(3170, 3172) | - sic.is_between(3190, 3199) | (sic == 3229) | (sic == 3260) | sic.is_between(3262, 3263) | (sic == 3269) | - sic.is_between(3230, 3231) | sic.is_between(3630, 3639) | sic.is_between(3750, 3751) | (sic == 3800) | - sic.is_between(3860, 3861) | sic.is_between(3870, 3873) | sic.is_between(3910, 3911) | (sic == 3914) | - (sic == 3915) | sic.is_between(3960, 3962) | (sic == 3991) | (sic == 3995)).then(9) - .when(sic.is_between(2300, 2390) | sic.is_between(3020, 3021) | sic.is_between(3100, 3111) | - (sic == 3130) | (sic == 3131) | sic.is_between(3140, 3151) | sic.is_between(3963, 3965)).then(10) - .when(sic.is_between(8000, 8099)).then(11) - .when((sic == 3693) | sic.is_between(3840, 3851)).then(12) - .when(sic.is_between(2830, 2836)).then(13) - .when(sic.is_between(2800, 2809) | sic.is_between(2810, 2819) | sic.is_between(2820, 2829) | - sic.is_between(2850, 2859) | sic.is_between(2860, 2869) | sic.is_between(2870, 2879) | - sic.is_between(2890, 2899)).then(14) - .when((sic == 3031) | (sic == 3041) | sic.is_between(3050, 3053) | sic.is_between(3060, 3069) | - sic.is_between(3070, 3079) | sic.is_between(3080, 3089) | sic.is_between(3090, 3099)).then(15) - .when(sic.is_between(2200, 2269) | sic.is_between(2270, 2279) | sic.is_between(2280, 2284) | - sic.is_between(2290, 2295) | (sic == 2297) | (sic == 2298) | (sic == 2299) | - sic.is_between(2393, 2395) | sic.is_between(2397, 2399)).then(16) - .when(sic.is_between(800, 899) | sic.is_between(2400, 2439) | sic.is_between(2450, 2459) | - sic.is_between(2490, 2499) | sic.is_between(2660, 2661) | sic.is_between(2950, 2952) | - (sic == 3200) | sic.is_between(3210, 3211) | sic.is_between(3240, 3241) | sic.is_between(3250, 3259) | - (sic == 3261) | (sic == 3264) | sic.is_between(3270, 3275) | sic.is_between(3280, 3281) | - sic.is_between(3290, 3293) | sic.is_between(3295, 3299) | sic.is_between(3420, 3429) | - sic.is_between(3430, 3433) | sic.is_between(3440, 3442) | (sic == 3446) | (sic == 3448) | - (sic == 3449) | sic.is_between(3450, 3452) | sic.is_between(3490, 3499) | (sic == 3996)).then(17) - .when(sic.is_between(1500, 1511) | sic.is_between(1520, 1529) | sic.is_between(1530, 1539) | - sic.is_between(1540, 1549) | sic.is_between(1600, 1699) | sic.is_between(1700, 1799)).then(18) - .when((sic == 3300) | sic.is_between(3310, 3317) | sic.is_between(3320, 3325) | sic.is_between(3330, 3339) | - sic.is_between(3340, 3341) | sic.is_between(3350, 3357) | sic.is_between(3360, 3369) | - sic.is_between(3370, 3379) | sic.is_between(3390, 3399)).then(19) - .when((sic == 3400) | sic.is_between(3443, 3444) | sic.is_between(3460, 3469) | sic.is_between(3470, 3479)).then(20) - .when(sic.is_between(3510, 3519) | sic.is_between(3520, 3529) | (sic == 3530) | (sic == 3531) | - (sic == 3532) | (sic == 3533) | (sic == 3534) | (sic == 3535) | (sic == 3536) | (sic == 3538) | - sic.is_between(3540, 3549) | sic.is_between(3550, 3559) | sic.is_between(3560, 3569) | - (sic == 3580) | (sic == 3581) | (sic == 3582) | (sic == 3585) | (sic == 3586) | - (sic == 3589) | sic.is_between(3590, 3599)).then(21) - .when((sic == 3600) | sic.is_between(3610, 3613) | sic.is_between(3620, 3621) | sic.is_between(3623, 3629) | - sic.is_between(3640, 3646) | sic.is_between(3648, 3649) | (sic == 3660) | (sic == 3690) | - sic.is_between(3691, 3692) | sic.is_between(3699, 3699)).then(22) - .when((sic == 2296) | (sic == 2396) | sic.is_between(3010, 3011) | (sic == 3537) | (sic == 3647) | - (sic == 3694) | (sic == 3700) | (sic == 3710) | (sic == 3711) | (sic == 3713) | - (sic == 3714) | (sic == 3715) | (sic == 3716) | (sic == 3792) | sic.is_between(3790, 3791) | - sic.is_between(3799, 3799)).then(23) - .when((sic == 3720) | (sic == 3721) | sic.is_between(3723, 3725) | sic.is_between(3728, 3729)).then(24) - .when(sic.is_between(3730, 3731) | sic.is_between(3740, 3743)).then(25) - .when(sic.is_between(3760, 3769) | (sic == 3795) | sic.is_between(3480, 3489)).then(26) - .when(sic.is_between(1040, 1049)).then(27) - .when(sic.is_between(1000, 1009) | sic.is_between(1010, 1019) | sic.is_between(1020, 1029) | - sic.is_between(1030, 1039) | sic.is_between(1050, 1059) | sic.is_between(1060, 1069) | - sic.is_between(1070, 1079) | sic.is_between(1080, 1089) | sic.is_between(1090, 1099) | - sic.is_between(1100, 1119) | sic.is_between(1400, 1499)).then(28) - .when(sic.is_between(1200, 1299)).then(29) - .when((sic == 1300) | sic.is_between(1310, 1319) | sic.is_between(1320, 1329) | sic.is_between(1330, 1339) | - sic.is_between(1370, 1379) | (sic == 1380) | (sic == 1381) | (sic == 1382) | (sic == 1389) | - sic.is_between(2900, 2912) | sic.is_between(2990, 2999)).then(30) - .when((sic == 4900) | sic.is_between(4910, 4911) | sic.is_between(4920, 4925) | (sic == 4930) | - (sic == 4931) | (sic == 4932) | (sic == 4939) | sic.is_between(4940, 4942)).then(31) - .when((sic == 4800) | sic.is_between(4810, 4813) | sic.is_between(4820, 4822) | sic.is_between(4830, 4839) | - sic.is_between(4840, 4841) | sic.is_between(4880, 4889) | (sic == 4890) | (sic == 4891) | - (sic == 4892) | sic.is_between(4899, 4899)).then(32) - .when(sic.is_between(7020, 7021) | sic.is_between(7030, 7033) | (sic == 7200) | sic.is_between(7210, 7212) | - (sic == 7214) | sic.is_between(7215, 7217) | sic.is_between(7219, 7221) | sic.is_between(7230, 7231) | - sic.is_between(7240, 7241) | sic.is_between(7250, 7251) | sic.is_between(7260, 7269) | - sic.is_between(7270, 7291) | sic.is_between(7292, 7299) | (sic == 7395) | (sic == 7500) | - sic.is_between(7520, 7529) | sic.is_between(7530, 7539) | sic.is_between(7540, 7549) | - (sic == 7600) | (sic == 7620) | (sic == 7622) | (sic == 7623) | (sic == 7629) | - sic.is_between(7630, 7631) | sic.is_between(7640, 7641) | sic.is_between(7690, 7699) | - sic.is_between(8100, 8199) | sic.is_between(8200, 8299) | sic.is_between(8300, 8399) | - sic.is_between(8400, 8499) | sic.is_between(8600, 8699) | sic.is_between(8800, 8899) | - sic.is_between(7510, 7515)).then(33) - .when(sic.is_between(2750, 2759) | (sic == 3993) | (sic == 7218) | (sic == 7300) | - sic.is_between(7310, 7319) | sic.is_between(7320, 7329) | sic.is_between(7330, 7339) | - sic.is_between(7340, 7342) | (sic == 7349) | sic.is_between(7350, 7353) | (sic == 7359) | - sic.is_between(7360, 7369) | (sic == 7374) | (sic == 7376) | (sic == 7377) | (sic == 7378) | - (sic == 7379) | (sic == 7380) | sic.is_between(7381, 7385) | sic.is_between(7389, 7394) | - (sic == 7396) | (sic == 7397) | sic.is_between(7399, 7399) | (sic == 7519) | (sic == 8700) | - sic.is_between(8710, 8713) | sic.is_between(8720, 8721) | sic.is_between(8730, 8734) | - sic.is_between(8740, 8748) | sic.is_between(8900, 8911) | sic.is_between(8920, 8999) | - sic.is_between(4220, 4229)).then(34) - .when(sic.is_between(3570, 3579) | (sic == 3680) | sic.is_between(3681, 3689) | (sic == 3695)).then(35) - .when(sic.is_between(7370, 7372) | (sic == 7375) | (sic == 7373)).then(36) - .when((sic == 3622) | sic.is_between(3661, 3666) | sic.is_between(3669, 3679) | (sic == 3810) | (sic == 3812)).then(37) - .when((sic == 3811) | (sic == 3820) | sic.is_between(3821, 3827) | sic.is_between(3829, 3839)).then(38) - .when(sic.is_between(2520, 2549) | sic.is_between(2600, 2639) | sic.is_between(2670, 2699) | - sic.is_between(2760, 2761) | sic.is_between(3950, 3955)).then(39) - .when(sic.is_between(2440, 2449) | sic.is_between(2640, 2659) | sic.is_between(3220, 3221) | - sic.is_between(3410, 3412)).then(40) - .when(sic.is_between(4000, 4013) | sic.is_between(4040, 4049) | (sic == 4100) | sic.is_between(4110, 4121) | - sic.is_between(4130, 4131) | sic.is_between(4140, 4142) | sic.is_between(4150, 4151) | - sic.is_between(4170, 4173) | sic.is_between(4190, 4199) | (sic == 4200) | sic.is_between(4210, 4219) | - sic.is_between(4230, 4231) | sic.is_between(4240, 4249) | sic.is_between(4400, 4499) | - sic.is_between(4500, 4599) | sic.is_between(4600, 4699) | (sic == 4700) | sic.is_between(4710, 4712) | - sic.is_between(4720, 4729) | sic.is_between(4730, 4739) | sic.is_between(4740, 4749) | - (sic == 4780) | (sic == 4782) | (sic == 4783) | (sic == 4784) | (sic == 4785) | sic.is_between(4789, 4789)).then(41) - .when((sic == 5000) | sic.is_between(5010, 5015) | sic.is_between(5020, 5023) | sic.is_between(5030, 5039) | - sic.is_between(5040, 5049) | sic.is_between(5050, 5059) | (sic == 5060) | (sic == 5063) | - (sic == 5064) | (sic == 5065) | sic.is_between(5070, 5078) | (sic == 5080) | sic.is_between(5081, 5088) | - (sic == 5090) | sic.is_between(5091, 5094) | (sic == 5099) | (sic == 5100) | sic.is_between(5110, 5113) | - sic.is_between(5120, 5122) | sic.is_between(5130, 5139) | sic.is_between(5140, 5149) | - sic.is_between(5150, 5159) | sic.is_between(5160, 5169) | sic.is_between(5170, 5172) | - sic.is_between(5180, 5182) | sic.is_between(5190, 5199)).then(42) - .when((sic == 5200) | sic.is_between(5210, 5219) | sic.is_between(5220, 5229) | sic.is_between(5230, 5231) | - sic.is_between(5250, 5251) | sic.is_between(5260, 5261) | sic.is_between(5270, 5271) | - (sic == 5300) | sic.is_between(5310, 5311) | (sic == 5320) | sic.is_between(5330, 5331) | - (sic == 5334) | sic.is_between(5340, 5349) | sic.is_between(5390, 5400) | sic.is_between(5410, 5412) | - sic.is_between(5420, 5429) | sic.is_between(5430, 5439) | sic.is_between(5440, 5449) | - sic.is_between(5450, 5459) | sic.is_between(5460, 5469) | sic.is_between(5490, 5500) | - sic.is_between(5510, 5529) | sic.is_between(5530, 5539) | sic.is_between(5540, 5549) | - sic.is_between(5550, 5559) | sic.is_between(5560, 5569) | sic.is_between(5570, 5579) | - sic.is_between(5590, 5599) | sic.is_between(5600, 5700) | sic.is_between(5710, 5722) | - sic.is_between(5730, 5736) | sic.is_between(5750, 5799) | (sic == 5900) | sic.is_between(5910, 5912) | - sic.is_between(5920, 5929) | sic.is_between(5930, 5932) | (sic == 5940) | sic.is_between(5941, 5949) | - sic.is_between(5950, 5959) | sic.is_between(5960, 5969) | sic.is_between(5970, 5979) | - sic.is_between(5980, 5990) | (sic == 5992) | (sic == 5993) | (sic == 5994) | (sic == 5995) | - sic.is_between(5999, 5999)).then(43) - .when(sic.is_between(5800, 5819) | sic.is_between(5820, 5829) | sic.is_between(5890, 5899) | - (sic == 7000) | sic.is_between(7010, 7019) | sic.is_between(7040, 7049) | (sic == 7213)).then(44) - .when((sic == 6000) | sic.is_between(6010, 6036) | sic.is_between(6040, 6062) | sic.is_between(6080, 6082) | - sic.is_between(6090, 6100) | sic.is_between(6110, 6113) | sic.is_between(6120, 6179) | - sic.is_between(6190, 6199)).then(45) - .when((sic == 6300) | sic.is_between(6310, 6331) | sic.is_between(6350, 6351) | sic.is_between(6360, 6361) | - sic.is_between(6370, 6379) | sic.is_between(6390, 6411)).then(46) - .when((sic == 6500) | (sic == 6510) | sic.is_between(6512, 6515) | sic.is_between(6517, 6519) | - sic.is_between(6520, 6532) | sic.is_between(6540, 6541) | sic.is_between(6550, 6553) | - sic.is_between(6590, 6599) | sic.is_between(6610, 6611)).then(47) - .when(sic.is_between(6200, 6299) | (sic == 6700) | sic.is_between(6710, 6726) | sic.is_between(6730, 6733) | - sic.is_between(6740, 6779) | sic.is_between(6790, 6795) | (sic == 6798) | sic.is_between(6799, 6799)).then(48) - .when(sic.is_between(4950, 4959) | sic.is_between(4960, 4961) | sic.is_between(4970, 4971) | - sic.is_between(4990, 4991)).then(49) - .otherwise(None) - ) - - -def ffi30(): - """ - Returns a Polars expression that classifies SIC codes into 30 Fama-French industries. - Usage: df.with_columns(ffi30_polars().alias('ffi30')) - """ - sic = pl.col('sic') - - return ( - pl.when(sic.is_between(100, 199) | sic.is_between(200, 299) | sic.is_between(700, 799) | - sic.is_between(910, 919) | sic.is_between(2000, 2099)).then(1) - .when(sic.is_between(2080, 2085)).then(2) - .when(sic.is_between(2100, 2199)).then(3) - .when(sic.is_between(920, 999) | sic.is_between(3650, 3652) | (sic == 3732) | - sic.is_between(3930, 3949) | sic.is_between(7800, 7999)).then(4) - .when(sic.is_between(2700, 2799) | (sic == 3993)).then(5) - .when((sic == 2047) | sic.is_between(2391, 2392) | sic.is_between(2510, 2519) | - sic.is_between(2590, 2599) | sic.is_between(2840, 2844) | sic.is_between(3160, 3172) | - sic.is_between(3190, 3199) | (sic == 3229) | (sic == 3260) | sic.is_between(3262, 3263) | - (sic == 3269) | sic.is_between(3230, 3231) | sic.is_between(3630, 3639) | - sic.is_between(3750, 3751) | (sic == 3800) | sic.is_between(3860, 3873) | - sic.is_between(3910, 3911) | (sic == 3914) | (sic == 3915) | sic.is_between(3960, 3962) | - (sic == 3991) | (sic == 3995)).then(6) - .when(sic.is_between(2300, 2390) | sic.is_between(3020, 3021) | sic.is_between(3100, 3111) | - (sic == 3130) | (sic == 3131) | sic.is_between(3140, 3151) | sic.is_between(3963, 3965)).then(7) - .when(sic.is_between(2830, 2836) | (sic == 3693) | sic.is_between(3840, 3851) | sic.is_between(8000, 8099)).then(8) - .when(sic.is_between(2800, 2829) | sic.is_between(2850, 2899)).then(9) - .when(sic.is_between(2200, 2284) | sic.is_between(2290, 2295) | (sic == 2297) | (sic == 2298) | - (sic == 2299) | sic.is_between(2393, 2399)).then(10) - .when(sic.is_between(800, 899) | sic.is_between(1500, 1799) | sic.is_between(2400, 2439) | - sic.is_between(2450, 2459) | sic.is_between(2490, 2499) | sic.is_between(2660, 2661) | - sic.is_between(2950, 2952) | (sic == 3200) | sic.is_between(3210, 3211) | - sic.is_between(3240, 3241) | sic.is_between(3250, 3259) | (sic == 3261) | (sic == 3264) | - sic.is_between(3270, 3275) | sic.is_between(3280, 3281) | sic.is_between(3290, 3299) | - sic.is_between(3420, 3442) | (sic == 3446) | (sic == 3448) | (sic == 3449) | - sic.is_between(3450, 3452) | sic.is_between(3490, 3499) | (sic == 3996)).then(11) - .when((sic == 3300) | sic.is_between(3310, 3317) | sic.is_between(3320, 3325) | - sic.is_between(3330, 3341) | sic.is_between(3350, 3357) | sic.is_between(3360, 3369) | - sic.is_between(3370, 3379) | sic.is_between(3390, 3399)).then(12) - .when((sic == 3400) | sic.is_between(3443, 3444) | sic.is_between(3460, 3479) | - sic.is_between(3510, 3599)).then(13) - .when((sic == 3600) | sic.is_between(3610, 3613) | sic.is_between(3620, 3621) | - sic.is_between(3623, 3629) | sic.is_between(3640, 3660) | (sic == 3690) | - sic.is_between(3691, 3692) | sic.is_between(3699, 3699)).then(14) - .when((sic == 2296) | (sic == 2396) | sic.is_between(3010, 3011) | (sic == 3537) | - (sic == 3647) | (sic == 3694) | (sic == 3700) | sic.is_between(3710, 3716) | - (sic == 3792) | sic.is_between(3790, 3791) | sic.is_between(3799, 3799)).then(15) - .when(sic.is_between(3720, 3721) | sic.is_between(3723, 3725) | sic.is_between(3728, 3731) | - sic.is_between(3740, 3743)).then(16) - .when(sic.is_between(1000, 1119) | sic.is_between(1400, 1499)).then(17) - .when(sic.is_between(1200, 1299)).then(18) - .when((sic == 1300) | sic.is_between(1310, 1389) | sic.is_between(2900, 2912) | sic.is_between(2990, 2999)).then(19) - .when((sic == 4900) | sic.is_between(4910, 4942)).then(20) - .when((sic == 4800) | sic.is_between(4810, 4899)).then(21) - .when(sic.is_between(7020, 7021) | sic.is_between(7030, 7033) | (sic == 7200) | - sic.is_between(7210, 7299) | (sic == 7395) | (sic == 7500) | sic.is_between(7510, 7549) | - (sic == 7600) | sic.is_between(7620, 7641) | sic.is_between(7690, 7699) | - sic.is_between(8100, 8199) | sic.is_between(8200, 8299) | sic.is_between(8300, 8399) | - sic.is_between(8400, 8499) | sic.is_between(8600, 8748) | sic.is_between(8800, 8999)).then(22) - .when(sic.is_between(3570, 3579) | (sic == 3622) | sic.is_between(3661, 3679) | - sic.is_between(3680, 3689) | (sic == 3695) | sic.is_between(3810, 3812) | - sic.is_between(3820, 3839) | (sic == 7373)).then(23) - .when(sic.is_between(2440, 2449) | sic.is_between(2520, 2549) | sic.is_between(2600, 2639) | - sic.is_between(2640, 2659) | sic.is_between(2670, 2699) | sic.is_between(2760, 2761) | - sic.is_between(3220, 3221) | sic.is_between(3410, 3412) | sic.is_between(3950, 3955)).then(24) - .when(sic.is_between(4000, 4013) | sic.is_between(4040, 4049) | (sic == 4100) | - sic.is_between(4110, 4173) | sic.is_between(4190, 4231) | sic.is_between(4240, 4249) | - sic.is_between(4400, 4499) | sic.is_between(4500, 4599) | sic.is_between(4600, 4699) | - (sic == 4700) | sic.is_between(4710, 4749) | (sic == 4780) | sic.is_between(4782, 4789)).then(25) - .when((sic == 5000) | sic.is_between(5010, 5199)).then(26) - .when((sic == 5200) | sic.is_between(5210, 5999)).then(27) - .when(sic.is_between(5800, 5829) | sic.is_between(5890, 5899) | (sic == 7000) | - sic.is_between(7010, 7019) | sic.is_between(7040, 7049) | (sic == 7213)).then(28) - .when((sic == 6000) | sic.is_between(6010, 6799)).then(29) - .when(sic.is_between(4950, 4961) | sic.is_between(4970, 4971) | sic.is_between(4990, 4991)).then(30) - .otherwise(None) - ) - - -def ffi12(): - """ - Returns a Polars expression that classifies SIC codes into 12 Fama-French industries. - Usage: df.with_columns(ffi12_polars().alias('ffi12')) - """ - sic = pl.col('sic') - - return ( - pl.when(sic.is_between(100, 999) | sic.is_between(2000, 2399) | sic.is_between(2700, 2749) | - sic.is_between(2770, 2799) | sic.is_between(3100, 3199) | sic.is_between(3940, 3989)).then(1) - .when(sic.is_between(2500, 2519) | sic.is_between(2590, 2599) | sic.is_between(3630, 3659) | - sic.is_between(3710, 3711) | (sic == 3714) | (sic == 3716) | sic.is_between(3750, 3751) | - (sic == 3792) | sic.is_between(3900, 3939) | sic.is_between(3990, 3999)).then(2) - .when(sic.is_between(2520, 2589) | sic.is_between(2600, 2699) | sic.is_between(2750, 2769) | - sic.is_between(3000, 3099) | sic.is_between(3200, 3569) | sic.is_between(3580, 3629) | - sic.is_between(3700, 3709) | sic.is_between(3712, 3713) | (sic == 3715) | - sic.is_between(3717, 3749) | sic.is_between(3752, 3791) | sic.is_between(3793, 3799) | - sic.is_between(3830, 3839) | sic.is_between(3860, 3899)).then(3) - .when(sic.is_between(1200, 1399) | sic.is_between(2900, 2999)).then(4) - .when(sic.is_between(2800, 2829) | sic.is_between(2840, 2899)).then(5) - .when(sic.is_between(3570, 3579) | sic.is_between(3660, 3692) | sic.is_between(3694, 3699) | - sic.is_between(3810, 3829) | sic.is_between(7370, 7379)).then(6) - .when(sic.is_between(4800, 4899)).then(7) - .when(sic.is_between(4900, 4949)).then(8) - .when(sic.is_between(5000, 5999) | sic.is_between(7200, 7299) | sic.is_between(7600, 7699)).then(9) - .when(sic.is_between(2830, 2839) | (sic == 3693) | sic.is_between(3840, 3859) | sic.is_between(8000, 8099)).then(10) - .when(sic.is_between(6000, 6999)).then(11) - .otherwise(12) - ) - -####################################################################################################################### -# TTM functions # -####################################################################################################################### - - -def ttm4(series, df): - """ - Calculate trailing 4-period sum (TTM4) using Polars. - - :param series: variables' name (string) - :param df: polars dataframe (used to compute the result as a Series) - :return: polars Expression that can be used in with_columns() - - Note: This function returns a Polars Expression for use in with_columns(). - Example: data_rawq.with_columns([ttm4('ibq', data_rawq).alias('ibq4')]) - """ - # Build expression for sum of current + 3 lags - return ( - pl.col(series) + - pl.col(series).shift(1).over('permno') + - pl.col(series).shift(2).over('permno') + - pl.col(series).shift(3).over('permno') - ) - - -def ttm12(series, df): - """ - Calculate trailing 12-period sum (TTM12) using Polars. - - :param series: variables' name (string) - :param df: polars dataframe (used to compute the result as a Series) - :return: polars Expression that can be used in with_columns() - - Note: This function returns a Polars Expression for use in with_columns(). - Example: crsp_mom.with_columns([(ttm12('mdivpay', crsp_mom) / pl.col('me')).alias('dy')]) - """ - # Build expression for sum of current + 11 lags - return ( - pl.col(series) + - pl.col(series).shift(1).over('permno') + - pl.col(series).shift(2).over('permno') + - pl.col(series).shift(3).over('permno') + - pl.col(series).shift(4).over('permno') + - pl.col(series).shift(5).over('permno') + - pl.col(series).shift(6).over('permno') + - pl.col(series).shift(7).over('permno') + - pl.col(series).shift(8).over('permno') + - pl.col(series).shift(9).over('permno') + - pl.col(series).shift(10).over('permno') + - pl.col(series).shift(11).over('permno') - ) - -def fillna_atq(df_q: pl.DataFrame, df_a: pl.DataFrame): - """ - Use annual chars to fill null values in quarterly chars. - Skips columns matching 'mom*' pattern. - """ - # find columns that are null in df_q AND exist in df_a - q_null_cols = [c for c in df_q.columns if df_q[c].is_null().any()] - a_cols = df_a.columns - candidates = list(set(q_null_cols) & set(a_cols)) - - # exclude mom* columns - na_columns_list = [c for c in candidates if not re.match(r'mom.', c)] - - if not na_columns_list: - return df_q - - # extract annual cols + keys, rename to '*_a' - df_temp = ( - df_a.select(['permno', 'date'] + na_columns_list) - .rename({c: f'{c}_a' for c in na_columns_list}) - ) - - # left join and coalesce - df_q = df_q.join(df_temp, on=['permno', 'date'], how='left') - df_q = df_q.with_columns([ - pl.coalesce([pl.col(c), pl.col(f'{c}_a')]).alias(c) - for c in na_columns_list - ]).drop([f'{c}_a' for c in na_columns_list]) - - return df_q - - -def fillna_ind( - df: pl.DataFrame, - method: str, - ffi: int, - not_fill_col: list -): - """ - Fill null values using industry-level mean or median grouped by date + ffi code. - """ - ffi_col = f'ffi{ffi}' - na_columns_list = [ - c for c in df.columns - if df[c].is_null().any() and c not in not_fill_col - ] - - if not na_columns_list: - return df - - if method == 'mean': - agg_exprs = [pl.col(c).mean().alias(f'{c}_fill') for c in na_columns_list] - elif method == 'median': - agg_exprs = [pl.col(c).median().alias(f'{c}_fill') for c in na_columns_list] - else: - raise ValueError(f"method must be 'mean' or 'median', got '{method}'") - - df_fill = ( - df.group_by(['date', ffi_col]) - .agg(agg_exprs) - ) - - df = df.join(df_fill, on=['date', ffi_col], how='left') - df = df.with_columns([ - pl.coalesce([pl.col(c), pl.col(f'{c}_fill')]).alias(c) - for c in na_columns_list - ]).drop([f'{c}_fill' for c in na_columns_list]) - - return df - - -def fillna_all( - df: pl.DataFrame, - method: str, - not_fill_col: list -): - """ - Fill null values using cross-sectional mean or median grouped by date. - """ - na_columns_list = [ - c for c in df.columns - if df[c].is_null().any() and c not in not_fill_col - ] - - if not na_columns_list: - return df - - if method == 'mean': - agg_exprs = [pl.col(c).mean().alias(f'{c}_fill') for c in na_columns_list] - elif method == 'median': - agg_exprs = [pl.col(c).median().alias(f'{c}_fill') for c in na_columns_list] - else: - raise ValueError(f"method must be 'mean' or 'median', got '{method}'") - - df_fill = ( - df.group_by('date') - .agg(agg_exprs) - ) - - df = df.join(df_fill, on='date', how='left') - df = df.with_columns([ - pl.coalesce([pl.col(c), pl.col(f'{c}_fill')]).alias(c) - for c in na_columns_list - ]).drop([f'{c}_fill' for c in na_columns_list]) - - return df - - -def standardize(df: pl.DataFrame): - """ - Cross-sectionally rank and standardize all char columns to [-1, 1]. - Excludes info columns. Null ranks are filled with 0. - """ - INFO_COLS = { - 'permno', 'date', 'datadate', 'gvkey', 'sic', 'count', - 'exchcd', 'shrcd', 'ffi49', 'ret', 'retadj', 'retx', - 'lag_me', 'ticker', 'conm', 'comnam', 'prc', 'shrout', - 'size_grp', - 'primaryexch', 'conditionaltype', 'tradingstatusflg', - 'sharetype', 'securitytype', 'securitysubtype', - 'usincflg', 'issuertype', - } - col_names = [c for c in df.columns if c not in INFO_COLS] - - for col_name in tqdm(col_names): - # dense rank within each date, then scale to [-1, 1] - df = df.with_columns([ - pl.col(col_name) - .rank(method='dense') - .over('date') - .alias(f'_rank_{col_name}') - ]) - # count non-null unique values per date - df = df.with_columns([ - pl.col(f'_rank_{col_name}') - .max() - .over('date') - .alias('_max_rank') - ]) - df = df.with_columns([ - pl.when(pl.col('_max_rank') > 1) - .then( - (pl.col(f'_rank_{col_name}') - 1) / - (pl.col('_max_rank') - 1) * 2 - 1 - ) - .otherwise(None) - .fill_null(0) - .alias(f'rank_{col_name}') - ]).drop([col_name, f'_rank_{col_name}', '_max_rank']) - - return df diff --git a/chars_ciz_monthly/iclink_ciz.sas b/chars_ciz_monthly/iclink_ciz.sas deleted file mode 100644 index a98ca68..0000000 --- a/chars_ciz_monthly/iclink_ciz.sas +++ /dev/null @@ -1,281 +0,0 @@ -/* ********************************************************************************* */ -/* ******************** W R D S R E S E A R C H M A C R O S ******************** */ -/* ********************************************************************************* */ -/* WRDS Macro: ICLINK_CIZ */ -/* Summary : Create IBES-CRSP Link Table */ -/* Author : Rabih Moussawi, WRDS */ -/* Date : September 25, 2006 */ -/* Update : November 2024 by Freda Drechsler for CRSP CIZ data format */ -/* Variables : - IBESID and CRSPID are IBES and CRSP Names Datasets */ -/* - OUTSET: IBES-CRSP link table output dataset */ -/* ********************************************************************************* */ - -%MACRO ICLINK_ciz (IBESID=IBES.ID,CRSPID=CRSP.STOCKNAMES_v2,OUTSET=WORK.ICLINK); - -/* ********************************************************************************* */ -/* FUNCTION: - Creates a link table between IBES TICKER and CRSP PERMNO */ -/* - Scores links from 0 (best link) to 6 */ -/* Possible IBES ID (names) file to use: */ -/* Detail History: ID File */ -/* Summary History: IDSUM File */ -/* Recommendation Detail and Summary Statistics: RECDID and RECDIDSUM Files */ -/* */ -/* INPUT: IBES and CRSP ID (or NAMES) Datasets, with historical identifiers list */ -/* - IBES: IBES.ID, IBES.IDSUM, IBES.RECID, or IBES.RECIDSUM files */ -/* - CRSP: CRSP.MSENAMES, CRSP.DSENAMES, or CRSP.STOCKNAMES files */ -/* */ -/* OUTPUT: ICLINK set stored in prespecified directory */ -/* - SCORE variable: lower scores are better and high scores may need further */ -/* checking before using them to link CRSP & IBES data. */ -/* In computing the score, a CUSIP match is considered better than a */ -/* TICKER match. The score also includes a penalty for differences in */ -/* company names-- CNAME in IBES and COMNAM in CRSP. Name penalty is */ -/* based upon SPEDIS, which is the spelling distance function in SAS. */ -/* SPEDIS=0 is a perfect score and SPEDIS<30 is usually good */ -/* enough to be considered a name match. */ -/* Note here that Exchange Ticker can also be used as Flag */ -/* "SCORE" levels: */ -/* - 0: BEST match: using (cusip, cusip dates and company names) */ -/* or (exchange ticker, company names and 6-digit cusip) */ -/* - 1: Cusips and cusip dates match but company names do not match */ -/* - 2: Cusips and company names match but cusip dates do not match */ -/* - 3: Cusips match but cusip dates and company names do not match */ -/* - 4: tickers and 6-digit cusips match but comp names do not match */ -/* - 5: tickers and names match but 6-digit cusips do not match */ -/* - 6: tickers match but names and 6-digit cusips do not match */ -/* ********************************************************************************* */ - -options nonotes; -/* Check Validity of Library Assignments */ -%if (%sysfunc(libref(crsp))) %then %do; - %let cs=/wrds/crsp/sasdata/; - libname crsp ("&cs/m_stock_v2","&cs/q_stock_v2","&cs/a_stock_v2"); -%end; -%if (%sysfunc(libref(ibes))) %then %do; libname ibes "/wrds/ibes/sasdata"; %end; - -/* Name End Dates variable in MSENAMES and DSENAMES is different than STOCKNAMES */ -%if %sysfunc(upcase(&CRSPID)) ne CRSP.STOCKNAMES_V2 %then %let condition = %str(RENAME=(nameendt=nameenddt)); -%else %let condition = ; - -%put ; -%put ### START. Creating IBES-CRSP Link Table: ICLINK ; -%put ## IBES NAMES (ID) Dataset Used: &IBESID; -%put ## CRSP NAMES (ID) Dataset Used: &CRSPID; -%put ## Step1: Linking using CUSIPs... ; - -/* Step 1: Link by CUSIP */ -/* IBES: Get the list of IBES TICKERS for US firms in IBES */ -proc sort data=&IBESID out=_IBES1 (keep=ticker cusip CNAME sdates); - where USFIRM=1 and not(missing(cusip)); - by ticker cusip sdates; -run; - -/* Create first and last 'start dates' for CUSIP link */ -proc sql; - create table _IBES2 - as select *, min(sdates) as fdate, max(sdates) as ldate - from _IBES1 - group by ticker, cusip - order by ticker, cusip, sdates; -quit; - -/* Label date range variables and keep only most recent company name for CUSIP link */ -data _IBES2; - set _IBES2; - by ticker cusip; - if last.cusip; - label fdate="First Start date of CUSIP record"; - label ldate="Last Start date of CUSIP record"; - format fdate ldate date9.; - drop sdates; -run; - -/* CRSP: Get all PERMNO-NCUSIP combinations */ -proc sort data=&CRSPID out=_CRSP1 (keep=PERMNO CUSIP IssuerNm name: &condition); - where not missing(CUSIP); - by PERMNO CUSIP namedt; -run; - -/* Arrange effective dates for CUSIP link */ -proc sql; - create table _CRSP2 - as select PERMNO,CUSIP,IssuerNm,min(namedt)as namedt,max(nameenddt) as nameenddt - from _CRSP1 - group by PERMNO, CUSIP - order by PERMNO, CUSIP, NAMEDT; -quit; - -/* Label date range variables and keep only most recent company name */ -data _CRSP2; - set _CRSP2; - by permno cusip; - if last.cusip; - label namedt="Start date of CUSIP record"; - label nameenddt="End date of CUSIP record"; - format namedt nameenddt date9.; -run; - -/* Create CUSIP Link Table */ -/* CUSIP date ranges are only used in scoring as CUSIPs are not reused for - different companies overtime */ -proc sql; - create table _LINK1_1 - as select * - from _IBES2 as a, _CRSP2 as b - where a.CUSIP = b.CUSIP - order by TICKER, PERMNO, ldate; -quit; - -/* Score links using CUSIP date range and company name spelling distance */ -/* Idea: date ranges the same cusip was used in CRSP and IBES should intersect */ -data _LINK1_2; - set _LINK1_1; - by TICKER PERMNO; - if last.permno; /* Keep link with most recent company name */ - name_dist = min(spedis(cname,IssuerNm),spedis(IssuerNm,cname)); - if (not ((ldate < namedt) or (fdate > nameenddt))) and name_dist < 30 then SCORE = 0; - else if (not ((ldate < namedt) or (fdate > nameenddt))) then score = 1; - else if name_dist < 30 then SCORE = 2; - else SCORE = 3; - keep TICKER PERMNO cname IssuerNm score; -run; - -%put ## Step2: Linking using TICKERs... ; -/* Step 2: Find links for the remaining unmatched cases using Exchange Ticker */ -/* Identify remaining unmatched cases */ -proc sql; - create table _NOMATCH1 - as select distinct a.* - from _IBES1 (keep=ticker) as a - where a.ticker NOT in (select ticker from _LINK1_2) - order by a.ticker; -quit; - -/* Drop Step1 Tables*/ -proc sql; drop table _IBES1,_IBES2,_CRSP1,_CRSP2; quit; - -/* Add IBES identifying information */ -proc sql; - create table _NOMATCH2 - as select b.ticker, b.CNAME, b.OFTIC, b.sdates, b.cusip - from _NOMATCH1 as a, &IBESID as b - where a.ticker = b.ticker and not (missing(b.OFTIC)) - order by ticker, oftic, sdates; -quit; - -/* Create first and last 'start dates' for Exchange Tickers */ -proc sql; - create table _NOMATCH3 - as select *, min(sdates) as fdate, max(sdates) as ldate - from _NOMATCH2 - group by ticker, oftic - order by ticker, oftic, sdates; -quit; - -/* Label date range variables and keep only most recent company name */ -data _NOMATCH3; - set _NOMATCH3; - by ticker oftic; - if last.oftic; - label fdate="First Start date of OFTIC record"; - label ldate="Last Start date of OFTIC record"; - format fdate ldate date9.; - drop sdates; -run; - -/* Get entire list of CRSP stocks with Exchange Ticker information */ -proc sort data=&CRSPID out=_CRSP1 (keep=ticker IssuerNm permno cusip name: &condition); - where not missing(ticker); - by permno ticker namedt; -run; - -/* Arrange effective dates for link by Exchange Ticker */ -proc sql; - create table _CRSP2 - as select permno,IssuerNm,ticker as crsp_ticker,cusip, - min(namedt)as namedt,max(nameenddt) as nameenddt - from _CRSP1 - group by permno, ticker - order by permno, crsp_ticker, namedt; -quit; -/* CRSP exchange ticker renamed to crsp_ticker to avoid confusion with IBES TICKER */ - -/* Label date range variables and keep only most recent company name */ -data _CRSP2; - set _CRSP2; - by permno crsp_ticker; - if last.crsp_ticker; - label namedt="Start date of exch. ticker record"; - label nameenddt="End date of exch. ticker record"; - format namedt nameenddt date9.; -run; - -/* Merge remaining unmatched cases using Exchange Ticker */ -/* Note: Use ticker date ranges as exchange tickers are reused overtime */ -proc sql; - create table _LINK2_1 - as select a.ticker,a.oftic, b.permno, a.cname, b.IssuerNm, a.cusip, b.cusip as crsp_cusip, a.ldate - from _NOMATCH3 as a, _CRSP2 as b - where a.oftic = b.crsp_ticker and - (ldate>=namedt) and (fdate<=nameenddt) - order by ticker, oftic, ldate; -quit; - -/* Score using company name using 6-digit CUSIP and company name spelling distance */ -data _LINK2_2; - set _LINK2_1; - name_dist = min(spedis(cname,IssuerNm),spedis(IssuerNm,cname)); - if substr(cusip,1,6)=substr(crsp_cusip,1,6) and name_dist < 30 then SCORE=0; - else if substr(cusip,1,6)=substr(crsp_cusip,1,6) then score = 4; - else if name_dist < 30 then SCORE = 5; - else SCORE = 6; -run; - -/* Some companies may have more than one TICKER-PERMNO link, */ -/* so re-sort and keep the case (PERMNO & Company name from CRSP) */ -/* that gives the lowest score for each IBES TICKER (first.ticker=1) */ -proc sort data=_LINK2_2; by ticker score; run; -data _LINK2_3; - set _LINK2_2; - by ticker score; - if first.ticker; - keep ticker permno cname IssuerNm permno score; -run; - -%put ## Step3: Finalizing Links and Scores... ; -/* Step 3: Add Exchange Ticker links to CUSIP links */ -/* Create Labels for ICLINK dataset and variables */ -/* Create final link table and save it in prespecified directory */ -data &OUTSET (label="IBES-CRSP Link Table"); - set _LINK1_2 _LINK2_3; -label CNAME = "Company Name in IBES"; -label IssuerNm= "Company Name in CRSP"; -label SCORE= "Link Score: 0(best) - 6"; -run; - -/* Final Sort */ -proc sort data=&OUTSET; by TICKER SCORE PERMNO; run; - -%put ## Step4: Link Table &OUTSET Ready... ; -/* House Cleaning */ -proc sql; -drop table _CRSP1,_CRSP2, - _LINK1_1,_LINK1_2,_LINK2_1,_LINK2_2,_LINK2_3, - _NOMATCH1,_NOMATCH2,_NOMATCH3; -quit; -%put ### DONE . ; %put ; -options notes; -%MEND ICLINK_ciz; - -%ICLINK_ciz; - -proc export data=WORK.ICLINK - outfile="iclink_ciz.csv" - dbms=csv - replace; -run; - -/* ********************************************************************************* */ -/* ************* Material Copyright Wharton Research Data Services *************** */ -/* ****************************** All Rights Reserved ****************************** */ -/* ********************************************************************************* */ \ No newline at end of file diff --git a/chars_ciz_monthly/impute_rank_output.py b/chars_ciz_monthly/impute_rank_output.py deleted file mode 100644 index 63fa6c8..0000000 --- a/chars_ciz_monthly/impute_rank_output.py +++ /dev/null @@ -1,271 +0,0 @@ -""" -Impute, rank, and output final characteristic datasets. - -Logic: -1. Load chars_a_raw.parquet and chars_q_raw.parquet. -2. For overlapping accounting variables (both annual & quarterly), pick the - value with the more recent datadate; fall back to whichever is non-null. -3. Shift return one period forward (t chars predict t+1 return). -4. Compute ffi49 industry codes for ALL outputs (before branching). -5. Produce four output files: - - chars_raw_no_impute.parquet - - chars_raw_imputed.parquet - - chars_rank_no_impute.parquet - - chars_rank_imputed.parquet -""" - -import polars as pl -from tqdm import tqdm - -from functions import ffi49, fillna_ind, fillna_all, standardize, INPUT_PATH, OUTPUT_PATH - -# ===================================================================== -# Variable lists -# ===================================================================== -OBS_VARS = [ - 'gvkey', 'permno', 'jdate', 'ticker', 'conm', 'comnam', - 'sic', 'ret', 'retx', 'retadj', - 'exchcd', 'shrcd', 'prc', 'shrout', 'size_grp', -] - -ACCOUNTING_VARS = [ - 'datadate', # must be here (not OBS_VARS) so it gets a_/q_ prefix for recency comparison - 'acc', 'bm', 'agr', 'alm', 'ato', 'cash', 'cashdebt', 'cfp', 'chcsho', - 'chtx', 'depr', 'ep', 'gma', 'grltnoa', 'lev', 'lgr', 'ni', 'noa', 'op', - 'pctacc', 'pm', 'rd_sale', 'rdm', 'rna', 'roa', 'roe', 'rsup', 'sgr', 'sp', - 'me_ia', 'bm_ia', - 'cashpr', 'cfp_ia', 'chatoia', 'egr', 'invest', 'chmom', 'rd', -] - -A_ONLY_VARS = [ - 'adm', 'herf', 'hire', - 'absacc', 'age', 'chempia', 'chinv', 'convind', 'currat', 'divi', 'divo', - 'grcapx', 'pchcapx_ia', 'pchcurrat', 'pchdepr', 'pchgm_pchsale', 'pchquick', - 'pchsale_pchinvt', 'pchsale_pchrect', 'pchsale_pchxsga', 'pchsaleinv', - 'quick', 'realestate', 'roic', 'salecash', 'salerec', 'saleinv', - 'secured', 'securedind', 'sin', 'tang', 'tb', 'chpmia', -] - -Q_ONLY_VARS = [ - 'abr', 'sue', 'cinvest', 'nincr', 'pscore', - 'roavol', -] - -M_VARS = [ - 'baspread', 'beta', 'ill', 'maxret', - 'mom12m', 'mom1m', 'mom36m', 'mom60m', 'mom6m', - 're', 'rvar_capm', 'rvar_ff3', 'rvar_mean', - 'seas1a', 'std_dolvol', 'std_turn', 'zerotrade', - 'me', 'dy', 'turn', 'dolvol', 'indmom', -] - - -# ===================================================================== -# Helpers -# ===================================================================== -def _available_cols(df, cols): - """Return only columns that exist in df, preserving order.""" - return [c for c in cols if c in df.columns] - - -def _replace_inf(df): - """Replace ±inf and NaN with null in all float columns.""" - float_cols = [c for c in df.columns if df[c].dtype in (pl.Float32, pl.Float64)] - if not float_cols: - return df - return df.with_columns([ - pl.when(pl.col(c).is_infinite() | pl.col(c).is_nan()) - .then(None) - .otherwise(pl.col(c)) - .alias(c) - for c in float_cols - ]) - - -def _reconcile(df_a, df_q): - """ - Merge annual and quarterly characteristics. - - For variables in ACCOUNTING_VARS: use the more recent datadate when both - frequencies are available, otherwise use whichever is non-null. - """ - a_prefix = {v: f'a_{v}' for v in ACCOUNTING_VARS} - q_prefix = {v: f'q_{v}' for v in ACCOUNTING_VARS} - - # --- annual side: obs + accounting + a_only + monthly --- - a_cols = _available_cols(df_a, OBS_VARS + ACCOUNTING_VARS + A_ONLY_VARS + M_VARS) - df_a_sel = df_a.select(a_cols).rename( - {v: a_prefix[v] for v in ACCOUNTING_VARS if v in a_cols} - ) - - # --- quarterly side: obs + accounting + q_only + quarterly-only monthly vars --- - q_cols = _available_cols(df_q, OBS_VARS + ACCOUNTING_VARS + Q_ONLY_VARS) - q_drop = [c for c in OBS_VARS if c not in ('gvkey', 'permno', 'jdate')] - df_q_sel = ( - df_q.select(q_cols) - .rename({v: q_prefix[v] for v in ACCOUNTING_VARS if v in q_cols}) - .drop([c for c in q_drop if c in q_cols]) - ) - - # merge - df = df_a_sel.join(df_q_sel, on=['gvkey', 'permno', 'jdate'], how='left') - - # reconcile each overlapping variable (skip datadate) - for var in tqdm(ACCOUNTING_VARS[1:], desc='Reconciling A/Q'): - a_col, q_col = f'a_{var}', f'q_{var}' - has_a = a_col in df.columns - has_q = q_col in df.columns - - if not has_a and not has_q: - df = df.with_columns(pl.lit(None).alias(var)) - elif not has_q: - df = df.with_columns(pl.col(a_col).alias(var)).drop(a_col) - elif not has_a: - df = df.with_columns(pl.col(q_col).alias(var)).drop(q_col) - else: - # Both exist: pick by recency, fall back to whichever is available - a_avail = pl.col(a_col).is_not_null() - q_avail = pl.col(q_col).is_not_null() - latest = ( - pl.when(pl.col('q_datadate') < pl.col('a_datadate')) - .then(pl.col(a_col)) - .otherwise(pl.col(q_col)) - ) - available = pl.when(a_avail).then(pl.col(a_col)).otherwise(pl.col(q_col)) - df = df.with_columns( - pl.when(a_avail & q_avail).then(latest).otherwise(available).alias(var) - ).drop([a_col, q_col]) - - # drop frequency-specific datadates - df = df.drop([c for c in ['a_datadate', 'q_datadate'] if c in df.columns]) - return df - - -def _shift_return(df): - """Shift return forward one month: t characteristics predict t+1 return.""" - df = df.sort(['permno', 'jdate']) - df = df.with_columns([ - pl.col('ret').shift(-1).over('permno').alias('ret_lead'), - pl.col('jdate').shift(-1).over('permno').alias('date'), - ]) - df = df.drop('ret').rename({'ret_lead': 'ret'}) - - df = df.filter(pl.col('ret').is_not_null()).drop('jdate') - return _replace_inf(df) - - -def _rank_df(df): - """Rank characteristics cross-sectionally, add log_me, fill nulls with 0.""" - out = df.clone() - out = out.with_columns(pl.col('me').alias('lag_me')) - # bm < 0 → null before ranking (GHZ convention) - if 'bm' in out.columns: - out = out.with_columns( - pl.when(pl.col('bm') < 0).then(None).otherwise(pl.col('bm')).alias('bm') - ) - out = standardize(out) - out = out.with_columns(pl.col('lag_me').log().alias('log_me')) - # (fixed-20260325) convert both NaN and ±inf to 0 after standardization so - # rank-stage float cleanup matches _replace_inf(). - float_cols = [c for c in out.columns if out[c].dtype in (pl.Float32, pl.Float64)] - out = out.with_columns([ - pl.when(pl.col(c).is_infinite() | pl.col(c).is_nan()) - .then(pl.lit(0)) - .otherwise(pl.col(c)) - .alias(c) - for c in float_cols - ]) - rank_cols = [c for c in out.columns if c.startswith('rank_')] - out = out.with_columns([pl.col(c).fill_null(0) for c in rank_cols]) - return out - - -# ===================================================================== -# Main -# ===================================================================== -if __name__ == '__main__': - # ------------------------------------------------------------------ - # Load - # ------------------------------------------------------------------ - print("Loading raw chars...", flush=True) - chars_a = pl.read_parquet(OUTPUT_PATH + 'chars_a_raw.parquet') - chars_q = pl.read_parquet(OUTPUT_PATH + 'chars_q_raw.parquet') - - chars_a = ( - chars_a.drop_nulls(subset=['permno']) - .with_columns(pl.col('permno').cast(pl.Int64), - pl.col('jdate').cast(pl.Date)) - .unique(subset=['permno', 'jdate']) - ) - chars_q = ( - chars_q.drop_nulls(subset=['permno']) - .with_columns(pl.col('permno').cast(pl.Int64), - pl.col('jdate').cast(pl.Date)) - .unique(subset=['permno', 'jdate']) - ) - - print(f" chars_a: {chars_a.shape}, chars_q: {chars_q.shape}", flush=True) - - # ------------------------------------------------------------------ - # Reconcile annual / quarterly - # ------------------------------------------------------------------ - print("Reconciling annual & quarterly...", flush=True) - df = _reconcile(chars_a, chars_q) - - # ------------------------------------------------------------------ - # Shift return forward - # ------------------------------------------------------------------ - print("Shifting return...", flush=True) - if 'retx' in df.columns: - df = df.drop('retx') - df = _shift_return(df) - df = df.sort(['permno', 'date']) - print(f" After shift: {df.shape}", flush=True) - - # ------------------------------------------------------------------ - # Fill SIC + compute ffi49 (shared by ALL outputs — Bug 5 fix) - # ------------------------------------------------------------------ - df = df.with_columns(pl.col('sic').forward_fill().over('permno')) - df = df.with_columns(pl.col('sic').fill_null(0).cast(pl.Int64)) - df = df.with_columns(ffi49().alias('ffi49')) - - # ------------------------------------------------------------------ - # Output 1: raw (no imputation) - # ------------------------------------------------------------------ - print("Saving chars_raw_no_impute.parquet ...", flush=True) - df.write_parquet(OUTPUT_PATH + 'chars_raw_no_impute.parquet') - - # ------------------------------------------------------------------ - # Output 2: imputed (industry-median → cross-sectional-median) - # ------------------------------------------------------------------ - print("Imputing...", flush=True) - df_impute = df.clone() - df_impute = df_impute.with_columns(pl.col('date').cast(pl.Date)) - df_impute = _replace_inf(df_impute) - - df_impute = fillna_ind(df_impute, method='median', ffi=49, not_fill_col=OBS_VARS) - df_impute = fillna_all(df_impute, method='median', not_fill_col=OBS_VARS) - - # IBES-based `re` has sparse coverage → fill remaining with 0 - if 're' in df_impute.columns: - df_impute = df_impute.with_columns(pl.col('re').fill_null(0)) - - print("Saving chars_raw_imputed.parquet ...", flush=True) - df_impute.write_parquet(OUTPUT_PATH + 'chars_raw_imputed.parquet') - - # ------------------------------------------------------------------ - # Output 3: ranked (no imputation) - # ------------------------------------------------------------------ - print("Ranking (no impute)...", flush=True) - df_rank = _rank_df(df) - print("Saving chars_rank_no_impute.parquet ...", flush=True) - df_rank.write_parquet(OUTPUT_PATH + 'chars_rank_no_impute.parquet') - del df_rank - - # ------------------------------------------------------------------ - # Output 4: ranked (imputed) — Bug 1 fix: rank the IMPUTED data - # ------------------------------------------------------------------ - print("Ranking (imputed)...", flush=True) - df_rank_imp = _rank_df(df_impute) - print("Saving chars_rank_imputed.parquet ...", flush=True) - df_rank_imp.write_parquet(OUTPUT_PATH + 'chars_rank_imputed.parquet') diff --git a/chars_ciz_monthly/merge_chars.py b/chars_ciz_monthly/merge_chars.py deleted file mode 100644 index 5f7e9d9..0000000 --- a/chars_ciz_monthly/merge_chars.py +++ /dev/null @@ -1,209 +0,0 @@ -""" -Merge accounting characteristics with satellite characteristics (rolling_chars, sue, abr, myre) -and CRSP return data. - -Inputs (all from OUTPUT_PATH = ../data/processed/): - - chars_a_accounting.parquet (from accounting.py) - - chars_q_accounting.parquet (from accounting.py) - - rolling_chars.parquet (from rolling_chars.py: beta, baspread, ill, maxret, rvar_capm, - rvar_ff3, rvar_mean, std_dolvol, std_turn, zerotrade) - - sue.parquet (from sue.py) - - myre.parquet (from myre.py) - - abr.parquet (from abr.py) - - crsp_msf.parquet (raw CRSP monthly, for delisting-adjusted returns / backfill) - -Outputs: - - chars_a_raw.parquet - - chars_q_raw.parquet -""" - -import polars as pl -from functions import INPUT_PATH, OUTPUT_PATH - -# ===================================================================== -# Satellite characteristic files to merge -# ===================================================================== -# Each entry: (filename, columns to keep besides permno/date, date_col) -_SATELLITE_FILES = [ - ('sue.parquet', ['sue'], 'date'), - ('myre.parquet', ['re'], 'date'), - ('abr.parquet', ['abr'], 'date'), -] - -# rolling_chars.parquet already has permno + date + many columns -_ROLLING_CHARS_FILE = 'rolling_chars.parquet' -_ROLLING_CHARS_COLS = [ - 'beta', 'baspread', 'ill', 'maxret', - 'rvar_capm', 'rvar_ff3', 'rvar_mean', - 'std_dolvol', 'std_turn', 'zerotrade', -] - -# CRSP CIZ -> accounting naming convention -_CRSP_SYNONYMS = { - 'mthcaldt': 'date', - 'mthprc': 'prc', - 'mthret': 'ret', - 'mthretx': 'retx', -} - - -def _load_satellite(filename, keep_cols, date_col): - """Load a satellite parquet, align date to month-end, deduplicate.""" - df = pl.read_parquet(OUTPUT_PATH + filename) - df = df.with_columns([ - pl.col('permno').cast(pl.Int64), - pl.col(date_col).cast(pl.Date).dt.month_end().alias('jdate'), - ]) - df = df.select(['permno', 'jdate'] + keep_cols) - df = df.unique(subset=['permno', 'jdate'], keep='last') - return df - - -def _load_rolling_chars(): - """Load rolling_chars parquet.""" - df = pl.read_parquet(OUTPUT_PATH + _ROLLING_CHARS_FILE) - df = df.with_columns([ - pl.col('permno').cast(pl.Int64), - pl.col('date').cast(pl.Date).dt.month_end().alias('jdate'), - ]) - available = [c for c in _ROLLING_CHARS_COLS if c in df.columns] - df = df.select(['permno', 'jdate'] + available) - df = df.unique(subset=['permno', 'jdate'], keep='last') - return df - - -def _build_crsp_backfill(): - """ - Build a CRSP backfill table for missing ret/retx/me. - CIZ (v2) mthret already includes delisting returns — no separate adjustment needed. - """ - crsp = pl.read_parquet(INPUT_PATH + 'crsp_msf.parquet') - crsp = crsp.rename(_CRSP_SYNONYMS, strict=False) - - # cast Decimal → Float64 - crsp = crsp.with_columns([ - pl.col(c).cast(pl.Float64) - for c in crsp.columns - if str(crsp[c].dtype).startswith('Decimal') - ]) - - crsp = crsp.filter( - pl.col('primaryexch').is_in(['N', 'A', 'Q']) & - (pl.col('conditionaltype') == 'RW') & - (pl.col('tradingstatusflg') == 'A') - ) - crsp = crsp.filter( - (pl.col('sharetype') == 'NS') & - (pl.col('securitytype') == 'EQTY') & - (pl.col('securitysubtype') == 'COM') & - (pl.col('usincflg') == 'Y') & - (pl.col('issuertype').is_in(['ACOR', 'CORP'])) - ) - - crsp = crsp.with_columns([ - pl.col('date').cast(pl.Date), - pl.col('permno').cast(pl.Int64), - pl.col('permco').cast(pl.Int64), - ]) - crsp = crsp.with_columns([ - pl.col('date').dt.month_end().alias('jdate'), - ]) - crsp = crsp.filter(pl.col('prc').is_not_null()) - crsp = crsp.with_columns([ - (pl.col('prc').abs() * pl.col('shrout')).alias('me'), - pl.col('ret').cast(pl.Float64).fill_null(0), - pl.col('retx').cast(pl.Float64).fill_null(0), - ]) - - # aggregate me: assign sum-of-permco-me to the permno with the largest me - crsp_summe = crsp.group_by(['jdate', 'permco']).agg(pl.col('me').sum()) - crsp_maxme = crsp.group_by(['jdate', 'permco']).agg(pl.col('me').max()) - crsp = crsp.join(crsp_maxme, on=['jdate', 'permco', 'me'], how='inner') - crsp = (crsp - .drop('me') - .join(crsp_summe, on=['jdate', 'permco'], how='inner') - .sort(['permno', 'jdate']) - .unique(subset=['permno', 'jdate'], keep='last') - ) - - crsp = crsp.select([ - 'permno', 'jdate', - pl.col('ret').alias('ret_fill'), - pl.col('retx').alias('retx_fill'), - pl.col('me').alias('me_fill'), - ]) - return crsp - - -def _merge_satellites(chars): - """Merge all satellite characteristics into chars.""" - # rolling chars - rolling = _load_rolling_chars() - chars = chars.join(rolling, on=['permno', 'jdate'], how='left') - - # other satellites - for filename, cols, date_col in _SATELLITE_FILES: - sat = _load_satellite(filename, cols, date_col) - chars = chars.join(sat, on=['permno', 'jdate'], how='left') - - return chars - - -def _backfill_crsp(chars, crsp_fill): - """Fill missing ret/retx/retadj/me from CRSP backfill table.""" - chars = chars.join(crsp_fill, on=['permno', 'jdate'], how='left') - for col_name in ['ret', 'retx', 'me']: - fill_col = f'{col_name}_fill' - if fill_col in chars.columns and col_name in chars.columns: - chars = chars.with_columns([ - pl.coalesce([pl.col(col_name), pl.col(fill_col)]).alias(col_name) - ]).drop(fill_col) - elif fill_col in chars.columns: - chars = chars.rename({fill_col: col_name}) - # drop rows without return - chars = chars.filter( - pl.col('ret').is_not_null() & - pl.col('retx').is_not_null() - ) - return chars - - -# ===================================================================== -# Main -# ===================================================================== -if __name__ == '__main__': - print("Loading accounting characteristics...", flush=True) - chars_a = pl.read_parquet(OUTPUT_PATH + 'chars_a_accounting.parquet') - chars_q = pl.read_parquet(OUTPUT_PATH + 'chars_q_accounting.parquet') - - # ensure types - for label, df in [('chars_a', chars_a), ('chars_q', chars_q)]: - df = df.with_columns([ - pl.col('permno').cast(pl.Int64), - pl.col('jdate').cast(pl.Date), - ]) - df = df.unique(subset=['permno', 'jdate'], keep='last') - if label == 'chars_a': - chars_a = df - else: - chars_q = df - - print("Loading satellite characteristics...", flush=True) - chars_a = _merge_satellites(chars_a) - chars_q = _merge_satellites(chars_q) - - print("Building CRSP backfill...", flush=True) - crsp_fill = _build_crsp_backfill() - - print("Backfilling CRSP data...", flush=True) - chars_a = _backfill_crsp(chars_a, crsp_fill) - chars_q = _backfill_crsp(chars_q, crsp_fill) - - # save - print(f"Saving chars_a_raw.parquet shape={chars_a.shape}", flush=True) - chars_a.write_parquet(OUTPUT_PATH + 'chars_a_raw.parquet') - - print(f"Saving chars_q_raw.parquet shape={chars_q.shape}", flush=True) - chars_q.write_parquet(OUTPUT_PATH + 'chars_q_raw.parquet') - - print("Done.", flush=True) diff --git a/chars_ciz_monthly/myre.py b/chars_ciz_monthly/myre.py deleted file mode 100644 index 57d8c51..0000000 --- a/chars_ciz_monthly/myre.py +++ /dev/null @@ -1,133 +0,0 @@ -# Calculate HSZ Replicating Anomalies -# RE: Revisions in analysts' earnings forecasts - -import polars as pl -from functions import INPUT_PATH, OUTPUT_PATH - -######################################################################### -# Merging IBES and CRSP by using ICLINK table. Merging last month price # -######################################################################### - -# Read ICLINK table from local file -iclink = pl.scan_csv(INPUT_PATH + 'iclink_ciz.csv') -# Convert all column names to lowercase first -iclink = iclink.select([pl.col(c).alias(c.lower()) for c in iclink.collect_schema().names()]) -# Then rename specific columns as needed -iclink = iclink.rename({'issuernm': 'comnam'}) - -# Read IBES data -ibes = pl.scan_parquet(INPUT_PATH + 'ibes.parquet') - -# Filtering IBES -ibes = ibes.filter( - pl.col('medest').is_not_null() & - pl.col('fpedats').is_not_null() -) - -# Add merge_date (end of month for statpers) -ibes = ibes.with_columns([ - pl.col('statpers').dt.month_end().alias('merge_date') -]) - -# Read CRSP monthly stock file from local -crsp_msf = pl.scan_parquet(INPUT_PATH + 'crsp_msf.parquet') -crsp_msf = crsp_msf.rename({'mthcaldt': 'date', 'mthprc': 'prc', 'mthcumfacpr': 'cfacpr'}) - -# Add merge_date (next month end) -crsp_msf = crsp_msf.with_columns([ - pl.col('date').dt.month_end().alias('date'), -]) -crsp_msf = crsp_msf.with_columns([ - (pl.col('date') + pl.duration(days=1)).dt.month_end().alias('merge_date') -]) - -# Merge IBES with ICLINK -ibes_iclink = ibes.join(iclink, on='ticker', how='left') - -# Merge with CRSP -ibes_crsp = ibes_iclink.join(crsp_msf, on=['permno', 'merge_date'], how='inner') -ibes_crsp = ibes_crsp.sort(['ticker', 'fpedats', 'statpers']) - -############################### -# Merging last month forecast # -############################### - -# Create last month columns using partitioned shift (no guard needed) -ibes_crsp = ibes_crsp.with_columns([ - pl.col('statpers').shift(1).over(['ticker', 'fpedats']).alias('statpers_last_month'), - pl.col('meanest').shift(1).over(['ticker', 'fpedats']).alias('meanest_last_month'), -]) - -# Re-sort -ibes_crsp = ibes_crsp.sort(['ticker', 'permno', 'fpedats', 'statpers']) - -########################### -# Drop empty "last month" # -# Calculate HXZ RE # -########################### - -ibes_crsp = ibes_crsp.filter(pl.col('statpers_last_month').is_not_null()) - -# Calculate adjusted price and monthly revision -ibes_crsp = ibes_crsp.with_columns([ - (pl.col('prc') / pl.col('cfacpr')).alias('prc_adj') -]) - -ibes_crsp = ibes_crsp.filter(pl.col('prc_adj') > 0) - -ibes_crsp = ibes_crsp.with_columns([ - ((pl.col('meanest') - pl.col('meanest_last_month')) / pl.col('prc_adj')).alias('monthly_revision') -]) - -# Create permno_fpedats identifier -ibes_crsp = ibes_crsp.with_columns([ - (pl.col('permno').cast(pl.Utf8) + '-' + pl.col('fpedats').cast(pl.Utf8)).alias('permno_fpedats') -]) - -# Drop duplicates and add count -ibes_crsp = ibes_crsp.unique(subset=['permno_fpedats', 'statpers'], keep='first') -ibes_crsp = ibes_crsp.with_columns([ - pl.col('permno_fpedats').cum_count().over('permno_fpedats').alias('count') -]) - -################## -# Calculate RE # -################## - -# Create lagged monthly_revision columns (partition by permno_fpedats to avoid mixing forecast periods) -ibes_crsp = ibes_crsp.with_columns([ - pl.col('monthly_revision').shift(1).over('permno_fpedats').alias('monthly_revision_l1'), - pl.col('monthly_revision').shift(2).over('permno_fpedats').alias('monthly_revision_l2'), - pl.col('monthly_revision').shift(3).over('permno_fpedats').alias('monthly_revision_l3'), - pl.col('monthly_revision').shift(4).over('permno_fpedats').alias('monthly_revision_l4'), - pl.col('monthly_revision').shift(5).over('permno_fpedats').alias('monthly_revision_l5'), - pl.col('monthly_revision').shift(6).over('permno_fpedats').alias('monthly_revision_l6') -]) - -# Calculate RE based on count -ibes_crsp = ibes_crsp.with_columns([ - pl.when(pl.col('count') == 4) - .then((pl.col('monthly_revision_l1') + pl.col('monthly_revision_l2') + pl.col('monthly_revision_l3')) / 3) - .when(pl.col('count') == 5) - .then((pl.col('monthly_revision_l1') + pl.col('monthly_revision_l2') + pl.col('monthly_revision_l3') + pl.col('monthly_revision_l4')) / 4) - .when(pl.col('count') == 6) - .then((pl.col('monthly_revision_l1') + pl.col('monthly_revision_l2') + pl.col('monthly_revision_l3') + pl.col('monthly_revision_l4') + pl.col('monthly_revision_l5')) / 5) - .when(pl.col('count') >= 7) - .then((pl.col('monthly_revision_l1') + pl.col('monthly_revision_l2') + pl.col('monthly_revision_l3') + pl.col('monthly_revision_l4') + pl.col('monthly_revision_l5') + pl.col('monthly_revision_l6')) / 6) - .otherwise(None) - .alias('re') -]) - -# Filter and finalize -ibes_crsp = ibes_crsp.filter(pl.col('count') >= 4) -ibes_crsp = ibes_crsp.sort(['ticker', 'statpers', 'fpedats']) -ibes_crsp = ibes_crsp.unique(subset=['ticker', 'statpers'], keep='first') - -# Select final columns and rename -ibes_crsp = ibes_crsp.select(['ticker', 'statpers', 'fpedats', 'anndats_act', 'curr_act', 'permno', 're']) -ibes_crsp = ibes_crsp.rename({'statpers': 'date'}) - -# Write output (collect LazyFrame to DataFrame first) -ibes_crsp.collect().write_parquet(OUTPUT_PATH + 'myre.parquet') - -print("RE calculation completed successfully!") \ No newline at end of file diff --git a/chars_ciz_monthly/rolling_chars.py b/chars_ciz_monthly/rolling_chars.py deleted file mode 100644 index bd2e451..0000000 --- a/chars_ciz_monthly/rolling_chars.py +++ /dev/null @@ -1,489 +0,0 @@ -""" -Rolling Window Characteristics using Polars -Rewritten from Pandas+multiprocessing to Polars for 10-100x speedup - -Characteristics included: -- beta: CAPM beta -- beta_ff5: Fama-French 5-factor market beta -- baspread: Bid-ask spread -- ill: Amihud (2002) illiquidity measure -- maxret: Maximum daily return -- rvar_capm: CAPM residual variance -- rvar_ff3: FF3 residual variance -- rvar_mean: Return variance -- std_dolvol: Std of log dollar volume -- std_turn: Std of turnover -- zerotrade: Zero trading days measure -""" - -import polars as pl -from polars import col -import polars_ols as pls -import time -from functions import INPUT_PATH, OUTPUT_PATH - -# ============================================================================= -# Utility Functions -# ============================================================================= - -def measure_time(func): - """Decorator to time function execution""" - def wrapper(*args, **kwargs): - start_time = time.time() - print(f"Function : {func.__name__.upper()}", flush=True) - print(f"Start : {time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(start_time))}", flush=True) - result = func(*args, **kwargs) - end_time = time.time() - print(f"End : {time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(end_time))}", flush=True) - total_seconds = end_time - start_time - minutes = int(total_seconds // 60) - seconds = total_seconds % 60 - print(f"Execution time : {minutes} minutes and {seconds:.2f} seconds", flush=True) - print() - return result - return wrapper - - -def gen_MMYY_column(date_col): - """Generate YYYYMM integer from date column for grouping""" - return (col(date_col).dt.year() * 100 + col(date_col).dt.month()).cast(pl.Int32) - - -def gen_consecutive_lists(input_list, k): - """Split a list into consecutive, non-overlapping sublists of length k""" - return [ - input_list[i : i + k] - for i in range(0, len(input_list), k) - if len(input_list[i : i + k]) == k - ] - - -def build_groups(input_list, k): - """Build k staggered groupings (offset windows) over a list""" - return [gen_consecutive_lists(input_list[offset:], k) for offset in range(k)] - - -def group_mapping_dfs(input_list, k): - """ - Create mapping DataFrames linking aux_date to group_number, - and group_number to new (max) aux_date - """ - groups = build_groups(input_list, k) - dfs = [ - pl.DataFrame({"aux_date": group}).with_columns( - group_number=pl.cum_count("aux_date"), - new_date=col("aux_date").list.max() - ) - for group in groups - ] - return [ - { - "group_map": df.explode("aux_date") - .select([col("aux_date").cast(pl.Int32), "group_number"]) - .lazy(), - "date_map": df.select([ - "group_number", - col("new_date").alias("aux_date") - ]) - .unique() - .sort(["group_number"]) - .lazy(), - } - for df in dfs - ] - - -# ============================================================================= -# Characteristic Calculation Functions -# ============================================================================= - -def capm_beta(df, min_obs=21): - """ - CAPM beta = cov(exret, mktrf) / var(mktrf) - Also computes variance of residuals from regression with intercept - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - (pl.cov("exret", "mktrf") / pl.var("mktrf")).alias("beta"), - # Residual variance from regression with intercept: - # residual = (exret - mean(exret)) - beta * (mktrf - mean(mktrf)) - ( - (col("exret") - pl.mean("exret")) - - (pl.cov("exret", "mktrf") / pl.var("mktrf")) * (col("mktrf") - pl.mean("mktrf")) - ).var().alias("rvar_capm"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def maxret(df, min_obs=21): - """ - Maximum daily return over the rolling window - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - pl.max("ret").alias("maxret"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def rvar_mean(df, min_obs=21): - """ - Variance of raw returns (no factor adjustment) - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - pl.var("ret").alias("rvar_mean"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def baspread(df, min_obs=21): - """ - Mean relative bid-ask spread: (askhi - bidlo) / ((askhi + bidlo) / 2) - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - ((col("askhi") - col("bidlo")) / ((col("askhi") + col("bidlo")) / 2)).mean().alias("baspread"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def std_dolvol(df, min_obs=21): - """ - Standard deviation of log dollar volume: std(log(|vol * prc|)) - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - ((col("vol") * col("prc").abs()).log()).std().alias("std_dolvol"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def std_turn(df, min_obs=21): - """ - Standard deviation of daily turnover: std(vol / shrout) - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - (col("vol") / col("shrout")).std().alias("std_turn"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def zerotrade(df, min_obs=21): - """ - Zero trading days measure: - zerotrade = (zero_count + (1/turnover_sum)/11000) * 63 / n_obs - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - (col("vol") == 0).sum().alias("zero_count"), - (col("vol") / col("shrout")).sum().alias("turn_sum"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .with_columns([ - ( - (col("zero_count") + (1.0 / col("turn_sum")) / 11000) - * 63.0 / col("n_obs") - ).alias("zerotrade") - ]) - .select(["permno", "group_number", "zerotrade"]) - ) - - -def ill(df, min_obs=21): - """ - Amihud (2002) illiquidity measure: - ill = mean(abs(ret) / (abs(prc) * vol)) - """ - return ( - df.group_by(["permno", "group_number"]) - .agg([ - (col("ret").abs() / (col("prc").abs() * col("vol"))).mean().alias("ill"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - - -def beta_ff5(df, min_obs=21): - """ - Fama-French 5-factor market beta using polars_ols - Returns the coefficient on mktrf from regression: - exret ~ mktrf + smb + hml + rmw + cma - """ - # Use polars_ols least_squares method with string column names - res_exp = pl.col("exret").least_squares.ols( - "mktrf", "smb", "hml", "rmw", "cma", - add_intercept=True, - mode="coefficients" - ) - result = ( - df.filter( - col("mktrf").is_not_null() & - col("smb").is_not_null() & - col("hml").is_not_null() & - col("rmw").is_not_null() & - col("cma").is_not_null() - ) - .group_by(["permno", "group_number"]) - .agg([ - res_exp.first().struct.field("mktrf").alias("beta_ff5"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - return result - - -def rvar_ff3(df, min_obs=21): - """ - Fama-French 3-factor residual variance using polars_ols - Computes var(residuals) from: exret ~ mktrf + smb + hml - """ - res_exp = pl.col("exret").least_squares.ols( - "mktrf", "smb", "hml", - add_intercept=True, - mode="residuals" - ) - result = ( - df.filter( - col("mktrf").is_not_null() & - col("smb").is_not_null() & - col("hml").is_not_null() - ) - .group_by(["permno", "group_number"]) - .agg([ - res_exp.var().alias("rvar_ff3"), - pl.len().alias("n_obs"), - ]) - .filter(col("n_obs") >= min_obs) - .drop("n_obs") - ) - return result - - -# ============================================================================= -# Processing Functions -# ============================================================================= - -def process_map_chunk(base_data, mapping, char_func, min_obs=21): - """ - Execute rolling computation for a mapping: - 1. Join base_data with mapping['group_map'] on aux_date - 2. Compute characteristic per (permno, group_number) - 3. Join mapping['date_map'] to remap to end date - """ - result = ( - base_data - .join(mapping["group_map"], on="aux_date", how="inner") - .pipe(char_func, min_obs=min_obs) - .join(mapping["date_map"], on="group_number", how="inner") - ) - return result - - -def compute_single_char(df, aux_maps, char_func, char_names, min_obs=21): - """ - Compute a single characteristic across all mapping chunks - - Parameters: - ----------- - char_names : str or list of str - Column name(s) to extract from the characteristic function - """ - results = [] - for mapping in aux_maps: - chunk_result = process_map_chunk(df, mapping, char_func, min_obs=min_obs) - results.append(chunk_result) - - # Handle both single column name (str) and multiple column names (list) - if isinstance(char_names, str): - char_names = [char_names] - - combined = pl.concat(results).select(["permno", "aux_date"] + char_names) - return combined - - -@measure_time -def compute_all_rolling_chars(input_path, output_path, n_months=3, min_obs=21): - """ - Main function to compute all rolling characteristics - - Parameters: - ----------- - input_path : str - Path to input parquet file with daily returns and factors - output_path : str - Path to output parquet file - n_months : int - Number of months in rolling window (default: 3) - min_obs : int - Minimum observations required per window (default: 21) - """ - - # 1. Load data - print(f"Loading data from {input_path}...", flush=True) - df = pl.scan_parquet(input_path) - - # 2. Prepare data - cast Decimal to Float64 for polars_ols compatibility - print("Preparing data...", flush=True) - df = ( - df - .rename({ - "dlycaldt": "date", - "dlyret": "ret", - "dlyvol": "vol", - "dlyprc": "prc", - "dlyhigh": "askhi", - "dlylow": "bidlo", - }) - .with_columns([ - col("date").cast(pl.Date), - col("permno").cast(pl.Int64), - (col("shrout") * 1000).cast(pl.Float64).alias("shrout"), - # Cast all Decimal columns to Float64 - col("ret").cast(pl.Float64), - col("vol").cast(pl.Float64), - col("prc").cast(pl.Float64), - col("askhi").cast(pl.Float64), - col("bidlo").cast(pl.Float64), - col("rf").cast(pl.Float64), - col("mktrf").cast(pl.Float64), - col("smb").cast(pl.Float64), - col("hml").cast(pl.Float64), - col("rmw").cast(pl.Float64), - col("cma").cast(pl.Float64), - ]) - .with_columns([ - (col("ret") - col("rf")).alias("exret"), - col("date").dt.month_end().alias("monthend"), - ]) - .with_columns([ - gen_MMYY_column("monthend").alias("aux_date"), - ]) - .filter( - col("ret").is_not_null() & - col("vol").is_not_null() - ) - .sort(["permno", "date"]) - ) - - # 3. Get unique months for window mapping - print("Generating rolling window mappings...", flush=True) - unique_months = ( - df.select("aux_date") - .unique() - .sort("aux_date") - .collect()["aux_date"] - .to_list() - ) - - aux_maps = group_mapping_dfs(unique_months, n_months) - - # 4. Compute each characteristic - print(f"Computing {n_months}-month rolling characteristics (min {min_obs} obs)...", flush=True) - - char_configs = [ - (capm_beta, ["beta", "rvar_capm"]), # Extract both columns from capm_beta - (maxret, "maxret"), - (rvar_mean, "rvar_mean"), - (baspread, "baspread"), - (std_dolvol, "std_dolvol"), - (std_turn, "std_turn"), - (zerotrade, "zerotrade"), - (ill, "ill"), - # (beta_ff5, "beta_ff5"), - (rvar_ff3, "rvar_ff3"), - ] - - all_results = None - for char_func, char_names in char_configs: - # Handle display name for logging - display_name = char_names if isinstance(char_names, str) else ", ".join(char_names) - print(f" Computing {display_name}...", flush=True) - char_result = compute_single_char(df, aux_maps, char_func, char_names, min_obs).collect() - - if all_results is None: - all_results = char_result - else: - all_results = all_results.join( - char_result, - on=["permno", "aux_date"], - how="full", - coalesce=True - ) - - # 5. Map back to actual dates - print("Mapping to final dates...", flush=True) - - date_map = ( - df.select(["aux_date", "date"]) - .group_by("aux_date") - .agg(col("date").max().alias("date")) - .collect() - ) - - output = ( - all_results - .join(date_map, on="aux_date", how="inner") - .drop("aux_date") - .sort(["permno", "date"]) - ) - - # 6. Write output - print(f"Writing output to {output_path}...", flush=True) - output.write_parquet(output_path) - - print(f"✓ Completed! Output shape: {output.shape}", flush=True) - print(f" Unique permnos: {output['permno'].n_unique()}", flush=True) - print(f" Date range: {output['date'].min()} to {output['date'].max()}", flush=True) - print(f" Columns: {output.columns}", flush=True) - - return output - - -if __name__ == "__main__": - # Configuration - INPUT_PATH = INPUT_PATH + "crsp_dsf.parquet" - OUTPUT_PATH = OUTPUT_PATH + "rolling_chars.parquet" - N_MONTHS = 3 - MIN_OBS = 21 - - result = compute_all_rolling_chars( - input_path=INPUT_PATH, - output_path=OUTPUT_PATH, - n_months=N_MONTHS, - min_obs=MIN_OBS - ) - - print("\nSample results:") - print(result.head(10)) \ No newline at end of file diff --git a/chars_ciz_monthly/sue.py b/chars_ciz_monthly/sue.py deleted file mode 100644 index fa8490b..0000000 --- a/chars_ciz_monthly/sue.py +++ /dev/null @@ -1,158 +0,0 @@ -# Calculate HSZ Replicating Anomalies -# SUE: Standardized Unexpected Earnings (Earnings surprise) - -import polars as pl -import duckdb -from pathlib import Path -from datetime import date -from functions import INPUT_PATH, OUTPUT_PATH - -################### -# Compustat Block # -################### -comp = ( - pl.scan_parquet(INPUT_PATH + "comp_fundq.parquet") - .select(["gvkey", "datadate", "fyearq", "fqtr", "epspxq", "ajexq"]) - .collect() -) - -################### -# CCM Block # -################### -ccm = ( - pl.scan_parquet(INPUT_PATH + "ccm.parquet") - .with_columns([ - pl.col("linkenddt").fill_null(date.today()) - ]) - .collect() -) - -# Merge comp with ccm -ccm1 = comp.join(ccm, on="gvkey", how="left") - -# Set link date bounds -ccm2 = ( - ccm1 - .filter( - (pl.col("datadate") >= pl.col("linkdt")) - & (pl.col("datadate") <= pl.col("linkenddt")) - ) - .select(["gvkey", "permno", "datadate", "fyearq", "fqtr", "epspxq", "ajexq"]) -) - -# Calculate EPS = epspxq / ajexq -ccm2 = ccm2.with_columns([ - (pl.col("epspxq") / pl.col("ajexq").replace(0, None)).alias("eps") -]) -ccm2 = ccm2.unique(subset=["permno", "datadate"]) - -# Filter out null eps and sort -ccm2 = ccm2.filter(pl.col("eps").is_not_null()) -ccm2 = ccm2.sort(["permno", "datadate"]) - -# Create count for each permno -ccm2 = ccm2.with_columns([ - pl.col("eps").cum_count().over("permno").alias("count") -]) - -# Create lag variables e1 to e12 -ccm2 = ccm2.with_columns([ - pl.col("eps").shift(1).over("permno").alias("e1"), - pl.col("eps").shift(2).over("permno").alias("e2"), - pl.col("eps").shift(3).over("permno").alias("e3"), - pl.col("eps").shift(4).over("permno").alias("e4"), - pl.col("eps").shift(5).over("permno").alias("e5"), - pl.col("eps").shift(6).over("permno").alias("e6"), - pl.col("eps").shift(7).over("permno").alias("e7"), - pl.col("eps").shift(8).over("permno").alias("e8"), - pl.col("eps").shift(9).over("permno").alias("e9"), - pl.col("eps").shift(10).over("permno").alias("e10"), - pl.col("eps").shift(11).over("permno").alias("e11"), - pl.col("eps").shift(12).over("permno").alias("e12"), -]) - -# Compute YoY differences: e{i}_diff = e{i} - e{i+4} -# Matches SIZ: std is computed on 4-quarter EPS differences, not raw levels -ccm2 = ccm2.with_columns([ - (pl.col(f"e{i}") - pl.col(f"e{i+4}")).alias(f"e{i}_diff") - for i in range(1, 9) -]) - -# Calculate sue_std based on count -# Using row-wise std on YoY differences, requires >= 6 non-null values -cols_6 = ["e6_diff", "e5_diff", "e4_diff", "e3_diff", "e2_diff", "e1_diff"] -cols_7 = ["e7_diff", "e6_diff", "e5_diff", "e4_diff", "e3_diff", "e2_diff", "e1_diff"] -cols_8 = ["e8_diff", "e7_diff", "e6_diff", "e5_diff", "e4_diff", "e3_diff", "e2_diff", "e1_diff"] - -def std_with_zero_handling(cols, min_valid=6): - """Row-wise std of YoY EPS differences; returns None if < min_valid non-null, 0 if all equal.""" - concat_expr = pl.concat_list(cols) - valid_count = concat_expr.list.eval(pl.element().is_not_null().cast(pl.Int32).sum()).list.first() - all_equal = concat_expr.list.max() == concat_expr.list.min() - std_val = concat_expr.list.eval(pl.element().std()).list.first() - return ( - pl.when(valid_count < min_valid).then(pl.lit(None, dtype=pl.Float64)) - .when(all_equal).then(pl.lit(0.0)) - .otherwise(std_val) - ) - -ccm2 = ccm2.with_columns([ - pl.when(pl.col("count") <= 6).then(None) - .when(pl.col("count") == 7).then(std_with_zero_handling(cols_6)) - .when(pl.col("count") == 8).then(std_with_zero_handling(cols_7)) - .otherwise(std_with_zero_handling(cols_8)) - .alias("sue_std") -]) - -# Calculate SUE -ccm2 = ccm2.with_columns([ - ((pl.col("eps") - pl.col("e4")) / pl.col("sue_std").replace(0, None)).alias("sue") -]) - -print(f"Calculated SUE: {ccm2.shape}") - -################### -# Monthly CRSP # -################### -crsp_msf = ( - pl.scan_parquet(INPUT_PATH + "crsp_msf.parquet") - .select(pl.col("mthcaldt").alias("date")) - .unique() - .collect() -) - -# Add plus12m column -ccm2 = ccm2.with_columns([ - (pl.col("datadate").dt.offset_by("12mo").dt.month_end()).alias("plus12m") -]) - -################### -# Populate to Monthly (using DuckDB for inequality join) -################### -con = duckdb.connect(":memory:") -con.register("ccm2", ccm2.to_arrow()) -con.register("crsp_msf", crsp_msf.to_arrow()) - -df = con.execute(""" - SELECT a.gvkey, a.permno, a.datadate, b.date, a.sue - FROM ccm2 a - LEFT JOIN crsp_msf b - ON a.datadate <= b.date - AND a.plus12m >= b.date - ORDER BY a.permno, b.date, a.datadate DESC -""").pl() - -con.close() - -# Keep first (most recent datadate) for each permno-date -df = df.unique(subset=["permno", "date"], keep="first") -df = df.filter(pl.col("date").is_not_null()) # Remove rows where no CRSP date matched -df = df.select(["gvkey", "permno", "datadate", "date", "sue"]) - -print(f"Final shape: {df.shape}") - -################### -# Save Output # -################### -df.write_parquet(OUTPUT_PATH + "sue.parquet") -print(f"Saved to {OUTPUT_PATH + 'sue.parquet'}") \ No newline at end of file