Patent Application
20190080051
Microarray techniques, and more recently, RNA sequencing, have fundamentally changed abilities to explore molecular mechanisms underlying complex diseases and are routinely used to identify disease-associated genome-wide changes in gene expression patterns. A common expectation is that identification of differentially expressed (DE) genes can help pinpoint the molecular processes perturbed in a disease, which, in turn, can be used as biomarkers for diagnosis and prognosis. However, conventional methods for differential expression analysis only allow for identification of average changes between two groups, not for identification of specific changes in a single subject.
Methods and systems disclosed herein address the foregoing problems and allow for a systematic quantification of the heterogeneity of disease states between different subjects on a molecular (e.g., gene or protein expression) level.
In one aspect, the invention features a method of generating a disease profile, the method involving detecting differential levels of one or more analytes (e.g., mRNA, a methylated nucleotide, protein, peptide, lipid, carbohydrate, or a metabolite) in one or more case subjects relative to the levels of the analytes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles for the case subjects; comparing the personalized perturbation profiles with a set of analytes whose differential presence is associated with said disease; and obtaining a set of overlapping analytes that defines the disease profile.
In another aspect, the invention features a method of generating a disease module, the method involving detecting differential expression of one or more genes in one or more case subjects relative to the expression levels of the genes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles for the case subjects; comparing the personalized perturbation profiles with a set of genes whose differential expression is associated with said disease; and compiling all genes that are perturbed in at least about 20% (e.g., 30%, 40%, 50%, 60%, 70%, 80%, 90% or 100%) of case subjects, thereby defining a disease module.
In another aspect, the invention features a method of classifying the disease state of a subject, the method involving detecting differential expression of one or more analytes in a subject to obtain a personalized perturbation profile for the subject; determining the fraction of analytes from a disease module that are differentially expressed in the personalized perturbation profile for the subject, thereby characterizing the disease state of the subject.
In another aspect, the invention features a method of classifying the disease state of a subject, the method involving detecting differential expression of one or more genes in a subject to obtain a personalized perturbation profile for the subject; determining the fraction of genes from a disease module that are differentially expressed in the personalized perturbation profile for the subject, thereby characterizing the disease state of the subject.
In another aspect, the invention features a method of determining whether the subject has the disease, the method involving detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the genes in one or more control subjects, thereby obtaining a personalized perturbation profile for the subject; compiling the personalized perturbation profile from each subject across a population of case subjects; comparing the compiled personalized perturbation profiles with a set of genes whose differential expression is associated with said disease and generating a statistical score; thereby determining whether the subject has the disease.
In another aspect, the invention features a computer-implemented method of generating a disease module, the method involving: (a) detecting differential levels of one or more analytes in one or more case subjects relative to the levels of the one or more analytes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more analytes for the case subjects; (b) comparing the personalized perturbation profiles with a set of one or more analytes whose differential presence is associated with said disease; and (c) obtaining a set of overlapping analytes that defines the disease module.
In one embodiment, step a. comprises: a. calculating a mean level of each analyte detected in the one or more control subjects; b. calculating, using the calculated mean level each analyte detected in the one or more control subjects, the deviation of the level of each analyte detected in the one or more test subjects; c. identifying, in the one or more test subjects, analytes with a deviation above or below a threshold deviation level from the calculated mean.
In another embodiment, the one or more analyte is a gene.
In another aspect, the invention features a computer-implemented method of classifying the disease state of a subject, the method involving: a. detecting differential levels of one or more genes in one or more case subjects relative to the levels of the one or more genes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more genes for the case subjects; b. comparing the personalized perturbation profiles with a set of one or more genes whose differential presence is associated with said disease; c. obtaining a set of overlapping genes that defines the disease module; and d. calculating a statistical score of the set of overlapping genes, and, based on the calculated score, classifying the disease state of the subject.
In one embodiment, step a. comprises: a. calculating a mean level of each gene detected in the one or more control subjects; b. calculating, using the calculated mean level each gene detected in the one or more control subjects, the deviation of the level of each gene detected in the one or more test subjects; c. identifying, in the one or more test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In another aspect, the invention features a computer-implemented method of generating a disease module, the method involving: a. detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the plurality genes in at least one control subject, thereby obtaining a personalized perturbation profile for the subject; b. compiling the personalized perturbation profile from the individual subject across a population of case subjects; c. comparing the compiled personalized perturbation profiles with a set of genes whose differential expression is associated with said disease; and d. obtaining a set of overlapping genes from the compiled perturbation profiles that defines the disease module.
In one embodiment, step a. involves: a. calculating a mean level of each gene detected in the at least one control subject; b. calculating, using the calculated mean level each gene detected in the at least one control subject, the deviation of the level of each gene detected in the test subjects; c. identifying, in the test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In one embodiment, the method further includes obtaining a set of partially overlapping genes and non-overlapping genes from the compiled perturbation profiles.
In one embodiment, the method further includes calculating, based on the overlapping genes, partially overlapping genes, and non-overlapping genes, an expression heterogeneity of the disease.
In another aspect, the invention provides specifically programmed computer system comprising:
a. at least one specialized computer machine comprising:
i. a non-transient memory, electronically storing particular computer executable program code; and
ii. at least one computer processor which, when executing the particular program code, becomes a specifically programmed computer processor configured to perform at least the following operations:
a. detecting differential levels of one or more analytes in one or more case subjects relative to the levels of the one or more analytes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more analytes for the case subjects;
b. comparing the personalized perturbation profiles with a set of one or more analytes whose differential presence is associated with said disease; and
c. obtaining a set of overlapping analytes that defines the disease module.
In one embodiment, step ii. a. comprises: a. calculating a mean level of each analyte detected in the one or more control subjects; b. calculating, using the calculated mean level each analyte detected in the one or more control subjects, the deviation of the level of each analyte detected in the one or more test subjects; c. identifying, in the one or more test subjects, analytes with a deviation above or below a threshold deviation level from the calculated mean.
In one embodiment, the one or more analyte is a gene.
In another aspect, the invention provides a specifically programmed computer system comprising:
a. at least one specialized computer machine comprising:
i. a non-transient memory, electronically storing particular computer executable program code; and
ii. at least one computer processor which, when executing the particular program code, becomes a specifically programmed computer processor configured to perform at least the following operations:
a. detecting differential levels of one or more genes in one or more case subjects relative to the levels of the one or more genes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more genes for the case subjects;
b. comparing the personalized perturbation profiles with a set of one or more genes whose differential presence is associated with said disease;
c. obtaining a set of overlapping genes that defines the disease module; and
d. calculating a statistical score of the set of overlapping genes, and, based on the calculated score, classifying the disease state of the subject.
In another embodiment, step ii. a. comprises: a. calculating a mean level of each gene detected in the one or more control subjects; b. calculating, using the calculated mean level each gene detected in the one or more control subjects, the deviation of the level of each gene detected in the one or more test subjects; c. identifying, in the one or more test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In another aspect, the invention features a specifically programmed computer system comprising:
a. at least one specialized computer machine comprising:
i. a non-transient memory, electronically storing particular computer executable program code; and
ii. at least one computer processor which, when executing the particular program code, becomes a specifically programmed computer processor configured to perform at least the following operations:
a. detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the plurality genes in at least one control subject, thereby obtaining a personalized perturbation profile for the subject;
b. compiling the personalized perturbation profile from the individual subject across a population of case subjects;
c. comparing the compiled personalized perturbation profiles with a set of genes whose differential expression is associated with said disease; and
d. obtaining a set of overlapping genes from the compiled perturbation profiles that defines the disease module.
In one embodiment, step ii. a. comprises: a. calculating a mean level of each gene detected in the at least one control subject; b. calculating, using the calculated mean level each gene detected in the at least one control subject, the deviation of the level of each gene detected in the test subjects; c. identifying, in the test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In one embodiment the specifically programmed computer system, further comprising obtaining a set of partially overlapping genes and non-overlapping genes from the compiled perturbation profiles.
In another embodiment the specifically programmed computer system, further comprising calculating, based on the overlapping genes, partially overlapping genes, and non-overlapping genes, am expression heterogeneity of the disease.
In various embodiments of the above aspect, a fraction greater than about 10%, 15%, 20%, 30%, 40%, 50%, 60%, 65%, 75%, 80%, 85%, 90%, or 95% indicates the presence of the disease in the subject. In one embodiment, the disease is a neurodegenerative disease (e.g., Parkinson's Disease or Huntington's Disease). In another embodiment, the disease is asthma. In another embodiment, the fraction defines a subset of patients within the disease module having similar personalized perturbation profiles.
Another example embodiment of the invention is a method for determining a disease state of a patient. The method includes generating personalized biomarker expression perturbation profiles for a plurality of individual subjects with a disease. The personalized biomarker expression perturbation profiles include representations of biomarkers that are perturbed beyond a threshold amount. The biomarker expression levels are associated with gene expression levels, and in some embodiments may be protein expression levels. The method further includes creating a disease module by combining representations of biomarkers from the personalized biomarker expression perturbation profiles. The disease module includes a network of representations of biomarkers having perturbations associated with the disease. The method further includes accessing biomarker data including representations of biomarker expressions for the patient from a sample obtained from the patient, and determining the disease state of the patient based on a comparison of the biomarker data and the disease module.
In some embodiments, the personalized biomarker expression perturbation profiles can be generated by comparing representations of biomarker expressions of the individual subjects with reference biomarker expression levels of a control group, and selecting for inclusion in the personalized biomarker expression perturbation profiles representations of biomarkers having expression levels exceeding corresponding biomarker expression levels of the control group by the threshold amount. Creating the disease module can include determining a number of random biomarker perturbations expected for the disease, and including a number of representations of biomarkers in the disease module that is greater than the expected number of random biomarker perturbations. Determining the disease state of the patient can include matching representations of perturbed biomarkers of the biomarker data with the representations of biomarkers of the disease module, and the method can determine that the patient has the disease if a number of representations of perturbed biomarkers of the biomarker data matching representations of biomarkers of the disease module exceeds a threshold level.
Another example embodiment of the invention is a system for determining a disease state of a patient. The system includes memory, a data source, a hardware processor in communication with the memory and the data source, and a control module in communication with the processor. The hardware processor is configured to perform a predefined set of operations in response to receiving a corresponding instruction selected from a predefined native instruction set of codes. The control module includes a first set of machine codes selected from the native instruction set for causing the hardware processor to obtain from the data source and store in the memory representations of biomarker expressions for a plurality of individual subjects with a disease. The biomarker expression levels are associated with gene expression levels, and in some embodiments may be protein expression levels. The control module further includes a second set of machine codes for causing the hardware processor to generate and store in the memory personalized biomarker expression perturbation profiles for the plurality of individual subjects. The personalized biomarker expression perturbation profiles include representations of biomarkers that are perturbed beyond a threshold amount. The control module further includes a third set of machine codes for causing the hardware processor to create and store in the memory a disease module by combining representations of biomarkers from the personalized biomarker expression perturbation profiles. The disease module includes a network of representations of biomarkers having perturbations associated with the disease. The control module further includes a fourth set of machine codes for causing the hardware processor to access from the data source biomarker data including representations of biomarker expressions for the patient from a sample obtained from the patient. The control module further includes a fifth set of machine codes for causing the hardware processor to determine the disease state of the patient based on a comparison of the biomarker data and the disease module.
The patent or application file contains at least one drawing executed in color. Copies of this patent or patent application publication with color drawing(s) will be provided by the Office upon request and payment of the necessary fee.
The foregoing will be apparent from the following more particular description of example embodiments of the invention, as illustrated in the accompanying drawings in which like reference characters refer to the same parts throughout the different views. The drawings are not necessarily to scale, emphasis instead being placed upon illustrating embodiments of the present invention.
By “alteration” is meant a change (increase or decrease) in the expression levels or activity of a gene or polypeptide as detected by standard art known methods such as those described herein. As used herein, an alteration includes a 10% change in expression levels, such as a 25% change, a 40% change, or a 50% or greater change in expression levels.
By “analyte” is meant a substance that is the subject of an analytical method. Exemplary analytes include proteins, polynucleotides (e.g., RNA, DNA, methylated DNA, and other modified polynucleotides), metabolites, carbohydrate, and lipids.
By “biologic sample” is meant any tissue, cell, fluid, or other material derived from an organism.
By “case subject” is meant a subject identified as having a disease.
By “control subject” is meant a healthy subject that does not have the disease.
In this disclosure, “comprises,” “comprising,” “containing” and “having” and the like can have the meaning ascribed to them in U.S. Patent law and can mean “includes,” “including,” and the like; “consisting essentially of” or “consists essentially” likewise has the meaning ascribed in U.S. Patent law and the term is open-ended, allowing for the presence of more than that which is recited so long as basic or novel characteristics of that which is recited is not changed by the presence of more than that which is recited, but excludes prior art embodiments.
“Detect” refers to identifying the presence, absence or amount of the object to be detected.
A “detectable” expression level, as used herein, means a level that is detectable by standard techniques currently known in the art or those that become standard at some future time, and include for example, differential display, RT (reverse transcriptase)-coupled polymerase chain reaction (PCR), Northern Blot, and/or RNase protection analyses. The degree of differences in expression levels need only be large enough to be visualized or measured via standard characterization techniques.
By “differential expression” is meant that expression is altered relative to a reference. In one embodiment, the alteration is significant when evaluated using a statistical method. In one embodiment, the alteration is increased or decreased relative to a threshold.
By “disease” is meant any condition or disorder that damages or interferes with the normal function of a cell, tissue, or organ. Examples of diseases include asthma, Parkinson's Disease, or Huntington's Disease.
By “disease state” is meant the presence, absence, or extent of disease in a subject.
By “disease module” is meant a pool of genes whose differential expression is associated with a disease.
By “disease profile” is meant a set of alterations in the level of an analyte that is associated with a disease state.
The term “expression” refers to the biosynthesis of a gene product. For example, in the case of a structural gene, expression involves transcription of the structural gene into mRNA and the translation of mRNA into one or more polypeptides.
In general, a “gene” is a region on the genome that is capable of being transcribed to an RNA that either has a regulatory function, a catalytic function, and/or encodes a protein. An eukaryotic gene typically has introns and exons, which may organize to produce different RNA splice variants that encode alternative versions of a mature protein.
By “nucleic acid molecule” is meant an oligomer or polymer of ribonucleic acid or deoxyribonucleic acid, or analog thereof.
By “polypeptide” is meant any chain of amino acids, regardless of length or post-translational modification (for example, glycosylation or phosphorylation).
By “reference” is meant a standard or control condition.
By “subject” is meant a mammal, including, but not limited to, a human or non-human mammal, such as a bovine, equine, canine, ovine, rodent, or feline.
By “marker” is meant any protein, polynucleotide or fragment thereof having an alteration in expression level or activity that is associated with a disease or disorder.
A description of example embodiments of the invention follows.
Gene expression data are routinely used to identify genes that on average exhibit different expression levels between a case and a control group. Yet, very few of such differentially expressed genes are detectably perturbed in individual patients. The disclosed methods and systems provide a framework to construct personalized perturbation profiles for individual subjects, identifying the set of genes that are significantly perturbed in each individual. This allows an analysis of the heterogeneity of the molecular manifestations of complex diseases by quantifying the expression-level similarities and differences among patients with the same phenotype. Despite the high heterogeneity of the individual perturbation profiles, patients with asthma, Parkinson's, and Huntington's disease, for example, share a broad pool of sporadically disease-associated genes. Individuals with considerable overlap with this pool have a 85%-100% chance of being diagnosed with the disease. The developed framework opens up the possibility to apply gene expression data in the context of precision medicine, with important implications for biomarker identification and diagnosis.
The disclosed methods and systems involve the identification of genes or proteins whose expression levels are perturbed in a single subject compared to a group of control subjects. The resulting personalized expression perturbation profiles (PEEPs) allow for a detailed investigation of the molecular roots of a disease state of a single subject, in contrast to conventional differential expression analysis methods that only yield average changes between two groups of subjects. The PEEPs can serve as a starting point to address various important challenges of personalized medicine, such as molecular-based diagnosis. The genes and/or proteins may be referred to herein as markers or biomarkers. Unlike conventional methods of differential expression analysis, which only allow for the identification of average changes between two groups, the disclosed methods and systems allow for a systematic quantification of the heterogeneity of disease states between different subjects on a molecular (e.g., gene or protein expression) level. The novel molecular signatures (disease modules) do not rely on a small set of marker genes, but on a larger set of genes that, by design, takes into account the heterogeneity of diseases.
Successful practice of the invention can be achieved with one or a combination of methods that can detect and/or quantify markers. These methods include, without limitation, sequencing methods (e.g., Sanger, Next Generation Sequencing, RNA-SEQ), hybridization-based methods, including those employed in biochip arrays, mass spectrometry (e.g., laser desorption/ionization mass spectrometry), fluorescence (e.g., sandwich immunoassay), surface plasmon resonance, ellipsometry and atomic force microscopy. Expression levels of markers (e.g., polynucleotides, polypeptides, or other analytes) are compared by procedures well known in the art, such as RT-PCR, Northern blotting, Western blotting, flow cytometry, immunocytochemistry, binding to magnetic and/or antibody-coated beads, in situ hybridization, fluorescence in situ hybridization (FISH), flow chamber adhesion assay, ELISA, microarray analysis, or colorimetric assays. Methods may further include one or more of electrospray ionization mass spectrometry (ESI-MS), ESI-MS/MS, ESI-MS/(MS)n, matrix-assisted laser desorption ionization time-of-flight mass spectrometry (MALDI-TOF-MS), surface-enhanced laser desorption/ionization time-of-flight mass spectrometry (SELDI-TOF-MS), desorption/ionization on silicon (DIOS), secondary ion mass spectrometry (SIMS), quadrupole time-of-flight (Q-TOF), atmospheric pressure chemical ionization mass spectrometry (APCI-MS), APCI-MS/MS, APCI-(MS)n, atmospheric pressure photoionization mass spectrometry (APPI-MS), APPI-MS/MS, and APPI-(MS), quadrupole mass spectrometry, fourier transform mass spectrometry (FTMS), and ion trap mass spectrometry, where n is an integer greater than zero.
Detection methods may include use of a biochip array. Biochip arrays useful in the invention include protein and polynucleotide arrays. One or more markers are captured on the biochip array and subjected to analysis to detect the level of the markers in a sample.
Markers may be captured with capture reagents that are immobilized to a solid support, such as a biochip, a multiwell microtiter plate, a resin, or a nitrocellulose membrane that is subsequently probed for the presence or level of a marker. Capture can be on a chromatographic surface or a biospecific surface. For example, a sample containing a protein marker may be used to contact the active surface of a biochip for a sufficient time to allow binding to a capture molecule. Unbound molecules are washed from the surface using a suitable eluant, such as phosphate buffered saline. In general, the more stringent the eluant, the more tightly the markers must be bound to be retained after the wash.
Upon capture on a biochip, a marker can be detected by a variety of detection methods selected from, for example, a gas phase ion spectrometry method, an optical method, an electrochemical method, atomic force microscopy and a radio frequency method. In one embodiment, mass spectrometry, and in particular, SELDI, is used. Optical methods include, for example, detection of fluorescence, luminescence, chemiluminescence, absorbance, reflectance, transmittance, birefringence or refractive index (e.g., surface plasmon resonance, ellipsometry, a resonant mirror method, a grating coupler waveguide method or interferometry). Optical methods include microscopy (both confocal and non-confocal), imaging methods and non-imaging methods. Immunoassays in various formats (e.g., ELISA) are popular methods for detection of analytes captured on a solid phase. Electrochemical methods include voltametry and amperometry methods. Radio frequency methods include multipolar resonance spectroscopy.
As used herein, the term “overlapping” refers to analytes that are present, and analytes that are absent (i.e. having a deviation above, or below the calculated mean level of the analyte, respectively) in the perturbation profiles of all case subjects. As used herein, the term “partially overlapping” refers to analytes that are present, and analytes that are absent (i.e. having a deviation above, or below the calculated mean level of the analyte, respectively) in the perturbation profiles of a portion of the case subjects. As used herein, the term “non-overlapping” refers to analytes that are present, and analytes that are absent (i.e. having a deviation above, or below the calculated mean level of the analyte, respectively) in the perturbation profiles of a one case subject.
In some embodiments, detecting differential levels of one or more analytes in one or more case subjects relative to the levels of the one or more analytes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more analytes for the case subjects comprises: a. calculating a mean level of each analyte detected in the one or more control subjects; b. calculating, using the calculated mean level each analyte detected in the one or more control subjects, the deviation of the level of each analyte detected in the one or more test subjects; and c. identifying, in the one or more test subjects, analytes with a deviation above or below a threshold deviation level from the calculated mean.
In some embodiments, detecting differential levels of one or more genes in one or more case subjects relative to the levels of the one or more genes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more genes for the case subjects comprises: a. calculating a mean level of each gene detected in the one or more control subjects; b. calculating, using the calculated mean level each gene detected in the one or more control subjects, the deviation of the level of each gene detected in the one or more test subjects; c. identifying, in the one or more test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
Referring to
Where, for each gene i of subject j the expression level lij is compared to the reference distribution of expression levels of that gene within the control group. The z-score captures how many standard deviations σcont(li) the individual expression level lij deviates from the mean value li
cont of the control group.
In some embodiments, the threshold deviation level is a global threshold zthresh that identifies the genes that are sufficiently perturbed in an individual subject.
In some embodiments, the resulting individual perturbation expression profile (PEEP) of a subject can be viewed as a “barcode,” representing the genes that are up- (zij>zthresh) or down-regulated (zij<zthresh) compared to the control group.
In some embodiments, zthresh is from 1.5 to 4. In some embodiments, is zthresh 2.5.
In some embodiments, the present invention is a computer-implemented method of generating a disease module, the method comprising: a. detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the plurality genes in at least one control subject, thereby obtaining a personalized perturbation profile for the subject; b. compiling the personalized perturbation profile from the individual subject across a population of case subjects; c. comparing the compiled personalized perturbation profiles with a set of genes whose differential expression is associated with said disease; and d. obtaining a set of overlapping genes from the compiled perturbation profiles that defines the disease module.
In some embodiments, detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the plurality genes in at least one control subject, thereby obtaining a personalized perturbation profile for the subject comprises: a. calculating a mean level of each gene detected in the at least one control subject; b. calculating, using the calculated mean level each gene detected in the at least one control subject, the deviation of the level of each gene detected in the test subjects; c. identifying, in the test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In some embodiments, the method further comprises obtaining a set of partially overlapping genes and non-overlapping genes from the compiled perturbation profiles. In some embodiments, the method further comprises calculating, based on the overlapping genes, partially overlapping genes, and non-overlapping genes, an expression heterogeneity of the disease.
In some embodiments, the expression heterogeneity of the disease is calculated by determining the mean pair-wise similarity of the data of the individuals in the case and control groups.
In some embodiments, the mean pair-wise similarity is determined by the distribution of Jaccard indicies
for all pairwise gene sets A and B of the individuals in the case and control groups.
In some embodiments, the statistical score is determined using the Fisher's exact test.
In some embodiments, the present invention provides a specifically programmed computer system comprising:
a. at least one specialized computer machine comprising:
In some embodiments, detecting differential levels of one or more genes in one or more case subjects relative to the levels of the one or more genes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more genes for the case subjects comprises: a. calculating a mean level of each analyte detected in the one or more control subjects; b. calculating, using the calculated mean level each analyte detected in the one or more control subjects, the deviation of the level of each analyte detected in the one or more test subjects; c. identifying, in the one or more test subjects, analytes with a deviation above or below a threshold deviation level from the calculated mean.
In some embodiments, the present invention provides a specifically programmed computer system comprising:
a. at least one specialized computer machine comprising:
In some embodiments, detecting differential levels of one or more genes in one or more case subjects relative to the levels of the one or more genes in one or more control subjects, thereby obtaining a set of personalized perturbation profiles of the one or more genes for the case subjects comprises: a. calculating a mean level of each gene detected in the one or more control subjects; b. calculating, using the calculated mean level each gene detected in the one or more control subjects, the deviation of the level of each gene detected in the one or more test subjects; and c. identifying, in the one or more test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In some embodiments, the present invention provides a specifically programmed computer system comprising: a. detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the plurality genes in at least one control subject, thereby obtaining a personalized perturbation profile for the subject; b. compiling the personalized perturbation profile from the individual subject across a population of case subjects; c. comparing the compiled personalized perturbation profiles with a set of genes whose differential expression is associated with said disease; and d. obtaining a set of overlapping genes from the compiled perturbation profiles that defines the disease module.
In some embodiments, detecting differential expression of a plurality of genes in an individual case subject relative to the expression levels of the plurality genes in at least one control subject, thereby obtaining a personalized perturbation profile for the subject comprises: a. calculating a mean level of each gene detected in the at least one control subject; b. calculating, using the calculated mean level each gene detected in the at least one control subject, the deviation of the level of each gene detected in the test subjects; c. identifying, in the test subjects, genes with a deviation above or below a threshold deviation level from the calculated mean.
In some embodiments, the method further comprises obtaining a set of partially overlapping genes and non-overlapping genes from the compiled perturbation profiles. In some embodiments, the method further comprises calculating, based on the overlapping genes, partially overlapping genes, and non-overlapping genes, an expression heterogeneity of the disease.
The following is a detailed explanation of concepts behind the disclosed methods and systems with specific examples.
The practice of many embodiments of the present invention employs, unless otherwise indicated, conventional techniques of molecular biology (including recombinant techniques), microbiology, cell biology, biochemistry and immunology, which are well within the purview of the skilled artisan. Such techniques are explained fully in the literature, such as, “Molecular Cloning: A Laboratory Manual”, second edition (Sambrook, 1989); “Oligonucleotide Synthesis” (Gait, 1984); “Animal Cell Culture” (Freshney, 1987); “Methods in Enzymology” “Handbook of Experimental Immunology” (Weir, 1996); “Gene Transfer Vectors for Mammalian Cells” (Miller and Calos, 1987); “Current Protocols in Molecular Biology” (Ausubel, 1987); “PCR: The Polymerase Chain Reaction”, (Mullis, 1994); “Current Protocols in Immunology” (Coligan, 1991). These techniques are applicable to the production of the polynucleotides and polypeptides of embodiments of the invention, and, as such, may be considered in making and practicing embodiments of the invention. Particularly useful techniques for particular embodiments will be discussed in the sections that follow.
The following examples are put forth so as to provide those of ordinary skill in the art with a complete disclosure and description of how to make and use the assay, screening, and therapeutic methods of the methods and system disclosed herein, and are not intended to limit the scope of what the inventors regard as their invention.
Microarray techniques, and more recently RNA sequencing, have fundamentally changed our ability to explore the molecular mechanisms underlying complex diseases, being routinely used to identify disease-associated genome-wide changes in gene expression patterns. An important goal of these studies is the identification of differentially expressed (DE) genes, whose expression level systematically differs between a case (disease) and a control (healthy) group. The expectation is that such DE genes will help pinpoint the molecular processes perturbed in a disease, which in turn can be used as biomarkers for diagnosis and prognosis (see Showe, M. K. et al. Gene expression profiles in peripheral blood mononuclear cells can distinguish patients with non-small cell lung cancer from patients with nonmalignant lung disease. Cancer research 69, 9202-9210 (2009) and Taylor, I. W. et al. Dynamic modularity in protein interaction networks predicts breast cancer outcome. Nature Biot. 27, 199-204 (2009)), patient classification and drug target identification. For example differential expression patterns of whole blood cells have long been considered promising candidates for cheap, easily accessible biomarkers for multiple diseases (see Zeller, T. & Blankenberg, S. Blood-based gene expression tests: promises and limitations. Circulation. Cardiovascular genetics 6, 139-140 (2013)).
Despite their extraordinary use in research and medicine, the interpretation and validation of gene expression patterns continues to offer major challenges. Indeed, results from similar studies are often inconsistent, the proposed biomarkers are often not reproduced, and the identified DE genes rarely point to a unique set of disease-associated genes (see Ein-Dor, L., Kela, I., Getz, G., Givol, D. & Domany, E. Outcome signature genes in breast cancer: is there a unique set? Bioinformatics 21 (2004)). For example, a meta study of multiple heart failure studies failed to identify any gene that is differentially expressed in all seven datasets, the most reproduced gene being differentially expressed only in four datasets (see Asakura, M. & Kitakaze, M. Global gene expression profiling in the failing myocardium. Circ. J. 73, 1568-1576 (2009)). Two main reasons are often listed as the source for these inconsistencies: (i) The comparison of different microarray-based measurement is hindered by important technical challenges, like the use of different platforms, dyes or statistical methods. (ii) There is intrinsic variability in gene expression levels, driven by both genetic factors, like the effect of single nucleotide polymorphisms (SNPs) and copy number variations (CNVs) on expression qualitative trait loci (eQTLs) (see Stranger, B. E. et al. Relative impact of nucleotide and copy number variation on gene expression phenotypes. Science 315, 848-853 (2007) and Stranger, B. E. et al. Population genomics of human gene expression. Nature Gen. (2007)), and non-genetic factors (see Spielman, R. S. et al. Common genetic variants account for differences in gene expression among ethnic groups. Nature Gen. (2007); Cheung, V. G. & Spielman, R. S. Genetics of human gene expression: mapping DNA variants that influence gene expression. Nature Rev. Gen. 10, 595-604 (2009); Wu, L. et al. Variation and genetic control of protein abundance in humans. Nature 499, 79-82 (2013); Alemu, E. Y, Carl, J. W., Bravo, H. C. & Hannenhalli, S. Determinants of expression variability. Nucleic Acids Res. (2014)), arising from epigenetic modifications (see Jaenisch, R. & Bird, A. Epigenetic regulation of gene expression: how the genome integrates intrinsic and environmental signals. Nature Gen. 33 (2003)) and the inherent stochasticity of biological processes (see Elowitz, M. B. Stochastic gene expression in a single cell. Science 297 (2002); Li, J., Min, R., Vizeacoumar, F. & Jin, K. Exploiting the determinants of stochastic gene expression in saccharomyces cerevisiae for genome-wide prediction of expression noise. Proc. Natl. Acad. Sci. U.S.A. (2010); and Balazsi, G., van Oudenaarden, A. & Collins, J. J. Cellular decision making and biological noise: from microbes to mammals. Cell 144, 910-925 (2011)). Here, we focus on a third important yet less explored factor: the heterogeneity of complex diseases, i.e., the possibility that multiple, only partially or non-overlapping molecular mechanisms can act in different patients with the same phenotype. For example, breast and colorectal tumors typically contain about 80 mutated genes (see Wood, L. D. et al. The genomic landscapes of human breast and colorectal cancers. Science 318, 1108-1113 (2007)). Yet, the mutations in different tumors have very little overlap, resulting in an astonishing number of more than 1,700 mutated genes identified in only 22 tumors. To date, about 140 “driver genes” have been identified, whose mutation promotes tumorigenesis in most cancer types, but only two to eight of these driver genes are mutated in any individual tumor (see Vogelstein, B. et al. Cancer genome landscapes. Science 339, 1546-1558 (2013)). A similar phenomenon is likely to occur at the gene expression level: many different perturbations may be associated with the same phenotype. We must therefore develop bottom-up methodologies that can interpret in a predictive fashion the inherent heterogeneity of individual perturbation profiles of both healthy and disease patients.
The disclosed methods and systems provide a framework to construct and integrate personalized perturbation profiles (PEEPs) from biomarker expression data, allowing us to systematically characterize the inherent heterogeneity of gene expression patterns. The approach is tested on asthma, a chronic inflammatory disease of the lung, Parkinson's disease (PD), a progressive disorder of the nervous system (see Scherzer, C. R. et al. Molecular markers of early Parkinson's disease based on gene expression in blood. Proc. Natl. Acad. Sci. U.S.A. 104, 955-960 (2007)), and Huntington's disease (HD), a neurodegenerative disorder caused by mutations in a single gene (HTT, Huntingtin) (see Borovecki, F. et al. Genome-wide expression profiling of human blood reveals biomarkers for Huntington's disease. Proc. Natl. Acad. Sci. U.S.A. 102, 11023-11028 (2005)). In all three diseases, there is a high heterogeneity between the PEEPs of individual patients. However, these heterogeneous patterns point to the existence of a single, highly predictive disease pool specific for each disease. Test results offer a conceptual change in the way disease-associated perturbations are interpret, in line with the emerging disease module hypothesis. Accordingly, disease-associated mutations perturb some cellular function that at the molecular level is encoded into a subnetwork of the underlying interactome. Therefore, multiple, often independent perturbations can impair the functional integrity of such a module, indicating that it is intrinsically impossible to associate a single gene or pathway to a specific pathophenotype.
Results
To illustrate the inherent limitations of group-based differential expression analysis, consider the POSTN gene, coding for the protein periostin. Periostin is an established biomarker for asthma (see Takayama, G. et al. Periostin: a novel component of subepithelial fibrosis of bronchial asthma downstream of IL-4 and IL-13 signals. The Journal of allergy and clinical immunology 118, 98-9104 (2006); Sidhu, S. S. et al. Roles of epithelial cell-derived periostin in TGF-activation, collagen production, and collagen gel elasticity in asthma. Proc. Natl. Acad. Sci. U.S.A. (2010); and Parulekar, A. D., Atik, M. A. & Hanania, N. A. Periostin, a novel biomarker of th2-driven asthma. Curr Opin. Pulm. Med. 20, 60-65 (2014)), its role in airway remodeling being exploited by an experimental asthma drug (see Corren, J. et al. Lebrikizumab treatment in adults with asthma. The New England Journal of Medicine (2011)). The strong differential expression pattern between asthmatic and healthy subjects confirms its asthma association (
To generalize the above observations, we inspected the expression levels of all genes that were differentially expressed according to a standard group-wise analysis in asthma, PD and HD (see Methods, below). As shown in
A framework for personalized gene-expression analysis—These limitations can be overcome and individual expression heterogeneity can be turned into a predictive information by constructing personalized perturbation profiles that reflect expression changes within a single subject, rather than mean changes between two groups (
capturing how many standard deviations σcont(li) the individual expression level lij deviates from the mean value li
cont of the control group. We then use a global threshold zthresh to identify the genes that are sufficiently perturbed in an individual subject. The resulting individual perturbation expression profile (PEEP) of a subject can be viewed as a “barcode,” representing the genes that are up- (zij>zthresh) or down-regulated (zij<−zthresh) compared to the control group (
To characterize the PEEPs, the group of genes perturbed in individuals is compared with DE genes obtained from a standard group-wise approach (see Methods, below). The first observation is that only for HD we find group-wise DE genes that are contained in all individual profiles. For asthma and PD, no single gene is perturbed in all case subjects.
Quantifying the heterogeneity among individual perturbation profiles—To quantify the underlying expression heterogeneity of a disease, we move beyond the group-wise DE genes and ask, instead, how similar are the PEEPs of two individuals with the same disease.
for all pairwise gene sets A and B of the individuals in the case and control groups of three diseases. For asthma, the mean pairwise similarity (J
=3×10−2) is three times higher in the case group than in the control group (
J
=1×10−2). While this difference is highly significant (p-value <10−77, Mann-Whitney U-test), in absolute numbers the overlap is small: While a typical asthmatic subject has on average 379 perturbed genes, the average number of shared perturbed genes between two asthmatic subjects is only 24 (
J
˜0.24 and
J
˜0.14 for HD and PD, respectively. The same analysis can also be performed on the full continuous z-score profiles using Pearson correlation as measure of similarity, yielding similar results (see
To quantify whether the observed overlap between the PEEPs of the case subjects could have emerged by chance, we calculated the statistical significance for each pair individually using Fisher's exact test. As expected, we find the overlap to be significant for all subject pairs, even after applying the most conservative Bonferroni correction (
This leads to our second main result: We observe highly significant similarities between the PEEPs of case subjects, similarities that are absent in healthy subjects. These similarities cannot be attributed to a few widely shared DE genes identified by the group-wise differential expression analysis, but arise from more complex patterns of pairwise overlaps.
Functional analysis of the perturbation profiles—The low overlap between the personalized profiles of case subjects prompts us to ask how the molecular level heterogeneity translates into relatively homogeneous disease phenotypes. To address this, we examine the extent to which the individual profiles reflect disruptions in common disease-specific pathways (
We next determined the Jaccard similarity of the respective individual perturbed pathway genes for each pair of subjects whose PEEPs are significantly enriched with genes of the pathway. The low similarity values (J˜0.1,
We repeated the analysis of the heterogeneity among perturbed pathways also for AD and PD, using general pathway annotations from MSigDB (see Subramanian, A. et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. U.S.A. 102, 15545-15550 (2005)) and functional Gene Ontology (GO) (see Ashburner, M. et al. Gene ontology: tool for the unification of biology. Nature Gen. 25, 25-29 (2000)). The results again indicate that the same biological function or pathway is perturbed in different ways in different patients (
Predicting diseases from PEEPs—Taken together, our results indicate that patients with the same disease exhibit highly heterogeneous perturbations that nevertheless point towards common functional disruptions. This suggests the existence of a broader group of genes, whose perturbations are associated with the specific disease. As demonstrated next, by compiling all genes that are perturbed in a significant fraction of the case subjects, we can accurately predict the disease state of each patient.
Given the relatively high number of genes perturbed in the individual profiles
Perturbations of these modules uniquely characterize the respective diseases. To show this, we used a repeated cross-validation approach and determined the different PEEP's overlap with the disease module (see Methods, below). We find that the fraction of genes from the disease module perturbed in an individual subject accurately predicts whether the subject has the disease. For asthma, the PEEPs of case subjects contain on average 21% of the asthma disease pool, compared to less than 7% for the control subjects. For PD and HD the overlap of the case subjects with the corresponding disease modules is much higher, obtaining 65% and 86% respectively, compared to 20% and 6% for the control subjects. This indicates that PD and HD are characterized by a more specific set of characteristic perturbations, while asthma displays a more heterogeneous range of associated perturbations. The receiver operating characteristics (ROC) in
Discussion
Group-wise expression analysis has two important limitations: (i) It can only identify genes that are consistently (i.e., in the same direction) perturbed in a large fraction of the patients. (ii) It does not yield patient specific information. Here, we introduced a simple, yet powerful method that overcomes these limitations and offers personalized perturbation profiles, or PEEPs. The method can be interpreted as a generalization of group-wise differential expression methods with PEEPs representing personalized differentially expressed genes. As a consequence, the PEEPs can be easily interpreted and further analyzed using established tools, such as the geneset enrichment analysis used above.
As illustrated in
Despite the high gene level variability, the commonalities at the functional and pathway level indicate that complex diseases arise from disruptions of certain biological processes or disease modules (see Menche, J. et al. Uncovering disease-disease relationships through the incomplete interactome. Science 347, 1257601 (2015)), hence the observed heterogeneity simply reflects the molecular diversity of such disruptions. We therefore expect considerable variability among the expression profiles of subjects with the same disease not despite, but because they all have the same disease. Recently, a number of studies proposed various strategies for dissecting disease heterogeneity, in particular in the field of cancer. The PARADIGM algorithm (see Vaske, C. J. et al. Inference of patient-specific pathway activities from multi-dimensional cancer genomics data using paradigm. Bioinformatics 26, i237-i245 (2010)), for example, infers patient-specific pathways using various omics-type information, such as expression and mutational data, together with curated pathway interactions. Another widespread algorithm, HotNet2 (see Leiserson, M. D. et al. Pan-cancer network analysis identifies combinations of rare somatic mutations across pathways and protein complexes. Nature Gen. 47, 106-114 (2015)) tackles the genetic heterogeneity of different cancer samples using the concept of information propagation starting from known mutations in order to identify cancer-related subnetworks in signaling networks. In this work, we document the existence of large disease module also on other disease areas and using transcriptional data only. We find that a sufficient level of random perturbations among these disease modules can accurately predict the presence/absence of a particular disease. We integrated the personalized profiles of all patients to reconstruct the respective disease module, finding that the fraction of genes in an individual's PEEP is a near perfect predictor for a patient's disease status. This suggests that personalized profiles could identify combinatorial biomarker signatures that go beyond single markers. Methods for constructing a disease module from given seed biomarkers are disclosed in PCT published application no. WO 2015/084461. With next-generation sequencing technology advancing at a fast pace, there is great potential for applying RNAseq technologies to identify transcriptional signatures also in a clinical setting (see Byron, S. A., Van Keuren-Jensen, K. R., Engelthaler, D. M., Carpten, J. D. & Craig, D. W. Translating RNA sequencing into clinical diagnostics: opportunities and challenges. Nature Rev. Gen. 17, 257-271 (2016)). Such signatures are of key importance for personalized medicine and could, for example, help diagnose previously unrecognized diseases. While the results presented here provide first evidence of the general feasibility of using our approach to obtain predictive biomarkers, a comprehensive reference base across all relevant diseases and more extensive tests concerning the robustness and reliability of the resulting disease pools will be required towards an actual clinical application. The observed heterogeneity among the individual perturbation profiles further indicates that single-target drugs may be effective only in a small number of patients. Instead, multi-target strategies may prove more promising for drug development (see Csermely, P., Korcsmáros, T., Kiss, H. J., London, G. & Nussinov, R. Structure and dynamics of molecular networks: a novel paradigm of drug discovery: a comprehensive review. Pharmacology & therapeutics 138, 333-408 (2013)). Our approach can be used to quantitatively assess the expected fraction of patients for which a drug is expected to be effective, helping guide the development of targets with maximal efficacy.
Materials & Methods
Gene expression data—We use data from an ongoing study by Janssen Research & Development for asthma (manuscript in preparation) and previously published expression profiling studies for Huntington's disease and idiopathic Parkinson's disease (see Lesnick, T. G. et al. A genomic pathway approach to a complex disease: axon guidance and parkinson disease. PLoS Genet. 3, e98 (2007)). The asthma dataset contains 55 case subjects with moderate or severe asthma and 25 gender and age matched healthy control subjects (see Silkoff, P. et al. Asthma characteristics and biomarkers from the airways disease endotyping for personalized therapeutics (adept) longitudinal profiling study. Respiratory research 16, 1 (2015) for a detailed description of the cohort). The asthma samples were collected from bronchoscopy (endobronchial biopsies and brushings), preserved immediately in RNAlater® solution and then maintained at −70° C. Qiagen miRNeasy kit (Qiagen; Germantown, Md.) and NuGen ovation pico WTA kit (NuGen Technologies; San Carlos, Calif.) were used to extract and amplify RNA. cDNA is profiled using Affymetrix HG-U133+PMchip (Affymetrix, Santa Clara, Calif.). CEL files were assessed using Almac Diagnostics Microarray Toolbox for quality control (chip image analysis, Affymetrix GeneChip QC, RNA degradation analysis, distribution analysis, principal components analysis, and correlation analysis) and technical outliers are excluded. Robust multi-array (RMA) method is used to re-normalize the profiles, followed by batch effect adjustment via linear modeling of batch (as random factor) and cohort. The Huntington's disease dataset19 (GEO accession number GSE1767) contains analysis of blood samples from 17 case subjects (5 presymptomatic and 12 symptomatic) and 14 control subjects. In HD, the gene expression is suggested to be altered in a variety of tissues including peripheral blood. Affymetrix U133A GeneChips and Amersham Biosciences CodeLink Uniset Human I and II bioarrays were used to analyze the gene expression in blood samples. The Parkinson's disease data (GSE7621) contains 16 case and 9 control subjects for which multiregional gene expression analysis was conducted in postmortem brain using Affymetrix HG U133 Plus 2.0 gene chips. For the PD and HD datasets, the details of the sample generation and expression profiling can be found in the original publications. We reprocessed the raw data set in GEO for Parkinson using RMA with quantile normalization as implemented in the R package ‘affy’. We verified the quality of the data sets by checking the gene expression distribution and sample clustering in PCA. All expression levels in the PD and HD data were log2 transformed to facilitate direct comparison of the three data, overall results do not depend on the transformation, however. Basic statistics of the used datasets are shown in
Group-wise differential expression analysis—Genes differentially expressed are identified between case and control subjects using the limma R Bioconductor package (see Smyth, G. K. limma: Linear models for microarray data. In Gentleman, R., Carey, V. J., Huber, W., Irizarry, R. A. & Dudoit, S. (eds.) Bioinformatics and computational biology solutions using R and Bioconductor, 397-420 (Springer, New York, 2005)). The difference between expression levels of case and control subjects are assessed by fitting the expression levels to a linear model using one coefficient for each group in the design matrix. The probesets were mapped to Entrez Gene IDs using the platform annotation files in each data set. In case there were multiple probesets corresponding to the same Gene ID, the probeset with the maximum expression was used in the analysis. The p-values were corrected for multiple hypothesis testing using the Benjamini-Hochberg method (see Benjamini, Y. & Hochberg, Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. Roy. Statist. Soc. Ser. B 289-300 (1995)). At a cut-off of FDR <0.2 we obtain 417, 524, and 7,419 DE genes for asthma, PD and HD, respectively.
Personalized perturbation analysis—To construct the personalized perturbation profile of a subject j we compare the expression level lij each of its genes i to the reference distribution of expression levels of the same gene within the control group. The extent to which gene i is perturbed in subject j is quantified by the z-score
that indicates how many standard deviations σcont(li) the individual expression level lij is away from the mean value li
cont of the control group. Note that if subject j itself is part of the control group, we do not consider it for the computation of the reference distribution but use only the remaining control subjects. We then use a global threshold zthresh to define the set of genes that are perturbed in an individual subject. Positive z-scores zji>zthresh indicate up-regulation, negative values zji<−zthresh indicate down-regulation. The role of zthresh is thus analogous to the one of cutoffs for calling differential expression in standard group-wise analyses. Higher zthresh values result in smaller, more stringent PEEPs that potentially miss out relevant genes with less pronounced perturbations, while lower zthresh values provide a more global picture that may, however, also contain an increased number of false positives. The precise choice of zthresh can be optimized for specific purposes, such as disease patient classification (compare with
We systematically evaluated the stability of the obtained z-scores against changes in the population size by removing increasing numbers of subjects from the control population: We first calculated the expression level that corresponds to z-score=2.5 compared to the average mean and standard deviation of all genes and all control subjects. We then observed how the z-score of this expression level changes while randomly removing an increasing number of subjects from the pool of control subjects (up to 50% of the original population). For each gene, we then calculated the ratio of the z-score obtained from the decreased population to the original z-score calculated from the original population.
Analytical comparison with randomly distributed genes—To determine the minimal number X of case subjects in which a gene must be perturbed in order to be collected into the global pool of disease associated genes we use a comparison with random expectation. We consider a null model where each subject has g perturbed genes that are drawn completely at random from all G genes. The probability for one gene to be perturbed in exactly k out of n subjects is then given by the binomial distribution
with p=g/G. Using the mean number of genes observed in the individual profiles for g, the histogram of the number of subjects per gene can now be obtained by simply multiplying G×f(k; n, p). We find excellent agreement between this formula and the distributions observed among the control subjects, but, as expected, not for case subjects (
which can be solved by simply testing increasing values. Finally, we choose the minimal value as X=Xrand+1, thereby ensuring a broad, yet high-quality pool of disease-associated genes. The calculated values are X=10 for asthma, X=9 for PD, and X=10 for HD.
Cross-validation analysis for disease state prediction—We performed a five-fold cross-validation analysis using the fraction of genes of the combinatorial pool of disease-associated genes that is contained in a subject's personal perturbation profile to predict the disease state of the subject. Note that we do not take the direction of the perturbation into account. If the fraction is larger than a given threshold that can be determined from the training data we classify the subject as “case,” otherwise as “control.” This threshold not only allows for patient classification, but can also be interpreted as a direct measure of the heterogeneity of a disease. For the cross-validation, we randomly split the subjects into five groups having similar proportions of cases and controls as in the full dataset. We then iteratively use each group as the validation set and the remaining four groups as training data to generate the PEEPs and the combinatorial disease pool. Next, we calculate the fraction of the combinatorial pool that is contained in the PEEP of each subject in the validation set. By using all identified fractions as putative thresholds for classification as “case” or “control” and comparing with the true labels we then construct the ROC (receiver operating characteristic) curve and calculate the AUC (area under the curve). Note that the classifier is completely blind to the information of the left-out validation subjects, thus avoiding overfitting due to the fact that the combinatorial pool itself is compiled from all genes that are perturbed in X or more case subjects. The entire procedure is repeated 100 times to get robust estimates of the ROC curve and the AUC.
We further compared the performance of the PEEP-based classification to a k-nearest-neighbor (kNN) based classification. For every sample in the test set we calculated the gene expression correlation with all samples in the training set and then ranked the training samples according to the strength of the correlation. The known disease states of the k most similar samples (i.e., highest correlation) is then used to score the test sample's likelihood to belong to the same class. After evaluating a range of values of k (=3, 5, 10, 15, 20) we found k=15 to offer the highest prediction accuracy. Note that while the kNN method allows for a high-quality classification, the subsequent interpretation of a classification result is less straightforward compared to the PEEP approach above, which is directly based on overlapping gene sets that can be immediately further investigated and potentially validated.
To estimate the influence of the sample size on the final accuracy of the classification analysis, we further repeated both the kNN and the PEEP-based classification using a two-fold validation scheme, such that only half of the case and control subjects are available for training. The results shown in Supplementary
Functional gene annotation data—To analyze the biological function of genes and gene sets, we can use Gene Ontology (GO) terms, general pathway annotations and asthma-specific pathways. GO annotations were downloaded from http://www.geneontology.org/. We only use high confidence annotations associated with the evidence codes EXP, IDA, IMP, IGI, IEP, ISS, ISA, ISM or ISO and further remove all associations with a non-empty “qualifier” column (see Berriz, G. F., Beaver, J. E., Cenik, C., Tasan, M. & Roth, F. P. Next generation software for functional trend analysis. Bioinformatics 25, 3043-3044 (2009)). Since the provided GO files only contain the most specific annotations explicitly, we add all implicit more general annotations by up-propagating the given annotations along the full GO tree.
The general pathway annotations were taken from the Molecular Signatures Database (MSigDB) published by the Broad Institute, Version 4.0. MSigDB integrates several pathway databases; we use those from KEGG, Biocarta, and Reactome. Asthma-specific pathways (
Gene set enrichment analysis—The enrichment analysis between a given gene set and a pathway or GO annotation (“term”) may be accomplished using Fisher's exact test. We considered a term to be significantly enriched if p-value <0.05 (Bonferroni correction for number of tested terms). For each bar in
R-package—We provide the R package “PePPeR” (Personalized Perturbation ProfileR) which includes functions to fetch expression data sets from the GEO database, identify group-wise DE genes and construct individual perturbation profiles. The R package along with its documentation is available at https://github.com/emreg00/pepper.
While this invention has been particularly shown and described with references to example embodiments thereof, it will be understood by those skilled in the art that various changes in form and details may be made therein without departing from the scope of the invention encompassed by the appended claims.
This application claims the benefit of U.S. Provisional Application No. 62/253,878, filed on Nov. 11, 2015. The entire teachings of the above application are incorporated herein by reference.
| Filing Document | Filing Date | Country | Kind |
|---|---|---|---|
| PCT/US2016/061401 | 11/10/2016 | WO | 00 |
| Number | Date | Country | |
|---|---|---|---|
| 62253878 | Nov 2015 | US |
