<jupyter_start><jupyter_text>Load Data<jupyter_code>!pip install shap
import pandas as pd
from patsy import dmatrix
import numpy as np
from sklearn import tree
from sklearn.ensemble import RandomForestRegressor
from sklearn.ensemble import RandomForestClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.linear_model import LinearRegression
from joblib import Parallel, delayed
import multiprocessing
#import shap<jupyter_output><empty_output><jupyter_text>Read SPXT: market indexDow Jones, S&P 500 and Nasdaq CompositeI only use spxt sheet, Vix sheet is ignored...<jupyter_code>spxt = pd.read_excel('/content/drive/My Drive/SWIB Quant Intern/SPXTVIX.xlsx',sheet_name='SPXT')<jupyter_output><empty_output><jupyter_text>One problem with SPXT, it is not strictly daily, some dates are missing.<jupyter_code>spxt[spxt.Dates >= '2001-08-31'].head(10)<jupyter_output><empty_output><jupyter_text>Calculate month return and daily return. we need daily return to calculate beta, and monthly return for calculate excess return.<jupyter_code>spxt['year_month'] = pd.to_datetime(spxt['Dates']).dt.to_period('M')
spxt['d_return'] = np.concatenate(   ([np.nan], [(spxt['SPXT'][t] - spxt['SPXT'][t-1])/spxt['SPXT'][t-1] for t in np.arange(spxt.shape[0]) if t != 0] )   )
spxt['m_return'] = [ ( spxt.loc[spxt['year_month'] == str(T), 'SPXT'].iloc[-1] -  spxt.loc[spxt['year_month'] == str(T), 'SPXT'].iloc[0] )/spxt.loc[spxt['year_month'] == str(T), 'SPXT'].iloc[0] for T in spxt['year_month']] 

spxt[(spxt.Dates >= '2000-01') & (spxt.Dates < '2000-02') ]
# -0.038310
(1925.41 - 2002.11) / 2002.11
# -0.058829
(6145.30 - 6529.42) / 6529.42

# -0.041092
(1919.84 - 2002.11)/2002.11<jupyter_output><empty_output><jupyter_text>Read individual stock index<jupyter_code>def load_basic_data(fileName):

        tempDF = pd.read_parquet(fileName)
        
        tempDF['date'] = pd.to_datetime(tempDF['date'])
                  
        return tempDF

tempDF = load_basic_data("/content/drive/My Drive/SWIB Quant Intern/SP500.pq")<jupyter_output><empty_output><jupyter_text>Build Signals For each stock_i where $i \in \{ 1,...1069\}$. we convert daily data to monthly data in order to remove noise. We group daily data by year_month. Suppose month T is the current month.stock_i_Month_T return = (Price_T.iloc[-1] - Price_T.iloc[0]) / Price_T.iloc[0]stock_i_month_T_excess_return = stock_i_month_T_return - market_month_T_returnNotes from Roger:1.   Build those signals based on price, return, volume, volatility, momentum2.   Use monthly data to remove noise. (first email)3.   the forecast is for **excess return**, violatility need to be calculated in **return space** not on price. performance to be measured by out of sample R24.  Y: next month excess return(T+1).  ?? Should use monthT  excess return??  Mom1m:  last month excess return (T-1)  Momentum: previous T-12:T-2 month excess return.  mom36: last 36 month excess return (T-36 month excess return?).    maxret1 : max monthly excess return in the last 12 month  maxret2 : max daily excess return in the last 1 month  turnover: dollar volume / market cap  beta: roll 252 days calculation against index. (return_i = beta * return_m + intercept + epsilon). Beta time-dependent: Beta_i_T I used the lastest available 252 days before last day of monthT.    beta^2  retvol: stdev(daily return)*sqrt(252). ? *sqrt(252)  idiovol: stdev(daily regression residual)*sqrt(252)<jupyter_code># Map 2001-09 - 2020-06 to T 0 - 225
T = pd.to_datetime(tempDF['date']).dt.to_period('M').unique()
T = dict(zip (T,np.arange(len(T))))
# for all stock_id
# stock id, total of 1069 stocks 

stock_id = list(tempDF.sedol.unique()) 

