Artificial neural network classification of copy number variation, part 2

Recap:

Welcome to the second part of this post series on building artificial neural network models for copy number classification.  In the first part, I described the problem with interpreting copy-ratio plots to find clinically-relevant CNV events.  The data from targeted capture deep sequencing are noisy and biased, and finding clinically-relevant genotypes in genes that have CNVs requires the analyst to visualize the CNV event and assign a classification on the basis of experience and expert knowledge.

The LASSO model

Once my training data were in place (see part 1), I used a multiple linear regression LASSO model as a machine-learning benchmark.  I did this to determine whether a more powerful neural network model would be warranted.   The LASSO model uses an “L1” prior to perform feature selection, setting some coefficients to zero as warranted by the data.   There is ample precedent for applying this type of model in bioinformatics settings where the goal is maximize predictive power without overfitting.

I fit the LASSO to the data, with 33% held out for validation.  The best fit was obtained with the alpha parameter set to 0.001.  k-fold cross validation (where k=10 and alpha=0.001) yielded an accuracy of 76%.   These results are surprisingly good, given the complexity of the CNV signals in the noisy data.  Unfortunately, 76% accuracy is simply not good enough for an automated method that will be used to predict genotypes in clinical data.

The ANN model

Next, I decided to construct an artificial neural network model.  My goal was to keep the model as simple as possible, while reaching a very high classification accuracy needed for clinical work.   To that end I constructed a one hidden-layer model with 19 input nodes corresponding to the 19 copy-ratio probes in the CNV data.  The output layer contained five nodes, corresponding to the five classes of defined CNV event or other event (for example, a very distinct sequencing artifact that kept appearing in the data):

 

In between the input and output layers I constructed a 10-node hidden layer.   A one hidden-layer neural network is the simplest form of the ANN model, and I tried to keep the number of hidden-layer nodes to a minimum as well.  Specific details about the model, hyper-parameter tuning, and the code will be available in the near future when I put a pre-print of this work on biorxiv.

Model training and cross-validation

I trained the model on the 175 sample dataset and on a 350 sample “synthetic” dataset created by adding gaussian “noise” to the real data.  The results are shown below, across 250 training epochs.

When the ANN model was tested with 10-fold cross validation, the accuracy reached a level of 96.5% (+/- 5.4%).   This is obviously a big improvement on the LASSO model, and reaches a level of accuracy that is good enough for clinical pipelines (with the caveat that low confidence predictions will still be checked “by hand.”)

Below, I’m showing a sample of the model output (left) and ground truth (right) from the test data.  The numbers (and colors) of the boxes correspond to the model’s probability in that classification.  You can see that most CNV events are called with high probability, but several (yellow boxes) are called correctly but with lower probability.   One event (red box) is called incorrectly with high probability.

Conclusions and caveats

Going into this project, I had no idea if the ANN model would be able to make predictions on the basis of so few examples in the training set.  The classic examples you see about ANN/CNN models rely on handwriting training sets with 10,000 or more images.   So I was surprised when the model did very well with extremely limited training data.   Since this method was developed for a clinical pipeline, it can be improved as the pipeline generates new training data with each new patient sample.  We would need many thousands of samples through our “legacy” pipeline to see enough examples of the rare star allele events in CYP2D6 that we could then classify them.  That is why I limited my CNV calling to three star alleles.

The low confidence, true positive predictions concern me less than the high confidence false negative.  Missing a real CNV that has impact on CYP2D6 function and therefore clinical relevance is very dangerous.  This can lead to incorrect prescribing and adverse drug reactions for the patient.   I really want to understand why the method makes predictions like this, and how to fix it.  Unfortunately, I believe it will require a lot more training data to solve this problem and that is something I lack.

My goals for this project now are 1) to publish a preprint on biorxiv describing this work and 2) to obtain some additional training/test datasets.   Because our pharmacogenomics test is not generating the kind of volume we expected, I may have to look around for another gene with clinically-relevant CNV events to test this method further.   For example, we do have an NGS-based test of hearing and deafness genes with thousands of validated patient samples.  One gene, STRC, has relevant CNVs that are complex and require analyst visualization to detect.  This may be a good system for follow up refinement of this type of model.

Artificial neural network classification of copy number variation, part 1

