<jupyter_start><jupyter_text>Predicting Stock Clossing Price Using News Headlines as Predictors This notebook explores the effect (if any) that news headlines have on forecasting stock market closing prices. I'll be using data from [this Kaggle competition](https://www.kaggle.com/aaron7sun/stocknews) and Apple stock data.  Import and Process News HeadlinesThe data was already downloaded to the `data` directory, so here we import it and show a few lines.<jupyter_code>import warnings
warnings.filterwarnings('ignore')
import pandas as pd
news_df = pd.read_csv('../data/Combined_News_DJIA.csv').drop(labels='Label', axis=1)
news_df.head(3)<jupyter_output><empty_output><jupyter_text>Remove leading `b`Some headlines were stored as bytes (preceded by a `b` but needed to be converted to strings. I tried [this solution](https://stackoverflow.com/questions/41918836/how-do-i-get-rid-of-the-b-prefix-in-a-string-in-python) but it didn't work, so I just created a function the removed the unnecessary characters.<jupyter_code>def byte_to_str(value):
    if isinstance(value, str) and (value.startswith("b'") or value.startswith('b"')):
        value = value[2:-1]
    return value

for column in news_df.iloc[:, 1:].columns:
    news_df[column] = news_df[column].apply(byte_to_str)
    
news_df.head(3)<jupyter_output><empty_output><jupyter_text>Create Latent Space RepresentationIn order to train a model that takes news headlines as predictors/features, we first have to convert the headlines to a numerical representation that contains syntactic and semantic information about the headlines. For that, we use [BERT](https://arxiv.org/abs/1810.04805), which is a transformer model pre-trained with all of Wikipedia and a large dataset of unpublished books. BERT produces a [latent state representation](https://towardsdatascience.com/understanding-latent-space-in-machine-learning-de5a7c687d8d) for each headline in the form of matrix where each row represents a word in a healine. The matrix is then transformed to a vector by calculating the mean of each column. Then, all the vector representations for all headlines of a given day are used to create a new matrix, which represents all articles for that day. This matrix is also transformed to a vector by calculating the mean of each column. The collection of vectors is then used as features to train the final model.<jupyter_code>import nltk
nltk.download('punkt')
from nltk.tokenize import word_tokenize
import torch
from transformers import BertConfig, BertModel, BertTokenizer
import numpy as np

model_name = "bert-base-cased"
device = "cuda:0" if torch.cuda.is_available() else "cpu"

config = BertConfig.from_pretrained(model_name)
bert_model = BertModel.from_pretrained(model_name, config=config)
bert_model.eval()
bert_tokenizer = BertTokenizer.from_pretrained(model_name, do_lower_case=True)
articles_matrix = np.empty((len(news_df.index), 768))

with torch.no_grad():
    for row in news_df.iterrows():
        day_index = row[0]
        headlines_series = row[1][1:]
        nan_count = headlines_series.isna().sum()
        day_article_matrix = np.empty((headlines_series.size - nan_count, 768))
        for index, headline in enumerate(headlines_series):
            if isinstance(headline, str):
                tokens = bert_tokenizer(headline, return_tensors='pt')
                bert_model.to(device)
                output = bert_model(tokens.input_ids)
                # outputs matrix of shape (num_tokens, 768)
                latent_matrix = output.last_hidden_state.to('cpu').detach().numpy()[0]
                mean_vector = np.mean(latent_matrix, axis=0).reshape((1, latent_matrix.shape[1]))
                day_article_matrix[index, :] = mean_vector
                
        mean_vector = np.mean(day_article_matrix, axis=0).reshape((1, day_article_matrix.shape[1]))
        articles_matrix[day_index, :] = mean_vector
print(articles_matrix.shape)
news_df['Date'] = pd.to_datetime(news_df['Date'], format='%Y-%m-%d')
news_numerical_df = pd.DataFrame(articles_matrix, index=news_df['Date'])
news_numerical_df<jupyter_output><empty_output><jupyter_text>Stock DataThe Apple stock data is downloaded using the [Yahoo Finance Package](https://pypi.org/project/yfinance/)<jupyter_code>import yfinance as yf
import datetime

tickers = ['AAPL']
start = datetime.datetime(2008,8,8)
end = datetime.datetime(2016,7,1)
stock_df = yf.download(tickers, start=start,end = end, interval='1d')
stock_df<jupyter_output>[*********************100%***********************]  1 of 1 completed<jupyter_text>Data DescriptionWe will use the daily stock data for Apple, encompasing the dates of 08/08/2008 - 07/01/2016. We can see below that the data is not stationary, but the deep learning architecture that we will be using to model the data (LSTM) is able to "learn non-linear and non-stationary nature of time series" [4](References), so we don't have to make the data stationary.<jupyter_code>import matplotlib.pyplot as plt

plt.figure(figsize=(10,6))
plt.title('Apple Clossing Values')
plt.ylabel('$')
print(stock_df['Close'].describe())
stock_df['Close'].plot()
plt.savefig('../plots/AAPL_20080808-20160701.png', dpi=300)
plt.show()<jupyter_output>count    1988.000000
mean       16.720844
std         8.410058
min         2.792857
25%         9.318839
50%        16.306428
75%        23.685357
max        33.250000
Name: Close, dtype: float64<jupyter_text>Data Pre-Processing  Scaling[Scaling the input data can reduce the complexity of the problem being solved by the deep learning algorithm](https://machinelearningmastery.com/how-to-improve-neural-network-stability-and-modeling-performance-with-data-scaling/), so it is a recommended step to try to reduce the scale of the predictors in order to improve the performance of the model. In this case, since the data from the BERT model is in a scale from -1 to 1, then we will scale the closing value of the stock to match that. And then we'll create a new dataset that includes the output from BERT model and the stock closing data.<jupyter_code>from sklearn.preprocessing import MinMaxScaler

scaler = MinMaxScaler(feature_range=(-1, 1))
scaler.fit(stock_df['Close'].dropna().values.reshape(-1, 1))
stock_df['Close_Normalized'] = scaler.transform(stock_df['Close'].values.reshape(-1, 1))

data_df = news_numerical_df.copy()
data_df = data_df.join(stock_df[['Close_Normalized', 'Close']])
data_df = data_df.dropna()
data_df<jupyter_output><empty_output><jupyter_text>Shifting Closing ValueSince we're training the model to predict the closing value of the Apple stock for any given day, taking the previous day's news headlines and the previous day's closing value, then we need to create a label column with the closing value from `T + 1`.We also remove the last row since the last series of news and closing value doesn't have a value for the following day.<jupyter_code>Y = np.empty((len(data_df.index),))
Y[0:-1] = data_df['Close'].values[1:]
data_df['Label'] = Y
data_df = data_df[:-1]
data_df<jupyter_output><empty_output><jupyter_text>Data SplitThe data is split in ~ 70% to train the model, ~ 15% to validate the model, and ~ 15% to test the model.<jupyter_code>num_train_items = round(len(data_df.index) * .7)
num_val_items = (len(data_df.index) - num_train_items) // 2
num_test_items = len(data_df.index) - num_train_items - num_val_items
print(num_train_items)
print(num_val_items)
print(num_test_items)
X_train = data_df.iloc[:num_train_items,:769].values
y_train = data_df.iloc[:num_train_items,-1].values
print(X_train.shape)
print(y_train.shape)
X_val = data_df.iloc[num_train_items:(num_train_items+num_val_items),:769].values
y_val = data_df.iloc[num_train_items:(num_train_items+num_val_items),-1].values
print(X_val.shape)
print(y_val.shape)
X_test = data_df.iloc[(num_train_items+num_val_items):,:769].values
y_test = data_df.iloc[(num_train_items+num_val_items):,-1].values
print(X_test.shape)
print(y_test.shape)
X_train = X_train.reshape((X_train.shape[0], 1, X_train.shape[1]))
X_val = X_val.reshape((X_val.shape[0], 1, X_val.shape[1]))
X_test = X_test.reshape((X_test.shape[0], 1, X_test.shape[1]))
print(X_train.shape)
print(X_val.shape)
print(X_test.shape)<jupyter_output>(1391, 1, 769)
(298, 1, 769)
(298, 1, 769)<jupyter_text>Model TrainingWe loosely follow [this tutorial](https://machinelearningmastery.com/multivariate-time-series-forecasting-lstms-keras/) to train the Long Short Term Memory (LSTM) model. After trying a few different values for the output dimmesion of the LSTM layer, 1700 was the one the provided the best accuracy.<jupyter_code>from tensorflow.keras.models import Sequential
from tensorflow.keras.layers import Dense
from tensorflow.keras.layers import LSTM

model = Sequential()
model.add(LSTM(1700, input_shape=(X_train.shape[1], X_train.shape[2])))
model.add(Dense(1))
model.compile(loss='mae', optimizer='adam')

history = model.fit(X_train, y_train, epochs=100, batch_size=72, validation_data=(X_val, y_val), verbose=2, shuffle=False)<jupyter_output>Epoch 1/100
20/20 - 3s - loss: 2.4849 - val_loss: 4.6264
Epoch 2/100
20/20 - 2s - loss: 8.8520 - val_loss: 14.3589
Epoch 3/100
20/20 - 2s - loss: 5.2638 - val_loss: 10.2710
Epoch 4/100
20/20 - 2s - loss: 6.4141 - val_loss: 10.9076
Epoch 5/100
20/20 - 2s - loss: 6.0652 - val_loss: 11.8424
Epoch 6/100
20/20 - 2s - loss: 5.3856 - val_loss: 10.6545
Epoch 7/100
20/20 - 2s - loss: 5.7609 - val_loss: 11.0129
Epoch 8/100
20/20 - 2s - loss: 5.4284 - val_loss: 10.8262
Epoch 9/100
20/20 - 2s - loss: 5.3332 - val_loss: 10.5421
Epoch 10/100
20/20 - 2s - loss: 5.2558 - val_loss: 10.2976
Epoch 11/100
20/20 - 2s - loss: 5.1756 - val_loss: 10.0783
Epoch 12/100
20/20 - 2s - loss: 5.0499 - val_loss: 9.7025
Epoch 13/100
20/20 - 2s - loss: 5.0142 - val_loss: 9.5289
Epoch 14/100
20/20 - 2s - loss: 4.8350 - val_loss: 9.0595
Epoch 15/100
20/20 - 2s - loss: 4.8029 - val_loss: 8.8277
Epoch 16/100
20/20 - 2s - loss: 4.5727 - val_loss: 8.1624
Epoch 17/100
20/20 - 2s - loss: 4.6399 - val_loss: 8.2102
Epoch 18/100
[...]<jupyter_text>The model converges at around 80 epochs.<jupyter_code>plt.figure(figsize=(10,6))
plt.plot(history.history['loss'], label='train')
plt.plot(history.history['val_loss'], label='test')
plt.title('Model Training Error per Epoch')
plt.ylabel('Error')
plt.xlabel('Epoch')
plt.legend()
plt.savefig('../plots/model_training_error_epochs.png', dpi=300)
plt.show()<jupyter_output><empty_output><jupyter_text>Model AccuracyWe can see below that the model is highly accurate on the test data with a root mean square error of 0.9. So this means that the the majority of the model's predicted values are within around -2 and +2 dollars of the actual closing values.<jupyter_code>y_hat = model.predict(X_test)
print(y_hat.shape)
from sklearn.metrics import mean_squared_error
from math import sqrt

rmse = sqrt(mean_squared_error(y_test, y_hat))
rmse

pd.DataFrame(
    {
        'Predicted Values': y_hat.reshape((y_hat.shape[0],)),
        'Actual Values': y_test.reshape((y_test.shape[0],))
    },
    index = data_df.iloc[(num_train_items+num_val_items):,:].index
).plot(figsize=(10,6), title=f'Model Accuracy {rmse} RMSE')
plt.savefig('../plots/model_accuracy.png', dpi=300)
plt.show()<jupyter_output><empty_output><jupyter_text>Univariate LSTM ModelThe model above is fairly accurate but in order to see how much is due to the news, we will also train a model with only the previous closing values, also using a LSTM architecture. We use some of the code from [this tutorial](https://machinelearningmastery.com/how-to-develop-lstm-models-for-time-series-forecasting/) to accomplish this.  Format and split the dataWe have to format the data such that we set the desired number of lags as predictors. So in the case where we select the number of lags to be 2, then each row of the training data would be T-1 and T-2, and the label for that row would be T. In the code below, we can see that the `lags` variable is set to 1 because that yielded better accuracy compared to lags of 2 and 3.<jupyter_code># split a univariate sequence
def split_sequence(sequence, lags):
    X, y = list(), list()
    for i in range(len(sequence)):
        # find the end of this pattern
        end_ix = i + lags
        # check if we are beyond the sequence
        if end_ix > len(sequence)-1:
            break
        # gather input and output parts of the pattern
        seq_x, seq_y = sequence[i:end_ix], sequence[end_ix]
        X.append(seq_x)
        y.append(seq_y)
    return np.array(X), np.array(y)
lags = 1
X, y = split_sequence(stock_df['Close'].dropna().values, lags)
X = scaler.transform(X)<jupyter_output><empty_output><jupyter_text>We split the data to match the proportions above.<jupyter_code>print(X.shape)
print(y.shape)
X_train_uni = X[:num_train_items]
y_train_uni = y[:num_train_items]
print(X_train_uni.shape)
print(y_train_uni.shape)
X_val_uni = X[num_train_items:(num_train_items+num_val_items)]
y_val_uni = y[num_train_items:(num_train_items+num_val_items)]
print(X_val_uni.shape)
print(y_val_uni.shape)
X_test_uni = X[(num_train_items+num_val_items):]
y_test_uni = y[(num_train_items+num_val_items):]
print(X_test_uni.shape)
print(y_test_uni.shape)
X_train_uni = X_train_uni.reshape((X_train_uni.shape[0], X_train_uni.shape[1], 1))
X_val_uni = X_val_uni.reshape((X_val_uni.shape[0], X_val_uni.shape[1], 1))
X_test_uni = X_test_uni.reshape((X_test_uni.shape[0], X_test_uni.shape[1], 1))<jupyter_output><empty_output><jupyter_text>Train Univariate ModelWe train the model using the hyperparameters listed below as they provided the best accuracy.<jupyter_code># define model
model_uni = Sequential()
model_uni.add(LSTM(400, input_shape=(lags, 1)))
model_uni.add(Dense(1))
model_uni.compile(optimizer='adam', loss='mae')
# fit model
history_uni = model_uni.fit(X_train_uni, y_train_uni, epochs=100, batch_size=72, validation_data=(X_val_uni, y_val_uni), verbose=2, shuffle=False)<jupyter_output>Epoch 1/100
20/20 - 1s - loss: 12.4037 - val_loss: 25.4204
Epoch 2/100
20/20 - 0s - loss: 12.1544 - val_loss: 25.2510
Epoch 3/100
20/20 - 0s - loss: 11.7523 - val_loss: 24.9645
Epoch 4/100
20/20 - 0s - loss: 11.0996 - val_loss: 24.5265
Epoch 5/100
20/20 - 0s - loss: 10.1176 - val_loss: 23.9100
Epoch 6/100
20/20 - 0s - loss: 9.0976 - val_loss: 23.1932
Epoch 7/100
20/20 - 0s - loss: 8.3325 - val_loss: 22.4183
Epoch 8/100
20/20 - 0s - loss: 7.7745 - val_loss: 21.6420
Epoch 9/100
20/20 - 0s - loss: 7.4402 - val_loss: 20.9233
Epoch 10/100
20/20 - 0s - loss: 7.1532 - val_loss: 20.1993
Epoch 11/100
20/20 - 0s - loss: 6.8723 - val_loss: 19.4618
Epoch 12/100
20/20 - 0s - loss: 6.6369 - val_loss: 18.7468
Epoch 13/100
20/20 - 0s - loss: 6.4400 - val_loss: 18.0660
Epoch 14/100
20/20 - 0s - loss: 6.2337 - val_loss: 17.3858
Epoch 15/100
20/20 - 0s - loss: 6.0156 - val_loss: 16.6856
Epoch 16/100
20/20 - 0s - loss: 5.7882 - val_loss: 15.9559
Epoch 17/100
20/20 - 0s - loss: 5.5564 - val_loss: 15.1991
E[...]<jupyter_text>The model converges at around 35 epochs, so we could go from 100 to 40 epochs in order to reduce computation.<jupyter_code>plt.figure(figsize=(10,6))
pyplot.plot(history_uni.history['loss'], label='train')
pyplot.plot(history_uni.history['val_loss'], label='val')
plt.title('Univariate Model Training Error per Epoch')
plt.ylabel('Error')
plt.xlabel('Epoch')
plt.legend()
plt.savefig('../plots/univariate_model_training_error_epochs.png', dpi=300)
pyplot.show()<jupyter_output><empty_output><jupyter_text>Univariate model accuracyThe univariate model is more accurate than the model that also takes news headlines as predictors. It has a RMSE of 0.7. This means that the majority of the model's predicted values are within around -1.4 and +1.4 dollars of the actual closing values.<jupyter_code>y_hat_uni = model_uni.predict(X_test_uni)
print(y_hat_uni.shape)
rmse_uni = sqrt(mean_squared_error(y_test_uni, y_hat_uni))
rmse_uni

pd.DataFrame(
    {
        'Predicted Values': y_hat_uni.reshape((y_hat.shape[0],)),
        'Actual Values': y_test.reshape((y_test.shape[0],))
    },
    index = data_df.iloc[(num_train_items+num_val_items):,:].index
).plot(figsize=(10,6), title=f'Univariate Model Accuracy {rmse_uni} RMSE')
plt.savefig('../plots/univariate_model_accuracy.png', dpi=300)
plt.show()<jupyter_output><empty_output><jupyter_text>Headlines-Only ModelBased on the two previous results, it seems that the headlines data provide noise that affects the accuracy of the model negatively. But in order confirm that, we will train a model that only takes the headlines as input and not previous closing prices.<jupyter_code>model_no_prev = Sequential()
model_no_prev.add(LSTM(1000, input_shape=(X_train[:,:,:-1].shape[1], X_train[:,:,:-1].shape[2])))
model_no_prev.add(Dense(1))
model_no_prev.compile(loss='mae', optimizer='adam')

history_no_prev = model_no_prev.fit(X_train[:,:,:-1], y_train, epochs=100, batch_size=72, validation_data=(X_val[:,:,:-1], y_val), verbose=2, shuffle=False)<jupyter_output>Epoch 1/100
20/20 - 4s - loss: 2.3778 - val_loss: 5.9805
Epoch 2/100
20/20 - 1s - loss: 8.0990 - val_loss: 14.2006
Epoch 3/100
20/20 - 1s - loss: 5.3513 - val_loss: 11.4627
Epoch 4/100
20/20 - 1s - loss: 5.9068 - val_loss: 11.2849
Epoch 5/100
20/20 - 1s - loss: 5.9516 - val_loss: 12.0395
Epoch 6/100
20/20 - 1s - loss: 5.5893 - val_loss: 11.9154
Epoch 7/100
20/20 - 1s - loss: 5.5885 - val_loss: 11.9473
Epoch 8/100
20/20 - 1s - loss: 5.5051 - val_loss: 11.9253
Epoch 9/100
20/20 - 1s - loss: 5.4408 - val_loss: 11.8427
Epoch 10/100
20/20 - 1s - loss: 5.3967 - val_loss: 11.7587
Epoch 11/100
20/20 - 1s - loss: 5.3478 - val_loss: 11.6483
Epoch 12/100
20/20 - 1s - loss: 5.3408 - val_loss: 11.6601
Epoch 13/100
20/20 - 1s - loss: 5.2619 - val_loss: 11.5342
Epoch 14/100
20/20 - 1s - loss: 5.2358 - val_loss: 11.4247
Epoch 15/100
20/20 - 1s - loss: 5.2092 - val_loss: 11.3486
Epoch 16/100
20/20 - 1s - loss: 5.1710 - val_loss: 11.2663
Epoch 17/100
20/20 - 1s - loss: 5.1205 - val_loss: 11.1291
Epoch 1[...]<jupyter_text>The model never converges, which means that model is overfitting on the training data and not able to extrapolate to the validation data. This happened even after trying a few different output dimensions for the LSTM layer.<jupyter_code>plt.figure(figsize=(10,6))
pyplot.plot(history_no_prev.history['loss'], label='train')
pyplot.plot(history_no_prev.history['val_loss'], label='val')
plt.title('Headlines-Only Model Training Error per Epoch')
plt.ylabel('Error')
plt.xlabel('Epoch')
plt.legend()
plt.savefig('../plots/headlines_only_model_training_error_epochs.png', dpi=300)
pyplot.show()

y_hat_no_prev = model.predict(X_test[:,:,:-1])
print(y_hat_no_prev.shape)<jupyter_output>(298, 1)<jupyter_text>Headlines-Only Model Accuracy<jupyter_code>rmse_no_prev = sqrt(mean_squared_error(y_test, y_hat_no_prev))
rmse_no_prev

pd.DataFrame(
    {
        'Predicted Values': y_hat_no_prev.reshape((y_hat_no_prev.shape[0],)),
        'Actual Values': y_test.reshape((y_test.shape[0],))
    },
    index = data_df.iloc[(num_train_items+num_val_items):,:].index
).plot(figsize=(10,6), title=f'Headlines-Only Model Accuracy {rmse_no_prev} RMSE')
plt.savefig('../plots/headlines_only_model_accuracy.png', dpi=300)
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12292013311148087, 'fraction_numerical': 0.07367928452579035, 'mean_word_length': 4.688553682342502, 'pattern_counts': {'":': 0, '<': 70, '<?xml version=': 0, '>': 71, 'https://': 8, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '19', 'dir': 'jupyter-structured-clean-dedup', 'id': '401792', 'license': 'MIT', 'n_tokens_mistral': 7948, 'n_tokens_neox': 6870, 'n_words': 1846, 'path': 'notebooks/Code.ipynb', 'repo_name': 'jorgeso/nlp_market_prediction', 'size': '19232'}