stock_month = pd.DataFrame({'year_month':[],'sedol':[],'sector':[],'m_volume':[],'spxt_m_return':[],
 'm_return':[],'excess_return':[],'mom1m':[],'next_month_excess_return':[],'Momentum1':[],'Momentum2':[],
 'Momentum3':[],'Momentum4':[],'Momentum5':[],'Momentum6':[],'Momentum7':[],'Momentum8':[],'Momentum9':[],'Momentum10':[],'Momentum11':[],
 'Momentum12':[],'Momentum36':[],'maxret1':[],'maxret2':[],'dollarVolume':[],'mktcap':[],'turnover':[],'beta':[],'beta_sq':[],'retvol':[],
 'idiovol':[] })


# for stock i, 0:400. 400:800 , 800:1049
for i in np.arange(400):
    print(i)
    # stock_i is daily, which is used for constructing stock_month
    stock_i = tempDF[tempDF['sedol'] == stock_id[i]]; stock_i = stock_i.sort_values(by = ['date']) # make sure the date is correct
    stock_i_sector = stock_i['sector'].unique() # sometimes na + sector
    stock_i['year_month'] = pd.to_datetime(stock_i['date']).dt.to_period('M'); 
    stock_i['date'] = pd.to_datetime(stock_i['date']).dt.to_period('D') # period is better, can + -.
    time_i = list(stock_i['year_month'].unique())

    # SP index in time_i
    #start_m = [str(T) for T in time_i][0] ;end_m = [str(T) for T in time_i][-1] # string. 
    #spxt_i = spxt.loc[start_m:end_m] # 

    # intitial new df, save stock_i monthly data
    stock_i_month = pd.DataFrame(data = {"year_month":time_i})

    # stock id and sector
    stock_i_month['sedol'] = stock_id[i]; stock_i_month['sector'] = stock_i_sector[-1]

    # Since it is VOLUME_AVG_30D, I use the .iloc[-1] of that month to represent that month's volume.
    stock_i_month['m_volume'] = [  stock_i[stock_i['year_month'] == T]['VOLUME_AVG_30D'].iloc[-1]   for T in time_i ]

    # spxt monthly return
    stock_i_month['spxt_m_return'] = [ spxt[spxt['year_month'] == T].m_return.iloc[0]   for T in time_i] # any iloc[.]

    # stock_i monthly return, (price.T.iloc[-1] - price.T.iloc[0] ) / price.T.iloc[0]
    stock_i_month['m_return'] = [ (stock_i[stock_i['year_month'] == T].LAST_PRICE.iloc[-1] - stock_i[stock_i['year_month'] == str(T)].LAST_PRICE.iloc[0])/ 
                                stock_i[stock_i['year_month'] == T].LAST_PRICE.iloc[0]  for T in time_i ]


    # stock_ i excess return
    stock_i_month['excess_return'] = stock_i_month['m_return'] - stock_i_month['spxt_m_return']

    # Mom1m: last month excess return (T-1)
    stock_i_month['mom1m'] = np.concatenate(( [np.nan],stock_i_month['excess_return'].to_list()[0:-1]))

    # next month excess return (T+1), Y 
    stock_i_month['next_month_excess_return'] = np.concatenate( ( stock_i_month['excess_return'].to_list()[1:] , [np.nan]) )

    # Momentum, previous T-12 :  T-2 month excess return. 
    # mom36: last 36 month excess return.
    def empty_series(ele):
        if ele.empty:
          return np.nan
        else:
          return ele.iloc[0]
    for j in np.concatenate((list(np.arange(1,13)),[36])):
        stock_i_month['Momentum'+ str(j)] =  [ empty_series(ele)  for ele in [ stock_i_month[stock_i_month['year_month'] == T -j].excess_return  for T in time_i]  ]

    # maxret1, max monthly excess return in the last 12 month, from T-12 : T
    stock_i_month['maxret1'] = [  stock_i_month[ (stock_i_month['year_month'] >= T -12 ) & (stock_i_month['year_month'] <= T  )].excess_return.max()  for T in time_i]
    # maxret2, max daily excess return in the last 1 month
    stock_i_month['maxret2'] = [ stock_i[stock_i['year_month'] == T].d_return.max()  for T in time_i]

    # turnover: dollarvolume/ marketcap
    # Dollarvolume,  dollarVolume_T.iloc[-1].  which equals to price_T.iloc[-1]* volume_T.iloc[-1]
    # marketcap mktcap_T.mean()
    stock_i_month['dollarVolume'] = [stock_i[stock_i['year_month'] == T].dollarVolume.iloc[-1]  for T in time_i]
    stock_i_month['mktcap'] = [ stock_i[stock_i['year_month'] == T].mktcap.mean()    for T in time_i]
    stock_i_month['turnover'] = stock_i_month['dollarVolume']/stock_i_month['mktcap']
    
    # beta, roll 252 days , r_i ~ beta r_m.
    beta = list()
    retvol = list()
    idiovol = list()
    for T in time_i:
      last_date = stock_i[stock_i['year_month'] == T].date.iloc[-1]
      start_date = last_date - 252
      # daily data for stock_i, for the most recent 252 days, however, maybe not available up to 252 days.
      stock_i_T = stock_i.loc[ (stock_i['date'] >= start_date) &  (stock_i['date'] <= last_date), ['date','sedol','d_return']]
      # def empty_series, if at date_t, spxt[ spxt['Dates'] == str(t) ].d_return maybe empty series, then iloc will get error
      stock_i_T['spxt_d_return'] = [empty_series(ele)  for ele in [ spxt[ spxt['Dates'] == str(t) ].d_return  for t in stock_i_T['date']  ] ]
      # retvol: std(d_return), for the most recent 252 days, however, maybe not available up to 252 days.
      retvol.append(stock_i_T['d_return'].std())
      # run linear regression r_i ~ r_m
      reg = LinearRegression().fit(stock_i_T[['d_return']],stock_i_T[['spxt_d_return']] )
      beta.append(reg.coef_[0][0])
      # residual
      predicted = reg.predict(stock_i_T[['spxt_d_return']]).reshape((-1,))
      residual = stock_i_T.d_return.to_numpy() - predicted
      idiovol.append(np.std(residual))

    stock_i_month['beta'] = beta
    stock_i_month['beta_sq'] = [ele**2 for ele in beta]
    stock_i_month['retvol'] = retvol
    stock_i_month['idiovol'] = idiovol

    # merge stock_i to previous DF
    stock_month = pd.concat([stock_month,stock_i_month])


