Classifying data with Python

Using scikit-learn to build a machine-learning pipeline that compares how well different classifiers predict the movement of a stock price.



This post walks through a solution to problem 10 from Chapter 4 of An Introduction to Statistical Learning. It demonstrates how to use Python features to build a practical machine learning pipeline from first principles.

The first part of this post specifies the problem, a binary classification problem, and introduces the fundamental differences between the classifiers considered. Subsequent sections introduce the dataset and perform exploratory data analysis to identify significant variables. I finish by introducing a simple machine learning pipeline built using standard Python libaries that applies combinations of preprocessing scripts and classifiers to the dataset.

The problem

I’m tackling Q10 from Chapter 4 of ISL that can be summarised as follows:

Compare the performance of logistic regression, linear discriminant analysis, quadratic discriminant analysis and k nearest neighbours on predicting whether a stock price will go up or down in the Weekly dataset.

The Weekly dataset records 1089 stock market returns for 21 years and can be downloaded from this location.

Maths

The Weekly dataset exemplifies a binary classification problem: given observations X, predict the binary response Y. We model the probability that Y belongs to a particular category. The equation below calculates the probability that X belongs to category 1.

p(X) = Pr(Y = 1 | X)

Lets describe the key differences between each classification method. Logistic regression models p(X) using the logistic function:

p(X) = \frac{e^{\beta_0 + \beta_1 X}}{1 + e^{\beta_0 + \beta_1 X}}

The coefficients \beta_0 and \beta_1 are fitted using the Maximum Likelihood method (see this explanation). This method selects \beta_0 and \beta_1 to make predicted probabilities match observed outcomes as closely as possible.

Linear discriminant analysis (LDA) models the distribution of predictors X given Y, using Bayes’ theorem to estimate Pr(Y=k|X=x). LDA outperforms logistic regression when classes are well-separated, sample sizes are small, or predictors follow a normal distribution. LDA also handles more than two response classes effectively.

LDA assumes observations within each class follow a multivariate Gaussian distribution with a class-specific mean and a covariance matrix shared across all classes. Quadratic Discriminant Analysis (QDA) relaxes this assumption: each class has its own covariance matrix.

The k Nearest Neighbour algorithm estimates the conditional distribution of Y given X and then classifies an observation to the class with the highest estimated probability. Given a positive integer K and test observation x_0 the k-NN first identifies the K points in the training data which are closest to x_0, represented by N_0. The classifier then estimates the conditional probability for class j as the fraction of points in N_0 whose response values equal j.

Pr (Y = j|X = x_0) = \frac{1}{K} \sum_{i \in N_0} I(y_i = j)

The dataset

The Weekly dataset contains the Direction column as the response variable (up or down) and uses the Volume column plus five lag variables as predictors. Here’s the code to load the data:

from pathlib import Path
import pyreadr
import pandas as pd

def get_dataset():

    # include a try/except clause to load data if it doesn't exist locally
    url = "https://github.com/cran/ISLR/blob/master/data/Weekly.rda?raw=true"
    dst_path = "Weekly.rda"
    if Path(dst_path).exists():
        Weekly = pyreadr.read_r(dst_path)
    else:
        Weekly = pyreadr.read_r(pyreadr.download_file(url, dst_path))
    df_weekly = Weekly['Weekly']
    df_weekly['Direction_text'] = df_weekly['Direction']
    # convert categorical Direction to numeric: 1 = Up, 0 = Down
    df_weekly['Direction'] = (df_weekly['Direction'] == 'Up').astype(int)
    return df_weekly

df_weekly = get_dataset()
print(df_weekly.describe().to_markdown())

Returns the table below. The lag variables share identical minimum and maximum values because they derive from the lagged Today column. Direction is coded as 1 for an up week and 0 for a down week, so its mean of 0.556 says that 605 of the 1089 weeks were up. Notice that columns are not normalised—values don’t fall between 0 and 1. This matters: distance-based classifiers like k-NN are sensitive to the original data scale.

YearLag1Lag2Lag3Lag4Lag5VolumeTodayDirection
count108910891089108910891089108910891089
mean2000.050.1505850.1510790.1472050.1458180.1398931.574620.1498990.555556
std6.033182.357012.357252.36052.360282.361281.686642.356930.497132
min1990-18.195-18.195-18.195-18.195-18.1950.087465-18.1950
25%1995-1.154-1.154-1.158-1.158-1.1660.332022-1.1540
50%20000.2410.2410.2410.2380.2341.002680.2411
75%20051.4051.4091.4091.4091.4052.053731.4051
max201012.02612.02612.02612.02612.0269.3282112.0261