In this series of posts, I want to describe some work I’ve been doing attempting to apply an artificial neural network model to the problem of classifying copy number variation in a clinically-important gene, CYP2D6.   I will be presenting a poster-paper on this topic at the 2018 ISMB/ISCB meeting in Chicago.

Here is the abstract from the poster:

Pharmacogenomics is a rapidly developing field that aims to deliver on the promise of personalized medicine by guiding pharmacological intervention using an understanding of a patient’s individual genotype for drug metabolism.  By avoiding ineffective or dangerous treatments, patient outcomes could be dramatically improved and hospital costs reduced. We have developed a sequencing-based pharmacogenomics screening panel that uses targeted capture to perform deep sequencing of 200+ critical drug metabolism genes.  Assessing copy-number variation is a critical part of correctly interpreting genotypes in key drug metabolism genes, such as CYP2D6.  Historically, this has been a time-consuming step in our clinical pipeline involving large amounts of expert analyst time and data visualization.  This study presents a novel application of an artificial neural network (ANN) machine learning algorithm to learn the complex patterns in CNV data.  The result is a trained network that can quickly and accurately classify copy number events according to known training categories in the CYP2D6 gene.  We show that a simple, one hidden-layer network is sufficient to achieve the extremely high accuracy and low false-positive rate required in a high-throughput clinical setting.

Motivating factors

Motivating this work is the fact that interpreting copy number variation data from targeted capture sequencing is difficult owing to several factors.  First, the data are noisy owing to biases in capture efficiency and GC content.  Second, the copy number variation events in a gene like CYP2D6 are complex and subtle, but have dramatic impacts on the functional status of the gene.  Third, most copy number variation detection methods are optimized for whole genome sequencing with smooth and even coverage across each chromosome.

The result is that, although most of the variant calling and interpretation is automated, a person still has to sit down with the copy number variation copy-ratio plots and make a manual determination of genotype for the CYP2D6 gene (and other genes if they contain clinically-actionable CNVs).   Below is one such copy-ratio plot that illustrates the problems we face:

copy ratio plot
A copy-ratio plot for 9 clinical samples from our targeted sequencing pipeline at the CYP2D6 gene and CYP2D7 pseudogene.  Note the noise in the data.

You can see from the plot that one sample (red) has what appears to be a deletion (copy ratio ~ 0.5).   We do look at these plots individually, but the problem with noise and complexity remains.  This means that the analyst must be, in effect, a domain expert on each gene with a clinically-actionable copy number variation.   This limits pipeline throughput and is impractical if the goal is to process hundreds or thousands of patient samples per month.

A role for machine learning?

It occurred to me that the copy ratio data fall into distinct patterns that a human (like myself) learns by eye and with experience.  These patterns could be classifiable by machine-learning methods.   I considered applying a convolution neural network (CNN) to the plot above.

However, I only had 175 samples in my test/train set drawn from real patients and CEPH/Coriell depositories.   Therefore, I thought it would be best to start with simple models and go from there.  To establish the ground truth CYP2D6 genotype, the samples had all be processed by our current pipeline methods, with many (but not all) confirmed by Taqman assay.

It began to dawn on me that I first needed a better representation of the training data.   I liked the output of the CNV-kit method (below) better than the method from our pipeline (above):

Visualization of CNV-kit output at CYP2D6. The gray and orange lines are not important to this discussion.

The important part of the plot above are the gray dots, representing copy ratio values across 19 different bins or segments along the CYP2D6/2D7 gene and pseudogene.   CNV-kit also had the advantage of being well-documented, fast to run, and available as open-source software on github.

Tidying the training set

The data plotted above looked like this after I did some wrangling and tidying in python:

Copy ratio data from CNV-kit output at the CYP2D6 and 2D7 gene/pseudogene.  There are five columns to the right (four not shown) that contain “one-hot” encoded ground-truth for model training.

With the data in “tidy” form (each column a variable, each row an observation, thank you Hadley Wickham), I was ready to train some machine learning models to see if I could classify common copy number variation events in CYP2D6 and bypass the time-consuming visualization step.

I had no idea if this would work, given that there are only 175 samples (total) to train on and some of the more rare copy number variation events only have a handful of examples.  In part 2, I will talk more about how I tried a LASSO regression model, which, while performing well, failed to yield the high accuracy needed for an automated clinical pipeline.  I then tried a simple one-layer neural network approach.  I’ll talk about the surprisingly good performance of neural network approach and also future directions for this project.