The direction of cross obesity after puberty in male but not female offspring

Background We investigated parent-of-origin and allele-specific expression effects on obesity and hepatic gene expression in reciprocal crosses between the Berlin Fat Mouse Inbred line (BFMI) and C57Bl/6NCrl (B6N). Results We found that F1-males with a BFMI mother developed 1.8 times more fat mass on a high fat diet at 10 weeks than F1-males of a BFMI father. The phenotype was detectable from six weeks on and was preserved after cross-fostering. RNA-seq data of liver provided evidence for higher biosynthesis and elongation of fatty acids (p = 0.00635) in obese male offspring of a BFMI mother versus lean offspring of a BFMI father. Furthermore, fatty acid degradation (p = 0.00198) and the peroxisome pathway were impaired (p = 0.00094). The circadian rhythm was affected as well (p = 0.00087). Among the highest up-regulated protein coding genes in obese males were Acot4 (1.82 fold, p = 0.022), Cyp4a10 (1.35 fold, p = 0.026) and Cyp4a14 (1.32 fold, p = 0.012), which hydroxylize fatty acids and which are known to be increased in liver steatosis. Obese males showed lower expression of the genetically imprinted and paternally expressed 3 (Peg3) gene (0.31 fold, p = 0.046) and higher expression of the androgen receptor (Ar) gene (2.38 fold, p = 0.068). Allelic imbalance was found for expression of ATP-binding cassette transporter gene Abca8b. Several of the differentially expressed genes contain estrogen response elements. Conclusions Parent-of-origin effects during gametogenesis and/or fetal development in an obese mother epigenetically modify the transcription of genes that lead to enhanced fatty acid synthesis and impair β-oxidation in the liver of male, but not female F1 offspring. Down-regulation of Peg3 could contribute to trigger this metabolic setting. At puberty, higher amounts of the androgen receptor and altered access to estrogen response elements in affected genes are likely responsible for male specific expression of genes that were epigenetically triggered. A suggestive lack of estrogen binding motifs was found for highly down-regulated genes in adult hepatocytes of obese F1 males (p = 0.074). Electronic supplementary material The online version of this article (doi:10.1186/s12864-015-2164-2) contains supplementary material, which is available to authorized users.


Background
Besides nutrition and sedentary lifestyle, genetic predisposition contributes to obesity. To assess the genetic determinants of obesity, many large-scale genome-wide association studies were performed in different human populations. These studies led to the identification of about 350 genomic loci contributing to obesity [1]. However, their effects are small and the overall contribution of all detected loci to the variance in a population was small (0.01 % to 0.34 %) [2][3][4][5][6]. Therefore, additional factors could contribute to obesity in an individual. Missing heritability can be hidden in genetic effects of rare alleles, epistatic interaction and gene-byenvironment interaction, but can be also contributed to by parent-of-origin specific effects. Parent-of-origin effects include the Y chromosome, mitochondria, parental-specific genetic imprinting, and maternal environmental effects during pregnancy and the suckling period.
In a mouse model of juvenile obesity, we tested, if parental obesity influences the transmission of predisposition for adiposity to F1 offspring. We used the Berlin Fat Mouse Inbred line (BFMI) which carries a recessive mutation responsible for high fat deposition in particular during early development between four and ten weeks [7]. At ten weeks BFMI males are 5.5 times as obese as C57BL/6Nrl mice (B6N).
In reciprocal crosses between BFMI and B6N mice we tested the development of obesity in F1 males and females and their response to the obese parent with respect to gene expression in the liver. The results provide evidence for regulation of genes by the influence of the direction-of-cross, which leads to obesity and fatty liver exclusively in male offspring of obese BFMI mothers.