Exploratory analysis

To identify which variables correlate significantly with Direction, I use logistic regression. I adapted the approach from the statsmodels documentation.

Here’s the logistic regression model. sm.Logit does not add an intercept on its own, so sm.add_constant supplies the \beta_0 term. The null hypothesis, H_0, posits no relationship between the predictors (Lag1 through Lag5, Volume) and stock price Direction. The chained methods format output directly as markdown, avoiding manual copying.

from io import StringIO
import statsmodels.api as sm

X = sm.add_constant(df_weekly[['Lag1', 'Lag2', 'Lag3', 'Lag4', 'Lag5', 'Volume']])
y = df_weekly['Direction']

# building the model and fitting the data
log_reg = sm.Logit(y, X).fit()
log_reg_summary = log_reg.summary()
results_as_html = log_reg_summary.tables[1].as_html()
print(pd.read_html(StringIO(results_as_html), header=0, index_col=0)[0].to_markdown())

Apart from the intercept, the logistic regression identifies only Lag2 as a significant predictor: its p-value of 0.030 is below 0.05, and its positive coefficient means a strong return two weeks ago slightly raises the probability of an up week. Volume, with a p-value of 0.54, is nowhere near significant 1. Lag2 is therefore the main classifier input; I also try Lag2 with Volume to see whether adding a non-significant predictor helps.

coefstd errzP>|z|[0.0250.975]
const0.26690.0863.1060.0020.0980.435
Lag1-0.04130.026-1.5630.118-0.0930.01
Lag20.05840.0272.1750.030.0060.111
Lag3-0.01610.027-0.6020.547-0.0680.036
Lag4-0.02780.026-1.050.294-0.080.024
Lag5-0.01450.026-0.5490.583-0.0660.037
Volume-0.02270.037-0.6160.538-0.0950.05

Applying classifiers

The task requires comparing four classifiers on the dataset:

  1. Logistic regression
  2. Linear discriminant analysis
  3. Quadratic discriminant analysis
  4. k Nearest neighbours

Though scikit-learn Pipeline exists, building a custom pipeline explores how core Python libraries work together.

I started by defining classes to describe the Preprocessing stage of the pipeline and the Classifier using the dataclasses module. The Preprocessing class holds the name and a function which performs a set of transforms on the dataset to create the test and training datasets. The Classifier holds the name, a scikit-learn classifier and optionally the classifier score and Preprocessing class. I’ve added optional type hints as these are a massive help when understanding code.

The run_classifier function takes a list of n Preprocessing and m Classifier classes and uses itertools.product to create an n \times m list of pipelines to be run on the Weekly dataset. itertools.product is similar to np.meshgrid but for Python objects rather than np.arrays.

Instead of filtering the results dataframe, I use attrgetter from the operator module to find the classifier with the highest score directly.

# create budget pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.neighbors import KNeighborsClassifier
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis, QuadraticDiscriminantAnalysis
from sklearn.metrics import accuracy_score
from sklearn.base import BaseEstimator
from typing import List, Callable, Tuple, Optional
from dataclasses import dataclass
from operator import attrgetter
import itertools

pd.options.display.float_format = '${:,.5f}'.format

@dataclass
class Preprocess:
    """
    Class to hold information on the preprocessing used on a dataset.

    """
    name: str
    transform: Callable


@dataclass
class Classifier:
    """
    Class to hold information on a Classifier.

    """
    name: str
    model: BaseEstimator
    preprocess: Optional[Preprocess] = None
    score: float = 0.0


def lag2_preprocess(df: pd.DataFrame) -> Tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame, pd.DataFrame]:
    """This is the default test training split required by problem."""
    idx_train = df['Year'] < 2009.0
    idx_test = df['Year'] > 2008.0
    all_cols = ['Lag1', 'Lag2', 'Lag3', 'Lag4', 'Lag5', 'Volume']
    X_train = df['Lag2'][idx_train].values[:, None] # need to pass a 2d array to .fix(X, y)
    y_train = df['Direction'][idx_train]
    X_test = df['Lag2'][idx_test].values[:, None]
    y_test = df['Direction'][idx_test]
    return X_train, y_train, X_test, y_test


def lag2_volume_preprocess(df: pd.DataFrame) -> Tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame, pd.DataFrame]:
    """This is the default test training split required by problem."""
    idx_train = df['Year'] < 2009.0
    idx_test = df['Year'] > 2008.0
    all_cols = ['Lag1', 'Lag2', 'Lag3', 'Lag4', 'Lag5', 'Volume']
    X_train = df[['Lag2', 'Volume']][idx_train]
    y_train = df['Direction'][idx_train]
    X_test = df[['Lag2', 'Volume']][idx_test]
    y_test = df['Direction'][idx_test]
    return X_train, y_train, X_test, y_test