# Save stock_month as .pkl
# It runs slowly So I partitioned the job and merge them
# stock_month.to_pickle('/content/drive/My Drive/SWIB Quant Intern/stock_month.pkl')


# Stock_month final data already saved, so directly load

stock_month = pd.read_pickle('/content/drive/My Drive/SWIB Quant Intern/stock_month.pkl')

# map T: 2001-06: 0,    2020-06 :226
stock_month['Month_T'] = stock_month['year_month'].map(T)
stock_month
stock_month.groupby('sedol').beta.agg(['mean','std','size']).sort_values(by = 'mean')<jupyter_output><empty_output><jupyter_text>ParalellI Tried to parallel the for loop, However, I didn't see improved performance.<jupyter_code># for all stock_id
# stock id, total of 1069 stocks 

stock_id = list(tempDF.sedol.unique()) 

def paralell_func(i,tempDF = tempDF,spxt = spxt,stock_id = stock_id):#,tempDF,spxt,stock_id):
    print(i)
    # stock_i is daily, which is used for constructing stock_month
    stock_i = tempDF[tempDF['sedol'] == stock_id[i]]; stock_i = stock_i.sort_values(by = ['date']) # make sure the date is correct
    stock_i_sector = stock_i['sector'].unique() # sometimes na + sector
    stock_i['year_month'] = pd.to_datetime(stock_i['date']).dt.to_period('M'); 
    stock_i['date'] = pd.to_datetime(stock_i['date']).dt.to_period('D') # period is better, can + -.
    time_i = list(stock_i['year_month'].unique())

    # SP index in time_i
    #start_m = [str(T) for T in time_i][0] ;end_m = [str(T) for T in time_i][-1] # string. 
    #spxt_i = spxt.loc[start_m:end_m] # 

    # intitial new df, save stock_i monthly data
    stock_i_month = pd.DataFrame(data = {"year_month":time_i})

    # stock id and sector
    stock_i_month['sedol'] = stock_id[i]; stock_i_month['sector'] = stock_i_sector[-1]

    # Since it is VOLUME_AVG_30D, I use the .iloc[-1] of that month to represent that month's volume.
    stock_i_month['m_volume'] = [  stock_i[stock_i['year_month'] == T]['VOLUME_AVG_30D'].iloc[-1]   for T in time_i ]

    # spxt monthly return
    stock_i_month['spxt_m_return'] = [ spxt[spxt['year_month'] == T].m_return.iloc[0]   for T in time_i] # any iloc[.]

    # stock_i monthly return, (price.T.iloc[-1] - price.T.iloc[0] ) / price.T.iloc[0]
    stock_i_month['m_return'] = [ (stock_i[stock_i['year_month'] == T].LAST_PRICE.iloc[-1] - stock_i[stock_i['year_month'] == str(T)].LAST_PRICE.iloc[0])/ 
                                stock_i[stock_i['year_month'] == T].LAST_PRICE.iloc[0]  for T in time_i ]


    # stock_ i excess return
    stock_i_month['excess_return'] = stock_i_month['m_return'] - stock_i_month['spxt_m_return']

    # Mom1m: last month excess return (T-1)
    stock_i_month['mom1m'] = np.concatenate(( [np.nan],stock_i_month['excess_return'].to_list()[0:-1]))

    # next month excess return (T+1), Y 
    stock_i_month['next_month_excess_return'] = np.concatenate( ( stock_i_month['excess_return'].to_list()[1:] , [np.nan]) )

    # Momentum, previous T-12 :  T-2 month excess return. 
    # mom36: last 36 month excess return.
    def empty_series(ele):
        if ele.empty:
          return np.nan
        else:
          return ele.iloc[0]
    for j in np.concatenate((list(np.arange(1,13)),[36])):
        stock_i_month['Momentum'+ str(j)] =  [ empty_series(ele)  for ele in [ stock_i_month[stock_i_month['year_month'] == T -j].excess_return  for T in time_i]  ]

    # maxret1, max monthly excess return in the last 12 month, from T-12 : T
    stock_i_month['maxret1'] = [  stock_i_month[ (stock_i_month['year_month'] >= T -12 ) & (stock_i_month['year_month'] <= T  )].excess_return.max()  for T in time_i]
    # maxret2, max daily excess return in the last 1 month
    stock_i_month['maxret2'] = [ stock_i[stock_i['year_month'] == T].d_return.max()  for T in time_i]

    # turnover: dollarvolume/ marketcap
    # Dollarvolume,  dollarVolume_T.iloc[-1].  which equals to price_T.iloc[-1]* volume_T.iloc[-1]
    # marketcap mktcap_T.mean()
    stock_i_month['dollarVolume'] = [stock_i[stock_i['year_month'] == T].dollarVolume.iloc[-1]  for T in time_i]
    stock_i_month['mktcap'] = [ stock_i[stock_i['year_month'] == T].mktcap.mean()    for T in time_i]
    stock_i_month['turnover'] = stock_i_month['dollarVolume']/stock_i_month['mktcap']
    
    # beta, roll 252 days , r_i ~ beta r_m.
    beta = list()
    retvol = list()
    idiovol = list()
    for T in time_i:
      last_date = stock_i[stock_i['year_month'] == T].date.iloc[-1]
      start_date = last_date - 252
      # daily data for stock_i, for the most recent 252 days, however, maybe not available up to 252 days.
      stock_i_T = stock_i.loc[ (stock_i['date'] >= start_date) &  (stock_i['date'] <= last_date), ['date','sedol','d_return']]
      # def empty_series, if at date_t, spxt[ spxt['Dates'] == str(t) ].d_return maybe empty series, then iloc will get error
      stock_i_T['spxt_d_return'] = [empty_series(ele)  for ele in [ spxt[ spxt['Dates'] == str(t) ].d_return  for t in stock_i_T['date']  ] ]
      # retvol: std(d_return), for the most recent 252 days, however, maybe not available up to 252 days.
      retvol.append(stock_i_T['d_return'].std())
      # run linear regression r_i ~ r_m
      reg = LinearRegression().fit(stock_i_T[['d_return']],stock_i_T[['spxt_d_return']] )
      beta.append(reg.coef_[0][0])
      # residual
      predicted = reg.predict(stock_i_T[['spxt_d_return']]).reshape((-1,))
      residual = stock_i_T.d_return.to_numpy() - predicted
      idiovol.append(np.std(residual))

    stock_i_month['beta'] = beta
    stock_i_month['beta_sq'] = [ele**2 for ele in beta]
    stock_i_month['retvol'] = retvol
    stock_i_month['idiovol'] = idiovol

    # merge stock_i to previous DF
    return stock_i_month