Animals
Reciprocal crosses were generated from the Berlin Fat Mouse Inbred line BFMI860-12/Hber (BFMI) and C57BL/6NCrl (B6N). Seven BFMI males were mated with seven B6N females and six B6N males with six BFMI860 females to generate 48 F1 animals (deemed patBFMI and matBFMI, respectively). Each reciprocal F1 offspring groups consisted of 12 males and 12 females. From these animals three F1 males (matBFMI) from three BFMI females and three F1 males (patBFMI) from three B6N females were used for RNA-seq transcriptome measurements. One high fat diet raised male and female from each inbred parental strain (BFMI and B6N) were used as parental strain control animals and measured via RNA-seq as well. For the cross-fostering experiment, six BFMIxB6N and six B6NxBFMI crosses generated two additional groups of reciprocal F1 animals with a total of 88 animals, 46 with a BFMI father (23 males, 23 females) and 42 with a BFMI mother (24 males, 18 females). For cross-fostering, litters were exchanged between B6N and BFMI mothers within 12 hours after birth. The litter size ranged between four and eleven offspring. Two animals per cross-fostering group were excluded due to outlier testing (R-package "outliers").

Phenotyping
Body weight and body composition (total lean mass, total fat mass) were recorded weekly from week three and four. Body composition was measured using EchoMRI whole-body composition analyzer (Echo Medical Systems, Houston, TX, USA). At ten weeks blood glucose levels were measured prior to dissection after 2 hours of fasting using a glucometer (Contour, Bayer, Leverkusen, Germany). Liver triglycerides and proteins were measured in homogenized liver samples by commercially available kits (DC Protein Assay, Bio-Rad, München, Germany; TRIGS Triglycerides Assay, Randox, Crumlin, United Kingdom). For test statistics between the matBFMI and patBFMI, we applied a linear model using anova in R. Since litter size had a significant effect on all traits, phenotypic values were adjusted for litter size before statistic tests were performed, except for blood glucose and liver triglycerides.

RNA preparation and sequencing
RNA sequencing was performed on 10 liver samples: 1 male and 1 female B6N mouse, 1 male and 1 female BFMI mouse, 3 F1 (BFMI x B6N) males (patBFMI), and 3 F1 (B6N x BFMI) males (matBFMI). The F1 males were of non-cross-fostered mice. Tissue samples were frozen in liquid nitrogen and stored at -80°C until RNA preparation. RNA was prepared using the RiboZero RNA library protocol. Paired end RNA sequencing was run on an Illumina HiSeq2000 machine.