def apply_classifier(model: Classifier, dataset: pd.DataFrame) -> Classifier:
    """
    Apply classifier to dataset given function to create training and test datasets.

    """
    X_train, y_train, X_test, y_test = model.preprocess.transform(dataset)
    model.model = model.model.fit(X_train, y_train)
    y_predict = model.model.predict(X_test)
    model.score = accuracy_score(y_predict, y_test)
    return model


def run_classifiers():
    # define preprocessing functions
    preprocessors = [
        Preprocess("Lag 2", lag2_preprocess),
        Preprocess("Lag 2 and Volume", lag2_volume_preprocess),
    ]

    # define classifier models
    classifiers = [
        Classifier("Logistic Regression", LogisticRegression(random_state=0)),
        Classifier("Linear Discriminant Analysis", LinearDiscriminantAnalysis()),
        Classifier("Quadratic Discriminant Analysis", QuadraticDiscriminantAnalysis()),
        Classifier("kNN, k=1", KNeighborsClassifier(n_neighbors=1)),
        Classifier("kNN, k=2", KNeighborsClassifier(n_neighbors=2))
    ]

    # associate preprocessor function with classifiers for analysis
    prepped_classifiers = []
    for preprocessor, classifier in itertools.product(preprocessors, classifiers):
        prepped_classifier = Classifier(classifier.name, classifier.model, preprocessor)
        prepped_classifiers.append(prepped_classifier)

    # load data
    df_weekly = get_dataset()

    # train classifiers on datasets
    trained_classifiers = [apply_classifier(classifier, df_weekly) for classifier in prepped_classifiers]

    # create results dataframe
    df = pd.DataFrame({"Classifier name": [classifier.name for classifier in trained_classifiers],
                       "Preprocess name": [classifier.preprocess.name for classifier in trained_classifiers],
                       "Score": [classifier.score for classifier in trained_classifiers]})

    with pd.option_context('display.float_format', '{:,.5f}'.format):
        print(df.to_html())

    # get classifiers with maximum scores
    # returns the first classifier with the highest score in the event of a tie
    best_classifier = max(trained_classifiers, key=attrgetter('score'))
    print("The best classifier is {0} using the {1} preprocessor with a score of {2:.4f}".format(best_classifier.name, best_classifier.preprocess.name, best_classifier.score))


run_classifiers()

Returns:

The best classifier is Logistic Regression using the Lag 2 preprocessor with a score of 0.6250

Logistic Regression and Linear Discriminant Analysis achieve the best performance overall, at 0.625 with Lag2 alone, and produce identical results with both preprocessing methods. Adding Volume makes logistic regression, LDA and QDA worse, as its p-value suggested it would. Only kNN with k=1 improves, and with both predictors it edges ahead of the linear methods, 0.558 to 0.538. The difference lies in how they fit coefficients: logistic regression uses maximum likelihood, while LDA derives coefficients from an estimated normal distribution. When data is approximately Gaussian, both methods converge. Linear Discriminant Analysis outperforms logistic regression when data is truly Gaussian; logistic regression wins otherwise.

Classifier name Preprocess name Score
0 Logistic Regression Lag 2 0.62500
1 Linear Discriminant Analysis Lag 2 0.62500
2 Quadratic Discriminant Analysis Lag 2 0.58654
3 kNN, k=1 Lag 2 0.49038
4 kNN, k=2 Lag 2 0.46154
5 Logistic Regression Lag 2 and Volume 0.53846
6 Linear Discriminant Analysis Lag 2 and Volume 0.53846
7 Quadratic Discriminant Analysis Lag 2 and Volume 0.47115
8 kNN, k=1 Lag 2 and Volume 0.55769
9 kNN, k=2 Lag 2 and Volume 0.46154

Review

This post demonstrates how core Python libraries combine to build a machine learning pipeline for binary classification. Four classifiers were applied to the Weekly dataset to predict stock price direction. The custom pipeline implementation showcased Python’s dataclasses, itertools, and operator modules. Extending this work could involve applying different scaling methods to improve classifier performance.

Footnotes

  1. Typical p-value cutoffs are 1\% or 5\%. Choosing which p-value to use and even whether they should be used in this way is an area of academic debate. I enjoyed this tweet on the subject. ↩