num_cores = multiprocessing.cpu_count()
     
results = Parallel(n_jobs=num_cores)(delayed(paralell_func)(i) for i in np.arange(len(stock_id[800:1068])) )<jupyter_output><empty_output><jupyter_text>Build Prediction Model let $T \in \{1,2,....,226 \}$ denote the month index, let $ i \in \{1,2,....,1069 \}$ denote the stock index. The response variable is Next_month_excess_return, denoted as $r_{i,T+1}$, of stock i and month T+1. The predictive characteristics are denoted as P-dimensional vector $z_{i,T}$, T is 1,2....226 following "citi-bank" paper, which measures the overal time trend for the macromarket.1.  we consider regression problem, which means that the objective function is minimization of MSE.    $r_{i,T+1} = g(z_{i,T}) + \epsilon_{i,T+1}$          2.  we consider classification problem, for example a binary classification of whether d_return is positive or not.       $P(r_{i,T+1} > 0 | z_{i,T}) = g(z_{i,T})$Following the paper(ML), The g function depends neither on i or T. By maintaining the same form over time and across different stocks, the model leverages information from the entire panel and every individual asset. Also g(.) depends on z only through $z_{i,T}$.We split sample into 2 disjoint time periods, 2001-2013 as training, 2014-2020 as testing.Regarding performance evaluation for regression, we use R square $R^2_{oos} = 1 - \frac{\sum_{(i,T)\in test}(r_{i,T+1}-\hat{r}_{i,T+1})^2}{\sum_{(i,T)\in test} r^2_{i,T+1}}$, and F1,accuracy for classification. How to handle missing data? I remove all rows contains NA, data size drop from 113594 to 78857, only 817 stocks left(1069).<jupyter_code>stock_month_dropna = stock_month.dropna()
stock_month_dropna
len(stock_month_dropna.sedol.unique())
T
# split data, remember after remove NAs, T starts from 36 to 225. 
# training 36: 180,   testing 180:225