Analysis of RNA-Seq data
After sequencing, bcl2fastq v1.8.4 [8] was used for base calling and demultiplexing. Paired end reads were trimmed using Trimmomatic [9] using standard settings. TruSeq3-PE adapter sequences provided by Illumina were used for adapter trimming. Reads were aligned against Mus musculus GRCm38 reference genome (Ensembl) using BWA [10], after which GATK [11] was applied to do base quality score recalibration, indel realignment, duplicate removal, SNP/INDEL discovery and genotyping across all 10 samples using standard hard filtering parameters or variant quality score recalibration according to GATK Best Practices recommendations [12].
Raw reads from RNA-Seq data were extracted from BAM files using the GenomicFeatures package using the Mus_musculus.GRCm38.81.gtf file obtained from ensembl, overlap between genes and raw reads was calculated in 'union' mode using the Rsamtools package [13]. Reads were quantile normalized between samples using the pre-processCore package. After sample quantile normalization, RPKM values for each gene were calculated, which were then log transformed ( ln(n + 1) ) to reduce the effect of extreme RPKM values in subsequent statistical tests. The bio-maRt R package was used to provide additional gene annotation such as MGI descriptions (using the M. musculus ensembl database GRCm38.p4).
We analyzed the gene expression of the reciprocal F1 populations with respect to differentially expressed genes, pathways, and allele specific expression. Testing for differences between the two groups of reciprocal crosses was done using two sided t-tests. Due to limited sample sizes of three animals per F1 group, adjustment for multiple testing was not applied for differential gene expression analysis. The original transcriptome sequencing data for all ten animals tested are available at the GEO database (http://www.ncbi.nlm.nih.gov/geo/).

Pathway analyses
Genes that were differentially expressed (p < 0.05) between the two groups of F1 males with RPKM values > 0.5 in at least one of the F1 groups and gene expression differences of > = 20 % were subjected to overrepresentation analyses in KEGG pathways using innateDB [14] to identify pathways associated with these genes [15]. We used the default hypergeometric algorithm and the Benjamini Hochberg algorithm for multiple testing corrections. Heatmaps were created for selected pathways using the mean gene expression per group.

Analysis of allele specific expression (ASE)
To assess ASE in the reciprocal F1 crosses, bcftools was used to perform SNPs calling on the entire population (all 10 animals: 2 BFMI, 2 B6N, 3 matBFMI, 3 patBFMI). Detected SNPs underwent a rigorous data quality checks to make sure detected SNPs showed sufficient coverage in the entire population. The individuals from the BFMI and B6N pure lines were analyzed to detect homozygous SNPs that showed different alleles between B6N and BFMI. We demanded that all SNPs found in the four parental strain animals had non-missing data (neither low confidence nor low coverage). SNPs were called using the samtools mpileup tool. Subsequently, the generated bcf file was processed using BCFtools [16] and filtered using the provided vcfutils.pl script for SNPs with a 100 reads minimum for the entire population (over all 10 animals). After using the entire population for SNP calling, we reanalyzed all SNPs per individual using the samtools mpileup tool to obtain DP4 measurements for each F1 individual. SNPs with low confidence (Phred probability threshold < 50) or less than 5 high quality reads in an individual were excluded from further analyses. Allele ratios alt/(ref + alt) were calculated for all remaining SNPs. Additionally, we filtered for consistency between the different individuals from the same F1 cross direction. If the standard deviation of the alt/(ref + alt) ratios across the 3 individuals per F1 group was > 20 %, the SNP was removed due to the fact that the data was not consistent across the individuals.
After filtering for high quality SNPs which are homozygous (but different) between the parental strains, χ 2 scores were calculated for each individual. χ 2 scores were adjusted to account for direction of the effect (< 50 % vs. > 50 %) and scores were square root transformed to minimize the impact of a single F1 individual (in the F1 groups) showing a highly significant deviation from the expected 50:50 ratio. Summed scores across all individuals in a F1 cross direction (matBFMI and patBFMI) were calculated to obtain a group based χ 2 score, which were then subtracted from each other to get a group difference score per SNP. To account for over dispersion in the read counts obtained from the DP4 values, we permuted our group based χ 2 scores to obtain a 5 % FDR threshold for SNPs that show significant ASE in either the patBFMI or the matBFMI group.

Analysis of transcription factor binding sites
Differentially expressed genes (p < 0.05) were submitted to transcription factor binding site (TFBS) analysis. For each gene we extracted the upstream (2000 bp) and downstream (1000 bp) region around the transcription start using the TxDb.Mmusculus.UCSC.mm10 package from bioconductor [17]. All known mouse TFBS motifs were downloaded from the Jaspar database into R using the MotifDb package [18,19].
Since only one weight matrices (jolma2013-Ar) was available in the Jasper database for androgen response elements (AREs) in mouse, we created position weight matrices for androgen response elements using the motifs from Shaffer et al. for ARE, ARE2, ADR3, IR3 [23]. We used a conservative threshold of 90 % for sensitivity and specificity. A list of the applied position weight matrices is given in Additional file 1: Table S3.
Additionally, we compared the occurrence of TFBS in both differentially and allele specific expressed genes with 1,000 permutations of random sets of genes of the respective same group size.

Quantitative PCR
RNA from ten male matBFMI and ten male patBFMI was extracted from liver tissue. Total RNA was isolated from liver tissue using the NucleoSpin RNA II kit (Macherey-Nagel, Düren, Germany). One μg of total RNA was reverse-transcribed into cDNA using Accu-Script High Fidelity Transcriptase (Stratagene Europe, Amsterdam, Netherlands) and Oligo (dT) primers (Ark Scientific, Sigma-Aldrich, München, Germany). The cDNA integrity was checked in a polymerase chain reaction (PCR). For semi-quantitative real-time PCR 20 ng cDNA was mixed with the qPCR mastermix of the Plus Kit for SYBR Green (Eurogentec, Cologne, Germany) and 300 nM of the primers at a final volume of 10 μl per reaction. The amplification protocol included a 10 min DNA denaturation step at 95°C, followed by 40 cycles consisting of 20 s at 60°C and 40 s at 72°C. Quantification of mRNA transcript amounts was performed using the Applied Biosystems® ViiA™ 7 Real-Time PCR System and ViiA™ Software v.1.2.1. Primers for target genes Ar, Peg3 and Srebf1 were designed using Primer3web [24]. Sequences for reference gene primers (Actb, Mrpl46, Srp72) were taken from Hruz et al. [25]. The geometric mean of the three reference genes was used to normalize target gene expression. Relative quantification of gene expression was performed according to ABI's user guide [26]. Further details are provided in Additional file 2: Table S4.

Phenotypes
Since BFMI animals carry a major recessive mutation leading to juvenile obesity in homozygous carriers [27], we expected lean F1-offspring in either cross direction. However, data provided evidence for the development of obesity in F1 males, if the mother originates from the obese BFMI line (matBFMI), while females of the same cross direction and both sexes of the reciprocal cross had normal body weight and fat mass (Table 1). Body mass did not differ significantly between offspring of the reciprocal crosses, but F1 males from BFMI mothers (matBFMI) had 1.83 times more adipose tissue mass at 10 weeks (p = 0.0037) than F1 males from BFMI fathers (patBFMI). Since the lean mass was only marginally lower in the obese matBFMI versus the lean patBFMI F1 males (0.98 fold, p = 0.054), the fat to lean mass ratio was also 1.86 times higher in matBFMI males (p = 0.0009). Moreover, matBFMI also showed ectopic fat storage in the liver with fat content being 1.30 times higher than in patBFMI males from the reciprocal cross (p = 0.084). The higher fat mass in matBFMI males began at puberty and was significant from seven weeks on (Fig. 1).
A cross-fostering experiment was performed to assess the influence of postnatal maternal effects such as milk energy content and behavior during the suckling period on later development of fat deposition in offspring. This shown are means ± standard deviation after correction for litter size, except for blood glucose and liver triglycerides, which were not affected by the litter size.
Asterisks additionally indicate significant effects of cross fostering on phenotypes comparing animals of the same direction of cross (patBFMI vs. cross-fostered patBFMI, matBFMI vs. cross-fostered matBFMI) (**p < 0.01, *p < 0.05,°p < 0.1) experiment validated the results obtained before, but fat deposition was a little lower ( Fig. 1). Cross-fostering led to a reduction of adipose tissue mass in male and female F1 offspring of BFMI mothers by approximately 22 % if they were raised by lean B6N foster mothers. Crossfostered F1 males of BFMI mothers still developed 1.36 times higher fat mass (p = 0.0142), had a 1.36 times higher fat to lean ratio (p = 0.0067), 1.40 times higher hepatic triglyceride levels (p = 0.03070) and higher blood glucose levels (1.10 fold, p = 0.03073) than cross-fostered F1 males of the reciprocal cross (Table 2). These data support an effect of the suckling period on later obesity and diabetes development that has been shown before [28][29][30].

Differential gene expression between F1 males of BFMI mothers and BFMI fathers
Since lipid accumulation of F1 males differed between the reciprocal crosses, we investigated the direction-of-cross effect on hepatic gene expression profiles in males. Overrepresentation analyses of genes that were differentially expressed between the two groups of reciprocal F1 males at p < 0.05 and a fold-change difference of at least 20 % provided evidence for alterations of pathways for lipid metabolism, peroxisome, circadian rhythm, cellular transport processes and the endocrine system as a results of the direction of cross ( Table 2, Table 3). Most genes with higher expression in the obese F1 males control the cellular utilization of fatty acids for energy production or fat synthesis including the hydrolysis of triglycerides into free fatty acids and glycerol, the conversion of free fatty acids into fatty acyl-CoA esters, the biosynthesis of saturated fatty acids or the degradation and conversion into acetyl-CoA by β-oxidation [31] (Fig. 2). But also genes for the transport between peroxisome and mitochondria via carnitine [32,33] and the electron carrier generation by the citrate cycle and the respiratory chain in mitochondria were differentially expressed. The metabolic changes in obese males were associated with a shift of the circadian rhythm, in which the gene Fig. 1 Development of phenotypes of males and females on a high fat diet of the two reciprocal F1-crosses until 10 weeks. The figures show means ± standard deviations of the total fat mass and the total fat mass-to-total-lean mass-ratios in F1-animals without (a and b, respectively) and with cross-fostering (c and d, respectively). Paternal and maternal BFMI (pat, mat) indicates direction of cross. Asterisks indicate significant differences between different reciprocal crosses for males (blue) and females (red) and significant differences within the same reciprocal cross between males and females of matBFMI (orange) and patBFMI (grey) (***p < 0.001, **p < 0.01, *p < 0.05)  Table 2). Disruptions of the normal circadian rhythm have been linked to obesity in Clock [34] and Arntl (alias Bmal1) knock-out mice and feeding experiments [35]. A shift of the circadian rhythm in offspring as a result of in utero exposure to maternal obesity has been shown in rats [36].
Indications for impaired β-oxidation in the obese males were found from up-regulated expression of cytoplasmic acetyl-CoA carboxylase (Acaca, 1.20 fold, p < 0.06) and citrate lyase (Acly, fold 1.11, p < 0.08), which usually inhibits the cytosol-mitochondria transport of the fatty acid derived substrates via allosteric inhibition of Crat (1.27 fold, p = 0.01) by malonyl-CoA [37][38][39]. Furthermore, acyl-CoA thioesterases (Acot) are transporters that regulate the amounts of intermediates from fatty acid oxidation inside mitochondria. The mitochondrial matrix bound Acot2 (1.62 fold, p = 0.05), for example, enhances hepatic fatty acid oxidation by its transporter activity to prevent substrate overload [40], which otherwise could impair β-oxidation [41]. In obese F1 males also Acot4 (1.82 fold, p < 0.02) and Acot3 (6.49 fold, p < 0.08) were upregulated, which could indicate high amounts of such intermediates inside the mitochondrial matrix, leading to impaired fatty acid oxidation and increased levels of free fatty acids.
In a situation, in which fatty acids are preferentially synthesized and accumulated, while their β-oxidation is impaired, the cell has to increase energy generation from alternative fuels. Previous measurements of the respiratory quotient of the BFMI inbred mice provided evidence for a higher use of carbohydrates as the main energy fuel [44]. Succinyl-CoA, for example, is a citrate cycle intermediate that can be produced from degradation of certain amino acids but also from the propanoate pathway. The moderate up-regulation of Mgll (1.17 fold, p = 0.01), Aldh3a2 (1.11 fold, p = 0.01) and Dak (1.09 fold, p = 0.02) from the glycerolipid pathway aims at the conversion of triglycerides to precursor substrates for subsequent glycolysis. Increased carbohydrate metabolism is further supported by up-regulation of the carbohydrate transporters Slc16a7 (1.20 fold, p = 0.0007) and Slc16a12 (1.13 fold, p = 0.02) in matBFMI (Table 3,  Additional file 3: Table S1).
In humans, Slc16a7 catalyzes the rapid transport of lactate, pyruvate, branched-chain oxo acids and other substrates across the plasma membrane of mitochondria and is regulated by estradiol [45,46]. This regulation could contribute to the fact that female F1 offspring of the BFMI mother are protected against elevated hepatic triglyceride storage.
Among the highest up-regulated genes in matBFMI were the cytochrome P450, family 4, subfamily A members Cyp4a10 (1.35 fold, p = 0.03) and Cyp4a14 (1.32 fold, p = 0.01), which together with an increase of their rate limiting regulator Por (1.16 fold, p = 0.06) is an indicator for increased ω-oxidation [47], a prerequisite for subsequent peroxisomal β-oxidation. These genes were found to play a part in the development of liver steatosis [48,49]. The cytochrome P450 proteins are involved in the synthesis of cholesterol, steroids, and other lipids. P-value for differential gene expression between reciprocal F1 groups of males was p < 0.05, minimum fold-change difference >20 %, The innateDB tool was used for the KEGG pathway analysis. matBFMI -BFMI mother, patBFMI -BFMI father 2 Steroid hormone transcription factor binding sites. Esrra/bestrogen related receptor a/b, ERestrogen receptor alpha, EREestrogen response element, IR3 -androgen inverted repeats of hexameric half-site separated by 3 bp of spacer The proteins encoded by Cyp4a10 and Cyp4a14 hydroxylize fatty acids [50,51]. Recently, fatty acid-hydroxy fatty acids (FAHFAs) have been discovered as a class of lipids with anti-diabetic and anti-inflammatory effects [52]. Higher concentrations of FAHFAs in the obese F1 males due to increased expression of the two Cyp4A genes would be consistent with similar blood glucose levels in male F1 offspring of both cross directions ( Table 1). The hepatic expression of the androgen receptor (Ar) gene was higher in the obese matBFMI F1 males compared to the lean counterparts (2.38 fold, p = 0.07). The androgen receptor functions as a steroid-hormone activated transcription factor [53]. Therefore, Ar could Bold gene symbols indicate up regulation in F1 males of a BFMI mother. A list of all differentially expressed genes (p < 0.1) is given in Additional file 1: Table S1. Chrchromosome, matBFMI -BFMI mother, patBFMI -BFMI father, Fchfold-changes 2 All values were log-transformed (ln) 3 TFBS -Transcription factor binding sites for steroid hormones. Esrrbestrogen related receptor b, EREestrogen response element, IR3androgen inverted repeats of hexameric half-site separated by 3 bp of spacer contribute to the sex specific obesity effect. The androgen receptor affects the liver metabolism differently in males and females [54]. Besides the liver effects, the androgen receptor protects females against high fat diet-induced obesity and dyslipidemia [55]. Additional studies showed that mice lacking the androgen receptor develop hepatic steatosis and insulin resistance in males but not females under high fat diet [56]. Among the genes that were highly differentially expressed between the two groups of F1 males (p < 0.05, fold-change difference min. 20 %, Table 3), six and two genes contain binding sites for estrogen (Acnat2, Cyp2c38, Cyp2c39, Cyp4a10, Fasn, Nlrp12) and androgen (Acnat2, Cyp4a10), respectively. Considering all differentially expressed genes under p < 0.05 estrogen (ER), androgen (AR) or glucocorticoid response (GR) elements occur in 220 up-(55 ER, 14 AR, 3 GR) and 189 downregulated (30 ER, 11 AR, 6 GR) genes in the obese versus lean F1 males (Additional file 3: Table S1).
We compared the occurrence of these transcription factor binding sites in genes being differentially expressed (p < 0.05, >20 % fold-change difference) between the obese and the non-obese F1 males with randomly selected genes. We found that binding motifs for estrogen seemed to be under-represented in downregulated genes (p = 0.074). However, when relaxing our restraints and define our down-regulated genes as p < 0.05 and do not filter for only the highest fold-changes (> 20 %), this under-representation of estrogen binding sites becomes highly significant (p = 0.004) while androgen binding sites tended to be over-represented in up-regulated genes (p = 0.056). Therefore, modifications of these genes as a result of epigenetic reprogramming during genetic imprinting or during pregnancy could modify their expression in response to sex-hormone production beginning at puberty.

