January 25, 2016.  We recently applied our networking and visualization tools to explore the biological underpinnings of genes linked to Intellectual Disability (ID), discovering subsets comprising molecularly related classes of ID genes.  Several of these mathematically derived classes fit with mechanisms known to play a role in ID such as metabolism, neuronal development, mitochondria, epigenetic regulation, and cohesin. This framework enabled the biological annotation of a number of uncharacterized ID genes such as RAI1 (Smith Magenis Syndrome, SMS), which exhibits tight linkage to epigenetic regulation networks.

    The largest subset, comprising approximately 16% of ID genes, exhibited strong linkage to the biology of histone modification (epigenetic regulation), a subset of cohesin genes, and the core regulators of circadian rhythm.  Recent work reported that small molecule HDAC inhibitors could reverse some of the deleterious phenotypes in the Histone MLL2/KMT2D heterozygous knock out mouse model of Kabuki syndrome (1).  This provides a rationale for exploring small molecule inhibitors of HDAC or Bromodomain proteins for their ability to modulate phenotypes resulting from genetic defects in ID histone modification linked genes such as RAI1.


    1. Genes linked to intellectual disability could be mathematically grouped into five biologically distinct subsets (Epigenetic Regulation, Mitochondrial function, Metabolism, Proliferation and Neuronal Function)

    2. Gene networks linked to Epigenetic Regulation, the enzymatic modification of chromatin to affect gene transcripton, comprised about 16% (63/380) of all ID genes examined.

    3. RAI1, the gene responsible for Smith Magenis Syndrome exhibited tight linkage to epigenetic regulation, suggesting that this PHD domain containing protein plays a role in regulation of gene expression through chromatin modulation.

    4. RAI1 and other epigenetic regulation linked genes also exhibited linkage to genes that control the circadian rhythm oscillator.  Several of these syndromes exhibit sleep disturbances (CDKL5, Williams syndrome linked GTF2I, MBD5 and RAI1), supporting hypotheses regarding a causal relationship between syndromic gene dysfunction and control of circadian rhythm.  

    Intellectual disability encompasses a wide range of human phenotypes that impact the development of intellect, language, social and motor skills (2).  Broadly defined, Intellectually disability has an estimated prevalence in Western populations of ~2%, with both environmental and genetic factors contributing to the condition. 

     Using a large set of transcriptome data, we conducted comprehensive gene expression correlation analysis to create mathematically based correlation networks for about 20,000 transcribed genes.  Many tightly linked genes were highly enriched in genes whose biological function/properties were known (see earlier Molquant work).  Figure 1 shows 36 such annotated networks, providing a biological framework for analyses.

    FIGURE 1

    Relationships among biologically curated gene networks

    Click for high resolution image

     Because the network tools rely only on mathematical correlation as opposed to knowledge based linkage, any gene can be placed into its relevant network.  This enables functional annotation of unknown or poorly characterized genes based on their expression correlation to sets of known genes.

    To apply this analytical framework to ID, we first collated 380 genes, each with a reported relationship to ID, from publications, reviews and meeting reports (Table 1).  Cross gene correlation analysis, clustering and plotting the genes revealed a complex set of networks among the 380 genes (figure 2).  Among the network clusters, five predominant clusters were observed, which together comprised approximately half of the ID genes using a network pearson correlation (mean of the 10 plotted network genes) of >0.4.  These mathematically derived clusters recapitulate several biological processes known to be associated with intellectual disability ((3))  The biological networks plotted in figure 1 were clustered together with the ID genes (figure 2), enabling annotation of the clusters.

    FIGURE 2

    Network relationships of 380 Genes linked to Intellectual Disability

    click for high resolution image

    The largest of the mathematically determined clusters contained predominantly genes known to be involved in epigenetic regulation, including ID genes CREBBP, ARID1A, ARID1B, MLL, MLL2, KDM5C, KDM6A, KDM6B, EHMT1, EP300, BRWD3, ATRX , MYST3, MYST4, and a number of others.  Also present in the cluster were several ID genes with limited biological characterization, including PRR12, PTPN23 and RAI1.  Using a “guilt by association” rationale, we hypothesize that these genes are also involved in epigenetic regulation.

    Focusing on RAI1, we then generated a network based on RAI1 itself .  The top 50 genome-wide expression correlates of RAI1 are clustered in figure 3.  Note that about 1/3 of the top genes in the RAI1 network are also ID genes (yellow), and that most of the genes in the network are linked to epigenetic regulation, further supporting the functional assignment of RAI1.

    FIGURE 3

    RAI1 anchored network plotted with biological networks (green) shows tight correlation to epigenetic regulation ID genes. Known ID genes in yellow, RAI1 in Red.

    click for high resolution image


    Bitter Taste Receptors Found on Neutrophils

    Our sense of taste might relate to the ability of our white blood cells to detect and respond to pathogens

    May 4, 2015.  As we continue to build more biology-centric networks, we ran into this one.  Humans have (at least) 28 Taste Receptors, G protein coupled receptors for sweet, umami and bitter taste.  TAS1Rs are thought to sense sweet and glutamate (umami), whereas the shorter TAS2Rs are responsible for detecting bitter taste.  While bitter taste receptors are thought to be specifically expressed on the taste papillae (buds) of the tongue, our analysis found an unmistakable signature for 22 of the TAS2R bitter taste receptors on blood neutrophils, as defined by neutrophil gene ELANE (FIGURE 1).  TAS1R3 was linked to a distinct set of  blood cells: CD8A postive T lymphocytes, NK cells (defined by KIR2DL3) and CD68 positive macrophages.  Consistent with this observation, a February 2015 publication reported the discovery of taste receptors on leukocytes, and a 2014 paper reported umami receptors on neutrophils.

    FIGURE 1

    The implications of this observation are uncertain.  However, recent work has identified bitter taste receptors as bacterial chemosensors in lung airway epithelium.  In addition, a variant in one bitter taste receptor was linked to susceptibility to respiratory infection. Perhaps it's not surprising that this mechanism might extend to our white blood cells, enabling them to use these chemical receptors to sense and seek out foreign invaders such as pathogenic bacteria.  We are curious to hear from the innate immunity & chemoreception communities on this.

    This observation also highlights one of the limitations of the Molquant tools.  Because our datasets do not contain any tongue samples, we are unable to identify the presence of these receptors on tongue papillae.  If these receptors represented previously uncharacterized biology, we would mistakenly conclude that these GPCRs were specifically expressed in neutrophils.  We look forward to the increasing size and diversity of available samples to improve the networks.


    April 2, 2015.  Some cancer susceptibility loci may be expected to share genes/mechanisms.  We recently analyzed breast cancer susceptibility loci, and were curious to see if other cancer risk GWAS showed overlapping genes/loci.  Because breast cancer risk genes BRCA1 and BRCA2 also confer increased risk for prostate cancer, we took a quick look at prostate cancer risk data. Twenty three loci were published last year in Nature Genetics "A meta-analysis of 87,040 individuals identifies 23 new susceptibility loci for prostate cancer". Of 56 breast cancer loci, we noted only a single overlapping locus on Chromosome 1 containing PEX14, APITD1, DFFA, PGD, CORT, suggesting the majority of the cancer susceptibility mechanisms were distinct between the two tumor types.

    To further explore candidate biology in the Prostate Cancer GWAS, 67 genes near the 23 loci were plotted on a bathymetry plot (Figure).

    Among several network clusters, the most prominent contained TMPRSS2, half of an oncogenic fusion partner (TMPRSS2-ERG) common in prostate cancer.  Other genes in the network MARVELD3, TTC9 (and neighbor MAP3K9), GOLPH3L, MYO6, HTR3A, SIX4 all exhibit prominent salivary and/or exocrine gland expression, including TMPRSS2 which is prostate specific.  While there is more to learn about the biology of this cluster, we present it here for the prostate cancer and GWAS community.


    March 23, 2015.  Blood gene expression represents a potentially exciting new approach to dynamically monitor health, first explored by 'omics pioneer Mike Snyder at Stanford. Last week's PLOS Genetics publication of blood RNAs linked to hypertension, "A Meta-analysis of Gene Expression Signatures of Blood Pressure and Hypertension" demonstrates the emerging use of RNA expression to characterize specific biological conditions. Led by Daniel Levy, MD, the authors identify 34 genes in blood RNA that correlate with blood pressure, mRNAs that in aggregate (top 55 SBP, 22 DBP) explain an estimated 5-9% of the variation in BP.

    To better understand how circulating lymphocytes can sense hypertension, we analyzed the 34 gene signature with Molquant analytical tools.  Using a large diverse set of expression data, networks are generated and visualized to identify co-expression relationships among genes--correlated networks are inferred to exhibit similar biology.  Note that here the GTEX tissue dataset (2900 RNAseq samples from human tissues) was not used to form the networks, such that the GTEX plots shown here represent an independent assessment of the networks identified.

    Figure 1 shows expression of the genes across the GTEX project tissue survey data.  Consistent with the lymphocytic origin of the RNA, a strong enrichment for blood expression was seen in the signature genes.


    To provide a more refined view of the 34 genes, Figure 2 shows the networks of each gene plotted (white) along with networks for several previously generated (1) hematopoietic lineages (red).

    FIGURE 2:

    Two dominant, and related biological networks were observed, one tightly linked to T cell and macrophage lineages, another comprising genes previously associated with inflammation (COX-2/PTGS2, GADD34/PPP1R15A, DUSP1, FOS).

    These networks bring a finer resolution to the findings of the authors who noted that GSEA analysis identified inflammation and apoptosis as top enriched pathways.  Together, these observations enable a hypothesis that the identified signature represents an activated hematopoietic signature associated with inflammation.  Inflammation has long been linked with cardiovascular conditions linked to hypertension (e.g. Inflammation in Hypertension).  

    The network analysis further suggests caution in linking signature genes to specific non-hematopoietic mechanisms of hypertension. For example, figure 2 identifies links between KCNJ2 and hematopoiesis (figure 2).  While KCNJ2 is expressed in heart, the GTEX tissue expression survey identifies blood as the top expressing tissue for KCNJ2, consistent with our analysis.  In addition, genes identified by the authors as representing a co-expression network in neutrophils exhibit limited linkage to neutropils in either figure 2 or the GTEX dataset (figure 3).  

    FIGURE 3:

    While this work represents an exciting addition to the use of genomics data to characterize the biology of hypertension, we would agree with the authors that one should be cautious about interpreting hypertension biology from hematopoietic cell transcriptomes.  It appears that blood based signatures can sense medically relevant conditions, but perhaps they do so by reflecting the impact of a condition on the hematopoietic system, not by revealing the genes driving the underlying mechanisms. Transcriptomic analysis from less accessible tissues such as endothelial, cardiac or neural tissue may provide a more direct window to mechanism.



    January 24, 2015.  We recently extended our network, cluster, and visualization tools to process thousands of genes into a single plot. To explore the utility of the tools, we generated large-scale network plots of disease-associated genes for five complex human disorders:  Rheumatoid Arthritis (173 genes), Schizophrenia (325 genes), Intellectual Disability (351 genes), Autism Spectrum Disorder (338 genes) and Cardiovascular Disease (213 genes).  In each case, the mathematically derived networks recapitulated sets of previously reported enriched genes.  However, the plots also identified additional gene networks that may point to distinct disease associated biology.  Stay tuned for detailed analysis of selected network clusters in future posts.

    The prolific generation of disease/gene association data provides an opportunity to rapidly advance our understanding of complex diseases, subsequently enabling the development of targeted therapeutic interventions.

    Of the thousands of studies now completed, we chose five disorders in which multiple studies have been conducted and analyzed, often in aggregate as meta-analyses.  Such approaches surface hundreds of genes or loci, where the link to the trait is based on a statistical linkage (either Genome Wide Association Study, GWAS; or an identified genomic alteration).  This approach typically provides little information on the biology of the association, especially in GWAS, where it can difficult to pinpoint the responsible gene in an identified locus.

    Using correlation statistics, a large dataset of diverse gene expression profiles and sequence variants, as well as curated biological networks, we looked for biologically relevant clusters within the reported lists of disease-linked genes.

    Biologically Defined Networks

    To demonstrate the network plots using known biology, figure 1 shows a plot representing 32 biologically annotated networks, each represented by 10 genes.  The degree of expression correlation among each of the 320 genes in the plot is indicated by color—note that each of the 10 genes of each network exhibit a high degree of correlation to one another, appearing as a rectangle along the diagonal.  Most networks exhibit little cross network correlation, whereas some (e.g. hypoxia response and endothelial cells) exhibit cross network linkage consistent with their related biological functions.  Of note, proteasome, SNRPB linked splicing and proliferation exhibit linkage—we are not aware of reports linking these biological processes.

    FIGURE 1

    Rheumatoid Arthritis

    A complex disease with relatively well-characterized biological underpinnings, many studies have noted strong enrichment of genes associated with immune function and immune cell activation.  Plotting 173 genes from a recent analysis is consistent with previously identified biology (figure 2).  

    FIGURE 2


    Recent GWAS and genetic analysis has identified strong enrichment for neuronal development and synaptic signaling among a large set of genes linked to schizophrenia.  Plotting 325 genes (figure 3) from a recent GWAS meta-analysis revealed several neuronal network clusters among the larger gene set.

    FIGURE 3


    A complex developmental disorder with difficult diagnostic criteria, a wealth genomic and genetic studies have identified strong enrichment for neuronal (and chromatin modifier) genes.  A plot of  338 genes curated from shows enrichment for neuronal development and function genes (figure 4).  

    FIGURE 4

    Cardiovascular Disease


    A leading cause of death in the developed world, cardiovascular disease has both genetic and environmental components.   Here, noting common loci among diverse cardiovascular phenotypes,  we examined a set of GWAS studies examining: electrocardiographic phenotypes (see blog entry August 2014), Sudden Cardiac death, Blood Pressure, Coronary Artery Disease GWAS, atrial fibrillation and QRS/QT duration.  The resulting 215 genes are listed here

    Plotting the 215 genes (figure 5) from the results of the diverse studies revealed several previously unreported network clusters.  Further exploration of the biology associated with the networks may provide new insight into the processes linked to cardiovascular disease.

    FIGURE 5


    September 24, 2014.  We recently developed a candidate assay for early detection of Ovarian Cancer.

    Early detection for cancer has the potential to dramatically alter outcomes in Oncology.  Like many topics in oncology, controversies remain over the magnitude of benefit of early detection (e.g. last week’s Salvo in the Wall Street Journal), however, when combined with molecular profiling (discussed below) a robust non-invasive cancer detection assay holds great promise for saving lives.

    Ovarian cancer is one tumor type for which early detection is likely to be most beneficial.  Figure 1 shows Ovarian cancer statistics from which collects epidemiological data on cancer detection and survival.  Due to its internal location and relative lack of symptoms, only 15% of ovarian cancer is detected when localized to the primary site, whereas more than 60% of ovarian cancer has already metastasized at diagnosis.  The difference in 5 year survival is striking:  92% survival if localized; 27% if metastatic.  The result is that ovarian cancer is one of the more lethal tumors with only 45% of patients surviving 5 years.  Compare that to prostate cancer at 99%, breast cancer at 89% or Colon cancer at 65%.  The most brutal is Pancreatic Adenocarcinoma at 7% (numbers according to SEER).  

    FIGURE 1

    Signatures may help address the “benign growth” critique.

    One concern regarding interpretation of the SEER data, and the source of most early detection controversy is the possibility of a selection bias in the detection of localized tumors.  Perhaps a larger proportion of early tumors would have never progressed; therefore, treatment and/or removal of all early lesions may be unnecessary and may actually increase risk to the patient.

    The recent availability of large molecular datasets in oncology provides the opportunity for robust molecular characterization of tumors, and has already led to RNA based prognostic signatures for several tumor types (e.g.  Genomic Health’s Oncotype DX breast, colon and prostate tests GHI prognostics test page).  The presence of prognostic signatures in localized tumors demonstrates the utility of RNA expression profiles to distinguish between benign and malevolent early lesions.  

    Twenty to Thirty RNAs can define high grade serous ovarian cancer

    We recently applied our analytical framework and tools to ovarian cancer and identified a concise set of genes that are highly enriched in most high grade serous ovarian cancer, relative to either other tumor types or a panel of normal tissues (Figure 2A, B). 

    FIGURE 2 

    The most studied ovarian tumor marker CA125, originally identified by Bob Bast and coworkers over thirty years ago (Reactivity of a monoclonal antibody with human ovarian carcinoma), exhibits relatively poorer RNA expression enrichment in either tumor panels (Figure 3).   Note that each individual gene in the ovarian signature exhibits greater ovarian cancer specific expression than the RNA for CA125.  Identification of CA125 was a remarkable feat in 1981, the new data should enable the creation of greatly improved marker sets as shown here.

    FIGURE 3 

    To demonstrate the power of the analytical tools we’ve created, we also generated signatures for Pancreatic Adenocarcinoma (figure 4A, B) and 6 other tumor types:  colorectal, endometrial, kidney (spanning clear cell, chromophobe and papillary), non small cell lung (spanning adenocarcinoma and squamous), melanoma and prostate (figure 4C).  

    This ovarian cancer-specific signature may represent a novel tool for the early blood based detection of malignant ovarian tumors.  The signature, however, represents only the first rate-limiting step towards the development of a blood based test for early ovarian cancer.

    Looking for partners to help test the ovarian (and potential other) signature(s)

    We are interested in identifying committed collaborators to help test the signature:

    1. Measuring signature genes in appropriate subsets of circulating nucleic acid compartments, e.g.  exosomes, cell free nucleic acids.
    2. Determining specificity and limit of detection in archived samples of both patients and unaffected individuals.
    3. Testing the signature in characterized archived frozen serum samples using a “prospective-retrospective” design used successfully for previous expression based diagnostic tests. 


    August 10, 2014.  We recently used our analysis tools to examine the genetics and genomics of cardiomyopathies and uncovered a new candidate pathway for common variants associated with cardiovascular phenotypes.  


    Network analysis of sixteen genes associated with Hypertrophic Cardiomyopathy (HCM) and thirty one genes associated with cardiovascular disease phenotypes in Genome Wide Association Studies (GWAS) identified three major classes of cardiovascular disease-associated genes.  Consistent with well established biology, one major group comprised sarcomere genes, however, the networks further grouped the known HCM genes into subsets by muscle type.  Networking generated a list of 250 candidate sarcomere associated genes for use in exploring sequence variants underlying HCM.  Analysis of gene lists from three GWA studies identified two additional classes, one comprising genes associated with endothelium/smooth muscle, and a more speculative brain expressed group.  Many candidate genes from GWAS loci did not fall into one of the three network classes.  These analyses demonstrate the utility of network analysis in prioritization and interpretation of GWAS data, provide candidate genes for assessing sequence variants in HCM, and identifies candidate biological processes linked to cardiovascular risk in the general population.



    Although recently surpassed by cancer as the leading cause of death in the US, Heart disease remains one of the major health issues faced by modern society (1).

    Decades of research have provided a good understanding of causes and mechanisms of cardiovascular disease, resulting in effective lifestyle recommendations and pharmacologic agents to reduce risks, especially regarding management of LDL and triglyceride levels. 

    While environmental contributors are well established, susceptibility to cardiovascular disease also exhibits a strong genetic component.  A better understanding of genes and pathways linked to cardiovascular phenotypes may provide new targets for intervention.  The recent flood of genomics data provides a great opportunity to develop new insight.

    Examples of the power of genomic information in cardiovascular drug discovery include PCSK9 (2, 3), first identified as a target in 2003 through positional cloning of familial cholesterol phenotypes (4, 5). ANGPTL3, another genomically validated target (6), awaits the development of targeted inhibitors; currently pursued by antisense and RNAi companies, but perhaps best targeted with an antibody as anticipated by this Regeneron patent application (7).  Recent functional analysis of 23 GWA studies on blood lipid levels may provide additional drug targets (8).

    Hypertrophic Cardiomyopathy

    In addition to lipid homeostasis, the genetics of other biological systems also impact cardiovascular disease.  Multisystem genetic disorders can often exhibit cardiomyopathies, such as some of the lysosomal disorders (LAMP2, PRKAG2, GLA/Fabry) (9), and Rasopathies, Ras/MAPK gene defects such as Noonan Syndrome (10).  Genetic Mitochondrial disorders also commonly present with cardiomyopathy (11). 

    The most common cardiovascular genetic disorder, Hypertrophic Cardiomyopathy (HCM) affects approximately 1:500 (0.2%) of the population, and typically results from dominant mutations in the structural components of contractile heart muscle (1213).  Symptoms in affected individuals are varied, many show no overt pathology.  This is the disorder responsible for the sudden death of young athletes during or just after vigorous exercise.  

    Muscle type-specific gene networks

    Drawing from the recently expanded GTEx human tissue data set ( we added a new set of networks generated using gene expression based correlation, this round focusing on organ type and tissue type specific networks.  Figure 1a shows a bathymetry plot of 12 networks, represented by 10 genes each, including: smooth muscle, skeletal muscle, cardiac muscle, brain, lung, pancreas, adipose tissue, liver, ribosomes, mitochondria, proliferation, endothelium and lysosomal biogenesis (seeded by TFEB).  Most networks exhibit little cross network linkage. To confirm tissue specific expression for the muscle networks, figure 1B shows GTEX portal generated expression plots for the seed genes used in the cardiac (NKX2-5), skeletal (CACNA1S)  and smooth muscle (MYH11) networks, demonstrating the extent of tissue specific expression of the seed genes in the GTEx RNAseq data. 

    Figure 1A  Biological networks bathymetry plot

    Figure 1B  Muscle subtype network seed gene expression in GTEx dataset

    Hypertrophic Cardiomyopathy gene networks

    Sixteen genes linked to HCM were selected from a recent review (14) and used as seeds for network generation.  Figure 2 shows a plot of networks generated using each gene as the seed, as well as networks for smooth, cardiac and skeletal muscle.  Consistent with the common sarcomere biology attributed to these genes, most networks exhibited considerable cross network linkage.  However, clear subsets emerged, where networks grouped based on what appears to be relative expression among the distinct muscle subtypes.  MYL2, MYH7, MYL3, MYOZ, MYH6, TCAP, TNNC1, ACTN2 and CSRP3 formed a group with dominant linkage to skeletal muscle, but significant cardiac muscle linkage as well.  ACTC1, NEXN, TPM1 and PLN formed a tight group characterized by strong linkage to smooth muscle.  MYBPC3, TNNT2 and TNNI3 formed a third tight group characterized by selective linkage to cardiac muscle.  We interpret these results to represent the normal, not pathologic networks for these genes as it is likely that most of the samples are from individuals without HCM associated gene mutations.

    Figure 2 Hypertrophic Cardiomyopathy gene networks form distinct subgroups linked to muscle type

    Sarcomere biology gene network

    There are 190 genes plotted in figure 2, each of which is likely associated with sarcomere biology.  Collating the top 50 network genes for each of the 19 seed genes in figure 2 resulted in a list of 250 unique genes.  A PDF of the gene list is available here.   We believe this represents a comprehensive list of sarcomere associated genes that includes known players, uncharacterized genes, and others with no reported link to muscle biology.   Such a list may support the interpretation of sequence variant analysis, considering the significant challenges in distinguishing between incidental and pathogenic variants (15). Five additional HCM associated genes that were not included in the network analysis, JPH2, LDB3, MYLK2, TTN and VCL are all present on the mathematically derived 250 gene list.

    Cardiovascular GWAS genes link to smooth muscle/endothelial and possible neuronal network

    In addition to the relatively rare genetic cardiomyopathies, genetic predisposition to cardiovascular disease across the general population is expected to be both common and multigenic.  Several large GWAS have been conducted linking a number of cardiovascular phenotypes to genetic loci.  The EchoGen consortium recently published 16 loci that exhibited significant association with at least one of five echocardiographic traits, including ventricular hypertrophy (16).  Five loci reached significance in the second phase.  To further explore all 16 loci for potential biological pathway links we first re-examined the candidate genes identified at each locus then generated networks for 18 genes associated with the loci.

    Most genes associated with the loci were either within the gene or the candidate gene was the only gene in the region.  For the 17p13 locus, each of the top three listed SNPs were most closely associated with SMG6, HN1L or SRR respectively, however, two of the implicated SNPs in 17p13 are correlated suggesting their may only one or two relevant genes.  We included all three genes in the networks.   Because PLN is a known gene in HCM, we didn’t pursue the other two candidates in the 6q22 locus (SLC35F1, C6ORF204/CEP85L).   SLC35F1 is an uncharacterized brain specific channel whose close network members are involved in neurogenesis (data not shown), CEP85L is more challenging to interpret, exhibiting broad expression with testis and brain highest; close network members involved in DNA damage response (ATM, TTBK2, TAOK1) and ataxias (ATM, SACS, QKI).   At the 5q23 locus (CCDC100/CEP120), the more recent genomic build 38 now positions the the top listed SNP adjacent to PRDM6 so we included it as well. 

    Figure 3 shows networks and linkages for the 18 selected genes and several relevant biological networks.  While half of the genes exhibit only modest network linkages, the other half formed two distinct biologically related groups.  PRDM6, PLN, PALMD, MEIS2 and PDE3A formed one meta-network that linked to both endothelium and smooth muscle networks.  PRDM6, a smooth muscle-specific transcriptional repressor (17 ) and PDE3A are both known to be involved in smooth muscle biology, and PLN, a HCM linked gene is a member of the smooth muscle sub-network in figure 2 above.  This enables a hypothesis that a subset of the genes associated with cardiac echocardiographic phenotypes function in smooth muscle/endothelial cell biology, and further that these genes are the relevant affected genes in their respective loci.  While PRDM6, PLN, PDE3A loci replicated in the study, PALMD (meta-analysis p value of 1.1 x 10-7) did not and the MEIS2 SNP is a rare allele.  This network supports continued exploration of PALMD and MEIS2 as well as PRDM6, PLN and PDE2A in future work.

    Figure 3 EchoGen GWAS candidate gene networks form subgroups that link to smooth muscle/endothelium and brain expressed networks

    A second meta-network shown in figure 3 comprises NOVA1, SRR, GRID1 and WWOX which exhibit neuronal expression and exhibit modest linkage to a brain specific network.   However, only the SRR locus replicated in the study.   In addition, many genes exhibit brain enriched expression patterns, so in the absence of other strong support caution is warranted.

    Endothelial/Smooth Muscle network seen in a second GWAS gene set

    An unpublished study on  Kaiser/UCSF GERA cohort was recently presented (18).  Ten genes from loci exhibiting significant association with Left Ventricular Hypertrophy were reported.  To further explore potential biological networks associated with cardiomyopathy phenotypes, we generated network plots for the reported genes:  CCDC141, CRIM1, CTNNA1, HAND1, NFIA, PROCR, SIPA1L1, TBX3, VGLL2, ZNF595 (figure 4).   One gene, VGLL2 exhibited good network linkage to skeletal muscle and was present in the sarcomere gene list generated above, aligning most closely with the skeletal dominant subgroup of HCM genes.   A known muscle differentiation gene, the association with hypertrophy suggests a role for VGLL2 in proper expression of muscle specific genes.  Of the other listed genes, CRIM1, CTNNA1 TBX3, NFIA, PROCR and SIPA1L1 formed a meta-network that also linked to endothelium and smooth muscle networks.  This second linkage of GWAS identified genes with endothelium and smooth muscle cell networks further supports these networks as relevant in cardiovascular disease.  

    Figure 4 GERA GWAS gene candidate networks form endothelium/smooth muscle linked subgroup  

    Sudden Cardiac Death Susceptibility locus 2q24.2 gene TANC1 is also a member of the endothelial/smooth muscle network

    A two stage GWAS for sudden cardiac death employed 4,400 cases and more than 30,000 controls to identify a single locus with candidate genes BAZ2B, TANC1 and WDSUB1 (22). BAZ2B is widely expressed (highest in cerebellum, ovary, testis) exhibiting a strong network association with chromatin modulators (data not shown).  WDSUB1, also widely expressed (highest in pituitary, ovary and adrenal gland) exhibits modest network associations to intracellular trafficking (data not shown).  TANC1 exhibits highest expression in skin, cervix, lung and artery, and exhibits strong network association with the GWAS candidate network identified above as well as endothelial/smooth muscle networks (figure 5). 

    Figure 5 TANC1 network joins GERA and EchoGen Endothelial/smooth muscle gene candidates


    June 3, 2014.  Last week, drafts of the Human Proteome were published in Nature (1).  These publications and associated data are an exciting addition to the rapidly expanding dataset on the human genome.  Curious to learn how these data fit with our work, we focused on one of the two publications (2).  The authors generated 25 million mass spectra yielding 293,000 peptides from 30 tissues and cells, and found peptides corresponding to 17,294 of the estimated 20,500 or so genes.  An accompanying website enables gene-based queries of the data (

    At Molquant, we are generating a comprehensive set of networks that organize the genome into biologically relevant groups using a diverse set of genomic data.  Using networks based on gene expression data, we used the web portal to explore how our networks look in the Proteome data, assessing coverage as well as consistency with our findings.

    The accessible data from the proteome website only allows heat map analysis of input genes, so we aren’t able to generate any numbers or correlation statistics.  Nevertheless, a visual representation of the protein levels of the queried genes provides a good sense of the degree of correlated expression, tissue specificity and extent of coverage of the queried networks in the proteome data.


    Molquant Networks

    Using correlation statistics, thousands of diverse transcriptome samples (tissues, cell lines, tumors) and biologically directed seed genes, we have been able to generate more than one hundred gene networks representing diverse biological pathways, subcellular structures, tissues and cell types.

    To assess the strength of the networks, and to explore inter-network relationships, we developed a heatmap plot (bathymetry plot) of the expression correlation value of each queried gene to every other gene in a set.  Note that the 10 selected genes representing each network exhibit a high degree of correlation; inter-network correlations can be “read” by examining the degree of correlation across networks.  

    One representative bathymetry plot is shown in figure 1, which plots 20 networks each represented by 10 genes.  Well-characterized genes were chosen to link the mathematically derived networks to known biology.   Pathways, subcellular structures and tissue specific networks are shown, labels above the 10 gene squares that represent each network.


    FIGURE 1

    We also visualize the expression of networks using expression heat maps and aggregate signature scores of the networks.  Shown in figure 2A is a heat map and score of 50 genes from a Molquant proliferation network across more than 1600 human tissue samples from the GTEx project.  In this case, GTEx represents a “validation set" since GTEx samples were excluded from network generation that created the gene list.  The highly correlated expression of the proliferation network in GTEx highlights the general applicability of the algorithms and samples used for network generation.


    FIGURE 2

    Molquant networks are also found in Proteome data.

    To compare the RNA based networks with the newly published proteome data, we first compared the proliferation network to a protein abundance heatmap generated from the proteome data (Figure 2B).  More than half of the genes from the Molquant network were prominently represented in the adult testis Proteome sample, consistent with both figure 2A and the known biology of the testis.  Fetal ovary also exhibited increased signal from many of the proliferation genes.  However, the proteome data available on the portal only reports a single point per tissue (presumably summed from the raw data) and does not reveal the variation and nuance observed in the GTEx RNAseq data, where each individual sample is represented.  For example, note the high proliferation signal in about 20% of the “blood” samples.  Nevertheless, the consistency of high testis expression of many of the proliferation network genes/proteins supported the quality and quantitative nature of the proteome data.  We recognize that it is an assumption that the Molquant proliferation network represents the gold standard, an assertion that is difficult to prove.  However, the preponderance of well-characterized proliferation genes in the mathematically generated network supports the validity of the gene list. 


    “Housekeeping” gene networks exhibit varied expression across tissues

    Many basic cellular functions such as translation, energy production and metabolism are common to virtually all cells.   Although we typically refer to the genes of these processes as “housekeeping” genes, one of the surprises from our network generation work was that these housekeeping networks exhibited significant variation across tissues.  Figure 3A shows an example for a mitochondria network where higher expression of the network is seen in heart and muscle, an observation consistent with the high energy requirement of contractile tissues.  The increased coordinate expression of the network is consistent between RNA expression and the Proteome data.  In Figure 3B, a Ribosome network exhibits coordinate but varied network expression across both RNAseq and proteome data, yet some discrepancies are seen between the two data sets.  For example, the proteome data shows high network expression in fetal liver and adult pancreas, whereas the RNAseq GTEx data exhibits highest network expression in prostate and subsets of blood.

    FIGURE 3

    Tissue specific networks compared across RNAseq and Proteome

    Continuing the comparison of Molquant derived networks between RNAseq and proteome data,  we next examined four networks generated to represent tissue specific/cell type specific functions.  Figure 4 shows comparisons of networks for Lung, Pancreas, Skeletal Muscle and Adipose cells.  For Lung and Pancreas (Figure 4A, B) the networks generally perform very similarly between RNAseq and Proteome data with nearly all of the queried genes exhibiting highest expression/abundance in the appropriate tissue.  Extending network comparisons for Skeletal Muscle and Adipose Cells (Figure 4C, D), tissues that were not directly collected in the Proteome set,  showed a much less consistent pattern in the Proteome data.  While the RNAseq data exhibited the expected tissue restriction, proteome data for muscle was enriched in fetal and adult heart, as well as in esophagus.  The latter tissue is perhaps enriched with smooth muscle constituents.   In contrast to the well conserved network patterns for the other gene sets, the adipose cell network did not exhibit coordinate expression of network genes in any of the proteome samples.   Adipose cells are often present within and “contaminate” a wide variety of tissues (note network expression enrichment in blood vessel, breast and nerve GTEx RNAseq samples),  a discernable signal is lacking in the the proteome data.  Among several possibilities, the precise process of sample collection in the proteome data may have excluded “contaminating” adipose tissue.


    FIGURE 4

    Molquant hematopoietic cell networks well defined in Proteome data

    Molquant tools can generate networks for most biology, we have also derived networks for many hematopoietic cell lineages (using well defined cell type-specific genes as seeds).  One feature of the Proteome data is the isolation and independent characterization of several hematopoietic lineages (2), whereas the GTEx dataset only collected “blood”, this dataset lacks the ability to discriminate distinct hematopoietic cell types.  Figure 5 shows a comparison of six hematopoietic networks (platelets, erythrocytes, neutrophils, natural killer cells, CD8 T cells, macrophages).  The GTEx dataset (figure 5A) exhibits strong network enrichment in “blood” as expected, however, it is difficult to recognize any sub-structure to the network expression.  The same networks examined in the Proteome data (figure 5A) exhibit enriched protein abundance signal accurately corresponding to the derived networks.  Erythrocytes exhibit enriched expression in fetal liver (the site of fetal erythropoiesis), and both neutrophil and macrophage networks show enriched network expression in “monocytes” as presumably no cell type selection was conducted on those samples.  

    FIGURE 5

    Robust performance of Molquant Networks across diverse human data types

    The gene expression based networks we've derived are drawn from thousands of diverse transcriptomes, yet are designed to accurately represent many aspects of human biology.  If successful, these networks become powerful tools to explore disease biology, gene function annotation and gene/phenotype linkages, all of which enable translational efforts in diagnosis, drug discovery and development.  However, RNA is an imprecise representation of protein (or post-translational modification...), so mapping RNA based networks onto the new Proteome data helps to assess the biological relevance of the tools we've built. The observation that most of the networks generated using transcriptomic data were also seen to exhibit correlated abundance in the proteomic data suggests that the algorithms and data used to generate the Molquant networks are robust and can be used for examining diverse genomic data types.


    1. Marx, V. Proteomics: An atlas of expression. Nature 509, 645-649 (2014).

    2. Kim, M.S., et al. A draft map of the human proteome. Nature 509, 575-581 (2014).


    April 10, 2014.  To demonstrate the power of the tools we’re building at Molquant, we recently took a look at Parkinson’s disease, a challenging and complex neurodegenerative disorder affecting nearly 10 million people worldwide. 

    Applying our correlation tools, we discovered previously unreported relationships among 16 known Parkinson’s disease genes.  Two major subgroups were found, which we could link to specific biological networks that we generated.  A picture emerged suggesting that the normal biology of these gene subgroups functioned to integrate nutrient and oxygen conditions with cellular processes that tend mitochondria and maintain appropriate levels of mitochondrial energy production in neurons. 

    This mathematically derived picture is consistent with the emerging view in the field, but further suggests that different Parkinson’s associated variant genes will impact different nodes of this biology.  Our analysis supports the recent hypothesis in the field that some diabetes drugs (thiazolidonediones) may be attractive therapeutic options for Parkinson’s disease, but also suggests that these agents may elicit different activity depending on a patient’s particular Parkinson’s associated gene variant.  In addition, the larger set of genes identified in this analysis represents a promising set of candidate genes to prioritize or expand the current Parkinson’s disease associated gene sets.

    The detailed manuscript can be downloaded PD Gene Networks Molquant.pdf with more figures, discussion, references for further reading.  An excerpt of the PDF follows:



    Biological Networks

    Using correlation algorithms applied to a large set of gene expression data, we developed a set of networks linked to biological processes.

    Figure 1 shows network and matrix “bathymetry” plots for 19 biologically defined networks.  In the network plot (1A), 10 genes for each network are plotted as nodes, correlation represent the length of the edges (edges not rendered here).  The bathymetry plot (1B) (a 2D matrix, so data are duplicated below the diagonal) shows correlation relationships among the 190 genes rendered as color and “height”. 

    Networks include:  vesicle transport, autophagy, axon motors, TFEB anchored lysosomal biogenesis, skeletal muscle, adipose tissue, ribosomes, PGC1a anchored mitochondrial biogenesis, mitochondria, glycolysis, mitosis, bone/cartilage, hypoxia/blood vessels, T cells, NK cells, CD68 macrophages, CD14 macrophages, neutrophils.  

    Figure 1A

    Figure 1B

    Parkinson’s Disease Gene Networks

    Sixteen genes in which genetic alterations were associated with risk of PD:  ATP13A2, GBA, FBXO7, LRRK2, MAPT, PARK2, PARK7, PINK1, SMPD1, SNCA, SYNJ1, UCHL1, and from an unpublished report:  NOVA2, OR56B4, PABPC1L, RPE65 (1).

    Each of the 16 genes served as the “seed” gene for network formation and the resulting network and bathymetry plots are shown in figure 2.  Two distinct network clusters emerged from the analysis:  The first, referred to as the PARK2 meta-network includes PARK2, LRRK2, PINK1 and the newly reported gene NOVA2.  The second, the SNCA meta-network includes SNCA, RPE65, SYNJ1 and MAPT.  Three other genes also exhibited correlation to these two clusters (FBXO7, SMPD1, and UCHL1), whereas PARK7, PABPC1L, OR56B4, GBA and ATP13A2 networks exhibited limited correlation linkage. Also, note that the two meta-networks exhibited correlation linkage to one another (figure 2B).  The data used to generate these networks are drawn from many sources; we interpret these networks to represent the normal biology of these genes, not the pathologic state associated with PD. 

    Figure 2A

    Figure 2B

    Links between Parkinson’s Gene Networks and Biological Networks

    To provide biological annotation of the Parkinson’s gene networks, biological networks and Parkinson’s gene networks were plotted together (Figure 3).

    The PARK2 meta-network exhibited correlation to four of the queried biological networks, a vesicle network seeded by VAMP2, a PGC1a (PPARGC1A) seeded network, a lysosomal biogenesis network seeded by TFEB, and a hypoxia network seeded by VEGFA. 

    These correlations enable a hypothesis that the non-pathogenic versions of LRRK2, PARK2, PINK1 and NOVA2 normally function in intracellular processes that control the integration of nutrient and oxygen sensing with the cellular machinery that controls mitochondrial homeostasis. The fact that the pathogenic forms of these genes are directly linked to disease suggests that the linked biological pathways are altered in PD.

    Figure 3

    The observations and hypothesis are consistent with several previous findings in the field. PARK2 and PINK1 have been shown previously to function in the control of mitochondrial homeostasis, including motility (2, 3), as has FBXO7 (4), which here exhibits modest correlation linkage to the PARK2 meta-network.  The transcriptional co-activator PGC1a/PPARGC1A is thought to be a central controller of this process (5, 6). In addition, meta-analysis of expression data from Parkinson’s affected tissues identified a PGC1a signature as one of the top altered pathways (7).


    Note that these biological processes are not directly related to the mitochondrion per se (figure 3 shows limited correlation linkage to integral mitochondrial genes), rather to components and regulators of intracellular organelle trafficking.  LRRK2 in particular, exhibits tight correlations to several cytoskeletal, vesicle/intracellular transport genes (COL4A3, COL4A4, SYNE1, ARHGAP24,31, GPRASP1, LMBRD1, SLC6A13, STX12, TRAK2/Milton), suggesting a role in modulation of intracellular transport.  Also consistent with the hypothesis, the Parkinson’s associated lysosomal gene SMPD1 also exhibits linkage to the PARK2 meta-networks.


    The SNCA meta-network also exhibited correlation to PPARGC1A (PGC1a) network, but the strongest linkages are to cytoskeletal motor protein networks: either the kinesin motor KIF1A, or the myosin motor gene myosin 5A (MYO5A).  The RPE65 network exhibited tightest linkage to the MYO5A motor network, whereas SYNJ1, a known synaptic vesicle protein exhibited tightest linkage to the KIF1A network; SNCA shared linkages with both motor networks.  Together, these observations further support a hypothesis that these Parkinson’s genes normally function as components of and regulators of axonal transport, which is intimately linked to mitochondrial homeostasis (8).


    The other queried Parkinson’s genes exhibit distinct linkage correlations.  PABPC1L, an RNA binding protein of unknown function, and to a lesser extent, ATP13A2 exhibit linkage to an autophagy network seeded by ATG4B.  PARK7 exhibits tight linkage to the integral mitochondrial protein network anchored by COX6A1 (as does the PD gene candidate HTRA2 (data not shown)).  Olfactory Receptor OR56B4 doesn’t exhibit significant linkage to any networks analyzed. Although not shown here, PD associated gene VPS35 exhibits weak linkage to SYNJ1, and PD gene PLA2G6 doesn’t exhibit significant linkage to any networks analyzed (data not shown).  Further analysis exploring a broader set of biological networks may reveal additional linkages.


    PGC1a, a potential protective network in some genetically defined forms of PD. 

    Consistent with the previously reported observation that a PGC1a associated signatures are altered in PD (7), this analysis hypothesizes a role for PGC1a in biological processes associated with some of the identified PD genes.  The PGC1a (PPARGC1A) network exhibits linkage to the both the PARK2 and SNCA meta-networks (figure 3).  Although several lines of evidence indicate that PGC1a is associated with mitochondrial biogenesis, networks here clearly distinguish between “intrinsic” mitochondrial networks (those comprising mitochondrial genes) and PGC1a associated networks.  The network analysis shown here postulates that PCG1a’s impact on mitochondrial biogenesis occurs through its modulation of these associated transport and sensor pathways. 


    Based on previous links between the thiozolidinedione drugs, PGC1a and PD, a phase 2 trial was initiated in 2011 (9).   Given the association of PGC1a pathways with specific subgroups of PD genes, assessment of the genotype of patients in such a study may provide additional information regarding genotype/efficacy relationships.

    Several lines of evidence have linked PGC1a with AMP kinase, a rheostat of nutrient availability that responds to AMP/ATP ratios in a cell.  The recently reported drug candidate R118, a potent AMPK activator (paper on related molecule 10) just entered the clinic, intended for treatment of peripheral artery disease.  If this molecule achieves adequate exposure and tolerability in human studies, it may represent an intriguing candidate for PD.


    Parkinson’s network gene lists

    In addition to providing tools for biological interpretation, the Parkinson’s gene networks identified here comprise candidate PD gene lists suitable for exploration of current and future GWAS analysis, or for targeted genotyping studies.

    Here, only examining genes of the two meta networks, figure 4 shows the overlap of the top 1000 network genes for each PD gene from either the PARK2 (figure 4A) or SNCA (figure 4B) meta-network.  The identified genes are those that correlate with at least three of the four PD genes in each network, and therefore hypothesized to participate in the same biological processes.  Here is a downloadable file of the gene set.

    Figure 4

    A                                                                                      B


     Intersection of these core network genes with the online PD Gene database resulted in 35 genes that are present both in the PARK2 and/or SNCA meta-networks and the PDgene database (figure 5). 

    Figure 5

    2014 chromosome 1 PD locus gene S1PR1 present in PARK2 meta-network

    A stratified analysis of PD GWAS data identified a previously unrecognized PD locus on Chr 1 that included candidate genes DPH5, OLFM3, S1PR1, SLC30A7, VCAM (11).  We note that the sphingosine phosphate receptor S1PR1 is included in the PARK2 meta-network.  S1PR1 expression also exhibits correlation to MAPT and SYNJ1 networks as well, providing a strong network support to hypothesize that  S1PR1 may represent the relevant affected gene for the locus.   S1PR1 is widely expressed, exhibiting multifunctional roles including lymphocyte trafficking (the target of the anti-MS drug fingolimod) and angiogenesis.  The network linkages observed here support a role for S1PR1 in PD gene associated function, including the previously mentioned hypoxia network.



     Only transcript level regulation is captured

    Although transcript level regulation is widespread, it represents only one of many ways to regulate biological processes.  Genes that are widely expressed and exhibit a high degree of post-translational regulation may not be captured in these networks.  Nevertheless, our observations have been that, given a large enough sample size with sufficient biological complexity, i.e. many diverse samples, we can obtain networks exhibiting tight correlation of genes comprising well known biological processes (glycolysis, mitochondria, ribosomes, skeletal muscle, adipose tissue, taste receptor cells and others).  However, several other attempts to generate correlation based biological networks were unsuccessful (circadian rhythm genes, Notch pathway genes, retina specific genes, nonsense mediated RNA decay and others), highlighting the challenges inherent in expression based networks.   Certainly, the identified networks will not comprehensively capture all of the component genes of a particular biology, but the observations presented here argue that this approach represents a useful tool to analyze and interpret complex biology.

    Untangling co-occurring biological processes

    Many biological processes, although distinct, occur in the same tissues or conditions.  If an uncharacterized gene correlates with two (or more) biological networks, it is difficult to assign it to a specific process.  For example, although we recognize proteasomal degradation, splicing, and DNA replication to be distinct biological processes, networks comprising core genes from each of these functions exhibit a high degree of cross correlation in our datasets (data not shown). Uncharacterized gene C20ORF24 correlates highly with networks seeded by proteasomal gene PSMA5, pre-replication complex gene MCM10, and splicing factor SNRPB, precluding assignment to any of those processes.  The co-regulation of these three networks probably tells an interesting story itself, but this analysis merely identifies the phenomenon, not its resolution.

    Data sources

    Although more than 7000 individual profiles were explored to derive the correlations shown here, the networks identified are subject to the composition of the samples included.  Many genes are known to exhibit distinct networks in different tissues, so the correlations observed here likely reflect what is occurring in the dominantly represented tissue(s) or cells.  For example, PGC1a may exhibit different correlations in brown adipose tissue or cardiac cells or neurons.  Furthermore, our examination of network correlations in large subsets of the overall dataset (e.g. 1600 GTEx RNAseq of human tissues ) demonstrates that while many networks are well replicated in the subsets, other networks are not well recapitulated.  We look forward to the rapidly increasing number of profiled samples to improve the networks.



    Molquant tools enable the exploration of additional Parkinson’s Disease candidates through the intersection of Disease Biology Networks and Parkinson’s Disease GWAS PDGene database (  Molquant, 2014

    Genomics research continues to expand ever more rapidly, with new studies linking genes to phenotypes coming out continuously. However, as more studies emerged, researchers began to recognize a problem: the “Missing Heritability.”

    We have long been able to estimate the genetic underpinnings (heritability) of a trait by assessing the correlation between the trait and relatedness. If a trait is genetically linked (i.e. if it runs in families), then the more closely related two individuals are, the more likely they will share that trait (Twins>Siblings>Parents>Grandparents>First Cousins and so on).  A trait like height has been estimated to have 80% heritability, that is 80% of your potential height is driven by your genes, 20% by diet and other factors.  One of the largest GWAS efforts yet conducted looked at height, and they found 180 gene variants that link to height.  However, when they add up the contribution of all those genes, they contribute only about ~10% of the known variation in height.  Where is the other 70%? (THAT’S the ‘Missing Heritability’!) We know height is inherited; the genes must be in there somewhere, but we haven’t been able to find them. (This explains the limited interpretation abilities of Illumina and 23andMe).  

    This mystery in genomics has generated many competing hypotheses as to why we haven’t done better at linking genes to traits.  It’s still controversial (e.g.1) , but the data are leaning toward the “additive model,” the notion that many genes contribute to any given phenotype, but each gene contributes only a little (2,3; a thoughtful discussion 4) (and a recent ASHG 2013 talk by MIT’s Yaniv Erlich in which his analysis of 44 million entries in public ancestry data support the additive model).  When the field of genomics matures, we will have access to genes and phenotypes from millions of people. This would help to surmount the missing heritability problem by sheer numbers (as suggested by this yeast model paper (5). Still, limitations in our analysis assumptions (e.g. how do we decide whether an association is significant), as well as complexity in population structure will additionally create analysis challenges. Can we make significant progress until we sequence and analyze millions of people?

    Molquant’s Approach to Genome Scale Data - A New Framework

    Taking a somewhat orthogonal approach, we believe that there are real opportunities today to better interpret genome scale data.  The dominant statistical approach of correlating millions of individual variables to a ‘yes’ or ‘no’ phenotype carries with it some fundamental limitations. Further, biology doesn’t appear to be organized in a one SNP/one trait manner.  Using a “machine” analogy, we know that each gene typically makes one “part” that assembles into a larger machine, where each part needs to be functioning properly for the machine to work optimally. Most biological systems that have been studied conform to this model, where each gene makes a protein that fits into a multi-protein complex (Ribosomes, transcriptional regulation, DNA replication/repair, mitochondrial complexes, etc.)

    Molquant plot assessing the relationships among a set of mathematically derived biological networks

    If we can develop reasonable algorithms and apply them to the rapidly growing tangle of data, previously unappreciated patterns emerge.  Our work so far has demonstrated this: using biologically-seeded statistical approaches Molquant has been able to organize large amounts of genome-scale data into biologically relevant groups or networks. Unlike curated annotation groupings like Gene Ontology (GO) terms, these networks are mathematically-derived and comprise both well-studied genes as well as uncharacterized genes.

    These networks provide powerful tools for a number of analytical activities.  They enable biological annotation of genes and disease processes, which is especially useful for providing biological hypotheses around disease genes identified through positional cloning or sequencing efforts.  Applied to gene::trait associations, these networks provide a novel framework for assessing the role of a particular SNP in a trait. Such a framework makes a biological assumption, that the majority of gene variation associated with a particular trait will not be randomly distributed across the genome, but will be found in the biological networks relevant to that trait.  Thus, at least some of the “missing heritability” of GWAS may be found in variants among the networked genes, which contribute to the phenotype, but wouldn’t meet the stringent statistical rigor imposed on single SNPs.  

    An example in which such a network emerged from straight GWAS analysis can be found in studies of autoimmunity.  Meta-analysis of SNPs identified across multiple autoimmune disease studies found that most of the risk genes functioned in control of the immune response (e.g. cytokines, T cell activation, antigen processing; a nice summary (6)) No surprise given our understanding of these conditions, but imagine if little were known about this disorder, biological networks comprising immune cell regulation genes would represent a highly enriched source for identifying relevant SNPs.

    As a proof-of-concept, we recently completed an analysis of Parkinson’s Disease, where there are several known susceptibility genes, but much less is known about the biology most of these genes. Using the Molquant tools, we identified networks among the various Parkinson’s associated genes and then further linked these networks to known biology.  We look forward to posting this analysis in the near future.

    Follow @molquant on Twitter to receive our news and updates.


    1 Charney, E., Still Chasing Ghosts: A New Genetic Methodology Will Not Find the “Missing Heritability” 

    2 Yang et al., “Genome partitioning of genetic variation for complex traits using common SNPs” Nature Genetics 2011 Jun; 43(6):519-25

    Greg Gibson, “Hints of Hidden Heritability in GWAS” Nature Genetics, Volume 42, Number 7. July 2010

    4 Greg Gibson, "Rare and common variants: twenty arguments"  Nature Reviews, 13, pp 135-145, 2012

    5 Bloom et al., “Finding the Sources of Missing Heritability in a Yeast Cross” Nature, 494, 234–237, 14 February 2013
    6 Visscher et al., “Five Years of GWAS Discovery” American Journal of Human Genetics, 2012 January 13; 90(1): 7–24

    Luke JostinsEstimating heritability using twins” Unzipped Genes Blog. Dec 13, 2010
    Brendan Maher, “Personal Genomes: The Case of the Missing Heritability” Nature. 456, 18-21 (2008)
    What’s Missing in Missing Heritability” A Molecular Matter Blog. Jan 14, 2013