Finding AD biomarkers in proteomics data

Feyn version: 1.5+

Can the QLattice deal with omics data that is noisy and contains thousands of features? It certainly can!

Omics data typically contains hundreds to thousands of features (proteins, transcripts, methylated DNA etc.) that are measured in samples derived from sources such as blood, tissue or cell culture. These types of approaches are often used for exploratory analysis e.g. in biomarker discovery or understanding the mechanism of action of a certain drug. It often resembles a bit of a "fishing exercise".

Thus, there is a need to quickly and reliably identify the most important features and their interactions that contribute to a certain signal (e.g. disease state, cell-type identity, cancer detection).

In this tutorial we present a brief workflow for building simple and interpretable models for proteomics data. This specific example is taken from a study by Bader & Geyer et al. 2020 (Mann group) and contains samples taken from the cerebrospinal fluid of Alzheimer Disease (AD) patients and non-AD patients. We will show you how to build QLattice model that can classify people into AD and non-AD according to their proteomic profiles.

The dataset contains over a thousand features (features in this example describe the intensity of different proteins measured by mass spectrometry).

import numpy as np
import pandas as pd
import feyn
from feyn.tools import split
from feyn.losses import binary_cross_entropy
from sklearn.model_selection import train_test_split
import warnings
warnings.filterwarnings('ignore')
import matplotlib.pyplot as plt

Load the data

Note, the data has been preprocessed and missing values have been imputed. It contains 1166 proteins and 88 non-AD and 49 AD subjects.

data = pd.read_csv("../data/ad_omics.csv")

# Let's record the categorical data types in our dataset (note features will be treated as numerical by default).
stypes = {}
for f in data.columns:
    if data[f].dtype =='object':
        stypes[f] = 'c'
        
# In addition to numerical features we also include two different categorical features
stypes

{'_collection site': 'c', '_gender': 'c'}

Split dataset into train and test set

# Set random seed for reproducibility
random_seed = 552

# Define the target variable
target = "_clinical AD diagnosis"

# Split
train, test = train_test_split(data, test_size=0.33, stratify=data[target], random_state=random_seed)

Train the QLattice

ql = feyn.QLattice() # Connecting
ql.reset(random_seed=random_seed) # Resetting

qgraph = ql.get_classifier(train.columns, output = target, stypes=stypes, max_depth=2) # Extract classifier models

# Train the QLattice
for _ in range(50):
    qgraph.fit(train, threads=8, loss_function=binary_cross_entropy)
    ql.update(qgraph.best())

Inspect the top graph

qgraph[0].plot_summary(test)

As expected, MAPT (i.e. Tau) seems to be driving most of the signal here. Let's investigate further.

Explore features

Let's look at how the different features play together.

qgraph[0].plot_partial(test, by = "MAPT", fixed = {"hCG_2039566" : 76000, "IGLV4-69" : 222999.55})

png

This partial derivate plot shows you how higher SELENBP1 levels shift the MAPT curve to the right. I.e. the higher your SELENBP1 levels, the higher your MAPT levels have to be for a positive AD prediction.

In the following 2d partial plots we show the same phenomoneon for the other 2 features as well.

qgraph[0].plot_partial2d(test, fixed = {"SELENBP1" : 285584, "IGLV4-69" : 222999.55})
qgraph[0].plot_partial2d(test, fixed = {"SELENBP1" : 285584, "hCG_2039566" : 32000})
qgraph[0].plot_partial2d(test, fixed = {"IGLV4-69" : 222999.55, "hCG_2039566" : 32000})

<AxesSubplot:xlabel='MAPT', ylabel='SELENBP1'>

png

How well are we doing? - ROC Curve

qgraph[0].plot_roc_curve(test)

png

Confusion Matrix

qgraph[0].plot_confusion_matrix(test)

png

There is much more to be explored. Watch this space and get in touch if you'd like to give it a shot.