Gene expression of hemizygous genetic elements
Differential gene expression on the X and Y chromosome and mitochondria was examined in particular because these homozygous genomic units can be polymorphic between the reciprocal crosses. Among them are the androgen receptor (Ar: 2.38 fold, p = 0.07), zinc finger protein X-linked (Zfx: 0.84 fold, p = 0.09) and polyglutamine binding protein 1 (Pqbp1: 1.07 fold, p < 0.1), which act as transcription factors. Differentially expressed X chromosomal genes that are involved in the lipid metabolism are ATP-binding cassette, sub-family D (ALD), member 1 (Abcd1: 1.22 fold, p = 0.04), HAUS augmin-like complex, subunit 7 (Haus7: 1.21 fold, p = 0.05), NAD(P) dependent steroid dehydrogenase-like (Nsdhl: 1.10 fold, p = 0.03), monocarboxylic acid transporter (Slc16a2) and glycerol kinase (Gyk) ( Table 4, Additional file 3: Table  S1). Zfx is well studies in cancer as a proliferation factor; it is also involved in animal size and fertility and lower gene expression was found to inhibit cell growth [57,58]. Dysfunctional Pqbp1 was associated to lipid metabolism in C. elegans. Repression of this gene led to reduction of lipid content in mammalian primary white adipocytes and also triglyceride synthesis from fatty acids was impaired [59]. Since no differences in expression were found for genes on the Y-chromosome, we provide further evidence that the Y chromosome can most likely ruled out as cause for differences in obesity of males. In the class of mitochondrial genes we found mitochondrial NADH dehydrogenase 1 (mt-Nd1) expressed at a lower level (0.95 fold, p = 0.02) in obese matBFMI F1 males (Table 4).