features = [ 'Month_T', 'mom1m', 'Momentum2', 'Momentum3', 'Momentum4', 'Momentum5',
       'Momentum6', 'Momentum7', 'Momentum8', 'Momentum9', 'Momentum10',
       'Momentum11', 'Momentum12', 'Momentum36', 'maxret1', 'maxret2',
        'turnover', 'beta', 'beta_sq', 'retvol',
       'idiovol']


X_train = stock_month_dropna.loc[ stock_month_dropna.Month_T.isin(np.arange(0,180)) , features] ## 60477 all 78857
Y_train = stock_month_dropna.loc[ stock_month_dropna.Month_T.isin(np.arange(0,180)) , 'next_month_excess_return']

X_test = stock_month_dropna.loc[ stock_month_dropna.Month_T.isin(np.arange(180,226)) , features] # 18380, match!
Y_test = stock_month_dropna.loc[ stock_month_dropna.Month_T.isin(np.arange(180,226)) , 'next_month_excess_return']

# build model
regr = RandomForestRegressor(max_depth=3,random_state=0)
regr.fit(X_train, Y_train)
regr.score(X_train,Y_train)
# Classification model
# Y to binary 
def pos_ind(a):
  return 1 if a > 0 else 0 
Y_train_bi = [pos_ind(ele) for ele in Y_train]
Y_test_bi = [pos_ind(ele) for ele in Y_test]

clf = RandomForestClassifier(max_depth= 3,random_state=0)
clf.fit(X_train, Y_train_bi)
#clf.score(X_test,Y_test)
print(clf.score(X_test,Y_test_bi))
clf.score(X_train,Y_train_bi)<jupyter_output>0.5145266594124048
{'fraction_non_alphanumeric': 0.1491971454058876, 'fraction_numerical': 0.030776092774308653, 'mean_word_length': 4.3816381638163815, 'pattern_counts': {'":': 2, '<': 37, '<?xml version=': 0, '>': 40, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '570862', 'license': 'MIT', 'n_tokens_mistral': 7031, 'n_tokens_neox': 6617, 'n_words': 1526, 'path': 'week2_build_signals.ipynb', 'repo_name': 'chengning-zhang/Intern-SWIB', 'size': '17936'}