|Home | About | Journals | Submit | Contact Us | Français|
Thyroid hormone (TH) is crucial for normal brain development. TH transporters control TH homeostasis in brain as evidenced by the complex endocrine and neurological phenotype of patients with mutations in monocarboxylate transporter 8 (MCT8). We investigated the mechanisms of disease by analyzing gene expression profiles in fibroblasts from patients with MCT8 mutations. Studying MCT8 and its transcriptional context in different comprehensive spatial and temporal human brain transcriptome data sets revealed distinct region-specific MCT8 expression. Furthermore, MCT8 demonstrated a clear age-dependent decrease, suggesting its importance in early brain development. Performing comparative transcriptome analysis, we linked the genes differentially expressed (DE) in patient fibroblasts to the human brain transcriptome. DE genes in patient fibroblasts were strongly over-represented among genes highly correlated with MCT8 expression in brain. Furthermore, using the same approach we identified which genes in the classical TH signaling pathway are affected in patients. Finally, we provide evidence that the TRα2 receptor variant is closely connected to MCT8. The present study provides a molecular basis for understanding which pathways are likely affected in the brains of patients with mutations in MCT8. Our data regarding a functional relationship between MCT8 and TRα2 suggest an unanticipated role for TRα2 in the (patho)physiology of TH signaling in the brain. This study demonstrates how genome-wide expression data from patient-derived non-neuronal tissue related to the human brain transcriptome may be successfully employed to improve our understanding of neurological disease.
Thyroid hormone (TH) is an important factor for brain development and function (1). The transcriptional actions of the bioactive hormone T3 are mediated by binding to its nuclear receptors (TRα1 and TRβ, encoded by THRA and THRB, respectively) (1). Cellular action of TH requires transport across the plasma membrane, which is facilitated by plasma membrane transporters (2,3).
Monocarboxylate transporter 8 (MCT8) is the most specific TH transporter known to date. Among different tissues, it is significantly expressed in brain (4). Its biological relevance was demonstrated in humans with psychomotor retardation and disturbed serum TH levels caused by mutations in MCT8 (5,6).
The molecular mechanisms underlying the neurological abnormalities in patients with mutations in MCT8 (for brevity hereafter called MCT8 patients) are unknown, because (i) knowledge of the normal spatial and temporal expression pattern of MCT8 in human brain is largely lacking; (ii) genes downstream of MCT8 have not been identified; (iii) Mct8 knockout (KO) mice replicate the abnormal thyroid state of MCT8 patients, but lack neurological abnormalities and are therefore not a suitable model to unravel the pathogenesis of the psychomotor retardation in MCT8 patients (7,8).
In an attempt to overcome these limitations, we first examined the gene expression profiles of human fibroblasts from MCT8 patients and controls. Next, we assessed MCT8 expression in the brain using several comprehensive human brain transcriptome data sets. Finally, comparative analysis of the transcriptomes of fibroblasts from MCT8 patients and the human brain specified the deranged molecular signature of the disease. Furthermore, the data suggest the role of the TRα2 receptor variant in the (patho)physiology of TH signaling.
The present study illustrates how genome-wide expression data of patient-derived non-neuronal tissue related to the human brain transcriptome may be successfully employed to improve our understanding of neurological diseases.
To assess differences in gene expression levels per se and in different thyroid states, we performed microarray analysis on fibroblasts from MCT8 patients and controls under three different conditions. We hypothesized that T3 addition would amplify differences between patients and controls, because intracellular T3 availability is decreased in fibroblasts lacking functional MCT8. Unexpectedly, PCA analysis indicated that differences between patients and controls are much larger than differences induced by T3. The genes that differed in expression after T3 treatment between patients and controls are presented in Supplementary Material, Table S1. They include AKR1C1-3 and RCAN2/ZAKI4α which have been previously identified as T3-responsive genes in fibroblasts. Validation of representative genes by qPCR is shown in Fig. 1A.
As only a limited number of genes responded differently to T3 treatment in patients versus controls, we focussed on the differentially expressed (DE) genes irrespective of thyroid state. Applying stringent selection criteria (see Supporting Information (SI) text for details), we identified 1617 unique DE genes (equaling 2159 probe sets; Supplementary Material, Table S2). Hierarchical clustering of DE genes reveals large differences between patients and controls per se and only subtle effects of T3 during 24 h incubation (Fig. 1B and C). Well-known T3 target genes in the brain such as RELN, BDNF and PTGDS were markedly decreased in MCT8 patient fibroblasts, and Ingenuity Pathway Analysis revealed an over-representation of genes in the TR/RXR pathway (data not shown). Together with our previous observations that TH uptake is impaired in fibroblasts from MCT8 patients (9), the present results suggest that fibroblasts are a suitable tool to study molecular derangements in MCT8 patients.
Among the DE genes, 305 unique genes were at least twice represented with different probe sets on the microarray, strengthening the present findings. To validate the microarray results, we confirmed significant changes in the expression of most of the selected genes by qPCR (Fig. 1D).
Gene ontology (GO)-enrichment analysis was used to functionally characterize the DE genes (10). Figure 1E shows the most prominent biological processes identified, including specific GO categories related to cell adhesion and actin cytoskeleton processes.
Next, we assembled different microarray data sets generated from in total 381 human brain samples to explore the spatial and age-dependent expression profiles of MCT8. Data set 1 was obtained using 169 tissue samples from 30 individuals and represents gene expression profiles in 53 anatomically distinct human brain regions from at least three subjects [average age: 39.0 (19–90) years]. Data set 2 consisted of 173 samples (hippocampus, entorhinal cortex, superior-frontal gyrus and postcentral gyrus) from 55 individuals (average age 64.4 (20–99) years]. Data set 3 consisted of dorsolateral prefrontal cortex samples from 39 individuals (age range: 0.1–83 years). Full details are given in Patients and Methods and SI text). Figure 2A shows the expression pattern of MCT8 in data set 1, revealing significant MCT8 expression in distinct brain areas. As MCT8 appears crucial for the brain development, we investigated whether MCT8 expression changed with advancing age. Using data set 3, demonstrating distinct age-related gene clustering (Supplementary Material, Fig. S2A), MCT8 is prominently ranked among the subset of genes with a clear age-dependent decrease in expression (Supplementary Material, Fig. S2B). MCT8 expression rapidly declines after the first years of life (Fig. 2B). To assess the functional context of genes among which MCT8 was expressed, we investigated which GO categories were over-represented in genes correlated (r≥ 0.4) with MCT8 in the different data sets. As might be anticipated, consistently enriched GO terms were related to neuron development, axonogenesis and cell adhesion (Supplementary Material, Table S3).
Then, we sought to integrate the data obtained in patients' versus control fibroblasts with the normal human brain transcriptome to define molecular pathways and specify candidate genes, which are possibly deranged in the brains of MCT8 patients.
Our strategy entailed the following steps. First, to avoid potential bias, we removed from the analysis genes which were not expressed in at least one of the brain samples. Second, for all remaining genes in the brain data sets, correlations with the expression of MCT8 (or in principle any gene of interest) were calculated. Third, we selected genes that were correlated (r≥ 0.4) with MCT8 expression. Based on the association-by-guilt principle, the resulting list was likely enriched in MCT8-related genes. Fourth, we intersected this list of selected genes in brain with the DE genes of fibroblasts from MCT8 patients. Fifth, the difference between the observed versus expected frequency by chance was calculated for the overlapping genes and a P-value was calculated for the fold enrichment. This approach was performed in all data sets and was extended to the TRs.
Figure 3 illustrates this approach for MCT8 in data set 1. Out of the 44,634 genes called present in the brain samples of this data set, we selected 1172 genes that highly correlated with MCT8 expression in brain (Fig. 3A and B). From the list of 1172 selected genes, 102 were identified among the DE genes in patients' versus control fibroblasts, representing a 2.1-fold enrichment (P = 1.9 × 10−12; Fig. 3C and D, Supplementary Material, Table S4). In contrast, the same approach with 1172 randomly selected genes did not result in an over-representation.
Thus, these findings indicate that genes that highly correlate with MCT8 expression in the human brain are more than twice as likely to be found among genes that differ between fibroblasts from MCT8 patients versus controls. This approach was validated using data set 2 (2.4-fold enrichment; P = 4.2 × 10−19; Supplementary Material, Table S4). Similarly, genes following the age-related pattern of MCT8 expression in brain were over-represented among the DE genes (data set 3; 1.5-fold enrichment; P = 1.7 × 10−4; Supplementary Material, Table S4).
Since MCT8 expression varied highly among the different brain regions (Fig. 2A), background noise might obscure the number of identified genes likely associated with MCT8. Therefore, we selected the brain regions with high versus low MCT8 expression from data set 1 and excluded samples from subjects >60 years of age to diminish age-related effects (data set 4). Theoretically, this data set should considerably increase the number and accuracy of MCT8-associated genes. Indeed, a total of 433 out of the 4152 selected transcripts were present among the DE genes, equaling a 3-fold enrichment (P = 2.9 × 10−67; Supplementary Material, Table S4). Again, categories involved in cell adhesion and actin cytoskeleton regulation were significantly over-represented.
Then, we aimed to identify genes which are consistently correlated with MCT8 among the different data sets and are DE in fibroblasts from MCT8 patients. Since the different data sets represent different brain areas, this approach precluded identification of brain region-specific genes. The search for MCT8-correlated genes in all four data sets resulted in 224 genes, of which 28 were DE in fibroblasts from MCT8 patients (2.7-fold enrichment; P = 6.5 × 10−7; Table 1). Interestingly, among these genes, two members of the plasticity-related gene family (LPPR2 and LPPR4) are highly expressed in neurons and have been implicated in axonal outgrowth (11).
Taken together, the approach of selecting genes in the brain which correlate with MCT8 expression and subsequent intersection with DE genes in MCT8 patient fibroblasts resulted in significant over-representations using all different data sets (Table 2).
We tested the current paradigm which assumes that MCT8 transports T3 across the plasma membrane, ultimately to ensure nuclear availability of T3 for the TRs. Thus, we selected genes which expression correlated highly with both MCT8 and one of the TR isoforms. A clear clustering of different probe sets for each isoform was observed in all human brain transcriptome data sets (Supplementary Material, Figs S3 and S4). Specific patterning of THRA1 and THRB was visible in line with the known expression of these genes in the brain, providing further confidence in this approach. For example, a much lower THRB expression was noticed in cerebellum compared with cortical regions (data not shown), which corresponds well with previous studies (12).
Since TRβ appears functional in fibroblasts (13), we started selecting genes that correlated with both MCT8 and THRB expression in the human brain. These lists were intersected with the DE genes in patient fibroblasts resulting in a ~2-fold enrichment (Table 3, Supplementary Material, Table S5). An exception was data set 3, in which no genes passed the selection criteria, which is explained by the observation that THRB does not exhibit an age-related expression pattern similar to MCT8 (Fig. S4). Examining those genes likely downstream of MCT8 and TRβ revealed several known T3-responsive genes (e.g. CCND2 and RELN). Furthermore, this approach also suggested many interesting candidates for further investigation such as the LIM homeobox gene LHX2 and THY1 and a hypothetical protein encoded by the locus FLJ90757 (14).
Next, we intersected genes correlated with both MCT8 and THRA1 expression in human brain with the DE genes in patient fibroblasts. As this query resulted in only a few genes in data sets 1–3, enrichment did not reach significance (data not shown). This could be explained if MCT8 and TRα1 co-expression is not brain-wide, but restricted to a few distinct brain areas. Indeed, in data set 4, genes that correlated with MCT8 and THRA1 expression were highly over-represented among the DE genes (3.7-fold enrichment; P = 1.6 × 10−11; Supplementary Material, Table S7). In addition, we performed a similar approach using a probe set which specifically detects THRA1 mRNA with an extended 3′-UTR (Supplementary Material, Fig. S3) also resulting in a strong enrichment (11-fold; P = 4.4 × 10−23; Supplementary Material, Table S6).
Taken together, our results specify genes in a functional pathway putatively downstream of MCT8 and each of the T3-binding TRs, which are dysregulated in MCT8 patients.
The THRA1 and THRA2 mRNAs encoding TRα1 and TRα2, respectively, result from differential splicing of the THRA gene (Supplementary Material, Fig. S3). Scrutinizing the different probe sets, we found that two probe sets specifically target the mRNA coding for TRα2 (designated as THRA2; Supplementary Material, Fig. S3). Blasting the sequence of these probes in an independent sequence database verified this finding. In addition, the intensities of the probe sets for THRA2 were 3–10-fold higher (P < 1 × 10−37 in data sets 1–3) than for THRA1, which is in perfect agreement with the previous findings of higher THRA2 than THRA1 expression in the brain (12,15). In contrast to TRα1 and TRβ, T3 does not bind to TRα2, leaving it an ‘orphan receptor’ (16).
The first striking observation was that THRA2 expression closely followed MCT8 expression in human brain, much more so than the other TR isoforms (Supplementary Material, Fig. S4, Table 4). To test the hypothesis that MCT8 and TRα2 are functionally linked, we selected genes correlating with both MCT8 and THRA2 expression in the brain. After intersection of these lists with the DE genes in patients' fibroblasts, strong and consistent enrichments (Table 5, Supplementary Material, Table S7) were found, which is illustrated in Fig. 4. We first separately identified 4152 and 4877 genes which highly correlated with MCT8 and THRA2 expression, respectively. Intersection of these lists with the DE genes in patients' fibroblasts yielded 433 and 552 overlapping genes for MCT8 and THRA2, respectively (Fig. 4A and B). Importantly, the large majority of genes were identical between these two lists (Fig. 4C, Supplementary Material, Table S7). Survey of this list uncovered not only genes important for cerebral function (e.g. NPDC1, PEA-15 and ZNF706) but also numerous genes whose function in the brain is as yet unknown (e.g. JOPS2 and CTXN1). The number of DE genes that overlapped with the genes in brain selected on the MCT8/THRA2 combination was much larger than on all other possible combinations (Fig. 4D).
As these observations point to an unanticipated role of TRa2 in the (patho)physiology of TH signaling, we aimed to substantiate the robustness of our strategy. Therefore, genes whose expression changed only after T3 treatment in control fibroblasts were retrieved (Supplementary Material, Table S1 and Fig. S5). This list of genes was searched for correlation with the TR isoforms in brain. As might be expected, the list of T3-responsive genes was not enriched for genes highly correlated with THRA2 in brain. In contrast, a significant enrichment was noted for genes correlated with THRA1 or THRB in multiple human brain data sets (Supplementary Material, Table S8). Furthermore, a comparative analysis of T3-responsive genes in fibroblasts identified previously (13) yielded a preferential enrichment in the THRB-correlated genes. In addition, several well-known T3-responsive genes identified in the mouse brain (HR, NGRN/RC3, BTEB and BDNF) all markedly correlated with THRA1 or THRB in the human brain. Thus, T3-responsive genes were mapped to the T3-binding TRs in human brain, substantiating the validity of our approach. Finally, we performed similar analyses with MCT8 and nuclear receptors unrelated to TH signaling to substantiate our method. In contrast with the TRs, only a few genes correlated with both MCT8 and either the glucocorticoid receptor or the estrogen receptor, none of which resulted in a significant enrichment among the DE genes. Taken together, these findings strongly suggest that MCT8 and TRa2 are closely linked and may act together in controlling specific gene sets.
Understanding the neurological phenotype in MCT8 patients is currently hampered as knowledge of MCT8 expression in the human brain is limited. A few studies reported MCT8 staining in some parts of the brain, but a systematic analysis of MCT8 expression in the human brain has not been performed yet (17–19). In the present study, we explored the MCT8 mRNA expression across multiple human brain regions, using genome-wide expression profiles. The present results indicate that MCT8 is significantly expressed in the limbic system, cortical regions involved in language function and areas important for motor control. Although confirmatory studies are needed before anatomical MCT8 expression may be linked to the disease phenotype, the affected cognitive functions in patients are comprehensive if MCT8 is normally expressed in these regions in unaffected subjects. Species-specific cerebral patterning of MCT8 may partially explain the differences between humans and mice lacking functional MCT8 (8,20). Furthermore, we noted that MCT8 has a strong age-dependent expression pattern in humans, at least in the dorsolateral prefrontal cortex, which provides a basis for the hypothesis that treatment during the first years of life may improve the neurological development of MCT8 patients (21).
In the present study, we analyzed and compared the transcriptome of fibroblasts from MCT8 patients with the transcriptome of normal human brain. Microarray studies have advanced the understanding of numerous processes in normal and diseased brain (22–26). Although conventional microarray analysis has yielded novel insights, important information may be missed. Therefore, integrating data from different sources may truly increase further our understanding (26).
Recently, several studies have improved the understanding of the complexity of human brain by the application of network analysis on gene expression profiles in different brain regions, demonstrating intrinsic organizational transcription patterns (27–30). An intuitive disadvantage of analyzing gene expression in brain is that cell-type-specific gene expression levels are missed. However, Oldham et al. (29) clearly demonstrated the existence of specific modules comprising co-expressed genes, which enabled to distinguish cell-specific gene expression signatures in bulk brain tissue. Correlative analysis of gene expression in the brain more precisely defined specific functional units than based on classical morphological markers (31).
In many neurological disorders, especially rare syndromes, diseased tissue is not available for analysis. Therefore, genetically engineered mice are frequently used as models for disease. However, this approach is often limited by species differences in brain function. Alternatively, peripheral cells obtained from patients with neurological disease have been used to seek for abnormalities, which is limited by the uncertainty to which extent differences reflect abnormalities in the brain (32,33). Obviously, the absence of neurological abnormalities in the Mct8 KO mice precludes their use as models for the human brain phenotype (7,8). We aimed to overcome these difficulties by integrating and comparing the DE genes in patient cells with the human brain transcriptome. To our knowledge, the present study is the first linking a human model of neurological disease (MCT8 patients) with different large human brain transcriptome data sets from as many as 381 samples.
The current working hypothesis assumes that in MCT8 patients transport of T3 across the blood–brain barrier and/or transport of T3 into central neurons is hampered. From this, it is commonly speculated that T3 target genes important for brain development are dysregulated. However, proof of this concept is currently lacking and the molecular mechanisms underlying the neurological phenotype are elusive. By using patient-derived fibroblasts, we now demonstrate that numerous genes are dysregulated in patients. We first focussed on the direct effects of T3 addition to the cultured fibroblasts. Different T3 exposure levels resulted in distinct clustering of responsive genes with modest differences between patients and controls.
The modest response to T3 is partially related to the strict cut-off values used to exclude false-positives. Approximately half of the genes responding to T3 treatment in the present study were reported previously as T3 responsive in fibroblasts or other cells and tissues, supporting the validity of our findings (9,13,34–37). Furthermore, 24 h T3 treatment may be insufficient to induce large changes in T3-responsive genes in fibroblasts (36). It is likely that fibroblasts from MCT8 patients are in a chronic hypothyroid state, which was also suggested by the enrichment of the TR/RXR pathway among DE genes irrespective of the cultured thyroid state. Interestingly, we observed reduced levels of several genes which are known to be T3-responsive in brain. Other classical genes which are regulated by T3 in brain (MBP, NGRN/RC3 and HR) were very low or not expressed in fibroblasts, precluding T3-dependent regulation. Nevertheless, the possibility remains that MCT8 transports an as-yet-unknown ligand the absence of which affects the gene expression in patient fibroblasts.
The comparative analysis with genome-wide expression data in the human brain enabled us to explore the relevance of DE genes in patient fibroblasts for brain. The list of down-regulated genes included genes that cause mental retardation when disrupted, such as the SLC6A8, GAMT, MECP2 and NLGN4X (38). Among the up-regulated genes in MCT8 patients, the Notch-ligand JAG1 and Notch-signaling targets HES1 and HEY1 are of particular interest, because these genes are important transcriptional repressors of neurogenesis (39). Interestingly, many of these genes (e.g. HES1, HEY1, JAG1, GAMT and NLGN4X) are strongly correlated with MCT8 expression in one or more of the studied brain transcriptome data sets. This specifies them as genes likely affected in brains of MCT8 patients.
An over-representation of genes related to cell adhesion and actin skeleton processes was noted among the DE genes. As these processes are intimately linked and act together in synapse formation, dysregulation of these molecules likely affects brain development (40). Furthermore, this approach identified a number of DE genes which are highly correlated with both MCT8 and the classic TRs in brain, providing evidence for the hypothesis that deficiency in MCT8 results in insufficient T3 supply to the TRs and consequent dysregulation of T3-responsive genes. Furthermore, this method suggests the concrete identities of putative T3 target genes (e.g. CTXN1 and DEAF1), which likely underlie the neurological abnormalities as examples for future studies.
The most unexpected observation was the strong link between MCT8 and the TRα2 receptor variant. In contrast to TRα1 and TRβ, it has been established that TRα2 does not bind T3 (15,41). TRα2 is able to inhibit effects mediated by either TRα1 or TRβ, but it is unknown whether this is a biologically relevant mechanism (41,42). There are several indications to believe that TRα2 has physiological relevance beyond this antagonistic mechanism of action. First, TRα2 has a ligand-binding domain. Although it does not bind T3, this does not preclude the binding of other ligands to TRα2 which could thus be regarded as an orphan receptor. Second, brain expression of TRα2 is roughly 10-fold higher than TRα1, whereas this ratio is lower in other tissues (12). The divergent expression pattern between cerebral TRα1 and TRα2 may hint at a specific role for TRα2 in brain.
The present results add several lines of evidence suggesting a direct role of TRα2 in TH signaling in brain. First, the expression patterns of THRA2 and MCT8 were highly similar in the human brain transcriptome. Second, the functional link between MCT8 and TRα2 was suggested by our comparative analysis which demonstrated that much more DE genes in patient fibroblasts were among the selection of genes in the human brain correlated with the combination MCT8/THRA2 than with the combinations MCT8/THRA1 or MCT8/THRB. The scenario in which TRα2 antagonizes TRα1 or TRβ function does not explain the present results. Rather less than more DE genes would be expected to overlap with the selected genes in the brain correlating with the combination MCT8/THRA2.
We like to emphasize that T3-responsive genes among the DE genes in patient fibroblasts clearly overlapped with the selected genes in brain correlating with the combination MCT8/THRA1 or MCT8/THRB, but not with MCT8/THRA2. These data are perfectly in line with the current dogma of T3 binding to the classic TRs and provide further confidence in our analysis method. What could be the functional relationship between MCT8 and TRα2? It may be speculated that MCT8 transports a ligand across the plasma membrane, which subsequently binds to TRα2, thereby controlling gene transcription. Therefore, future studies may address the possible MCT8/THRA2 relationship by extending co-expression and functional studies, for example, using animal models such as the TRα2−/− mice (43). Taken together, our findings indicate a novel role for TRα2 in TH signaling in the brain, independent of the actions mediated by the T3 receptors.
In conclusion, the present study demonstrates for the first time which genes are DE in cells derived from patients with mutations in MCT8. Furthermore, we analyzed the expression of MCT8 and its transcriptional context in several comprehensive data sets of the human brain transcriptome. Comparative analysis of the transcriptome of cells from MCT8 patients with the human brain transcriptome enabled us to identify genes that are putatively dysregulated in the brains of MCT8 patients. Our approach indicated a functional relationship between MCT8 and TRα2. This may have not only implications for the understanding of the disease, but also suggests an important role for the ‘orphan receptor’ TRα2 and its relevance for TH signaling in the brain. Finally, this study exemplifies how genome-wide expression data from patient-derived non-neuronal tissue related to the human brain transcriptome may be successfully employed to improve our understanding of neurological diseases.
Microarray analysis was performed on skin fibroblasts from patients with MCT8 mutations, which have been described previously (9). Complete data sets are available online as GEO entry GSE20538. The microarray data of data set 1 were obtained from Gene Logic Inc. (Gaithersburg, MD, USA). The microarrays used to assemble data sets 2 and 3 were collected from two published studies (44,45). Raw data from these studies are available in GSE11512 (data set 2) and GSE11882 (data set 3). For additional information on microarray processing, see SI text. qPCR was performed with SYBR Green I (Eurogentec, Maastricht, The Netherlands). For details, see the SI text.
Functional annotation was done using Ingenuity Pathway Analysis (Ingenuity, Mountain View, CA, USA) and the web-based DAVID program (david.abcc.ncifcrf.gov) (10).
A Pearson-based correlation metric was used to compare the seed genes (MCT8 and the TRs) with all other probe sets present on the chip. Selected genes were intersected with DE genes in fibroblasts from MCT8 patients versus controls. Likelihood and significance of the overlap were tested using χ2-test (or the Fisher exact test if appropriate).
W.E.V. and T.J.V. designed research; W.E.V., Z.Ő., W.F.J.v.I. performed research; R.S. and F.W.V. contributed new reagents; W.E.V., S.M.A.S., P.J.v.d.S. and T.J.V. analyzed data and W.E.V. and T.J.V. wrote the paper.
This work was supported by grants from the Netherlands Organization for Scientific Research (W.E.V.) and the Netherlands Genomics Initiative (S.M.A.S.). Funding to pay the Open Access Charge was provided by Netherlands Organisation for Scientific Research.
We thank Fernando Rivadeneira (Erasmus MC) for excellent help with statistical analysis.
Conflict of Interest statement. None declared.