Quantitative PCR validation
We verified our gene expression fold-change results from the RNA-seq analyses performing qPCR for Ar, Peg3 and Srebf1. The qPCR results confirmed the direction of fold-changes in the obese matBFMI animals for  Table S4). Srebf1 was found to be 1.66 fold upregulated, although with lower significance (p = 0.065) as compared to the RNA-seq results (p = 0.027). It is possible that the upregulation of the androgen receptor is a reaction to low androgen levels that are known to be caused by obesity, though a study in cell culture found androgens to inhibit Peg3 expression [60].

Allele specific expression
Since allele-specific expression could result from parental imprinting, we exploited the sequence information of expressed genes. After rigorous quality check, we found  Table S2). These SNPs were used to search for parent-of origin-specific expression in the liver of 10 week old males. Deviations from balanced biallelic expression can cause different phenotypes by causing total gene expression differences but also independently since allele functions can be altered due to sequence variants, altering the function of gene products. We found evidence for allelic imbalances for protein coding genes Abca8b, Commd1, H13, Harbi1, Mcts2, Peg3, Plg, Prkd3 and Zrsr1 (Table 5). Paternal allele expression was identified for paternally expressed 3 (Peg3) and zinc finger (CCCH type), RNA binding motif and serine/arginine rich 1 (Zrsr1). Maternal allele expression was found for the gene H13. Previously, all three genes have been described to be genetically imprinted during gametogenesis and differentially expressed in offspring [61][62][63][64]. However, even if parent-of-origin specific expression is proven, the expression level of the paternal allele could differ between offspring due to cis or trans acting regulation that cannot be derived from RNA-seq data [64]. This was the case for Peg3 in our study. Interestingly, the obese matBFMI males showed only 0.31 fold expression of Peg3 (p = 0.046), which was one of the genes with the most extreme decrease in expression. Peg3 is a Krüppel C2H2-type zinc finger protein, which plays a role in the olfactory system related sexual learning, cell proliferation and p53-mediated apoptosis [65,66]. In zebrafish, Peg3 inhibits the Wnt pathway [67]. In our data, lower expression of Peg3 in in obese males was accompanied by lower expressions of the Wnt target genes cadherin 1 (Cdh1: 0.75 fold, p = 0.03), epidermal growth factor receptor (Egfr: 0.85 fold, p = 0.05) and neuropilin 1 (Nrp1: 0.90 fold, p = 0.02). High Peg3 expression was recently suggested to protect from high-fat diet induced obesity in F1 offspring of reciprocal crosses between PWK and B6 [68]. The 'paternal transmission' of the Peg3 effect was deduced from the observation that B6 develop Very little is known about the function of the imprinted gene Zrsr1 although it was associated with the pluripotency state in induced pluripotent stem cells [69]. Furthermore, it could be shown that trans-acting DNA binding proteins, such as methyl-CpG binding domain proteins 2 (Mecp2), methyl-CpG binding domain protein Mbd1 and Mbd3 bind at maternal and paternal Zrsr1 alleles differently [70]. In our F1 male mice, Zrsr1 was 0.80 fold expressed in the obese F1 males (p = 0.08). The predominant expression of the maternal H13 allele was confirmed by our data and the gene was only 0.94 fold expressed (p < 0.05) in the matBFMI males.
Differences in expression of imprinted genes could be causal for regulation of down-stream genes contributing to obesity. Recently, paternal intergenerational metabolic reprogramming was found in drosophila, where paternal fat-diet led to changes in the germ cell chromatin state, transcription patterns and modified offspring metabolism. The genetic factors identified in drosophila have been verified in mouse and human studies [71].
We also found a lower proportion of the BFMI Abca8b allele expressed by the obese matBFMI males. While the patBFMI expressed between 58-68 %, the matBFMI expressed only 31-46 % of the BFMI allele. This corresponds to a mean difference for the alternative allele proportion of 24 % under a fold-change of 0.87 (p = 0.06) in total gene expression. ATP-binding cassette (ABC) transporter genes represent the largest family of transmembrane proteins, which arrange the transport of various molecules across all cell membranes utilizing energy from bound ATP. Most ABC genes move compounds from the cytoplasm to the outside of the cell or into an extracellular compartment (endoplasmic reticulum, mitochondria, peroxisome). ABC transporters shuttle hydrophobic compounds within the cell as part of a metabolic process or outside the cell for transport to other organs or secretion from the body. The Abca subfamily includes the largest Abc genes (several >2000 amino acids) and its members are involved in transport of vitamin A (retinol, etc.) derivatives as well as in disorders of cholesterol transport and high-density lipoproteins (HDL) biosynthesis [72]. Our data from pathway analyzes and ASE indicate that altered functions of ABC transporters (mainly subfamily Abca and Abcd) contribute to the accumulation of body fat in the obese matBFMI males.

