Patent Application
20220027388
The technology disclosed relates to the use of machine learning techniques, which may be referred to as artificial intelligence, implemented on computers and digital data processing systems for the purpose of assessing the pathogenicity of biological sequence variants and deriving other pathogenicity related data using the pathogenicity assessments. These approaches may include or utilize corresponding data processing methods and products for emulation of intelligence (i.e., knowledge based systems, reasoning systems, and knowledge acquisition systems) and/or systems for reasoning with uncertainty (e.g., fuzzy logic systems), adaptive systems, machine learning systems, and artificial neural networks. In particular, the technology disclosed relates to using deep learning-based techniques for training deep convolutional neural networks for pathogenicity assessment as well as the use or refinement of such pathogenicity information.
The subject matter discussed in this section should not be assumed to be prior art merely as a result of its mention in this section. Similarly, a problem mentioned in this section or associated with the subject matter provided as background should not be assumed to have been previously recognized in the prior art. The subject matter in this section merely represents different approaches, which in and of themselves can also correspond to implementations of the claimed technology.
Genetic variations can help explain many diseases. Every human being has a unique genetic code and there are many genetic variants within a group of individuals. Many or most genetic variants that are deleterious have been depleted from genomes by natural selection. However, it is still difficult to identify which genetic variations are likely to be of clinical interest.
Further, modeling the properties and functional effects (e.g., pathogenicity) of variants is a challenging task in the field of genomics. Despite the rapid advancement of functional genomic sequencing technologies, interpretation of the functional consequences of variants remains a great challenge due to the complexity of cell type-specific transcription regulation systems.
Systems, methods, and articles of manufacture are described for constructing a variant pathogenicity classifier and for using or refining such pathogenicity classifier information. Such implementations may include or utilize non-transitory computer readable storage medium storing instructions executable by a processor to perform actions of the system and methodology described herein. One or more features of an implementation can be combined with the base implementation or other implementations, even if not explicitly listed or described. Further, implementations that are not mutually exclusive are taught to be combinable such that one or more features of an implementation can be combined with other implementations. This disclosure periodically may reminds the user of these options. However, omission from some implementations of recitations that repeat these options should not be taken as limiting the potential combinations taught in the following sections. Instead, these recitations are hereby incorporated forward by reference into each of the following implementations.
This system implementation and other systems disclosed optionally include some or all of the features as discussed herein. The system can also include features described in connection with methods disclosed. In the interest of conciseness, alternative combinations of system features are not individually enumerated. Further, features applicable to systems, methods, and articles of manufacture are not repeated for each statutory class set of base features. The reader will understand how features identified can readily be combined with base features in other statutory classes.
In one aspect of the discussed subject matter, methodologies and systems are described that train a convolutional neural network-based variant pathogenicity classifier, which runs on numerous processors coupled to memory. Alternatively, in other system implementations trained or suitably parameterized statistical models or techniques and/or other machine learning approaches may be employed in addition to or in the alternative of neural network-based classifiers. The system uses benign training examples and pathogenic training examples of protein sequence pairs generated from benign variants and pathogenic variants. The benign variants include common human missense variants and non-human primate missense variants occurring on alternative non-human primate codon sequences that share matching reference codon sequences with humans. The sampled humans may belong to different human subpopulations which may include or be characterized as: African/African American (abbreviated AFR), American (abbreviated AMR), Ashkenazi Jewish (abbreviated ASJ), East Asian (abbreviated EAS), Finnish (abbreviated FIN), Non-Finnish European (abbreviated NFE), South Asian (abbreviated SAS), and Others (abbreviated OTH). The non-human primate missense variants include missense variants from a plurality of non-human primate species, including, but not necessarily limited to: Chimpanzee, Bonobo, Gorilla, B. Orangutan, S. Orangutan, Rhesus, and Marmoset.
As discussed herein, a deep convolutional neural network, running on the numerous processors, may be trained to classify a variant amino acid sequence as benign or pathogenic. Thus, an output of such a deep convolutional neural network may include, but is not limited to, a pathogenicity score or classification for the variant amino acid sequence. As may be appreciated, in certain implementations suitably parameterized statistical models or techniques and/or other machine learning approaches may be employed in addition to or in the alternative of neural network-based approaches.
In certain embodiments discussed herein, the pathogenicity processing and/or scoring operations may include further features or aspects. By way of example, various pathogenicity scoring thresholds may be employed as part of the evaluation or assessment process, such as to assess or score a variant as benign or pathogenic. By way of example, in certain implementations, a suitable percentile of a pathogenicity score per gene for use as a threshold for likely pathogenic variants may be in the range from 51% to 99%, such as, but not limited to, the 51st, 55th, 65th, 70th, 75th, 80th, 85th, 90th, 95th, or 99th percentile. Conversely, a suitable percentile of the pathogenicity score per gene for use as a threshold for likely benign variants may be in the range from 1% to 49%, such as, but not limited to, the 1st, 5th, 10th, 15th, 20th, 25th, 30th, 35th, 40th, or 45th percentile.
In further embodiments, the pathogenicity processing and/or scoring operations may include further features or aspects that allow selection effects to be estimated. In such embodiments, forward time simulation of allele frequencies within a given population, using suitable inputs characterizing mutation rates and/or selection, may be employed to generate an allele frequency spectrum at a gene of interest. Depletion metrics may then be calculated for variants of interest, such as by comparing allele frequency spectra with and without selection, and a corresponding selection-depletion function fitted or characterized. Based on a given pathogenicity score and this selection-depletion function, a selection coefficient can be determined for a given variant based on the pathogenicity score generated for the variant.
In additional aspects, the pathogenicity processing and/or scoring operations may include further features or aspects that allow genetic disease prevalence to be estimated using pathogenicity scores. With respect to calculating a genetic disease prevalence metric for each gene, in a first methodology the trinucleotide context configurations of a set of deleterious variants is initially obtained. For each trinucleotide context in this set, a forward-time simulation assuming certain selection coefficients (e.g. 0.01) is performed to generate the expected allele frequency spectrum (AFS) for that trinucleotide context. Summing up AFS across the trinucleotides weighted by frequencies of trinucleotides in a gene produces the expected AFS for the gene. The genetic disease prevalence metric in accordance with this approach may be defined as the expected cumulative allele frequencies of variants having pathogenicity scores exceeding the threshold for that gene.
In a further aspect, the pathogenicity processing and/or scoring operations may include features or methodology to recalibrate the pathogenicity scoring. With respect to such as recalibration, in one example embodiment, a recalibration approach may focus on the percentiles of pathogenicity scores of variants, as these may be more robust and less affected by selection pressure exerted on the whole genes. In accordance with one implementation, a survival probability for each percentile of pathogenicity scores is calculated, which constitutes a survival probability correction factor that implies that the higher the percentile of pathogenicity scores, the less chance the variant survives purifying selection. The survival probability correction factor may be employed to perform a recalibration so as to help mitigate the effects of noise on the estimation of selection coefficients in missense variants.
The preceding description is presented to enable the making and use of the technology disclosed. Various modifications to the disclosed implementations will be apparent, and the general principles defined herein may be applied to other implementations and applications without departing from the spirit and scope of the technology disclosed. Thus, the technology disclosed is not intended to be limited to the implementations shown, but is to be accorded the widest scope consistent with the principles and features disclosed herein. The scope of the technology disclosed is defined by the appended claims.
These and other features, aspects, and advantages of the present invention will become better understood when the following detailed description is read with reference to the accompanying drawings in which like characters represent like parts throughout the drawings, wherein:
The following discussion is presented to enable any person skilled in the art to make and use the technology disclosed, and is provided in the context of a particular application and its requirements. Various modifications to the disclosed implementations will be readily apparent to those skilled in the art, and the general principles defined herein may be applied to other implementations and applications without departing from the spirit and scope of the technology disclosed. Thus, the technology disclosed is not intended to be limited to the implementations shown, but is to be accorded the widest scope consistent with the principles and features disclosed herein.
The following discussion includes aspects related to the training and use of neural networks, including convolutional neural networks, which may be used to implement certain of the analytics discussed below, such as the generation of variant pathogenicity scores or classifications and the derivation of useful clinical analytics or metrics based on such pathogenicity scores or classifications. With this in mind, certain aspects and features of such neural networks may be mentioned or referenced in describing the present techniques. To streamline discussion, a baseline knowledge of such neural networks is presumed in describing the present techniques. However, additional information and description of relevant neural network concepts is provided toward the end of the description for those wanting additional explanation of the relevant neural network concepts. Further, it should also be appreciated that, though neural networks are primarily discussed herein to provide a useful example and to facilitate explanation, other implementations may be employed in place of or in addition to neural network approaches, including but not limited to trained or suitably parameterized statistical models or techniques and/or other machine learning approaches.
In particular, the following discussion may utilize certain concepts related to neural networks (e.g., convolutional neural networks) in implementations used to analyze certain genomic data of interest. With that in mind, certain aspects of the underlying biological and genetic problems of interest are outlined here to provide a useful context to the problems being addressed herein, for which the neural network techniques discussed herein may be utilized.
Genetic variations can help explain many diseases. Every human being has a unique genetic code and there are many genetic variants within a group of individuals. Many or most genetic variants that are deleterious have been depleted from genomes by natural selection. However, it is still desirable to identify which genetics variations are likely to be pathogenic or deleterious. In particular, such knowledge may help researchers focus on the likely pathogenic genetic variants and accelerate the pace of diagnosis and cure of many diseases.
Modeling the properties and functional effects (e.g., pathogenicity) of variants is an important but challenging task in the field of genomics. Despite the rapid advancement of functional genomic sequencing technologies, interpretation of the functional consequences of variants remains a great challenge due to the complexity of cell type-specific transcription regulation systems. Therefore, a powerful computational model for predicting the pathogenicity of variants can have substantial benefits for both basic science and translational research.
Further, advances in biochemical technologies over the past decades have given rise to next generation sequencing (NGS) platforms that quickly produce genomic data at much lower costs than ever before, resulting in ever increasing amounts of genomic data being generated. Such overwhelmingly large volumes of sequenced DNA remain difficult to annotate. Supervised machine learning algorithms typically perform well when large amounts of labeled data are available. However, in bioinformatics and many other data-rich disciplines, the process of labeling instances is costly. Conversely, unlabeled instances are inexpensive and readily available. For a scenario in which the amount of labeled data is relatively small and the amount of unlabeled data is substantially larger, semi-supervised learning represents a cost-effective alternative to manual labeling. This therefore presents an opportunity to use semi-supervised algorithms to construct deep learning-based pathogenicity classifiers that accurately predict pathogenicity of variants. Databases of pathogenic variants that are free from human ascertainment bias may result.
Regarding machine-learning based pathogenicity classifiers, deep neural networks are a type of artificial neural networks that use multiple nonlinear and complex transforming layers to successively model high-level features. Deep neural networks provide feedback via backpropagation which carries the difference between observed and predicted output to adjust parameters. Deep neural networks have evolved with the availability of large training datasets, the power of parallel and distributed computing, and sophisticated training algorithms.
Convolutional neural networks (CNNs) and recurrent neural networks (RNNs) are components of deep neural networks. Convolutional neural networks may have an architecture that comprises convolution layers, nonlinear layers, and pooling layers. Recurrent neural networks are designed to utilize sequential information of input data with cyclic connections among building blocks like perceptrons, long short-term memory units, and gated recurrent units. In addition, many other emergent deep neural networks have been proposed for limited contexts, such as deep spatio-temporal neural networks, multi-dimensional recurrent neural networks, and convolutional auto-encoders.
Given that sequence data are multi- and high-dimensional, deep neural networks have promise for bioinformatics research because of their broad applicability and enhanced prediction power. Convolutional neural networks have been adapted to solve sequence-based problems in genomics such as motif discovery, pathogenic variant identification, and gene expression inference. Convolutional neural networks use a weight-sharing strategy that is useful for studying DNA because it can capture sequence motifs, which are short, recurring local patterns in DNA that are presumed to have significant biological functions. A hallmark of convolutional neural networks is the use of convolution filters. Unlike traditional classification approaches that are based on elaborately-designed and manually-crafted features, convolution filters perform adaptive learning of features, analogous to a process of mapping raw input data to the informative representation of knowledge. In this sense, the convolution filters serve as a series of motif scanners, since a set of such filters is capable of recognizing relevant patterns in the input and updating themselves during the training procedure. Recurrent neural networks can capture long-range dependencies in sequential data of varying lengths, such as protein or DNA sequences.
Training deep neural network, as depicted schematically in
Another element in the training of deep neural networks is regularization, which refers to strategies intended to avoid overfitting and thus achieve good generalization performance. For example, weight decay adds a penalty term to the objective loss function so that weight parameters converge to smaller absolute values. Dropout randomly removes hidden units from neural networks during training and can be considered an ensemble of possible subnetworks. Furthermore, batch normalization provides a new regularization method through normalization of scalar features for each activation within a mini-batch and learning each mean and variance as parameters.
With the preceding high-level overview in mind with respect to the presently described techniques, the presently described techniques differ from prior pathogenicity classification models, which employ a large number of human-engineered features and meta-classifiers. In contrast, in certain embodiments of the techniques described herein a simple deep learning residual network may be employed which takes as input only the amino acid sequence flanking a variant of interest and the orthologous sequence alignments in other species. In certain implementations, to provide the network with information about protein structure two separate networks may be trained to learn, respectively, secondary structure and solvent accessibility from sequence alone. These may be incorporated as sub-networks in the larger deep learning network to predict effects on protein structure. Using sequence as a starting point avoids potential biases in protein structure and functional domain annotation, which may be incompletely ascertained or inconsistently applied.
The accuracy of the deep learning classifier scales with the size of the training dataset, and variation data from each of six primate species independently contributes to boosting the accuracy of the classifier. The large number and diversity of extant non-human primate species, along with evidence showing that the selective pressures on protein-altering variants are largely concordant within the primate lineage, suggests systematic primate population sequencing as an effective strategy to classify the millions of human variants of unknown significance that currently limit clinical genome interpretation.
Further, common primate variation provides a clean validation dataset for evaluating existing methods that is completely independent of previously used training data, which has been hard to evaluate objectively because of the proliferation of meta-classifiers. The performance of the present model described herein was evaluated, along with four other popular classification algorithms (Sift, Polyphen2, CADD, M-CAP), using 10,000 held-out primate common variants. Because roughly 50% of all human missense variants would be removed by natural selection at common allele frequencies, the 50th-percentile score was calculated for each classifier on a set of randomly picked missense variants that were matched to the 10,000 held-out primate common variants by mutational rate, and that threshold was used to evaluate the held-out primate common variants. The accuracy of the presently disclosed deep learning model was significantly better than the other classifiers on this independent validation dataset, using either deep learning networks that were trained only on human common variants or using both human common variants and primate variants.
With the preceding in mind, and by way of summary, the methodology described herein differs from the existing methods for predicting pathogenicity of variants in various ways. First, the presently described approach adopts a novel architecture of semi-supervised deep convolutional neural networks. Second, reliable benign variants are obtained from human common variants (e.g., from gnomAD) and primate variants, while the highly confident pathogenic training set is generated through iterative balanced sampling and training, to avoid circular training and testing of models using the identical human curated variant databases. Third, deep learning models for secondary structure and solvent accessibility are integrated into the architecture of the pathogenicity model. The information obtained from the structure and solvent models are not limited to label prediction for specific amino acid residues. Rather, the readout layer is removed from the structure and solvent models, and the pre-trained models are merged with the pathogenicity model. While training the pathogenicity model, the structure and solvent pre-trained layers also backpropagate to minimize the error. This helps the pre-trained structure and solvent model to focus on the pathogenicity prediction problem.
As also discussed herein, outputs (e.g., pathogenicity scores and/or classifications) of a model trained and used as described herein, may be used to generate additional data or diagnoses of value, such as estimation of selection effects against a range of clinically significant variants and estimation of genetic disease prevalence. Other relevant concepts, such as recalibration of the model outputs and the generation and use of threshold values for characterizing pathogenic and benign variants are also described.
As used herein:
A base refers to a nucleotide base or nucleotide, A (adenine), C (cytosine), T (thymine), or G (guanine).
The terms “protein” and “translated sequence” may be used interchangeably.
The terms “codon” and “base triplet” may be used interchangeably.
The terms “amino acid” and “translated unit” may be used interchangeably.
The phrases “variant pathogenicity classifier”, “convolutional neural network-based classifier for variant classification”, and “deep convolutional neural network-based classifier for variant classification” may be used interchangeably.
A. Training and Inputs
Turning to example implementations, a deep learning network is described herein that may be used for variant pathogenicity classification (e.g., pathogenic or benign) and/or generation of a quantitative metric (e.g., a pathogenicity score) that numerically characterizes pathogenicity or the lack of pathogenicity. In one context of training a classifier using only variants with benign labels, the prediction problem was framed as whether a given mutation is likely to be observed as a common variant in the population. Several factors influence the probability of observing a variant at high allele frequencies, though deleteriousness is the primary focus of the present discussion and description. Other factors include, but are not limited to: mutation rate, technical artifacts such as sequencing coverage, and factors impacting neutral genetic drift (such as gene conversion).
With respect to training the deep learning network, the importance of variant classification for clinical applications has inspired numerous attempts to use supervised machine learning to address the problem, but these efforts have been hindered by the lack of an adequately sized truth dataset containing confidently labeled benign and pathogenic variants for training.
Existing databases of human expert curated variants do not represent the entire genome, with ˜50% of the variants in the ClinVar database coming from only 200 genes (˜1% of human protein-coding genes). Moreover, systematic studies identify that many human expert annotations have questionable supporting evidence, underscoring the difficulty of interpreting rare variants that may be observed in only a single patient. Although human expert interpretation has become increasingly rigorous, classification guidelines are largely formulated around consensus practices and are at risk of reinforcing existing tendencies. To reduce human interpretation biases, recent classifiers have been trained on common human polymorphisms or fixed human-chimpanzee substitutions, but these classifiers also use as their input the prediction scores of earlier classifiers that were trained on human curated databases. Objective benchmarking of the performance of these various methods has been elusive in the absence of an independent, bias-free truth dataset.
To address this issue, the presently described techniques leverage variation from non-human primates (e.g., chimpanzee, bonobo, gorilla, orangutan, rhesus, and marmoset) that contributes over 300,000 unique missense variants that are non-overlapping with common human variation and largely represent common variants of benign consequence that have been through the sieve of purifying selection. This greatly enlarges the training dataset available for machine learning approaches. On average, each primate species contributes more variants than the whole of the ClinVar database (˜42,000 missense variants as of November 2017, after excluding variants of uncertain significance and those with conflicting annotations). Additionally, this content is free from biases in human interpretation.
With respect to generating a benign training dataset for use in accordance with the present techniques, one such dataset was constructed of largely common benign missense variants from human and non-human primates for machine learning. The dataset comprised common human variants (>0.1% allele frequency; 83,546 variants), and variants from chimpanzee, bonobo, gorilla, and orangutan, rhesus, and marmoset (301,690 unique primate variants).
Using one such dataset that comprised common human variants (allele frequency (AF)>0.1%) and primate variation, a deep residual network, referred to herein as PrimateAI or pAI, was trained. The network was trained to receive as an input the amino acid sequence flanking the variant of interest and the orthologous sequence alignments in other species. Unlike existing classifiers that employ human-engineered features, the presently described deep learning network was trained to extract features directly from the primary sequence. In certain implementations, to incorporate information about protein structure, and as described in greater detail below, separate networks were trained to predict the secondary structure and the solvent accessibility from the sequence alone, and these were included as subnetworks in the full model. Given the limited number of human proteins that have been successfully crystallized, inferring structure from the primary sequence has the advantage of avoiding biases due to incomplete protein structure and functional domain annotation. The total depth of one implementation of the network, with protein structure included, was 36 layers of convolutions, comprising roughly 400,000 trainable parameters.
B. Protein Structure 2nd Structure and Solvent Accessibility Sub-Networks
In one example of an implementation, the deep learning network for pathogenicity prediction contains 36 total convolutional layers, including 19 convolutional layers for secondary structure and solvent accessibility prediction sub-networks, and 17 for the main pathogenicity prediction network which takes as input the results of the secondary structure and solvent accessibility sub-networks. In particular, because the crystal structures of most human proteins are unknown, the secondary structure network and the solvent accessibility prediction network were trained to enable the network to learn protein structure from primary sequence.
The secondary structure and solvent accessibility prediction networks in one such implementation have identical architecture and input data, but differ on the prediction states. For example, in one such implementation, the input to the secondary structure and solvent accessibility networks is an amino acid position frequency matrix (PFM) of suitable dimensions (e.g., a 51 length×20 amino acid PFM) encoding conservation information from the multiple sequence alignment of human with 99 other vertebrates.
In one embodiment, and with reference to
Following the residual blocks 140, the softmax layer 154 computes probabilities of the three states for each amino acid, among which the largest softmax probability determines the state of the amino acid. The model in one such implementation is trained with accumulated categorical cross entropy loss function for the whole protein sequence using the ADAM optimizer. In one implementation, which is illustrated, once the networks were pre-trained on secondary structure and solvent accessibility, instead of directly taking the output of the networks as inputs for the pathogenicity prediction network 160, the layer before the softmax layer 154 was instead taken so that more information would pass through to the pathogenicity prediction network 160. In one example, the output of the layer before the softmax layer 154 is an amino acid sequence of suitable length (e.g., 51 amino acids in length) and becomes the input for the deep learning network for pathogenicity classification.
With the preceding in mind, the secondary structure network is trained to predict a 3-state secondary structure: (1) alpha helix (H), (2) beta sheet (B), or (3) coils (C). The solvent accessibility network is trained to predict a 3-state solvent accessibility: (1) buried (B), (2) intermediate (I), or (3) exposed (E). As noted above, both networks only take primary sequence as their inputs 150, and were trained using labels from known crystal structures in the Protein DataBank. Each model models predict one respective state for each amino acid residue.
With the preceding in mind, and by way of further illustration of an example implementation, for each amino acid position in the input dataset 150, a window from the position frequency matrix was taken corresponding to the flanking amino acids (e.g., flanking 51 amino acids) and this was used to predict the label for either secondary structure or solvent accessibility for the amino acid at the center of the length amino acid sequence. The labels for secondary structure and relative solvent accessibility were directly obtained from the known 3D crystal structure of the protein using the DSSP software and did not require prediction from primary sequence. To incorporate the secondary structure network and solvent accessibility networks as part of the pathogenicity prediction network 160, position frequency matrices were computed from the human-based 99 vertebrate multiple sequence alignments. Although the conservation matrices generated from these two methods are generally similar, backpropagation was enabled through the secondary structure and solvent accessibility models during training for pathogenicity prediction to allow fine-tuning of the parameter weights.
By way of example, Table 1 shows example model architecture details for a 3-state secondary structure prediction deep learning (DL) model. The shape specifies the shape of the output tensor at each layer of the model and the activation is the activation given to the neurons of the layer. Inputs to the model were the position-specific frequency matrices of suitable dimensions (e.g., 51 amino acid length, 20 depth) for the flanking amino acid sequence around the variant.
Similarly, the model architecture illustrated in Table 2 shows example model architecture details for the 3-state solvent accessibility prediction deep learning model, which as noted herein may be identical in architecture to the secondary structure prediction DL model. The shape specifies the shape of the output tensor at each layer of the model and the activation is the activation given to the neurons of the layer. Inputs to the model were the position-specific frequency matrices of suitable dimensions (e.g., 51 amino acid length, 20 depth) for the flanking amino acid sequence around the variant.
The best testing accuracy for the 3-state secondary structure prediction model was 80.32%, similar to the state-of-the-art accuracy predicted by DeepCNF model on a similar training dataset. The best testing accuracy for the 3-state solvent accessibility prediction model was 64.83%, similar to the current best accuracy predicted by RaptorX on a similar training dataset.
With reference to the Tables 1 and 2 (reproduced below) and
Each of the input channel was passed through a 1D convolution layer (layer 1a and 1b) with 40 kernels and linear activations. This layer was used to upsample the input dimensions from 20 to 40. All of the other layers throughout the model used 40 kernels. The two layers (1a and 1b) activations were merged together by summing the values across each of the 40 dimensions (i.e., merge mode=‘sum’). The output of the merge node was passed through a single layer of 1D convolution (layer 2) followed by linear activation.
The activations from layer 2 was passed through a series of 9 residual blocks (layers 3 to 11). The activations of layer 3 was fed to layer 4 and layer 4's activation was fed to layer 5, and so on. There were also skip connections that directly summed the output of every 3rd residual block (layers 5, 8, and 11). The merged activations were then fed to two 1D convolutions (layers 12 and 13) with ReLU activations. The activations from layer 13 were given to the softmax readout layer. The softmax computed the probabilities of the three-class outputs for the given input.
As a further note, in one implementation of the secondary structure model, the 1D convolutions had an atrous rate of 1. In an implementation of the solvent accessibility model, the last 3 residual blocks (layers 9, 10 and 11) had an atrous rate of 2 to increase the coverage of the kernels. With respect to these aspects, an atrous/dilated convolution is a convolution where the kernel is applied over an area larger than its length by skipping input values with a certain step, also called atrous convolution rate or dilation factor. Atrous/dilated convolutions add spacing between the elements of a convolution filter/kernel so that neighboring input entries (e.g., nucleotides, amino acids) at larger intervals are considered when a convolution operation is performed. This enables incorporation of long-range contextual dependencies in the input. The atrous convolutions conserve partial convolution calculations for reuse as adjacent nucleotides are processed. Atrous/dilated convolutions allow for large receptive fields with few trainable parameters. The secondary structure of a protein is strongly dependent on the interactions of amino acids in close proximity. Thus models with higher kernel coverage improved the performance slightly. Conversely, solvent accessibility is affected by the long-range interactions between amino acids. Therefore, for the model with high coverage of kernels using atrous convolutions, its accuracy was more than 2% higher than that of the short-coverage models.
C. Pathogenicity Prediction Network Architecture
With respect to the pathogenicity prediction model, semi-supervised deep convolutional neural network (CNN) models were developed to predict pathogenicity of variants. The input features to the models include protein sequences and conservation profiles flanking variants and depletion of missense variants in specific gene regions. The changes caused by variants to secondary structure and solvent accessibility were also predicted by deep learning models and this was integrated into the pathogenicity prediction model. Predicted pathogenicity in one such implementation is on a scale from 0 (benign) to 1 (pathogenic).
The architecture for one such pathogenicity classification neural network (e.g., PrimateAI) is described schematically in
In certain embodiments, some or all layers of the deep residual network use a ReLU activation function, which greatly accelerates the convergence of stochastic gradient descent compared to saturating nonlinearities such as sigmoid or hyperbolic tangent. Other examples of activation functions that can be used by the technology disclosed include parametric ReLU, leaky ReLU, and exponential linear unit (ELU).
As described herein, some or all layers may also employ batch normalization by which distribution of each layer in a convolution neural network (CNN) changes during training and it varies from one layer to another. This reduces the convergence speed of the optimization algorithm.
With the preceding in mind, and with reference to
For indirect inputs 1f and 1g, the pre-trained layers of the secondary structure prediction model are loaded, excluding the softmax layer. For input 1f, the pre-trained layers are based on the human reference sequence for the variants along with the PSSM generated by PSI-BLAST for the variant. Likewise, for input 1g, the pre-trained layers of the secondary structure prediction models are based on the human alternative sequence as the input along with the PSSM matrix. Inputs 1h and 1i correspond to the similar pre-trained channels containing the solvent accessibility information for reference and alternative sequences of the variant, respectively.
In this example, the five direct input channels are passed through an upsampling convolution layer of 40 kernels with linear activations. The layers 1a, 1c, 1d and 1e are merged with values summed across the 40 feature dimensions to produce layer 2a. In other words, the feature map of the reference sequence is merged with the three types of conservation feature maps. Similarly, 1b, 1c, 1d, and 1e are merged with values summed across the 40 feature dimensions to generate layer 2b, i.e., the features of the alternative sequence is merged with the three types of conservation features.
The layers 2a and 2b are batch normalized with the activation of ReLU and each passed through a 1D convolution layer of filter size 40 (3a and 3b). The outputs of layers 3a and 3b are merged with 1f, 1g, 1h, and 1i with the feature maps concatenated to each other. In other words, the feature maps of reference sequence with conservation profile, and alternative sequence with the conservation profile are merged with the secondary structure feature maps of the reference and alternative sequence and the solvent accessibility feature maps of the reference and alternative sequence (layer 4).
The outputs of layer 4 are passed through six residual blocks (layers 5, 6, 7, 8, 9, 10). The last three residual blocks have an atrous rate of 2 for the 1D convolutions, to give higher coverage for the kernels. The output of layer 10 is passed through a 1D convolution of filter size 1 and activation sigmoid (layer 11). The output of layer 11 is passed through a global maxpool that picks a single value for the variant. This value represents the pathogenicity of the variant. The details of one implementation of the pathogenicity prediction model are shown in Table 3.
D. Training (Semi-Supervised) and Data Distribution
With respect to semi-supervised learning approaches, such techniques allow utilization of both labeled and un-labeled data for training the network(s). The motivation for choosing semi-supervised learning is that human-curated variant databases are unreliable and noisy, and in particular, lack reliable pathogenic variants. Because semi-supervised learning algorithms use both labeled and unlabeled instances in the training process, they can produce classifiers that achieve better performance than completely supervised learning algorithms that have only a small amount of labeled data available for training. The principle behind semi-supervised learning is that intrinsic knowledge within unlabeled data can be leveraged in order to strengthen the prediction capability of a supervised model that only uses labeled instances, thereby providing a potential advantage for semi-supervised learning. Model parameters learned by a supervised classifier from a small amount of labeled data may be steered towards a more realistic distribution (which more closely resembles the distribution of the test data) by the unlabeled data.
Another challenge that is prevalent in bioinformatics is the data imbalance problem. The data imbalance phenomenon arises when one of the classes to be predicted is underrepresented in the data because instances belonging to that class are rare (noteworthy cases) or hard to obtain. Minority classes are typically the most important to learn because they may be associated with special cases.
An algorithmic approach to handle imbalanced data distributions is based on ensembles of classifiers. Limited amounts of labeled data naturally lead to weaker classifiers, but ensembles of weak classifiers tend to surpass the performance of any single constituent classifier. Moreover, ensembles typically improve the prediction accuracy obtained from a single classifier by a factor that validates the effort and cost associated with learning multiple models. Intuitively, aggregating several classifiers leads to better overfitting control, since averaging the high variability of individual classifiers also averages the classifiers' overfitting.
While the preceding relates to the training and validation of a pathogenicity classifier implemented as a neural network, the following sections relate various implementation-specific and use-case scenarios for further refining and/or utilizing pathogenicity classification using such a network. In a first aspect, discussion of scoring thresholds and the use of such score thresholds is described.
As noted herein, a presently disclosed pathogenicity classification network, such as (but not limited to) the PrimateAI or pAI classifier described herein, may be used to generate a pathogenicity score useful for distinguishing or screening pathogenic variants from benign variants within genes. Because pathogenicity scoring as described herein is based on the extent of purifying selection in humans and non-human primates, pathogenicity scores associated with pathogenic and benign variants are expected to be higher in the genes that are under strong purifying selection. On the other hand, for genes under neutral evolution or weak selection, the pathogenicity scores for pathogenic variants tend to be lower. This concept is illustrated visually in
To assess possible thresholds that might be useful in assessing pathogenicity scores, potential score thresholds were studied using 84 genes which contained at least 10 benign/likely benign variants and at least 10 pathogenic and likely pathogenic variants in ClinVar. These genes were used to help evaluate suitable pathogenicity score thresholds for each gene. For each of these genes, the mean pathogenicity scores were measured for benign and pathogenic variants in ClinVar.
In one implementation, it was observed that mean pathogenicity scores for pathogenic and benign ClinVar variants in each gene correlated well with the 75th and 25th percentiles of pathogenicity scores (here PrimateAI or pAI scores) in that gene, as shown graphically in
In view of the present approach, a suitable percentile of pathogenicity score per gene for use as a cutoff for likely pathogenic variants may be in the range defined by and including the 51st percentile to the 99th percentile (e.g., the 65th, 70th, 75th, 80th, or 85th percentile). Conversely, a suitable percentile of pathogenicity score per gene for use as a cutoff for likely benign variants may be in the range defined by and including the 1st percentile to the 49th percentile (e.g., the 15th, 20th, 25th, 30th, or 35th percentile).
With respect to use of these thresholds,
Clinical studies and patient care are example use-case scenarios in which a pathogenicity classification network, such as PrimateAI, may be employed to classify and/or separate pathogenic variants from benign variants within genes. In particular, clinical genome sequencing has become standard-of-care for patients with rare genetic diseases. Rare genetic diseases are often, if not primarily, caused by highly deleterious rare mutations, which are generally easier to detect due to their severity. However, those rare mutations that underlie common genetic diseases remain largely uncharacterized due to their weak effects and large numbers.
With this in mind, it may be desirable to understand the mechanism between rare mutations and common diseases and to study the evolutionary dynamics of human mutations, particularly in the context of pathogenicity scoring of variants as discussed herein. During the evolution of human population, new variants have been constantly generated by de novo mutation while some of them have also been removed due to natural selection. If the human population size was constant, the allele frequencies of variants affected by the two forces would ultimately reach equilibrium. With this in mind, it may be desirable to use observed allele frequencies to determine the severity of nature selection on any variants.
However, the human population is not in the steady state at any moment, and has been instead exponentially expanding since the advent of agriculture. Therefore, in accordance with certain approaches discussed herein, forward time simulation may be used as a tool to investigate the effects of the two forces on the distribution of allele frequencies of variants. Aspects of this approach are described with respect to the steps illustrated in
With this in mind, a forward time simulation of neutral evolution using de novo mutation rates 280 may be employed as part of modeling (step 282) the distribution of allele frequencies of variants over time. As a baseline, a forward-time population model may be simulated assuming neutral evolution. Model parameters 300 were derived by fitting (step 302) simulated allele frequency spectra (AFS) 304 to that of the observed synonymous mutations (synonymous AFS 308) in the human genome. A simulated AFS 304 generated using the set of optimal model parameters 300 (i.e., those parameters corresponding to the best fit) may then be used in conjunction with other concepts as discussed herein, such as variant pathogenicity scoring, to derive useful clinical information.
As the distribution of rare variants is of primary interest, the evolutionary history of the human population in this example implementation is simplified into four stages of exponential expansion with different growth rates in this simulation, as shown in
In this example, a long burn-in period (about 3,500 generations) with a small change in effective population size was employed in the first phase. The population size change may be denoted as n. As the time after the burn-in is unknown, this time may be denoted as T1 and the effective population size at T1 as 10,000*n. The growth rate 284 during the burn-in is g1=n{circumflex over ( )}(1/3,500).
At 1400 years AD, the census population size across the world is estimated to have been about 360 million. At 1700 years AD, the census population size increased to about 620 million, and at 2000 years AD, it is 6.2 billion. Based on these estimates, the growth rates 284 at each stage can be derived as shown in Table 4:
For generation j 286, Nj chromosomes are randomly sampled from the previous generation to form the new-generation population, where Nj=gj*Nj-1, where gj is the growth rate 284 at generation j. Most of mutations are inherited from the previous generation during chromosome sampling. De novo mutations are then applied to these chromosomes according to de novo mutation rates (μ) 280.
With respect to the de novo mutation rates 280, in accordance with certain implementations, these may be derived in accordance with the following approach or an equivalent approach. In particular, in one such implementation three large parental-offspring trio datasets totaling 8,071 trios were obtained with whole genome sequencing from literature sources (Halldorsson set (2976 trios), Goldmann set (1291 trios), and Sanders set (3804 trios)). Merging these 8,071 trios, the de novo mutations were obtained that are mapped to the intergenic regions and the de novo mutation rates 280 derived for each of 192 trinucleotide context configurations.
These estimates of mutation rates were compared with other literature mutation rates (Kaitlin's mutation rates derived from intergenic regions of 1,000 genome projects), as shown in
With respect to the mutation rates at CpG islands, methylation levels at CpG sites affect the mutation rates substantially. To precisely compute the CpGTi mutation rates, the methylation levels at those sites should be taken into account. With this in mind, in an example implementation the mutation rate and CpG islands may be calculated in accordance with the following approach.
First, the impact of methylation levels on CpG mutation rates was evaluated by using whole genome bisulfite sequencing data (obtained from Roadmap Epigenomics project). The methylation data for each CpG island was extracted and averaged across 10 embryonic stem cell (ESC) samples. Those CpG islands were then separated into 10 bins based on 10 defined methylation levels, as shown in
CpG sites were classified into two types: (1) high methylation (if averaged methylation levels >0.5); and (2) low methylation (if averaged methylation levels ≤0.5). The de novo mutation rates for each of the 8 CpGTi tri-nucleotide contexts were computed for high and low methylation levels, respectively. Averaging across the 8 CpGTi tri-nucleotide contexts, the CpGTi mutation rates were obtained: 1.01e-07 for high methylation and 2.264e-08 for low methylation, as shown in Table 6.
The allele frequency spectrum (AFS) of exome sequencing data was then fitted. In one such sample implementation, simulations were performed assuming 100,000 independent sites and using various combinations of parameters, T1, r and n, where T1ϵ(330, 350, 370, 400, 430, 450, 470, 500, 530, 550), rϵ(20, 25, 30, . . . , 100, 105, 110), and nϵ(1.0, 2.0, 3.0, 4.0, 5.0) were considered.
Each of the three major classes of mutations were simulated separately using the different de novo mutation rates, i.e. CpGTi, non-CpGTi and Tv (as shown in Table 6). For CpGTi, high and low methylation levels were separately simulated and the two AFSs merged by applying the proportions of high or low methylation sites as weights.
For each combination of parameters as well as each mutation rate, a human population was simulated through present day. One thousand sets of 246,000 chromosomes were then randomly sampled (step 288) (e.g., from a target or final generation 290), corresponding to the sample size of gnomAD exomes. A simulated AFS 304 was then generated by averaging (step 292) across the one thousand respective sampled sets 294.
In terms of validation, human exome polymorphism data from genome Aggregation Database (gnomAD) v2.1.1 was acquired, which collected the whole-exome sequencing (WES) data of 123,136 individuals from eight subpopulations around the world (http://gnomad.broadinstitute.org/). Variants were excluded that did not pass filters, have median coverage <15, or fall within low complexity regions or segmental duplication regions, where the boundary of regions were defined in the file downloaded from https://storage.googleapis.com/gnomad-public/release/2.0.2/README.txt.). Variants were retained that were mapped to the canonical coding sequences defined by the UCSC genome browser for the hg19 build.
The synonymous allele frequency spectrum 308 of gnomAD was generated (step 306) by counting the number of synonymous variants in seven allele frequency categories, including singletons, doubletons, 3≤allele count (AC)≤4, . . . , and 33≤AC≤64. Variants with AC >64 were discarded as the focus was on rare variants.
A likelihood ratio test was applied to evaluate (step 302) the fit of the simulated AFS 304 to the gnomAD AFS (i.e., synonymous AFS 308) of rare synonymous variants across the three mutational classes. The heatmaps of Pearson's chi-squared statistic (corresponding to −2*log likelihood ratio), as shown in
In one example implementation, and with reference to
To improve the precision of the simulation, separate simulation was applied (step 282) for each of the 192 tri-nucleotide contexts using their specific de novo mutation rates 230 derived from 8071 trios (i.e., parent-offspring trios). Under each value of selection coefficient 320 and each mutation rate 280, a human population with the initial size of ˜20,000 chromosomes was simulated expanding to present day. One thousand sets 294 were randomly sampled (step 288) from the resulting population (i.e., target or final generation 290). Each set contained ˜500,000 chromosomes, corresponding to the sample size of gnomAD+Topmed+UK biobank. For each of the 8 CpGTi tri-nucleotide contexts, high and low methylation levels were simulated separately. The two AFSs were merged by applying the proportions of high or low methylation sites as weights.
After the AFSs for 192 tri-nucleotide contexts were obtained, these AFS were weighted by the frequency of the 192 tri-nucleotide contexts in the exome to generate the AFS of exome. This procedure was repeated for each of the 36 selection coefficients (step 306).
A selection-depletion curve 312 was then derived. In particular, as the selective pressure (s) against mutations rises, variants are expected to be progressively depleted. With the simulated AFS (304) under various levels of selection, a metric characterized as “depletion” was defined to measure the proportion of variants wiped out by purifying selection compared to the scenario under neutral evolution (i.e., without selection). In this example, depletion may be characterized as:
Depletion values 316 were generated for each of the 36 selection coefficients (step 314) to draw the selection-depletion curve 312 shown in
With the preceding discussion regarding selection and depletion characterization using forward time simulation in mind, these factors may be used to estimate selection coefficients for missense variants based on pathogenicity (e.g., PrimateAI or pAI) scores.
By way of example, and with reference to
Variants were merged among the three datasets to form the final dataset 352 to compute depletion metrics. Allele frequencies were recomputed to average across the three datasets in accordance with:
where i represents the dataset index. If a variant didn't appear in one dataset, a zero (0) was assigned to AC for that dataset.
With respect to depletion in nonsense mutations and splice variants, the fraction of PTVs depleted by purifying selection compared with the expected number in each gene may be calculated. Due to the difficulty in computing the expected number of frameshifts in a gene, focus was instead placed on nonsense mutations and splice variants, denoted as loss-of-function mutations (LOF).
The number of nonsense mutations and splice variants were counted (step 356) in each gene of the merged dataset to get the observed number of LOFs 360 (the numerator of the formula below). To calculate (step 362) the expected number of LOFs 364, the file containing constraint metrics (https://storage.googleapis.com/gnomad-public/release/2.1.1/constraint/gnomad.v2.1.1.1of_metrics.by_gene.text.bgz) was downloaded from the gnomAD database website. The observed synonymous variants in the merged dataset was used as baseline and was converted to the expect number of LOFs 364 by multiplying the ratio of the expected number of LOFs and the expected number of synonymous variants from gnomAD. The depletion metric 380 was then calculated (step 378) and verified to be within [0,1]. If less than 0, a 0 is assigned and vice versa. The above may be represented as:
Based on the depletion metric 380 of LOFs, the estimates of selection coefficients of LOFs 390 can be derived (step 388) for each gene using the respective selection-depletion curve 312.
With respect to the calculation (step 378) of depletion 380 in missense variants, in one example (illustrated in
For each gene, the pathogenicity score percentiles 420 (pAI percentiles in the present example) were divided (step 424) into 10 bins (e.g., (0.0, 0.1], (0.1, 0.2], . . . , (0.9, 1.0]), and the number 428 of the observed missense variants falling in each bin were counted (step 426). The depletion metric 380 is calculated 430 similar to that of LOFs, except that a respective depletion metric 380 was computed for each of the 10 bins. The correction factor of missense/synonymous variants from gnomAD was applied, similar as to that used in the calculation of LOF depletion described herein, to the expected number of missense variants in each bin. The above may be represented as:
Based on the depletion of 10 bins within each gene, the relationship 436 between percentiles 420 of pathogenicity scores and the depletion metric(s) 380 may be derived (step 434). In one example, the median percentile of each bin was determined and smooth splines fitted to the median points of the 10 bins. An example of this is shown in
Based on this methodology, for every possible missense variant, its depletion metric 380 can be predicted based on its pathogenicity percentile score 420 using a gene-specific fitted spline. The selection coefficient 320 of this missense variant can then be estimated using the selection-depletion relationship (e.g., selection-depletion curve 312 or other fitted function).
In addition, the expected number of deleterious rare missense variants and PTVs can be estimated individually. For example, it may be of interest to estimate, on average, the number of deleterious rare variants a normal individual might carry in their coding genome. For this example of an implementation, focus was on rare variants with AF<0.01%. To calculate the expected number of deleterious rare PTVs per individual, it is equivalent to summing up the allele frequencies of PTVs with selection coefficients 320 exceeding a certain threshold, as shown by:
E[#PTV|s>k]=ΣAF|s>k (5)
As PTVs include nonsense mutations, splice variants, and frameshifts, computations were made separately for each category. From the results shown in Table 6 below, it may be observed that each individual has about 1.9 rare PTVs with s>0.01 (equal to or worse than a BRCA1 mutation).
In addition, the expected number of deleterious rare missense variants were computed by summing up allele frequencies of rare missenses with selection coefficients 320 exceeding different thresholds:
E[#missense|s>k]=ΣAF|s>k (6)
From the results shown in Table 7 below, there are about 4 times as many missense variants as PTVs with s>0.01.
To promote the adoption and usage of pathogenicity scores of missense variants in clinical settings, the relationship between pathogenicity scores and clinical disease prevalence among those genes of clinical interest was investigated. In particular, methodologies were developed for estimating genetic disease prevalence using various metrics based on pathogenicity scores. Non-limiting examples of two such methodologies based on pathogenicity scores to predict genetic disease prevalence are described herein.
By way of preliminary context in terms of the data employed in this study, the DiscovEHR data referenced herein is a collaboration between the Regeneron Genetics Center and Geisinger Health System that aims to catalyze precision medicine by integrating whole-exome sequencing with clinical phenotypes from longitudinal electronic health records (EHRs) of 50,726 adult participants in Geisinger's MyCode Community Health Initiative. With this in mind, and with reference to
With this context in mind, two examples of approaches based on pathogenicity scores 206 (e.g., PrimateAI or pAI scores) for predicting genetic disease prevalence were developed. In these methodologies, and with reference to
With respect to calculating this genetic disease prevalence metric for each gene, depicted as step 478 of
For each trinucleotide context 500, the forward-time simulation, as described herein, is performed (step 502) to generate the expected (i.e., simulated) allele frequency spectrum (AFS) 304, assuming a selection coefficient 320 equal to 0.01 and using the de novo mutation rates 280 for that trinucleotide context. In one implementation of this methodology the simulation simulated 100,000 independent sites among ˜400K chromosomes (˜200K samples). Thus, the expected AFS 304 for a specific trinucleotide context 500 in such a context is the simulated AFS/100,000*the occurrence of the trinucleotide in the list of deleterious variants. Summing up across the 192 trinucleotides produce the expected AFS 304 for the gene. The genetic disease prevalence metric of a specific gene (i.e., expected CAF 480) in accordance with this approach is defined as the sum (step 506) of simulated rare allele frequencies (i.e., AF≤0.001) in the expected AFS 304 for the gene.
The genetic disease prevalence metric as derived in accordance with the second methodology is similar to that derived using the first methodology, but differs in the definition of the list of deleterious missense variants. In accordance with the second methodology, and with reference to
With respect to the genetic disease prevalence metric as derived using this second methodology,
As described herein, pathogenicity scores generated in accordance with the present teachings are derived using neural networks trained primarily based on DNA flanking sequences around variants, conservation among species, and protein secondary structures. However, the variation associated with pathogenicity scores (e.g., PrimateAI scores) can be large (e.g., about 0.15). In addition, certain implementations of the generalized model discussed herein for calculating pathogenicity scores do not utilize the information of observed allele frequencies in human population during training. In certain circumstances, some variants with high pathogenicity scores might appear to have allele counts >1, implying that there is a need to penalize these pathogenicity scores based on allele counts. With this in mind, it may be useful to recalibrate pathogenicity scores to address such circumstances. In one example embodiment discussed herein, a recalibration approach may focus on the percentiles of pathogenicity scores of variants, as these may be more robust and less affected by selection pressure exerted on the whole genes.
With this in mind, and with reference to
obsAI˜trueAI+e,e˜N(0,sd=0.15) (7)
The distribution 554 of observed pathogenicity score percentiles 550 in this context is the discrete uniform distribution of true pathogenicity score percentiles overlaid with Gaussian noise, as shown in
To visualize the probability that a variant with a true pathogenicity score percentile (x-axis in
Turning to
This depletion metric 522 measures the probability that a variant falling within each bin might be removed by purifying selection. An example of this is shown in
In accordance with one implementation, a smooth spline was fitted (step 584) to the survival probability vs. the median pathogenicity score percentile of each bin across the bins (e.g., 10 bins) and produced the survival probability for each percentile of pathogenicity scores. In accordance with this approach, this constitutes a survival probability correction factor 590, implying that the higher the percentile of pathogenicity scores 206, the less chance the variant survives purifying selection. In other implementations, other techniques may be employed, such as interpolation, in place of fitting a smooth spline. High pathogenicity scores 206 of those observed variants may then be penalized or corrected in accordance with this correction factor 590.
With the preceding in mind, and turning to
The true pathogenicity score percentiles were divided (step 604) into bins (i.e., merge 1%-10% (the first 10 columns of the recalibrated heatmap 596) into a first bin), 11%-20% (the next 10 columns of the recalibrated heatmap 596) into second bin, . . . , etc.), which represent the probability a variant might come from each of the true pathogenicity score percentile bins. For variants of that gene having an observed pathogenicity score percentile (e.g. x %, corresponding to the x-th row of the recalibrated heatmap 596), the probability that this variant might fall within each of the true pathogenicity score percentile bins (e.g., 10 bins) can be obtained (step 608). This may be denoted as the variant contribution 612 to each bin.
In this example, the expected number of missense variants 620 within each of the bins (e.g., 10 bins) (derived at step 624) is the sum of the variant contributions to that bin across all the observed missense variants in the respective gene. Based on this expected number of missense variants 620 falling within each of the bins for a gene, the depletion formula for missense variants discussed herein can be used to calculate (step 630) the corrected depletion metric 634 for each missense bin. This is illustrated in
In this manner of recalibrating the pathogenicity scores 206, the noise distribution in percentiles of pathogenicity scores 206 is modeled and the noise in predicted depletion metrics 522 reduced. This may help mitigate the effects of noise on the estimation of selection coefficients 320 in missense variants as discussed herein.
In the preceding discussion, various aspects of neural network architecture and use are referenced in the context of a pathogenicity classification or scoring network. While extensive knowledge of these various aspects of neural network design and use are not believed to be necessary for those wishing to understand and employ a pathogenicity classification network as discussed herein, for the benefit of those wishing additional detail, the following neural network primer is provided by way of additional reference.
With this in mind, the term “neural network” in a general sense may be understood to be a computational construct that is trained to receive a respective output and, based upon its training generate an output, such as a pathogenicity score, in which the input is modified, classified, or otherwise processed. Such a construct may be referred to as a neural network due to being modeled after a biological brain, with different nodes of the construct being equated to “neurons”, which may be interconnected with a wide range of other nodes so as to allow for complex potential interconnections between nodes. In general, a neural network may be considered a form of machine learning, as the pathways and associated nodes are typically trained by way of example (e.g., using sample data where the inputs and outputs are known or where a cost function can be optimized) and may learn or evolve over time as the neural network is used and its performance or outputs are modified or retrained.
With this in mind and by way of further illustration,
Neural networks 700 can be categorized into different types based on their mode of operation. By way of example, a convolutional neural network is one type of neural network that employs or incorporates one or more convolution layers, as opposed to dense or densely connected layers. In particular, densely connected layers learn global patterns in their input feature space. Conversely, convolution layers learn local patterns. By way of example, in the case of images, convolution layers may learn patterns found in small windows or subsets of the inputs. This focus on local patterns or features gives convolutional neural networks two useful properties: (1) the patterns they learn are translation invariant, and (2) they can learn spatial hierarchies of patterns.
With regard to this first property, after learning a certain pattern in one portion or subset of a data set, a convolution layer can recognize the pattern in other portions of the same or different data sets. In contrast, a densely connected network would have to learn the pattern anew if it were present elsewhere (e.g., at a new location). This property makes convolutional neural networks data efficient because they need fewer training samples to learn representations that can then be generalized so as to be identified in other contexts and locations.
Regarding the second property, a first convolution layer can learn small local patterns while a second convolution layer learns larger patterns made of the features of the first layers, and so on. This allows convolutional neural networks to efficiently learn increasingly complex and abstract visual concepts.
With this in mind, a convolutional neural network is capable of learning highly non-linear mappings by interconnecting layers of artificial neurons 704 arranged in many different layers 702 with activation functions that make the layers dependent. It includes one or more convolutional layers, interspersed with one or more sub-sampling layers and non-linear layers, which are typically followed by one or more fully connected layers. Each element of the convolutional neural network receives inputs from a set of features in the previous layer. The convolutional neural network learns concurrently because the neurons in the same feature map have identical weights. These local shared weights reduce the complexity of the network such that when multi-dimensional input data enters the network, the convolutional neural network avoids the complexity of data reconstruction in feature extraction and regression or classification process.
Convolutions operate over 3D tensors, called feature maps, with two spatial axes (height and width) as well as a depth axis (also called the channels axis). The convolution operation extracts patches from its input feature map and applies the same transformation to all of these patches, producing an output feature map. This output feature map is still a 3D tensor: it has a width and a height. Its depth can be arbitrary because the output depth is a parameter of the layer, and the different channels in that depth axis stand for filters. Filters encode specific aspects of the input data.
For instance in an example where a first convolution layer takes a feature map of a given size (28, 28, 1) and outputs a feature map of size (26, 26, 32): it computes 32 filters over its input. Each of these 32 output channels contains a 26×26 grid of values, which is a response map of the filter over the input, indicating the response of that filter pattern at different locations in the input. That is what the term feature map means in this context: every dimension in the depth axis is a feature (or filter), and the 2D tensor output [:, :, n] is the 2D spatial map of the response of this filter over the input.
With the preceding in mind, convolutions are defined by two key parameters: (1) the size of the patches extracted from the inputs, and (2) the depth of the output feature map (i.e., the number of filters computed by the convolution). In typical implementations these start with a depth of 32, continue to a depth of 64, and terminate with a depth of 128 or 256, though certain implementations may vary from this progression.
Turning to
With the preceding in mind, a convolutional neural network comprises convolution layers which perform the convolution operation between the input values and convolution filters (matrix of weights) that are learned over multiple gradient update iterations during a training process. Where (m, n) is the filter size and W is the matrix of weights, a convolution layer performs a convolution of W with the input X by calculating the dot product W·x+b, where x is an instance of X and b is the bias. The step size by which the convolution filters slide across the input is called the stride, and the filter area (m×n) is called the receptive field. A same convolution filter is applied across different positions of the input, which reduces the number of weights learned. It also allows location invariant learning, i.e., if an important pattern exists in the input, the convolution filters learn it no matter where it is in the sequence.
As may be understood from the preceding discussion, the training of a convolutional neural network is an important aspect of the network performing a given task of interest. The convolutional neural network is adjusted or trained so that the input data leads to a specific output estimate. The convolutional neural network is adjusted using back propagation based on a comparison of the output estimate and the ground truth until the output estimate progressively matches or approaches the ground truth.
The convolutional neural network is trained by adjusting the weights between the neurons based on the difference (i.e., error, δ) between the ground truth and the actual output. The intermediary step in the training process includes generating a feature vector from the input data using the convolution layers, as described herein. The gradient with respect to the weights in each layer, starting at the output, is calculated. This is referred to as the backward pass or going backwards. The weights in the network are updated using a combination of the negative gradient and previous weights.
In one implementation, the convolutional neural network 150 uses a stochastic gradient update algorithm (such as ADAM) that performs backward propagation of errors by means of gradient descent. The algorithm includes computing the activation of all neurons in the network, yielding an output for the forward pass. The error and the correct weights are then calculated per layer. In one implementation, the convolutional neural network uses a gradient descent optimization to compute the error across all the layers.
In one implementation, the convolutional neural network uses a stochastic gradient descent (SGD) to calculate the cost function. A SGD approximates the gradient with respect to the weights in the loss function by computing it from only one, randomized, data pair. In other implementations, the convolutional neural network uses different loss functions such as Euclidean loss and softmax loss. In a further implementation, an Adam stochastic optimizer is used by the convolutional neural network.
The convolution layers of the convolutional neural network serve as feature extractors. In particular, convolution layers act as adaptive feature extractors capable of learning and decomposing the input data into hierarchical features. Convolution operations typically involve a “kernel”, applied as a filter on the input data, producing an output data.
The convolution operation includes sliding (e.g., incrementally moving) the kernel over the input data. For each position of the kernel, the overlapping values of the kernel and the input data are multiplied and the results are added. The sum of products is the value of the output data at the point in the input data where the kernel is centered. The resulting different outputs from many kernels are called feature maps.
Once the convolutional layers are trained, they are applied to perform recognition tasks on new inference data. Since the convolutional layers learn from the training data, they avoid explicit feature extraction and implicitly learn from the training data. Convolution layers use convolution filter kernel weights, which are determined and updated as part of the training process. The convolution layers extract different features of the input, which are combined at higher layers. The convolutional neural network uses a various number of convolution layers, each with different convolving parameters such as kernel size, strides, padding, number of feature maps and weights.
Further aspects of implementation of a convolution neural network may involve sub-sampling of layers. In this context, sub-sampling layers reduces the resolution of the features extracted by the convolution layers to make the extracted features or feature maps robust against noise and distortion. In one implementation, sub-sampling layers employs two types of pooling operations, average pooling and max pooling. The pooling operations divide the input into non-overlapping spaces or regions. For average pooling, the average of the values in the region is calculated. For max pooling, the maximum value of the values is selected.
In one implementation, the sub-sampling layers include pooling operations on a set of neurons in the previous layer by mapping its output to only one of the inputs in max pooling and by mapping its output to the average of the input in average pooling. In max pooling, the output of the pooling neuron is the maximum value that resides within the input. In average pooling, the output of the pooling neuron is the average value of the input values that reside with the input neuron set.
A further aspect of a neural network implementation relevant to the present concepts is the use of non-linear layers. Non-linear layers use different non-linear trigger functions to signal distinct identification of likely features on each hidden layer. Non-linear layers use a variety of specific functions to implement the non-linear triggering, including, but not limited to, the rectified linear units (ReLUs), hyperbolic tangent, absolute of hyperbolic tangent, sigmoid and continuous trigger (non-linear) functions. In one implementation, a ReLU activation implements the function y=max(x, 0) and keeps the input and output sizes of a layer the same. One potential advantage of using ReLU is that the convolutional neural network may be trained many times faster. ReLU is a non-continuous, non-saturating activation function that is linear with respect to the input if the input values are larger than zero and zero otherwise.
In other implementations, the convolutional neural network may use a power unit activation function, which is a continuous, non-saturating function. The power activation function is able to yield x and y-antisymmetric activation if c is odd and y-symmetric activation if c is even. In some implementations, the unit yields a non-rectified linear activation.
In yet other implementations, the convolutional neural network may use a sigmoid unit activation function, which is a continuous, saturating function. The sigmoid unit activation function does not yield negative activation and is only antisymmetric with respect to the y-axis.
A further feature of convolutional neural networks is the use of residual connections that reinject prior information downstream via feature-map addition, as illustrated in
Another concept present in convolutional neural networks relevant to the present techniques and approaches is the use of skip connections. The principle behind residual learning is that the residual mapping is easier to be learned than the original mapping. Residual networks stack a number of residual units to alleviate the degradation of training accuracy. Residual blocks make use of special additive skip connections to combat vanishing gradients in deep neural networks. At the beginning of a residual block, the data flow is separated into two streams: (1) the first carries the unchanged input of the block, (2) the second applies weights and non-linearities. At the end of the block, the two streams are merged using an element-wise sum. One advantage of such constructs is to allow the gradient to flow through the network more easily.
With the benefit of such residual networks, deep convolutional neural networks (CNNs) can be easily trained and improved accuracy can be achieved for data classification, object detection, and so forth. Convolutional feed-forward networks connect the output of the lth layer as input to the (l+1)th layer. Residual blocks add a skip-connection that bypasses the non-linear transformations with an identify function. An advantage of residual blocks is that the gradient can flow directly through the identity function from later layers to the earlier layers.
An additional aspect related to implementations of convolutional neural networks that may be applicable to the present pathogenicity classification approaches is batch normalization, which is a method for accelerating deep network training by making data standardization an integral part of the network architecture. Batch normalization can adaptively normalize data even as the mean and variance change over time during training and works by internally maintaining an exponential moving average of the batch-wise mean and variance of the data seen during training. One effect of batch normalization is that it helps with gradient propagation, much like residual connections, and thus facilitates the use of deep networks.
Batch normalization can therefore be seen as yet another layer that can be inserted into the model architecture, just like the fully connected or convolutional layer. The batch normalization layer may typically be used after a convolutional or densely connected layer, though it can also be used before a convolutional or densely connected layer.
Batch normalization provides a definition for feed-forwarding the input and computing the gradients with respect to the parameters and its own input via a backward pass. In practice, batch normalization layers are typically inserted after a convolutional or fully connected layer, but before the outputs are fed into an activation function. For convolutional layers, the different elements of the same feature map (i.e. the activations) at different locations are normalized in the same way in order to obey the convolutional property. Thus, all activations in a mini-batch are normalized over all locations, rather than per activation.
1D Convolution
A further technique used in implementations of convolutional neural networks that may be applicable to the present approach relates to the use of 1D convolutions to extract local 1D patches or subsequences from sequences. A 1D convolution approach obtains each output step from a window or patch in the input sequence. 1D convolution layers recognize local patterns in a sequence. Because the same input transformation is performed on every patch, a pattern learned at a certain position in the input sequences can be later recognized at a different position, making 1D convolution layers translation invariant for translations. For instance, a 1D convolution layer processing sequences of bases using convolution windows of size 5 should be able to learn bases or base sequences of length 5 or less, and it should be able to recognize the base motifs in any context in an input sequence. A base-level 1D convolution is thus able to learn about base morphology.
Another aspect of convolutional neural networks that may be useful or utilized in the present contexts relates to global average pooling. In particular, global average pooling can be use used to replace fully connected (FC) layers for classification by taking the spatial average of features in the last layer for scoring. This reduces the training load and bypasses overfitting issues. Global average pooling applies a structural prior to the model and it is equivalent to linear transformation with predefined weights. Global average pooling reduces the number of parameters and eliminates the fully connected layer. Fully connected layers are typically the most parameter and connection intensive layers, and global average pooling provides much lower-cost approach to achieve similar results. The main idea of global average pooling is to generate the average value from each last layer feature map as the confidence factor for scoring, feeding directly into the softmax layer.
Global average pooling may provide certain benefits including, but not limited to: (1) there are no extra parameters in global average pooling layers thus overfitting is avoided at global average pooling layers; (2) since the output of global average pooling is the average of the whole feature map, global average pooling is robust to spatial translations; and (3) because of the huge number of parameters in fully connected layers which usually take over 50% in all the parameters of the whole network, replacing them by global average pooling layers can significantly reduce the size of the model, and this makes global average pooling very useful in model compression.
Global average pooling makes sense, since stronger features in the last layer are expected to have a higher average value. In some implementations, global average pooling can be used as a proxy for the classification score. The feature maps under global average pooling can be interpreted as confidence maps, and force correspondence between the feature maps and the categories. Global average pooling can be particularly effective if the last layer features are at a sufficient abstraction for direct classification; however, global average pooling alone may not be sufficient or suitable if multilevel features should be combined into groups like parts models, which may be more suitably addressed by adding a simple fully connected layer or other classifier after the global average pooling.
As may be appreciated, the neural network aspects of the present discussion, as well as the analytics and processing performed on the pathogenicity classifier output by the described neural network may be implemented on a computer system or systems. With this in mind, and by way of further context,
With respect to the production environment, and as shown in
With respect to the actual hardware architecture that may be employed to run or support the training servers 802, production servers 804, administrative interface 816, and/or client interface(s) 826, such hardware may be physically embodied as one or more computer systems (e.g., servers, workstations, and so forth). Examples of components which may be found in such a computer system 850 are illustrated in
With this in mind
In one implementation in which the computer system 850 is used to implement or train a pathogenicity classifier, the neural networks 102 such as benign dataset generator 812, variant pathogenicity classifier 160, secondary structure classifier 130, solvent accessibility classifier 132, and semi-supervised learner 814 are communicably linked to the storage subsystem 862 and user interface input devices 882.
In the depicted example, and in the context in which the computer system 850 is used to implement or train a neural network as discussed herein, one or more deep learning processors 894 may be present as part of the computer system 850 or otherwise in communication with the computer system 850. In such an embodiment, the deep learning processors can be GPUs or FPGAs and can be hosted by a deep learning cloud platforms such as Google Cloud Platform, Xilinx, and Cirrascale. Examples of deep learning processors include Google's Tensor Processing Unit (TPU), rackmount solutions like GX4 Rackmount Series, GX8 Rackmount Series, NVIDIA DGX-1, Microsoft' Stratix V FPGA, Graphcore's Intelligent Processor Unit (IPU), Qualcomm's Zeroth platform with Snapdragon processors, NVIDIA's Volta, NVIDIA's DRIVE PX, NVIDIA's JETSON TX1/TX2 MODULE, Intel's Nirvana, Movidius VPU, Fujitsu DPI, ARM's DynamicIQ, IBM TrueNorth, and others.
In the context of the computer system 850, the user interface input devices 882 can include a keyboard; pointing devices such as a mouse, trackball, touchpad, or graphics tablet; a scanner; a touch screen incorporated into the display; audio input devices such as voice recognition systems and microphones; and other types of input devices. In general, use of the term “input device” may be construed as encompassing all possible types of devices and ways to input information into computer system 850.
User interface output devices 886 can include a display subsystem, a printer, a fax machine, or non-visual displays such as audio output devices. The display subsystem can include a cathode ray tube (CRT), a flat-panel device such as a liquid crystal display (LCD), a projection device, or some other mechanism for creating a visible image. The display subsystem can also provide a non-visual display such as audio output devices. In general, use of the term “output device” may be construed as encompassing all possible types of devices and ways to output information from computer system 850 to the user or to another machine or computer system.
Storage subsystem 862 stores programming and data constructs that provide the functionality of some or all of the modules and methods described herein. These software modules are generally executed by a processor 854 alone or in combination with other processors 854.
Memory 866 used in the storage subsystem 862 can include a number of memories including a main random access memory (RAM) 878 for storage of instructions and data during program execution and a read only memory (ROM) 874 in which fixed instructions are stored. A file storage subsystem 870 can provide persistent storage for program and data files, and can include a hard disk drive, a floppy disk drive along with associated removable media, a CD-ROM drive, an optical drive, or removable media cartridges. The modules implementing the functionality of certain implementations can be stored by file storage subsystem 870 in the storage subsystem 862, or in other machines accessible by the processor 854.
Bus subsystem 858 provides a mechanism for letting the various components and subsystems of computer system 850 communicate with each other as intended. Although bus subsystem 858 is shown schematically as a single bus, alternative implementations of the bus subsystem 858 can use multiple busses.
The computer system 850 itself can be of varying types including a personal computer, a portable computer, a workstation, a computer terminal, a network computer, a television, a mainframe, a stand-alone server, a server farm, a widely-distributed set of loosely networked computers, or any other data processing system or user device. Due to the ever-changing nature of computers and networks, the description of computer system 850 depicted in
This written description uses examples to disclose the invention, including the best mode, and also to enable any person skilled in the art to practice the invention, including making and using any devices or systems and performing any incorporated methods. The patentable scope of the invention is defined by the claims, and may include other examples that occur to those skilled in the art. Such other examples are intended to be within the scope of the claims if they have structural elements that do not differ from the literal language of the claims, or if they include equivalent structural elements with insubstantial differences from the literal languages of the claims.
This application claims priority to U.S. Provisional Patent Application No. 63/055,724, filed on Jul. 23, 2020, which is herein incorporated by reference in its entirety for all purposes.
| Number | Date | Country | |
|---|---|---|---|
| 63055724 | Jul 2020 | US |
