This article has Open Peer Review reports available.
Identification of co-expression gene networks, regulatory genes and pathways for obesity based on adipose tissue RNA Sequencing in a porcine model
© Kogelman et al.; licensee BioMed Central Ltd. 2014
Received: 25 June 2014
Accepted: 24 September 2014
Published: 30 September 2014
Obesity is a complex metabolic condition in strong association with various diseases, like type 2 diabetes, resulting in major public health and economic implications. Obesity is the result of environmental and genetic factors and their interactions, including genome-wide genetic interactions. Identification of co-expressed and regulatory genes in RNA extracted from relevant tissues representing lean and obese individuals provides an entry point for the identification of genes and pathways of importance to the development of obesity. The pig, an omnivorous animal, is an excellent model for human obesity, offering the possibility to study in-depth organ-level transcriptomic regulations of obesity, unfeasible in humans. Our aim was to reveal adipose tissue co-expression networks, pathways and transcriptional regulations of obesity using RNA Sequencing based systems biology approaches in a porcine model.
We selected 36 animals for RNA Sequencing from a previously created F2 pig population representing three extreme groups based on their predicted genetic risks for obesity. We applied Weighted Gene Co-expression Network Analysis (WGCNA) to detect clusters of highly co-expressed genes (modules). Additionally, regulator genes were detected using Lemon-Tree algorithms.
WGCNA revealed five modules which were strongly correlated with at least one obesity-related phenotype (correlations ranging from -0.54 to 0.72, P < 0.001). Functional annotation identified pathways enlightening the association between obesity and other diseases, like osteoporosis (osteoclast differentiation, P = 1.4E-7), and immune-related complications (e.g. Natural killer cell mediated cytotoxity, P = 3.8E-5; B cell receptor signaling pathway, P = 7.2E-5). Lemon-Tree identified three potential regulator genes, using confident scores, for the WGCNA module which was associated with osteoclast differentiation: CCR1, MSR1 and SI1 (probability scores respectively 95.30, 62.28, and 34.58). Moreover, detection of differentially connected genes identified various genes previously identified to be associated with obesity in humans and rodents, e.g. CSF1R and MARC2.
To our knowledge, this is the first study to apply systems biology approaches using porcine adipose tissue RNA-Sequencing data in a genetically characterized porcine model for obesity. We revealed complex networks, pathways, candidate and regulatory genes related to obesity, confirming the complexity of obesity and its association with immune-related disorders and osteoporosis.
Obesity is a complex common health problem, and because of its exponential growth in prevalence in the last decades, the World Health Organization (WHO) has recognized it as a global epidemic since 1997. Obesity is an excess of adipose tissue, resulting from an imbalance between energy intake and energy expenditure . As adipose tissue plays a key role in the control of energy balance by secreting e.g. hormones, cytokines and growth factors, a disturbance in this tissue is strongly associated with other severe diseases such as type 2 diabetes, cardiovascular diseases and various types of cancer .
In this study, we use subcutaneous adipose tissue from a previously established porcine model  for studying the genetics of obesity. The value of the pig as a model for obesity has been investigated and proven by diverse studies, as the pig (omnivorous like humans) shares metabolic, digestive and cardiovascular features with humans [4, 5]. This was also obvious from a previous study of our F2 pig population which showed a strong genetic divergence for diverse obesity-related traits extensively phenotyped [3, 6].
Obesity is highly heritable, with heritability estimates between 40 and 70 percent , and although many studies have focused on finding loci involved in obesity (e.g. using genome-wide association studies (GWAS) ), the genetic risk variants identified to date explain a limited proportion of the heritability. In addition to genetic approaches, transcriptomic analyses have shown to be useful in studying complex diseases. Transcriptomic analyses give insight into the intermediate step between genes and its function, providing the opportunity to better understanding the biological mechanisms . Microarray technologies have been the main platform in recent years  and microarray expression data gave the opportunity of both assessing thousands of genes at the same time while providing a better understanding of the underlying biological processes of many complex diseases [11–13]. Expression data are extensively used to detect differentially expressed genes, but network approaches have also gained ground to reveal more of the complex transcriptional regulation by detecting sets of highly co-regulated genes. Clusters of co-regulated genes that share a common function, called ‘modules’, are thought to work together in a network and correspond to, for example, biological pathways. Several network approaches are available, such as the Weighted Gene Co-expression Network Analysis (WGCNA) method (based on correlation patterns between expression profiles)  which has proven its superiority over Partial Correlation and Information Theory (PCIT) methods , and the Lemon-Tree project (based on probabilistic graphical models) . We have previously applied WGCNA approaches to complex traits in animal populations [17, 18], and the method has also proved to be reliable in various human diseases, e.g. different types of cancer [19, 20].
A relatively novel platform for gene expression data is RNA Sequencing (RNA-Seq), which allows us to study the complete transcriptome in more detail and with more precise measurements in comparison with microarray platforms, also facilitating the discovery of novel genes . This has already been proven by Iancu et al.  who compared the WGCNA approach applied to microarray expression data and RNA-Seq data. Their results showed that the greater sensitivity and dynamic range results in a better estimation of network properties, such as network density and centralization. The huge advantage of detecting novel genes for complex traits will also be a major opportunity in network approaches, where the biological information included in the network helps in revealing the function of the novel genes. To date, a limited number of studies have analyzed RNA-Seq data using the WGCNA approach to investigate complex diseases and traits, e.g. psoriasis , heart failure in a mouse model  and alcohol use disorders in a mouse model . Moreover, the WGCNA approach has also been used on RNA-Seq data of human subcutaneous adipose tissue, revealing a cluster of genes associated with serum triglyceride regulation . Also the Lemon-Tree algorithms, previously published as LeMoNe, have showed their value in different biological traits [27–29]. The combined use of high-throughput omics data from clearly characterized groups of individuals for complex diseases and mathematical models to build gene co-expression and regulatory networks is at the core of systems biology methods [13, 17, 30].
In the present study we applied the WGCNA method and Lemon-Tree algorithms, using RNA-Seq data, with the aim to elucidate transcriptomic regulations of obesity by detecting pathways, novel and regulator genes involved in its pathogenesis. To our knowledge, this is the first study using systems biology and network approaches to study the overall complex transcriptional regulation of obesity using RNA-Seq data in a genetically characterized porcine model.
Results and discussion
The genetic background of selected animals
Descriptive statistics (mean and standard deviation) and test of difference of means for a selection of obesity-related traits for the three subgroups
Weighted Gene Co-expression Network Analysis (WGCNA)
Functional enrichment of modules
Overview of the most significantly overrepresented KEGG pathways and GO terms associated with the modules detected using WGCNA
KEGG pathway or GO-term
Frequency in module
Number of genes in term
KEGG: Osteoclast differentiation
KEGG: Natural Killer cel mediated cytotoxity
KEGG: B cell receptor signaling pathway
GO: Immune system response
GO: Cell activation
GO: Regulation of immune system response
GO: Extracellular region
GO: Extracellular matrix
GO: Proteinaceous extracellular matrix
GO: Ribosome biogenesis
GO: Nuclear par
These expression profiles are as expected, as various studies have shown that an increase of adipose tissue is correlated with an increase of several types of immune cells, particularly with macrophages [32, 33] which are key players in the initiation of a chronic inflammatory state in obesity . As expected for the 69 genes, the significant GO terms and KEGG pathways were still nearly all immune related, except that the most significant KEGG pathway was osteoclast differentiation (Padj = 1.4E-7) (Figure 3C & 3D). Osteoclasts are derived from macrophages, one of the most up-regulated immune cells in adipose tissue of obese individuals, and are therefore also closely linked to many immune diseases . Bone marrow houses two kinds of stem cells: the mesenchymal stromal cells which are precursors for osteoblasts and adipocytes and the hematologic stem cells originating from osteoclasts. Moreover, there is an important communication between adipose tissue and skeleton where factors secreted by adipocytes affect bone remodeling, i.e. leptin, adiponectin, pro-inflammatory cytokines as Interleukin 6 (IL-6) [36, 37]. IL-6 is known to be an important regulator of the immune and hematopoietic systems and it has been associated with osteoporosis disease and rheumatoid arthritis [38, 39]. Osteoporosis is a polygenic trait , whereby increased bone fragility results from increased adipocytes and osteoclastogenesis and insufficient osteoblastogenesis . When looking at the functions of the different genes present in the Blue module, we find many genes which have a clear function in the immune system and also have been associated with osteoclast differentiation, e.g. BCL2A1, DOK2, PTPN6 and several genes encoding cell surface molecules (e.g. CD45 and CD68). Another gene in this module is SPI1. SPI1 is encoding the transcription factor PU.1 protein which activates gene expression during myeloid and B-lymphoid cell development. A study of Wang et al.  has shown that PU.1 is expressed in white adipose tissue and plays a role in adipogenesis. Moreover, variations in SPI1 play a role in osteoclastogenesis as for example, PU.1 deficient mice develop osteoporosis , and it increases the risk of fracture by its effect on ALOX15.
Other highly significant associated pathways in this module are all immune-related, of which the most significant is Natural killer cell mediated cytotoxity (P-value = 3.8E-5). In fact, obesity causes morphological changes in adipose tissue, resulting in a state of chronic low-grade inflammation . Furthermore, natural killer (NK) cells are critical in the innate immune response, less examined in association with obesity, but it has been shown that diet-induced obese mice show a reduced NK cytotoxity after infection . Another study showed an increased level of NK cells in healthy obese compared with unhealthy obese individuals, suggesting its importance in metabolic processes . Several studies have shown and investigated the link between the immune system and metabolism [48, 49], also in combination with obesity [50, 51]. This also explains the significant association of the other KEGG pathways and GO terms in this module.
The Black module (MTROI = 0.35) shows a strong reverse correlation (-0.42) with fasting glucose levels (FGL). The KEGG pathways are not significant after BH correction, but before BH correction the most significant pathway is ECM-receptor interaction (P = 0.001). Several GO terms related to this extracellular matrix (ECM) are found to be significantly overrepresented, also after BH correction, e.g., extracellular region (Padj = 5.5E-6), extracellular matrix (Padj = 3.6E-5) and proteinaceous extracellular matrix (Padj = 3.6E-5). As we are interested in the genes which are involved in the pathways representing the high positive correlation with fatness, but with a high negative correlation with glucose levels, we examined the association of the genes between the two traits. We selected leaf fat at slaughter (SLfat) and FGL as traits of interest because of their high correlations. The correlations of the expression profiles with these traits show that there is a wide variation in their correlations with both traits, and that there is a weak negative correlation (-0.23) between the Gene-Trait correlations of SLfat and FGL. Next, we only selected genes having a correlation >0.4 with both SLfat and FGL, resulting in a selection of 36 genes, of which 24 were assigned a gene name, for further functional annotation. Of these genes we will only comment on the most relevant in relation to obesity. ADAMTS-12 is a metalloprotease necessary for normal immunological response . The PFPK gene (phosphofructokinase, platelet) is a key regulatory enzyme in glycolysis. In the first GWAS presented on obesity, this gene was found to be associated, but did not get validated in the replication stage . However, differential gene expression in the visceral adipose tissue shows differential expression of the PFPK gene in obese vs. lean individuals . Two of the genes in the modules are proprotein convertase subtilisin/kexins (PCSK): the PCSK5 and PCSK6 gene. Several studies have shown their relevance in metabolic-associated processes, for example high-density lipoprotein (HDL) metabolism [55, 56]. The BDKRB2 (BK type 2 receptor) gene has been suggested to be acting as a genetic modulator of glucose homeostasis, potentially increasing susceptibility to diseases like diabetes . The Probable G-protein coupled receptor 133 (encoded by the GPR133 gene) has been associated using GWAS to body weight control in mice  and height control in human . The MFAP5 gene is also present in the Black module, encoding the microfibrillar-associated protein 5. This is an extracellular matrix glycoprotein, shown to be expressed in adipose tissue, and its expression is positively correlated with BMI, and change in body fat mass . Moreover, Vaittinen et al.  also found correlations between MAFP5 expression and insulin resistance markers.
The third investigated module is the Lightyellow Module (MTROI = -0.22), showing a high MTR with omental fat at slaughter (-0.54). This module does not show any significant KEGG pathways after BH correction, but several significant overrepresented GO terms related to transcription were detected, e.g. nucleolus (Padj = 0.013), ribosome biogenesis (Padj = 0.034), nuclear part (Padj = 0.034), DNA-dependent transcription, elongation (Padj = 0.038). All these processes are involved in gene expression, which is altered under obesity. The most highly interconnected genes in this module are AMD1, NOL9, EIF4A1, POLR1C and ABCE1. AMD1, the hub-gene in this module, plays a key role in polyamine biosynthesis. Polyamines have shown to affect growth and development of adipose tissue, and increased levels of polyamines have been associated with childhood obesity . Moreover, in an Indian childhood cohort AMD1 has been associated with obesity and plasma leptin levels, speculating that AMD1 influences the susceptibility to obesity by modulating the polyamine metabolism or DNA methylation . NOL9 plays a role in ribosomal RNA (rRNA) processing  which takes place in the nucleolus. It has previously been shown that the EIF4A1 gene plays a role in the protein synthesis pathway which was altered in insulin-resistant obese individuals . Furthermore, EIF4A1 was negatively correlated with the insulin receptor gene (InsR). POLR1C is important in RNA polymerization, which is necessary for transcription. ABCE1 is involved in viral assembly by inhibiting the action of RNase L, a regulator of innate immunity. It has been shown that RNase L plays an important role in adipogenesis regulation  and is able to restore insulin response in muscle cells of obese individuals, suggesting a role for RNase L in “healthy obese” subjects . This finding is in agreement with a study of Mahdi et al.  that showed a correlation between the expression levels of ABCE1 and insulin secretion, following the same WGCNA approach as this study.
Neither the Greenyellow module (MTROI = 0.45) nor the Brown module (MTROI = 0.47) showed any significant KEGG pathways or GO terms after BH correction.
Differential network analysis
Previous studies have shown the potential of studying the co-expression patterns of sub-networks, thereby comparing two different states of the phenotype of interest . Therefore, we created and compared the lean sub-network and obese sub-network using the expression data of the 12 low OI animals and the 12 high OI animals. For both sub-networks the gene dendrogram was created and modules were subsequently defined to colors (results not shown). Both sub-networks revealed clearly different modules; however, the obese modules were not preserved in the lean sub-network, and vice versa. Biologically, this shows different active pathways in the lean versus obese animals, as expected.
We found eight genes (IFI35, SNX19, MARC2, MORN2, PLA2G16, CCL5, CACNA2D1, and REV1) that were hub genes in the obese sub-network, but lowly interconnected in the lean sub-network. All except IFI35 and REV1 are up-regulated in the obese animals. Many of the obese hub-genes are, as expected, active in obesity-related pathways, where we can see a subdivision in fatness related pathways and immune-related pathways. Both IFI35 and CCL5 are immune-related genes. IFI35 (interferon-induced protein 35) is a cytokine important for communication between cells in the innate immune response. A recent study showed that IFI35 negatively regulates the RIG-I antiviral signaling , which in turn activates the innate immune response. These patterns are in concordance with the down-regulation of IFI35 in our obese pigs, as it has been accepted that the native immune response is triggered in obese individuals . CCL5 (chemokine (C-C motif) ligand 5) is secreted by bone marrow stromal cells and encodes a protein which is also known as RANTES (regulated on activation, normal T cell expressed and secreted). It plays a key role in recruiting leukocytes during immune response and plays a role in inducing the activation of Natural Killer cells . It may also play a role in the inflammation of obese human white adipose tissue , as it shows up-regulation in obese subjects . SNX19 (sortin nexin 19) encodes a membrane-associated protein complex which is associated with coronary heart disease and myocardial infarction. It can also interact with IA-2 which is a major auto-antigen in type 1 diabetes and a regulator of insulin secretion . Several of the obese hub-genes are also related to fat and body weight. The first one is MARC2 (Mitochondrial Amidoxime Reducing Component 2), which has been associated with decreased total body fat and increased circulating glucose levels in mice . The PLA2G16 (adipocyte phospholipase A2 (AdPLA) gene, also present in the Greenyellow module, encodes an enzyme that is mainly found in adipose tissue and catalyzes the release of fatty acids from phospholipids in adipose tissue . Genetic variants of this gene have been associated with neutral lipid storage disease  and growth traits in cattle . AdPLA-deficient ob/ob mice show an increased lipolysis, a reduced adipose tissue mass and triglyceride content, but show a normal adipogenesis, showing its role in the control of adipocyte lipolysis and thereby its important role in the development of obesity . The CACNA2D1 gene (Voltage-dependent calcium channel subunit alpha-2/delta-1) encodes a protein in the voltage-dependent calcium channel complex. In a meta-analysis on body fat in mice it was shown to be associated with body fat percentage  and has been associated with carcass and meat quality traits in cattle [80, 81]. The exact function of the MORN2 gene is unknown; this gene was also present in the Greenyellow module. The last obese hub-gene is REV1, which encodes a DNA repair protein. The fact that this gene is down-regulated in obese animals may be explained by the fact that the DNA repair process is perturbed under obesity. It has been seen, for example, that mutations in genes encoding for DNA repair proteins are associated with developmental inhibition, immunodeficiency and increased risk of cancer in humans .
In addition, we found 11 genes (CSF1R, EVI2B, SAMHD1, PCD1A, CD68, FAM105A, P2Y12R, NCEH1, SPI1, FRMPD4 and PigE-108A11.6) which were hub-genes in the lean network but weakly interconnected in the obese sub-network. All these genes are up-regulated in the obese animals. All of them, except FRMPD4, were also found in one of the WGCNA modules. Except for the SAMHD1gene, they were all found in the Blue module; SAMDHD1 was present in the Black module. Several of the 11 genes are immune-related, e.g. CSF1R, SAMHD1, PCD1A, CD68, SPI1 and PigE-108A11.6. CSF1R (Colony stimulating factor 1 receptor) controls the production, differentiation and function of macrophages  and its expression has been correlated with body fat percentage in humans . SAMHD1 (SAM domain and HD domain-containing protein 1) is an enzyme exhibiting phosphohydrolase activity, converting nucleotide triphosphates to a nucleoside and triphosphate. PCD1A is the porcine cd1 antigen involved in the presentation of lipid antigens to T cells, but its precise function is unknown. CD68 encodes a glycoprotein that binds to low-density lipoprotein, expressed on monocytes/macrophages. Its expression was significantly up-regulated in acquired obesity in monozygotic twin pairs and correlated to liver fat and insulin resistance . SPI1 has been found in the Blue module in WGCNA, and its relation with adipogenesis and ostoclastogenesis has been discussed. PigE-108A11.6 is the orthologous of LILRB5 (leukocyte immunoglobulin-like receptor subfamily B member 5), which is expressed in immune cells where they bind to MHC class 1 molecules on antigen-presenting cells and inhibit stimulation of an immune response. It has been found to be up-regulated in omental adipose tissue in obese individuals . Also various other genes were detected among the hub-genes in the lean animals. The P2Y12R gene encodes a protein which is an important regulator in blood clotting and consequently associated with heart attacks . Moreover, insulin inhibits blood platelet aggregation by suppressing the P2Y12 pathway, and therefore type 2 diabetes leads to up-regulation of the P2Y12 pathway, resulting in increased platelet reactivity . The NCEH1 (Neutral cholesterol ester hydrolase 1) gene encodes an enzyme which is located in the endoplasmic reticulum which plays a role in the regulation of the levels of platelet activating factor and lysophospholipids. NCEH1 has a key role in the reverse cholesterol transport in macrophages, thereby playing a critical role in human atherosclerosis .
Detection of regulator genes
Several of the detected modules showed an overrepresentation of obesity-related KEGG pathways or GO terms, but the WGCNA approach did not reveal potential gene regulators of detected modules, as WGCNA is an undirected network. The challenge of finding potential regulators of significant gene modules related to obesity was addressed by using Lemon-Tree algorithms (available at https://code.google.com/p/lemon-tree/). Lemon-Tree created a set of potential regulators consisting of transcription factors and signal transducers using the GO categories ‘transcription factor activity’ and ‘signal transducer activity. Those potential regulators were then assigned to nodes with corresponding probabilistic scores for being a “regulator”. Genes will have a high probabilistic score in case there is a differential expression pattern on each side of the particular tree node. After generation of multiple statistically equal partitions of conditions for each cluster of genes, the algorithm uses an ensemble approach to sum the strength with which a regulator participates in each regulatory tree. The overall statistical confidence was calculated and used for prioritizing regulators, represented by the global probabilistic score (Prob. score). Mathematical details of the Lemon-Tree algorithms can be found in Joshi et al. . In this way, we built a regulatory module network using the RNA-Seq expression data and assigned high-scoring regulator genes to the detected clusters of co-expressed genes. We detected in total 43 tight clusters of no less than 10 genes, which resulted in a total of 1417 genes in the clusters. A pre-defined group of genes was tested to see whether they may regulate the expression levels of the genes in a particular module, represented by the global probabilistic score (Prob. Score). We tested the significance of the assigned regulators using a t-test, showing that the assigned regulators are significantly different from randomly assigned regulators (P-value = 0.000).
We also examined the genes present in the Black module in the Lemon-Tree results. Many of the genes in the Black module clustered together in cluster10. The Lemon-Tree algorithms point to ACKR2 as regulator gene (Figure 5B), again encoding a chemokine receptor, but in this case its exact function is not known.
Examining the obesity-related measurements of the 36 RNA-sequenced animals (divided among three genetically extreme groups) showed that there are clear differences between the three groups. Further investigation of their expression profiles using systems biology and network approaches revealed pathways that potentially play a key role in the development of obesity and obesity-related diseases (e.g. inflammation and osteoporosis). Specifically, our network approach WGCNA revealed several clusters of highly interconnected genes of which some were associated to immune-related pathways, which have previously been implicated in the pathogenesis of obesity. A subset of these pathways was also related to osteoclast differentiation (Padj = 1.4E-7), which has been associated with the immune system and osteoporosis. The co-expression of immune-related genes and their association with osteoclast differentiation shows the close association between obesity, the immune system and bone remodeling diseases like osteoporosis – these results have important implications for obese human patients. Furthermore, Lemon-Tree algorithms detected three regulator genes: CCR1, MSR1 and SPI1. Those were previously associated with the immune system and/or osteoclast differentiation, but here we show that those genes may have a key role in the link between obesity and osteoporosis. Another module contained genes previously associated with lipid metabolism and glucose metabolism, which reveals the known close association between obesity and other metabolic diseases like type 2 diabetes. Differential network analysis revealed several genes which showed a different level of activity in the pathways they represented. Many of those genes were previously associated with obesity or the immune system, but not detected using single-gene association studies. This stresses the potential of using systems biology or network biology approaches on complex diseases in order to detect genes which may be important in disease development, and consequently could be potential drug markers.
In conclusion, this study shows the advantage of using systems or network biology approaches in complex diseases to unravel their genetic architecture and transcriptional regulation. We revealed several new pathways and, more importantly, key regulator genes present in those pathways which give a better understanding of the complex transcriptional regulation of obesity. Furthermore, we have confirmed, using our genetically characterized omnivorous porcine model, the known association of obesity with other bone-related and immunity-related diseases at the transcriptomic level, and revealed several genes actively present in those associations.
The F2 pig population and subcutaneous adipose tissue samples
The samples used in this study were derived from a pig resource population established with the purpose of identifying the molecular background for obesity and obesity related traits. For a detailed description of the population and the phenotypes see (Kogelman et al.) . In this study we have only used the Duroc *Göttingen Minipig intercross that comprise a total of 279 animals generated by intercrossing Göttingen minipig boars (Ellegaard A/S) and Duroc sows (DanBred breeding herd). The Göttingen minipig breed is characterized for their small size, and its genetic predisposition for obesity, and for sharing several metabolic impairments seen in humans . The Duroc pig is a production breed, intensively selected for leanness and growth during the last 60 years. This population has been extensively phenotyped for OOR traits, e.g. weight, conformation, dual energy x-ray absorptiometry (DXA) scanning and slaughter measurements. The F2 pigs were produced at the research farm, University of Copenhagen Tåstrup, Denmark. Animal care and maintenance have been conducted according to the Danish “Animal Maintenance Act” (Act 432 dated 09/06/2004) and biological samples were collected according to the Danish “Veterinary Procedures Act” (Act 433 dated 09/06/2004). After slaughtering, a biobank of tissues was created by sampling several tissues for each animal, including subcutaneous adipose tissue. Following, we created one aggregate predicted genetic value based on the selection index theory  to represent the degree of obesity in each animal: the Obesity Index (OI) . The OI followed a normal distribution, and based on this distribution animals were categorized in two extreme groups: obese and lean. A third group was categorized as being around the mean of the distribution. We selected 36 animals from those three groups for RNA-sequencing of subcutaneous adipose tissue: 12 low OI (lean), 12 intermediate OI, and 12 high OI (obese). For animal selection, the family structure and the genders of the animals was taken into account, by maximizing the number of directly related animals (full sibs) to three per group and taking care of equal gender distribution within the groups.
Phenotypic characterization of the selected animals
The phenotypic measurements of the animals were investigated to compare the phenotypic distributions between the three selected groups. Data that were not normally distributed was transformed to approach normality. Differences in means were tested using an analysis of variance model using the R statistical package . All phenotypic measurements were corrected for gender, and in case of a significant effect corrected for batch and age effect, as described before in Kogelman et al. . Results were noted as highly significant with a P-value below 0.001, significant with a P-value below 0.01 and tend to be significant with a P-value below 0.05.
RNA isolation and RNA-Sequencing
Total RNA was isolated from porcine subcutaneous adipose tissue using the RNeasy lipid Mini kit (Qiagen, Germany) following manufacturer’s recommendations. Briefly, 100 mg of tissue were homogenized in Qiazol buffer using the gentleMACS™ Octo Dissociator system (Milteny Biotec, GmbH, Germany) with M tubes (Milteny Biotec, GmbH, Germany) and the RNA_02 program recommended by the manufacturer. The RNA was DNase treated in order to degrade the remaining genomic DNA. At the end of the protocol the RNA was eluated in 30 μl of RNAse free water. Quantity and quality were assessed by Nanodrop ND-1000 spectrophotometer. Furthermore integrity of the isolated RNA was inspected by electrophoresis in a 1.4% agarose gel and by measuring the RQI value on an Experion™ system (BioRad) using Eukaryote Total RNA StdSens kit (BioRad). The RQI average for all the 36 extracted samples was 8.56.
For the RNA-seq, libraries were constructed using 400 ng of total RNA and TruSeq RNA Sample Prep (Illumina) with Poly-A pull down rRNA depletion following manufacturer’s recommendations. Samples were sequenced on the HiSeq2500 platform, dividing the 36 samples over 4 lanes, using a read length of 100 bp paired-end reads . Before alignment, reads with a low quality and adapters were detected using FastQC and removed. Remaining reads were mapped to the SScrofa10.2.72 genome using default parameters in STAR aligner. This resulted in an average of 30,557,234 uniquely mapped reads per sample, of which on average 81.60 percent was mapped in the intragenic region (within introns or exons). On average 20,390 transcripts were detected among the mapped reads. Read counts were estimated at gene-level using HTSeq .
RNA-Sequencing data normalization and gene selection
Previous studies have been shown that genes with extreme low expression levels are less reliable . Therefore, genes with expression levels equal or lower than 5 were removed from the dataset, resulting in 12,253 genes per sample. The between-sample bias was removed by estimating the library size factor using the estimateSizeFactor() function in DESeq . Normalization was then performed using the voom() variance-stabilization function in the R-package Limma , whereby samples were corrected for gender.
Weighted Gene Co-expression Network Analysis
Because of computational limitations the dataset had to be further reduced for network construction. As it is believed that non-changing genes provide limited information in a co-expression network setting, genes were selected based on their variation (SD > 0.25), resulting in 8,745 genes. Furthermore, based on the connectivity (sum of connection strengths with all other genes) genes were selected for the lean, intermediate and obese sub dataset, as genes with a high connectivity (hub genes) are thought to be more biologically important . The 1500 highest connected genes in each of the lean, intermediate and obese sub-dataset resulted in a joint dataset of 3532 unique genes. Consequently, we built an unsigned co-expression network using the WGCNA R-package , using 3532 genes. The adjacency matrix was created by calculating the Pearson’s correlations between all genes, and raised to a power β of 7. The power β was chosen based on the scale-free topology criterion , resulting in a scale-free topology index (R2) of 0.92. Next, the Topological Overlap Measure (TOM), representing the overlap in shared neighbors, was calculated using the adjacency matrix. The dissimilarity TOM was used as input for the dendrogram, and modules (clusters of highly interconnected genes) were detected as branches of the dendrogram using the DynamicTreeCut algorithm . All modules were assigned to a color. The module eigengene was used to represent each module, which was calculated by the first principal component, thereby capturing the maximal amount of variation of the module. Using the module eigengene, the Module-Trait relationships were estimated by calculating the Pearson’s correlations between the module eigengene and the traits of interest. Those Module-Trait relationships were used to select potential biologically interesting modules for downstream analysis. Modules were selected when they had a correlation >0.5 with at least one of the selected traits. Genes in the module were selected when their intra-modular connectivity with that particular module was >0.6, the intra-modular connectivity with all other modules <0.6. The intra-modular connectivity is calculated as the correlation between the gene’s expression profile and the expression profile of the module eigengene. Another gene characteristic is the Gene Trait correlation: the correlation between the gene’s expression profile and the phenotype of interest.
After selection of modules and genes retained in those modules, we investigated those modules using several different module characteristics. An important intra-modular analysis was detection of the correlation between the Module Membership and the Gene-Trait relationships. The Module Membership is based on measuring the intra-modular connectivity: the correlation of a particular gene in the module with the module eigengene. The Gene-Trait relationship is the correlation between the expression profiles of a particular gene with the phenotype of interest. We selected a subset of the modules based on functional annotation. Using BioMart  the associated gene names were detected. Gene length bias has been shown to be an important bias in RNA-Seq, as longer genes will be sequenced deeper than shorter genes . Therefore, we used GOSeq  to detect overrepresented gene ontology (GO) terms and KEGG pathways, which is able to correct for gene length bias. The Probability Weighting Fuction (PWF) is obtained depending on the gene length, which is used in the Wallenius approximation to calculate overexpressed GO terms and KEGG pathways among the selected modules. To correct for the multiple testing problem, the p-values are adjusted using the Benjamini-Hochberg (BH) correction. GO terms and KEGG pathways were thought to be significant when the adjusted p-value was below 0.05.
The differential connectivity is a measure of the differences in gene interactions between the lean and obese animals, potentially identifying genes which underlie a difference in transcriptional regulation . Therefore, two sub-networks were created: one with only the 12 lean animals (lean sub-network) and one with only the 12 obese animals (obese sub-network). For network construction, we used the 12,253 genes which passed the quality control measures as described above. Again, we selected the most varying genes by using the variance threshold of S.D. > 0.25, resulting in 8,745 genes for network construction. The adjacency matrix was constructed in both sub-networks, and was created by calculating the Pearson’s correlations and raising them to the power β = 9, based on the scale-free topology assumption. The connectivity of each gene was calculated in both the lean and obese sub-network, as the sum of connection strengths of a particular gene with all other genes.
The differential connectivity (k_diff) was calculated by subtracting the connectivity of the genes in the lean sub-network from the connectivity in the obese sub-network. This resulted in a normal distribution of values between -1 and 1. When the k_diff was positive, it meant that the genes were highly connected in the obese sub-network than in the lean sub-network, while a negative k_diff meant that the genes were highly connected in the lean sub-network than in the obese sub-network. Furthermore, we performed a standard t-test, comparing the gene expression of the obese vs. the lean animals, to assign the difference in expression between the two sub-networks. Genes were labeled as differential connected when the absolute k_diff was above 0.6, resulting in the detection of genes which were hubgenes in one sub-network, but not in the other sub-network. Differentially connected genes were selected for functional annotation when their absolute t-test statistic was above 1.96 indicating that the genes are differentially expressed between the obese and lean animals. Associated genes were detected again using BioMart .
Detection of regulator genes
Availability of supporting data
The data discussed in this publication have been deposited in NCBI's Gene Expression Omnibus  and are accessible through GEO Series accession number GSE61271 (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE61271) .
The project is supported by a grant from the Danish Council for Strategic Research (BioChild Project: http://www.biochild.ku.dk) and from the EU-FP7 Marie Curie Actions – Career Integration Grant (CIG-293511), both granted to Haja N. Kadarmideen, from the University of Copenhagen Ph.D. stipend awarded to Lisette J.A. Kogelman and from the Danish Ministry of Science and Technology (the “UNIK Project for Food Fitness and Pharma for Health”) to Merete Fredholm.
- Bray GA: An Atlas of Obesity and Weight control. 2002, London, UK: CRC PressGoogle Scholar
- Kershaw EE, Flier JS: Adipose Tissue as an Endocrine Organ. J Clin Endocrinol Metab. 2004, 89 (6): 2548-2556. 10.1210/jc.2004-0395.PubMedGoogle Scholar
- Kogelman LJA, Kadarmideen HN, Mark T, Karlskov-Mortensen P, Bruun CS, Cirera S, Jacobsen MJ, Jørgensen CB, Fredholm M: An F2 pig resource population as a model for genetic studies of obesity and obesity-related diseases in humans: Design and genetic parameters. Front Genet. 2013, 4: 29.PubMedPubMed CentralGoogle Scholar
- Spurlock ME, Gabler NK: The Development of Porcine Models of Obesity and the Metabolic Syndrome. J Nutr. 2008, 138 (2): 397-402.PubMedGoogle Scholar
- Rocha D, Plastow G: Using commercial pigs in the search for genes behind human obesity. Nat Rev Genet. 2005, 6 (3).Google Scholar
- Kogelman LJA, Pant SD, Fredholm M, Kadarmideen HN: Systems genetics of obesity in an F2 pig model by genome-wide association, genetic network and pathway analyses. Front Genet. 2014, 5: 214.PubMedPubMed CentralGoogle Scholar
- Speliotes EK, Willer CJ, Berndt SI, Monda KL, Thorleifsson G, Jackson AU, Allen HL, Lindgren CM, Luan J, Mägi R, Randall JC, Vedantam S, Winkler TW, Qi L, Workalemahu T, Heid IM, Steinthorsdottir V, Stringham HM, Weedon MN, Wheeler E, Wood AR, Ferreira T, Weyant RJ, Segrè AV, Estrada K, Liang L, Nemesh J, Park JH, Gustafsson S, et al: Association analyses of 249,796 individuals reveal 18 new loci associated with body mass index. Nat Genet. 2010, 42 (11): 937-948. 10.1038/ng.686.PubMedPubMed CentralGoogle Scholar
- Day FR, Loos RJF: Developments in Obesity Genetics in the Era of Genome-Wide Association Studies. J Nutrigenet Nutrigenomics. 2011, 4 (4): 222-238. 10.1159/000332158.PubMedGoogle Scholar
- Adams JU: Transcriptome: Connecting the Genome to Gene Function. Nat Educ. 2008, 1 (1): 195.Google Scholar
- Heller MJ: DNA microarray technology: devices, systems, and applications. Annu Rev Biomed Eng. 2002, 4: 129-153. 10.1146/annurev.bioeng.4.020702.153438.PubMedGoogle Scholar
- Slonim DK: From patterns to pathways: gene expression data analysis comes of age. Nat Genet. 2002, 32: 502-508. 10.1038/ng1033.PubMedGoogle Scholar
- Nica AC, Dermitzakis ET: Using gene expression to investigate the genetic basis of complex disorders. Hum Mol Genet. 2008, 17 (R2): R129-R134. 10.1093/hmg/ddn285.PubMedPubMed CentralGoogle Scholar
- Kadarmideen HN: Genomics to systems biology in animal and veterinary sciences: progress, lessons and opportunities. Livest Sci. 2014, 116: 232-248.Google Scholar
- Zhao W, Langfelder P, Fuller T, Dong J, Li A, Hovarth S: Weighted Gene Coexpression Network Analysis: State of the Art. J Biopharm Stat. 2010, 20 (2): 281-300. 10.1080/10543400903572753.PubMedGoogle Scholar
- Kadarmideen HN, Watson-Haigh N: Building gene co-expression networks using transcriptomics data for systems biology investigations: Comparison of methods using microarray data. Bioinformation. 2012, 8 (18): 855-861. 10.6026/97320630008855.PubMedPubMed CentralGoogle Scholar
- Joshi A, De Smet R, Marchal K, Van de Peer Y, Michoel T: Module networks revisited: computational assessment and prioritization of model predictions. Bioinformatics. 2009, 25 (4): 490-496. 10.1093/bioinformatics/btn658.PubMedGoogle Scholar
- Kadarmideen HN, Watson-Haigh NS, Andronicos NM: Systems biology of ovine intestinal parasite resistance: disease gene modules and biomarkers. Mol Biosyst. 2011, 7 (1): 235-246. 10.1039/c0mb00190b.PubMedGoogle Scholar
- Kogelman LJA, Byrne K, Vuocolo T, Watson-Haigh N, Kadarmideen HN, Kijas J, Oddy H, Gardner G, Gondro C, Tellam R: Genetic architecture of gene expression in ovine skeletal muscle. BMC Genomics. 2011, 12 (1): 607-10.1186/1471-2164-12-607.PubMedPubMed CentralGoogle Scholar
- Clarke C, Madden SF, Doolan P, Aherne ST, Joyce H, O’Driscoll L, Gallagher WM, Hennessy BT, Moriarty M, Crown J, Kennedy S, Clynes M: Correlating transcriptional networks to breast cancer survival: a large-scale coexpression analysis. Carcinogenesis. 2013, 34 (10): 2300-2308. 10.1093/carcin/bgt208.PubMedGoogle Scholar
- Wang L, Tang H, Thayanithy V, Subramanian S, Oberg AL, Cunningham JM, Cerhan JR, Steer CJ, Thibodeau SN: Gene Networks and microRNAs Implicated in Aggressive Prostate Cancer. Cancer Res. 2009, 69 (24): 9490-9497. 10.1158/0008-5472.CAN-09-2183.PubMedPubMed CentralGoogle Scholar
- Wang Z, Gerstein M, Snyder M: RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009, 10 (1): 57-63. 10.1038/nrg2484.PubMedPubMed CentralGoogle Scholar
- Iancu OD, Kawane S, Bottomly D, Searles R, Hitzemann R, McWeeney S: Utilizing RNA-Seq Data for De-Novo Coexpression Network Inference. Bioinformatics. 2012, 28 (12): 1592-1597. 10.1093/bioinformatics/bts245.PubMedPubMed CentralGoogle Scholar
- Li B, Tsoi LC, Swindell WR, Gudjonsson JE, Tejasvi T, Johnston A, Ding J, Stuart PE, Xing X, Kochkodan JJ, Voorhees JJ, Kang HM, Nair RP, Abecasis GR, Elder JT: Transcriptome Analysis of Psoriasis in a Large Case–control Sample: RNA-Seq Provides Insights into Disease Mechanisms. J Invest Dermatol. 2014, 134 (7): 1828-1838. 10.1038/jid.2014.28.PubMedPubMed CentralGoogle Scholar
- Lee J-H, Gao C, Peng G, Greer C, Ren S, Wang Y, Xiao X: Analysis of Transcriptome Complexity Through RNA Sequencing in Normal and Failing Murine Hearts. Circ Res. 2011, 109 (12): 1332-1341. 10.1161/CIRCRESAHA.111.249433.PubMedPubMed CentralGoogle Scholar
- Darlington TM, Ehringer MA, Larson C, Phang TL, Radcliffe RA: Transcriptome analysis of Inbred Long Sleep and Inbred Short Sleep mice. Genes Brain Behav. 2013, 12 (2): 263-274. 10.1111/gbb.12018.PubMedPubMed CentralGoogle Scholar
- Haas B, Horvath S, Pietilainen K, Cantor R, Nikkola E, Weissglas-Volkov D, Rissanen A, Civelek M, Cruz-Bautista I, Riba L, Kuusisto J, Kaprio J, Tusie-Luna T, Laakso M, Aguilar-Salinas CA, Pajukanta P: Adipose Co-expression networks across Finns and Mexicans identify novel triglyceride-associated genes. BMC Med Genomics. 2012, 5 (1): 61-10.1186/1755-8794-5-61.PubMedPubMed CentralGoogle Scholar
- Michoel T, De Smet R, Joshi A, Van de Peer Y, Marchal K: Comparative analysis of module-based versus direct methods for reverse-engineering transcriptional regulatory networks. BMC Syst Biol. 2009, 3: 49-10.1186/1752-0509-3-49.PubMedPubMed CentralGoogle Scholar
- Vermeirssen V, Joshi A, Michoel T, Bonnet E, Casneuf T, Van de Peer Y: Transcription regulatory networks in Caenorhabditis elegans inferred through reverse-engineering of gene expression profiles constitute biological hypotheses for metazoan development. Mol Biosyst. 2009, 5 (12): 1817-1830. 10.1039/b908108a.PubMedGoogle Scholar
- Bonnet E, Tatari M, Joshi A, Michoel T, Marchal K, Berx G, Van de Peer Y: Module Network Inference from a Cancer Gene Expression Data Set Identifies MicroRNA Regulated Modules. PLoS One. 2010, 5 (4): e10162-10.1371/journal.pone.0010162.PubMedPubMed CentralGoogle Scholar
- Chuang HY, Hofree M, Ideker T: A decade of systems biology. Annu Rev Cell Dev Biol. 2010, 26: 721-744. 10.1146/annurev-cellbio-100109-104122.PubMedPubMed CentralGoogle Scholar
- Oshlack A, Wakefield M: Transcript length bias in RNA-seq data confounds systems biology. Biol Direct. 2009, 4 (1): 14-10.1186/1745-6150-4-14.PubMedPubMed CentralGoogle Scholar
- Weisberg SP, McCann D, Desai M, Rosenbaum M, Leibel RL, Ferrante AW: Obesity is associated with macrophage accumulation in adipose tissue. J Clin Invest. 2003, 112 (12): 1796-1808. 10.1172/JCI200319246.PubMedPubMed CentralGoogle Scholar
- Ferrante AW: The immune cells in adipose tissue. Diabetes Obes Metab. 2013, 15 (s3): 34-38. 10.1111/dom.12154.PubMedPubMed CentralGoogle Scholar
- Patel PS, Buras ED, Balasubramanyam A: The Role of the Immune System in Obesity and Insulin Resistance. J Obes. 2013, 2013: 9.Google Scholar
- Ginaldi L, Di Benedetto M, De Martinis M: Osteoporosis, inflammation and ageing. Immun Ageing. 2005, 2 (1): 14-10.1186/1742-4933-2-14.PubMedPubMed CentralGoogle Scholar
- Gimble JM, Robinson CE, Wu X, Kelly KA: The function of adipocytes in the bone marrow stroma: an update. Bone. 1996, 19 (5): 421-428. 10.1016/S8756-3282(96)00258-X.PubMedGoogle Scholar
- Devlin MJ, Rosen CJ: The bone-fat interface: basic and clinical implications of marrow adiposity. Lancet Diabetes Endocrinol. in pressGoogle Scholar
- Naka T, Nishimoto N, Kishimoto T: The paradigm of IL-6: from basic science to medicine. Arthritis Res. 2002, 4 (Suppl 3): S233-S242. 10.1186/ar565.PubMedPubMed CentralGoogle Scholar
- Ferrari SL, Karasik D, Liu J, Karamohamed S, Herbert AG, Cupples LA, Kiel DP: Interactions of Interleukin-6 Promoter Polymorphisms With Dietary and Lifestyle Factors and Their Association With Bone Mass in Men and Women From the Framingham Osteoporosis Study. J Bone Miner Res. 2004, 19 (4): 552-559. 10.1359/JBMR.040103.PubMedGoogle Scholar
- Grant SFA, Ralston SH: Genes and osteoporosis. Trends Endocrinol Metab. 1997, 8 (6): 232-236. 10.1016/S1043-2760(97)00058-1.PubMedGoogle Scholar
- Pino AM, Rosen CJ, Rodríguez JP: In Osteoporosis, differentiation of mesenchymal stem cells (MSCs) improves bone marrow adipogenesis. Biol Res. 2012, 45: 279-287. 10.4067/S0716-97602012000300009.PubMedGoogle Scholar
- Wang F, Tong Q: Transcription factor PU.1 is expressed in white adipose and inhibits adipocyte differentiation. Am J Physiol Cell Physiol. 2008, 295 (1): C213-C220. 10.1152/ajpcell.00422.2007.PubMedGoogle Scholar
- Tondravi MM, McKercher SR, Anderson K, Erdmann JM, Quiroz M, Maki R, Teitelbaum SL: Osteopetrosis in mice lacking haematopoietic transcription factor PU.1. Nature. 1997, 386 (6620): 81-84. 10.1038/386081a0.PubMedGoogle Scholar
- Tranah GJ, Taylor BC, Lui LY, Zmuda JM, Cauley JA, Ensrud KE, Hillier TA, Hochberg MC, Li J, Rhees BK, Erlich HA, Sternlicht MD, Peltz G, Cummings SR, Study of Osteoporotic Fractures (SOF) Research Group: Genetic variation in candidate osteoporosis genes, bone mineral density, and fracture risk: the study of osteoporotic fractures. Calcif Tissue Int. 2008, 83 (3): 155-166. 10.1007/s00223-008-9165-y.PubMedPubMed CentralGoogle Scholar
- Stienstra R, Duval C, Müller M, Kersten S: PPARs, Obesity, and Inflammation. PPAR Res. 2007, 2007: 95974.PubMedGoogle Scholar
- Smith AG, Sheridan PA, Harp JB, Beck MA: Diet-Induced Obese Mice Have Increased Mortality and Altered Immune Responses When Infected with Influenza Virus. J Nutr. 2007, 137 (5): 1236-1243.PubMedGoogle Scholar
- Lynch LA, O'Connell JM, Kwasnik AK, Cawood TJ, O'Farrelly C, O'Shea DB: Are Natural Killer Cells Protecting the Metabolically Healthy Obese Patient?. Obesity. 2009, 17 (3): 601-605. 10.1038/oby.2008.565.PubMedGoogle Scholar
- Osborn O, Olefsky JM: The cellular and signaling networks linking the immune system and metabolism in disease. Nat Med. 2012, 18 (3): 363-374. 10.1038/nm.2627.PubMedGoogle Scholar
- Matarese G, La Cava A: The intricate interface between immune system and metabolism. Trends Immunol. 2004, 25 (4): 193-200. 10.1016/j.it.2004.02.009.PubMedGoogle Scholar
- Martí A, Marcos A, Martínez JA: Obesity and immune function relationships. Obes Rev. 2001, 2 (2): 131-140. 10.1046/j.1467-789x.2001.00025.x.PubMedGoogle Scholar
- Wellen KE, Hotamisligil G: Obesity-induced inflammatory changes in adipose tissue. J Clin Invest. 2003, 112 (12): 1785-1788. 10.1172/JCI20514.PubMedPubMed CentralGoogle Scholar
- Moncada-Pazos A, Obaya AJ, Llamazares M, Heljasvaara R, Suárez MF, Colado E, Noël A, Cal S, López-Otín C: ADAMTS-12 Metalloprotease Is Necessary for Normal Inflammatory Response. J Biol Chem. 2012, 287 (47): 39554-39563. 10.1074/jbc.M112.408625.PubMedPubMed CentralGoogle Scholar
- Scuteri A, Sanna S, Chen W-M, Uda M, Albai G, Strait J, Najjar S, Nagaraja R, Orrú M, Usala G, Dei M, Lai S, Maschio A, Busonero F, Mulas A, Ehret GB, Fink AA, Weder AB, Cooper RS, Galan P, Chakravarti A, Schlessinger D, Cao A, Lakatta E, Abecasis GR: Genome-Wide Association Scan Shows Genetic Variants in the FTO Gene Are Associated with Obesity-Related Traits. PLoS Genet. 2007, 3 (7): e115-10.1371/journal.pgen.0030115.PubMedPubMed CentralGoogle Scholar
- Baranova A, Collantes R, Gowder S, Elariny H, Schlauch K, Younoszai A, King S, Randhawa M, Pusulury S, Alsheddi T, Ong JP, Martin LM, Chandhoke V, Younossi ZM: Obesity-related Differential Gene Expression in the Visceral Adipose Tissue. OBES SURG. 2005, 15 (6): 758-765. 10.1381/0960892054222876.PubMedGoogle Scholar
- Iatan I, Dastani Z, Do R, Weissglas-Volkov D, Ruel I, Lee JC, Huertas-Vazquez A, Taskinen M-R, Prat A, Seidah NG, Pajukanta P, Engert JC, Genest J: Genetic Variation at the Proprotein Convertase Subtilisin/Kexin Type 5 Gene Modulates High-Density Lipoprotein Cholesterol Levels. Circ Cardiovasc Genet. 2009, 2 (5): 467-475. 10.1161/CIRCGENETICS.109.877811.PubMedPubMed CentralGoogle Scholar
- Choi S, Korstanje R: Proprotein convertases in high-density lipoprotein metabolism. Biomark Res. 2013, 1 (1): 1-8. 10.1186/2050-7771-1-1.Google Scholar
- Alvim RO, Santos PCJL, Nascimento RM, Coelho GLLM, Mill JG, Krieger JE, Pereira AC: BDKRB2 + 9/-9 Polymorphism Is Associated with Higher Risk for Diabetes Mellitus in the Brazilian General Population. Exp Diabetes Res. 2012, 2012: 4.Google Scholar
- Chan Yingguang F, Jones Felicity C, McConnell E, Bryk J, Bünger L, Tautz D: Parallel Selection Mapping Using Artificially Selected Mice Reveals Body Weight Control Loci. Curr Biol. 2012, 22 (9): 794-800. 10.1016/j.cub.2012.03.011.PubMedGoogle Scholar
- Tönjes A, Koriath M, Schleinitz D, Dietrich K, Böttcher Y, Rayner NW, Almgren P, Enigk B, Richter O, Rohm S, Fischer-Rosinsky A, Pfeiffer A, Hoffmann K, Krohn K, Aust G, Spranger J, Groop L, Blüher M, Kovacs P, Stumvoll M: Genetic variation in GPR133 is associated with height: genome wide association study in the self-contained population of Sorbs. Hum Mol Genet. 2009, 18 (23): 4662-4668. 10.1093/hmg/ddp423.PubMedPubMed CentralGoogle Scholar
- Vaittinen M, Kolehmainen M, Schwab U, Uusitupa M, Pulkkinen L: Microfibrillar-associated protein 5 is linked with markers of obesity-related extracellular matrix remodeling and inflammation. Nutr Diabetes. 2011, 1 (8): e15-10.1038/nutd.2011.10.PubMedPubMed CentralGoogle Scholar
- Codoñer-Franch P, Tavárez-Alonso S, Murria-Estal R, Herrera-Martín G, Alonso-Iglesias E: Polyamines Are Increased in Obese Children and Are Related to Markers of Oxidative/Nitrosative Stress and Angiogenesis. J Clin Endocrinol Metab. 2011, 96 (9): 2821-2825. 10.1210/jc.2011-0531.PubMedGoogle Scholar
- Tabassum R, Jaiswal A, Chauhan G, Dwivedi OP, Ghosh S, Marwaha RK, Tandon N, Bharadwaj D: Genetic Variant of AMD1 is Associated with Obesity in Urban Indian Children. PLoS One. 2012, 7 (4): e33162-10.1371/journal.pone.0033162.PubMedPubMed CentralGoogle Scholar
- Heindl K, Martinez J: Nol9 is a novel polynucleotide 5′‒kinase involved in ribosomal RNA processing. EMBO J. 2010, 29 (24): 4161-4171. 10.1038/emboj.2010.275.PubMedPubMed CentralGoogle Scholar
- MacLaren R, Cui W, Simard S, Cianflone K: Influence of obesity and insulin sensitivity on insulin signaling genes in human omental and subcutaneous adipose tissue. J Lipid Res. 2008, 49 (2): 308-323.PubMedGoogle Scholar
- Fabre O, Breuker C, Amouzou C, Salehzada T, Kitzmann M, Mercier J, Bisbal C: Defects in TLR3 expression and RNase L activation lead to decreased MnSOD expression and insulin resistance in muscle cells of obese people. Cell Death Dis. 2014, 5: e1136-10.1038/cddis.2014.104.PubMedPubMed CentralGoogle Scholar
- Mahdi T, Hänzelmann S, Salehi A, Muhammed Sarheed J, Reinbothe Thomas M, Tang Y, Axelsson Annika S, Zhou Y, Jing X, Almgren P, Krus U, Taneera J, Blom AM, Lyssenko V, Esguerra JL, Hansson O, Eliasson L, Derry J, Zhang E, Wollheim CB, Groop L, Renström E, Rosengren AH: Secreted Frizzled-Related Protein 4 Reduces Insulin Secretion and Is Overexpressed in Type 2 Diabetes. Cell Metab. 2012, 16 (5): 625-633. 10.1016/j.cmet.2012.10.009.PubMedGoogle Scholar
- Fuller T, Ghazalpour A, Aten J, Drake T, Lusis A, Horvath S: Weighted gene coexpression network analysis strategies applied to mouse weight. Mamm Genome. 2007, 18 (6–7): 463-472.PubMedPubMed CentralGoogle Scholar
- Das A, Dinh PX, Panda D, Pattnaik AK: Interferon-Inducible Protein IFI35 Negatively Regulates RIG-I Antiviral Signaling and Supports Vesicular Stomatitis Virus Replication. J Virol. 2014, 88 (6): 3103-3113. 10.1128/JVI.03202-13.PubMedPubMed CentralGoogle Scholar
- Lumeng CN: Innate immune activation in obesity. Mol Aspects Med. 2013, 34 (1): 12-29. 10.1016/j.mam.2012.10.002.PubMedGoogle Scholar
- Robertson MJ: Role of chemokines in the biology of natural killer cells. J Leukoc Biol. 2002, 71 (2): 173-183.PubMedGoogle Scholar
- Keophiphath M, Rouault C, Divoux A, Clement K, Lacasa D: CCL5 promotes macrophage recruitment and survival in human adipose tissue. Arterioscler Thromb Vasc Biol. 2010, 30 (1): 39-45. 10.1161/ATVBAHA.109.197442.PubMedGoogle Scholar
- Matter CM, Handschin C: RANTES (Regulated on Activation, Normal T Cell Expressed and Secreted), Inflammation, Obesity, and the Metabolic Syndrome. Circulation. 2007, 115 (8): 946-948. 10.1161/CIRCULATIONAHA.106.685230.PubMedGoogle Scholar
- Leslie RDG, Atkinson MA, Notkins AL: Autoantigens IA-2 and GAD in Type I (insulin-dependent) diabetes. Diabetologia. 1999, 42 (1): 3-14. 10.1007/s001250051105.PubMedGoogle Scholar
- Morgan H, Beck T, Blake A, Gates H, Adams N, Debouzy G, Leblanc S, Lengger C, Maier H, Melvin D, Meziane H, Richardson D, Wells S, White J, Wood J, de Angelis MH, Brown SD, Hancock JM, Mallon AM, EUMODIC Consortium: EuroPhenome: a repository for high-throughput mouse phenotyping data. Nucleic Acids Res. 2010, 38 (Database issue): D577-D585.PubMedGoogle Scholar
- Ahmadian M, Duncan RE, Sul HS: The skinny on fat: lipolysis and fatty acid utilization in adipocytes. Trends Endocrinol Metab. 2009, 20 (9): 424-428. 10.1016/j.tem.2009.06.002.PubMedPubMed CentralGoogle Scholar
- Ahmadian M, Wang Y, Sul HS: Lipolysis in adipocytes. Int J Biochem Cell Biol. 2010, 42 (5): 555-559. 10.1016/j.biocel.2009.12.009.PubMedGoogle Scholar
- Sun J, Zhu J, Xue J, Zhang C, Lan X, Lei C, Chen H: Haplotype combinations of AdPLA gene polymorphisms associate with growth traits in Chinese cattle. Mol Biol Rep. 2012, 39 (6): 7069-7076. 10.1007/s11033-012-1538-7.PubMedGoogle Scholar
- Jaworski K, Ahmadian M, Duncan RE, Sarkadi-Nagy E, Varady KA, Hellerstein MK, Lee H-Y, Samuel VT, Shulman GI, Kim KH, de Val S, Kang C, Sul HS: AdPLA ablation increases lipolysis and prevents obesity induced by high-fat feeding or leptin deficiency. Nat Med. 2009, 15 (2): 159-168. 10.1038/nm.1904.PubMedPubMed CentralGoogle Scholar
- Wuschke S, Dahm S, Schmidt C, Joost HG, Al-Hasani H: A meta-analysis of quantitative trait loci associated with body weight and adiposity in mice. Int J Obes. 2006, 31 (5): 829-841.Google Scholar
- Hou G-Y, Yuan Z-R, Gao X, Li J-Y, Gao H-J, Chen J-B, Xu S-Z: Genetic Polymorphisms of the CACNA2D1 Gene and Their Association with Carcass and Meat Quality Traits in Cattle. Biochem Genet. 2010, 48 (9–10): 751-759.PubMedGoogle Scholar
- Yuan Z, Xu S: Novel SNPs of the bovine CACNA2D1 gene and their association with carcass and meat quality traits. Mol Biol Rep. 2011, 38 (1): 365-370. 10.1007/s11033-010-0117-z.PubMedGoogle Scholar
- Shao N, Jiang WY, Qiao D, Zhang SG, Wu Y, Zhang XX, Hua LX, Ding Y, Feng NH: An updated meta-analysis of XRCC4 polymorphisms and cancer risk based on 31 case–control studies. Cancer Biomark. 2012, 12 (1): 37-47.PubMedGoogle Scholar
- Sweet MJ, Hume DA: CSF-1 as a regulator of macrophage activation and immune responses. Arch Immunol Ther Exp (Warsz). 2003, 51 (3): 169-177.Google Scholar
- Ortega Martinez de Victoria E, Xu X, Koska J, Francisco AM, Scalise M, Ferrante AW, Krakoff J: Macrophage Content in Subcutaneous Adipose Tissue: Associations With Adiposity, Age, Inflammatory Markers, and Whole-Body Insulin Action in Healthy Pima Indians. Diabetes. 2009, 58 (2): 385-393.PubMedPubMed CentralGoogle Scholar
- Pietilainen KH, Kannisto K, Korsheninnikova E, Rissanen A, Kaprio J, Ehrenborg E, Hamsten A, Yki-Jarvinen H: Acquired obesity increases CD68 and tumor necrosis factor-alpha and decreases adiponectin gene expression in adipose tissue: a study in monozygotic twins. J Clin Endocrinol Metab. 2006, 91 (7): 2776-2781. 10.1210/jc.2005-2848.PubMedGoogle Scholar
- Gomez-Ambrosi J, Catalan V, Diez-Caballero A, Martinez-Cruz LA, Gil MJ, Garcia-Foncillas J, Cienfuegos JA, Salvador J, Mato JM, Fruhbeck G: Gene expression profile of omental adipose tissue in human obesity. FASEB J. 2004, 18 (1): 215-217.PubMedGoogle Scholar
- Cattaneo M: P2Y12 receptor antagonists: a rapidly expanding group of antiplatelet agents. Eur Heart J. 2006, 27 (9): 1010-1012.PubMedGoogle Scholar
- Ferreira IA, Mocking AI, Feijge MA, Gorter G, van Haeften TW, Heemskerk JW, Akkerman JW: Platelet inhibition by insulin is absent in type 2 diabetes mellitus. Arterioscler Thromb Vasc Biol. 2006, 26 (2): 417-422.PubMedGoogle Scholar
- Igarashi M, Osuga J, Uozaki H, Sekiya M, Nagashima S, Takahashi M, Takase S, Takanashi M, Li Y, Ohta K, Kumagai M, Nishi M, Hosokawa M, Fledelius C, Jacobsen P, Yagyu H, Fukayama M, Nagai R, Kadowaki T, Ohashi K, Ishibashi S: The critical role of neutral cholesterol ester hydrolase 1 in cholesterol removal from human macrophages. Circ Res. 2010, 107 (11): 1387-1395. 10.1161/CIRCRESAHA.110.226613.PubMedGoogle Scholar
- Kauts M-L, Pihelgas S, Orro K, Neuman T, Piirsoo A: CCL5/CCR1 axis regulates multipotency of human adipose tissue derived stromal cells. Stem Cell Res. 2013, 10 (2): 166-178. 10.1016/j.scr.2012.11.004.PubMedGoogle Scholar
- Hoshino A, Iimura T, Ueha S, Hanada S, Maruoka Y, Mayahara M, Suzuki K, Imai T, Ito M, Manome Y, Yasuhara M, Kirino T, Yamaguchi A, Matsushima K, Yamamoto K: Deficiency of chemokine receptor CCR1 Causes osteopenia due to impaired functions of osteoclasts and osteoblasts. J Biol Chem. 2010, 285 (37): 28826-28837. 10.1074/jbc.M109.099424.PubMedPubMed CentralGoogle Scholar
- Takemura K, Sakashita N, Fujiwara Y, Komohara Y, Lei X, Ohnishi K, Suzuki H, Kodama T, Mizuta H, Takeya M: Class A scavenger receptor promotes osteoclast differentiation via the enhanced expression of receptor activator of NF-kappaB (RANK). Biochem Biophys Res Commun. 2010, 391 (4): 1675-1680. 10.1016/j.bbrc.2009.12.126.PubMedGoogle Scholar
- Rosen CJ, Bouxsein ML: Mechanisms of Disease: is osteoporosis the obesity of bone?. Nat Clin Pract Rheum. 2006, 2 (1): 35-43. 10.1038/ncprheum0070.Google Scholar
- Clowes JA, Riggs BL, Khosla S: The role of the immune system in the pathophysiology of osteoporosis. Immunol Rev. 2005, 208: 207-227. 10.1111/j.0105-2896.2005.00334.x.PubMedGoogle Scholar
- Johansen T, Hansen HS, Richelsen B, Malmlöf K: The Obese Gottingen Minipig as a Model of the Metabolic Syndrome: Dietary Effects on Obesity, Insulin Sensitivity, and Growth Hormone Profile. Comp Med. 2001, 51 (2): 150-155.PubMedGoogle Scholar
- Cameron ND: Selection Indices andn Prediction of Genetic Merit in Animal Breeding. 1997, Madison: CABI; the University of WisconsinGoogle Scholar
- R-Core-Team: R: A language and environment for statistical computing. R Foundation for Statistical Computing. 2013, Vienna, Austria, http://www.R-project.org/.Google Scholar
- Kadarmideen HN, Kogelman LJA: RNA-Sequencing of lean, intermediate, and obese pigs. NCBI's Gene Expression Omnibus. 2014, http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE61271, 1.Google Scholar
- Anders S, Pyl PT, Huber W: HTSeq - A Python framework to work with high-throughput sequencing data. BioRxiv Preprint. 2014Google Scholar
- Łabaj PP, Leparc GG, Linggi BE, Markillie LM, Wiley HS, Kreil DP: Characterization and improvement of RNA-Seq precision in quantitative transcript expression profiling. Bioinformatics. 2011, 27 (13): i383-i391. 10.1093/bioinformatics/btr247.PubMedPubMed CentralGoogle Scholar
- Analysing RNA-Seq data with the DESeq Package. [http://www.bioconductor.org/help/course-materials/2011/BioC2011/LabStuff/DESeq.pdf].
- Smyth GK: limma: Linear Models for Microarray Data. Bioinformatics and Computational Biology Solutions Using R and Bioconductor. Edited by: Gentleman R, Carey V, Huber W, Irizarry R, Dudoit S. 2005, New York: Springer, 397-420.Google Scholar
- Jeong H, Mason SP, Barabasi AL, Oltvai ZN: Lethality and centrality in protein networks. Nature. 2001, 411 (6833): 41-42. 10.1038/35075138.PubMedGoogle Scholar
- Langfelder P, Horvath S: WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008, 9 (1): 559-10.1186/1471-2105-9-559.PubMedPubMed CentralGoogle Scholar
- Zhang B, Horvath S: A General Framework for Weighted Gene Co-Expression Network Analysis. Stat Appl Genet Mol Biol. 2005, 4 (1; Article 17): 1-43.Google Scholar
- Langfelder P, Zhang B, Horvath S: Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R. Bioinformatics. 2008, 24 (5): 719-720. 10.1093/bioinformatics/btm563.PubMedGoogle Scholar
- Haider S, Ballester B, Smedley D, Zhang J, Rice P, Kasprzyk A: BioMart Central Portal—unified access to biological data. Nucleic Acids Res. 2009, 37 (suppl 2): W23-W27.PubMedPubMed CentralGoogle Scholar
- Young M, Wakefield M, Smyth G, Oshlack A: Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010, 11 (2): R14-10.1186/gb-2010-11-2-r14.PubMedPubMed CentralGoogle Scholar
- Segal E, Shapira M, Regev A, Pe'er D, Botstein D, Koller D, Friedman N: Module networks: identifying regulatory modules and their condition-specific regulators from gene expression data. Nat Genet. 2003, 34 (2): 166-176. 10.1038/ng1165.PubMedGoogle Scholar
- Joshi A, Van de Peer Y, Michoel T: Analysis of a Gibbs sampler method for model-based clustering of gene expression data. Bioinformatics. 2008, 24 (2): 176-183. 10.1093/bioinformatics/btm562.PubMedGoogle Scholar
- Edgar R, Domrachev M, Lash AE: Gene Expression Omnibus: NCBI gene expression and hybridization array data repository. Nucleic Acids Res. 2002, 30 (1): 207-210. 10.1093/nar/30.1.207.PubMedPubMed CentralGoogle Scholar
- The pre-publication history for this paper can be accessed here:http://www.biomedcentral.com/1755-8794/7/57/prepub
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly credited. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.