Conclusions
In reciprocal crosses between the obese line BFMI and the lean line B6N, we observed high fat deposition with and without cross-fostering beginning at puberty only in male F1 offspring of obese BFMI mothers, but not in females of the same cross and males and females of the reciprocal cross. Gene expression analysis of liver RNA provided evidence for higher expression of genes controlling fat deposition and impaired β-oxidation of lipids as well as effects on the circadian rhythm in the obese F1 males. The obese F1 males likely use more carbohydrates and amino acids for energy production than lean males. Among the extremely up-regulated genes were Cytochrome A4 genes leading to hydroxylation of fatty acids. Furthermore, typical imprinted genes were differentially expressed between the two groups of males in reciprocal crosses.
The observations led us to the conclusion that parent-of-origin effects during gametogenesis and pregnancy play a significant role for the development of obesity in later life. Postnatal maternal effects significantly influence fat deposition, but could not entirely overwrite perinatal parent-of-origin-effects. Reprogramming of gametes and metabolic imprints during pregnancy set epigenetic marks that lead to differential expression of genes in response to sexhormone activation beginning at puberty. We suggest that the obese F1 males store more fat due to impaired β-oxidation and cellular compound transport. They compensate the lower energy production from fatty acids by an increased use of alternative fuels for the mitochondrial and non-mitochondrial energy production. The hydroxylases Cyt4a10 and Cyp4a14 are suggested to contribute to increased fatty liver, but the hydroxylated fatty acids are likely nonhazardous. The imprinted genes Peg3, Zrsr1, und H13 and the androgen receptor might play important sex-specific roles for affected down-stream signaling and metabolic pathways. We further suggest that the Y and X chromosomes as well as mitochondria play a minor role in the identified phenotypic differences.
Since the higher adipose tissue mass in F1 males of the BFMI mother was maintained during cross-fostering, parent-of-origin effects of genomic reprogramming during gametogenesis and maternal effects during pregnancy could be causal for the elevated fat accumulation.
Although, we found indicators for substrate inhibition circuits in mitochondria, the hypothesis of impaired substrate selection or mitochondrial substrate overload by fatty acids in BFMI mitochondria needs further research. Furthermore, we like to point out that allele-specificexpression analysis with RNA-seq data requires high control standards. Allele-specific expression could only be tested, if sequence variation was identified between the examined mouse strains. Since the analysis performed in this study was very conservative, the data we present here are highly reliable. Nonetheless, we miss