Spontaneous liver disease in wild-type C57BL/6JOlaHsd mice fed semisynthetic diet

Mouse models are frequently used to study mechanisms of human diseases. Recently, we observed a spontaneous bimodal variation in liver weight in C57BL/6JOlaHsd mice fed a semisynthetic diet. We now characterized the spontaneous variation in liver weight and its relationship with parameters of hepatic lipid and bile acid (BA) metabolism. In male C57BL/6JOlaHsd mice fed AIN-93G from birth to postnatal day (PN)70, we measured plasma BA, lipids, Very low-density lipoprotein (VLDL)-triglyceride (TG) secretion, and hepatic mRNA expression patterns. Mice were sacrificed at PN21, PN42, PN63 and PN70. Liver weight distribution was bimodal at PN70. Mice could be subdivided into two nonoverlapping groups based on liver weight: 0.6 SD 0.1 g (approximately one-third of mice, small liver; SL), and 1.0 SD 0.1 g (normal liver; NL; p<0.05). Liver histology showed a higher steatosis grade, inflammation score, more mitotic figures and more fibrosis in the SL versus the NL group. Plasma BA concentration was 14-fold higher in SL (p<0.001). VLDL-TG secretion rate was lower in SL mice, both absolutely (-66%, p<0.001) and upon correction for liver weight (-44%, p<0.001). Mice that would later have the SL-phenotype showed lower food efficiency ratios during PN21-28, suggesting the cause of the SL phenotype is present at weaning (PN21). Our data show that approximately one-third of C57BL/6JOlaHsd mice fed semisynthetic diet develop spontaneous liver disease with aberrant histology and parameters of hepatic lipid, bile acid and lipoprotein metabolism. Study designs involving this mouse strain on semisynthetic diets need to take the SL phenotype into account. Plasma lipids may serve as markers for the identification of the SL phenotype.


Introduction
Inbred mouse strains are frequently used in biomedical research, often in combination with standardized (semisynthetic) diets. The rationale behind these choices is to minimize genetic a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 animal supplier. The various cohorts were subsequently reared, maintained and sacrificed. Virgin C57BL/6JOlaHsd breeders (12 weeks of age at delivery, Envigo, The Netherlands) were acclimatized for at least 2 weeks at our facility and were fed AIN-93G from arrival onward. C57BL/6JOlaHsd breeders were fed AIN-93G during breeding, pregnancy and lactation. Breeders were housed in groups (females), or individual (males) prior to breeding, and were mated in 2-3F+1M groups. The breeding paradigm used in this study was similar to that applied in earlier studies by us [8,9]. Males were removed after 2 days. Pregnancy was confirmed by a >2 g increase in body weight after 1 week. Upon confirmation of pregnancy, females were housed individually and were not bred again. Only the F 1 offspring were used in this study. Nonpregnant females were mated again the following week, for a maximum of 4 times or until pregnancy was confirmed. Each female breeder was used a maximum of one pregnancy. Delivery day was recorded as postnatal day (PN) 0. Pups were randomized between dams, and litters were culled to 4M+2F at PN2, weaned at PN21 and males were housed either in pairs or solitarily (solo). Breeders and female offspring, not further used in this study, were terminated (CO 2 ) at weaning. A C57BL/6J control cohort was bred from the colony maintained at our local animal facility (Central Animal Facility, University of Groningen).

Study design
The term 'cohorts' refer to a group of mice sacrificed on the same postnatal day. At weaning (PN21), C57BL/6JOlaHsd mice were pair-housed (entire PN42 cohort, entire PN63 cohort, pair-housed PN70 cohort) or solo-housed (solo-housed PN70 cohort only). Mice were sacrificed at PN21, PN42, PN63 or PN70. In a subset of the PN70 cohort, VLDL-TG secretion was determined (see below). S7 Fig shows a visual representation of the study design. Reported 'n = ' represent total number of mice within a cohort, unless specified to be 'SL' or 'NL'. Pair housed C57BL/6J control mice, bred and kept at the Central Animal Facility of the University of Groningen, were fed low-fat semisynthetic control diet (D12450J) during PN56-140. The calculated macro-and micronutrient contents of the diets are given in S4 Table. The fatty acid composition of the AIN-93G diet is provided in S3 Table.

Body composition
Lean and fat mass were quantified by time-domain nuclear magnetic resonance (LF90II, Bruker Optics, Billerica, MA). Mice were placed into a semi-transparent plastic tube, which was subsequently placed into the bore of the NMR machine. The tube minimized body movements. Body composition analysis took approximately 90 seconds per mouse. Body composition analysis did not require fasting or anesthesia [9]. Simultaneously, body weights and food weights were recorded. Measurements were performed from PN21 once per week until termination at PN70. Germany) [10] (1 g/kg BW in~200 μl sterile PBS), retro-orbital blood was drawn at 0, 1, 2, 3 & 5 h in 9 h fasted (midnight-9AM, food taken away at midnight) mice. Blood was drawn during the light-phase. In total, less than 0.2 ml blood was withdrawn per mouse. After the last time point, mice were anaesthetized (isoflurane & O 2 ) and sacrificed by heart puncture and cervical dislocation. A terminal blood sample was collected. Triglycerides (TG) were measured enzymatically (see below). The VLDL-TG secretion rate was calculated from the plasma TG slope over 5 h (mM.h -1 ) multiplied by the estimated plasma volume (28.17 ml/kg BW) [10,11], giving the TG secretion rate (μmol.h -1 .kg BW -1 ). It was assumed that all mice had approximately the same plasma volume. Poloxamer 407 was assumed to have a strong effect on gene expression. Thus, tissues obtained from the VLDL-TG secretion cohort were not used for further postmortem analyses.

Termination
At PN21 unfasted mice were anaesthetized at 9AM and sacrificed by heart puncture (n = 12). At PN42, fasted mice (9AM-1PM) were anaesthetized (isoflurane & O 2 ) and sacrificed by heart puncture (n = 14). At PN63 (n = 14) and PN70 (n = 14 pair, n = 8 solo-housed), fasted (midnight-9AM) mice were anaesthetized (isoflurane & O 2 ) and sacrificed by heart puncture. A terminal blood sample was collected. It was assumed that our primary readout parameter (liver weight) would not be largely affected by a duration of fasting up to 9 h. At PN21, pups were housed with their dam and may have been breastfed prior to termination. To minimize animal suffering or discomfort, we chose not to wean the pups prior to PN21 to allow for fasting [12]. At PN42, given the still relatively young age of the mice, we chose a fasting period of 4 h to minimize animal suffering or discomfort. In adult mice (at PN63 and PN70), we chose a standard 9 h fasting period as it was not expected that this would cause unacceptable levels of animal suffering or discomfort. Blood samples were centrifuged at 2,000×g for 15 minutes at 4˚C in a tabletop centrifuge (Eppendorf, Nijmegen, the Netherlands). Liver weights were recorded. Whole tissues were snap-frozen in liquid N 2 and later cryogenically crushed to powder using a mortar and pestle.

Plasma bile acids
Plasma bile acid species were quantified by liquid chromatography-mass spectrometry. To 25 μl of plasma, we added a mixture of internal standards (isotopically labeled bile acids). Samples were centrifuged at 16,000×g for 10 min in a tabletop centrifuge (Eppendorf, Nijmegen, the Netherlands) and the supernatant was transferred and evaporated at 40˚C under a stream of N 2 . Samples were reconstituted in 200 μl methanol:water (1:1), mixed and centrifuged at 1,800×g for 3 min. The supernatant was filtered using a 0.2 μm spin-filter at 2,000×g for 10 min. Filtrates were transferred to vials and 10 μl was injected into the LC-MS/MS system. The LC-MS/MS system consisted of a Nexera X2 Ultra High Performance Liquid Chromatography system (SHIMADZU, Kyoto, Japan), coupled to a Sciex Qtrap 4500 MD triple quadrupole mass spectrometer (SCIEX, Framingham, MA, USA). Data were analyzed with Analyst MD 1.6.2 software.

Fatty-acyl chain profiling
Fatty acid methyl esters (FAMEs) were quantified using gas chromatography (GC) [13]. Cryogenically crushed tissues were homogenized in Potter-Elvehjem tubes (Sigma, St. Louis, MO, USA) in ice-cold phosphate buffered saline (PBS, Gibco, Fisher Scientific, Landsmeer, the Netherlands) solution. A known quantity of homogenized tissue, food, or plasma was transferred to Sofirell tubes, and capped with silicone-PTFE septum screw caps. An internal standard (heptadecanoic acid, C17, Sigma, St. Louis, MO, USA) was added. Lipids were transmethylated at 90˚C for 4 h in 6 M HCl:methanol (ratio 1:5), liquid-liquid extracted twice using hexane, transferred to a clean tube, dried at 45˚C under a stream of N 2 , reconstituted in hexane (n-hexane PA, Merck) and transferred to GC vials (Aluglas, cat. no. 1013679, APG Europe, Uithoorn, the Netherlands) with inserts (Aluglas, cat. no. 1013586, APG Europe) and caps (VWR, cat. no. 548-0085, Amsterdam, the Netherlands). Samples were analyzed by gas chromatography [13]. The GC system consisted of 6890N network gas chromatograph (Agilent, Middelburg, the Netherlands) and was equipped with a HP-ULTRA 1 dimethylpolysiloxane, nonpolar column (50 m length x 0.2 mm diameter, 0.11 μm film thickness; Agilent, Middelberg, the Netherlands). Samples were injected at 275˚C using a 7683 ALS autosampler. The initial column temperature was 160˚C and was ramped up at 2˚C/min to 240˚C, followed by a second-stage ramp-up at 10˚C/min to 290˚C. The carrier gas was helium at 34 kPa prepressure. The flame ionization detector (FID) was set to 300˚C and recorded each sample for 1 hour at 1 Hz. Data were integrated and analyzed using Atlas Chromatography Data System software.

Liver lipids and total protein content
Total lipids were extracted from liver homogenates using the Bligh & Dyer method [14]. In brief, 50 μl liver homogenate was mixed with 750 μl water, 3 ml chloroform:methanol (1:2), 1.2 ml water and 1 ml chloroform (vortexed at every step) and centrifuged 10 min at 1,000×g in a swinging bucket centrifuge (Roto Silenta, Hettich, Geldermalsen, the Netherlands). The bottom fraction was transferred to a clean tube, dried at 50˚C under a stream of N 2 , and reconstituted in 1 ml water with 2% triton X-100 (Fisher Scientific, Landsmeer, the Netherlands). From there, using the abovementioned commercially available kits, liver triglycerides, cholesterol and NEFA were quantified. Hepatic protein content was determined from diluted liver homogenates using a commercial BCA protein assay (Pierce, Thermo Fisher, Landsmeer, the Netherlands).

Acylcarnitine profiling
Using liquid chromatography-tandem mass spectrometry (LC-MS/MS), acylcarnitine species were quantified [15]. To 10 μl plasma, or 50 μl liver homogenate, a mixture of internal standards (isotopically labeled acylcarnitine species) and acetonitrile was added. Samples were mixed and centrifuged (15,000×g, 10 min) to precipitate proteins. Supernatant was transferred to GC vials. Samples were analyzed using LC-MS/MS [15]. The LC-MS/MS system consisted of an API 3000 LC-MS/MS system equipped with a Turbo ion spray source (Applied Biosystems/MDS Sciex, Ontario, Canada). Data were analyzed with Analyst and Chemoview software (Applied Biosystems/MSDSciex, Ontario, Canada).

RNA isolation and gene expression analysis
Using TRI-Reagent (Sigma, St. Louis, MO, USA), total RNA was extracted from cryogenically crushed whole livers. RNA was quantified by NanoDrop (NanoDrop Technologies, Wilmington, DE, USA). RNA integrity was confirmed by observing the 18S and 28S ribosomal RNA bands on 1% agarose gel (Ultra pure agarose, Thermo Fisher) in Tris-acetate-EDTA buffer (Thermo Fisher, Landsmeer, the Netherlands). The cDNA was synthesized using M-MLV (Invitrogen, Breda, the Netherlands) and random nonamers (Sigma, Darmstadt, Germany). cDNA was quantified by relative standard curve method using quantitative real-time PCR [16]. In brief, cDNA was pooled and serially diluted (10x, 20x, 40x, 80x, 160x). The cDNA samples were diluted to a working concentration (20x) prior to relative quantification. Diluted cDNA samples and the serially diluted relative standard curve were pipetted into MicroAmp plates (Thermo Fisher, Landsmeer, the Netherlands). To each well we added either TaqMan fast advanced master mix (Thermo Fisher) or Sybr green master mix (Thermo Fisher). Primer and TaqMan probe (Eurogentec, Luik, Belgium) sequences are given in S2 Table. Plates were run in a Quantstudio 5 real-time PCR system (Thermo Fisher).

Histological analysis
Liver was formalin-fixed and paraffin-embedded and sectioned. Paraffin sections were stained with haematoxylin and eosin (H&E) for routine histological analysis, with Picosirius Red for detection of fibrillar collagen, and with Ki-67 for mitotic figures. H&E-stained specimens were scored blindly for steatosis, non-alcoholic steatosis (NAS) [17], ballooning [18] and findings were reviewed by a certified veterinary pathologist (AdB). Ki-67-positive hepatocyte nuclei were counted in 5 separate x40 fields by a single assessor. Fibrosis was assessed using Pico Sirius-red and quantified using ImageJ. Small liver sections were stained with Pico-Sirius red. Portal and parenchymal fibrosis was assessed using digital image analysis software (ImageJ). Briefly, two x20 digital photomicrographs, representing between 50 and 100% of total sample surface area (where possible avoiding large, longitudinal vascular structures), were obtained from comparable regions in every section using polarized light. The Images were analyzed using an in-house developed macro that quantified the area of positively stained collagen fibers (fibrosis) in every image.

Statistical analysis
Statistics were performed using SPSS 23 (SPSS Inc., USA) and R (R Core Team, 2018) [20]. Repeated measures were plotted as median and interquartile range. Single-time data are plotted as Tukey boxplots and scatter plots. No data were excluded. Data were not assumed to be normally distributed, thus tested non-parametrically. Groups were compared using the exact two-sided Mann-Whitney U test. A p<0.05 was considered to indicate rejection of the null hypothesis. A parameter was considered bimodal when the null hypothesis of Hartigan's dip test was rejected and the scatter plot showed 2 distinct clusters. Principal component analysis (PCA, correlation matrix) was restricted to 2 factors and varimax rotated. PCA variables were exported by method regression and plotted. Classical hierarchical cluster analysis was computed using the unweighted pair group method with arithmetic mean (UPGMA) on the Gower's similarity coefficient for mixed data [21].

Liver weight and VLDL secretion of C57BL/6JOlaHsd mice fed semisynthetic diet
We noted a bimodal distribution in absolute and relative liver weight distributions irrespective of whether the mice were pair-or solo-housed at postnatal day (PN) 70 (Hartigan's dip test, p<0.05, Fig 1A). Mice could be subdivided into two, nonoverlapping groups based on liver weight: 0.6 SD 0.1 g (small liver; SL, 2.3 SD 0.1%BW) and 1.0 SD 0.1 g (normal liver; NL, 3.4 SD 0.2%BW). Liver weights at PN70 were statistically significant between SL and NL both in wet weights as well as relative to BW (both p<0.001). Values for bodyweight, as well as lean and fat mass at PN70 were unimodal ( Fig 1B). Offspring labeled SL and NL appeared randomly distributed across (surrogate) dams and cages. To determine whether the observations on the bimodality of liver weights were consistent we repeated the experimental procedures in a separately bred and reared cohort (the "PN63" cohort). We noted a bimodal distribution in absolute and relative liver weight distributions at PN63 (p<0.05, Fig 1A). VLDL-TG secretion rate was, irrespective of whether the mice were pair-or solo-housed, lower in SL mice both absolutely (-66%, p<0.001) and after correcting for liver weight (-44%, p<0.001), indicating that the lower VLDL-TG secretion rate possibly is, at least in part, a direct consequence of the differences in liver weight in these mice ( Fig 1C). We noted that plasma was visually more yellow in SL versus NL mice at PN70 (S8 Fig). Based on liver weights, the SL phenotype occurred at similar rates in pair-and solo-housing (5/14, 36% and 3/8, 38%, respectively, the 'PN70' cohort, Fig 1A). The SL/NL phenotype also occurred in the independent PN63 cohort (8/15, 53%). In subsequent analyses we analyzed parameters from pair-housed mice only, unless specifically mentioned.

Liver morphology and histological analyses for steatosis, inflammation, proliferation, and fibrosis
Diffuse macroscopic hepatic pallor was exclusively seen in SL mice (Fig 2). Histological analyses revealed notably more severe centrilobular hepatocellular hypertrophy and more frequent karyocytomegaly in SL, as compared with NL mice (p<0.001, Table 1). Ballooning was not seen in either SL or NL livers. SL livers had higher steatosis grades (p<0.05), more necrotic cells and mitotic figures per microscopic field (p<0.001), and more prominent mixed inflammatory cell lobular inflammation (p<0.01) The presence of pigmented macrophages and moderate bile duct hyperplasia was a consistent feature of the SL phenotype (Table 1). Hepatocyte proliferation was higher in SL mice ( Fig 2C, Table 1, p<0.01). Moderate liver fibrosis was present in SL mice (Fig 2D, p<0.01). Fibrous depositions were seen around portal triads and central vasculature, and associated with bile duct proliferation. Bridging fibrosis was seen extending from portal to central regions along the sinusoids (particularly evident subcapsularly).

Hepatic triglyceride content, plasma liver enzymes and fasting plasma lipids in SL and NL mice
In accordance with the histological steatosis, we found a higher hepatic triglyceride content in SL livers compared with NL (+105%, Fig 3A). Hepatic free fatty acid content (29 SD 11 versus 16 SD 5 μmol/g, p<0.05) was higher in SL versus NL. Plasma liver enzymes ASAT and ALAT (aspartate/alanine aminotransferase) were higher in SL ( Fig 3B). Fasting plasma TG and cholesterol levels were lower. Plasma NEFA was subtly higher in SL (Fig 3C), whereas total plasma protein (52 SD 4 versus 50 SD 4 mg/ml) concentrations were similar. Hepatic expression of genes related to lipogenesis (Fasn, Srebp1f, and Dgat1/2) were similar. Expression of Pparg1, the master regulator of the adipogenic program to store fats in lipid droplets, was higher in SL ( Fig 3D).

Essential fatty acids and long-chain polyunsaturated fatty acids in SL and NL mice
Essential fatty acid deficiency (EFAD) may cause or contribute to liver disease [22]. We assessed EFA status/concentrations by determining the hepatic fatty acyl-chain (FA) profile, primarily defined by moieties of plasma membranes, cholesteryl esters, free fatty acids and triglycerides. The hepatic FA profile at PN70 showed large differences between SL and NL. All ω6 species, apart from 18:2ω6, were higher (p<0.05) in the livers of SL mice ( Fig 3E). The hepatic osbond:DHA ratio, a marker of poor DHA status, was higher in SL ( Fig 3F). These observations were also found in the separate groups of mice at PN63 (S1H & S1I Fig, respectively). With respect to ω3/6 moieties, we detected 18:3ω3, 18:2ω6 and trace amounts 20:2ω6 in the diet, (S3 Table). The ω3/6 levels were typical for a soybean oil-based diet with no other FA-bearing ingredients, such as the semisynthetic AIN-93G diet. Osbond acid (22:5ω6)

PLOS ONE
Spontaneous liver disease in wild type C57BL/6JOlaHsd mice concentration was higher in plasma, heart and skeletal muscle, resulting in a higher osbond: DHA ratio (S1A-S1F Fig). Analyzing the hepatic fatty acyl-chain profiles using principal component analysis (PCA), NL and SL mice from PN63 and PN70 appeared to form discrete clusters based on their NL/SL status, whereas their age (PN63 versus PN70) seemingly had no effect on these principal components (S1G Fig). To approximate when the SL phenotype develops, we repeated the experimental procedures in a separately bred and reared cohorts (the "PN21" and "PN42" cohort). Mice from the PN42 cohort also appeared to form two discrete clusters, near the SL and NL cluster, whereas the PN21 cohort clustered separately (S1G The molar plasma triene/tetraene ratio, a surrogate EFAD marker, was similar in SL and NL mice at PN70 and 63 (S1B & S1J Fig, respectively). Essential (dietary) fatty acids can be elongated, desaturated and used to synthesize (among others) ligands involved in inflammatory processes. Hepatic expression of elongases Elovl5/6 and of cyclooxygenase 1 (Ptgs1), were similar between SL and NL mice. Expression of desaturase enzymes Fads1/2 (Δ5/6 desaturase, D5/6D) and of lipoxygenase 15 (Alox15) was higher in SL mice (Fig 3G). Lipids can enter the liver as free fatty acids via transporters including Cd36, which was higher expressed in SL ( Fig  3H). The expression of markers of β-oxidation was ambiguous; expression of the transporter of long-chain FAs (Fabp1) and the master regulator of lipid metabolism (Ppara) were lower in SL, whereas that of a member of the carnitine shuttle (Cpt1a) and the master regulator of mitochondrial/peroxisomal biogenesis (Pgc1a, Fig 3H) was higher in SL.

PLOS ONE
Spontaneous liver disease in wild type C57BL/6JOlaHsd mice

Hepatic mRNA expression of BA synthesis and transport genes in SL and NL mice
Hepatic Cyp7a1 expression, encoding the rate limiting enzyme in BA synthesis, was higher in SL mice at PN70. Other bile acid synthesis enzymes Cyp27a1 and Cyp2c70 were lower in SL, whereas Cyp8b1 was similar. The bile acid receptor Fxr and its downstream target Shp were lower in SL, with similar levels of Lxr and Rxr (Fig 4C). Expression of hepatic BA transporters (Bsep, Ntcp and Mrp3) was lower in SL, suggesting lower transhepatic BA fluxes. Expression of Mdr2, the biliary phospholipid transporter, was similar between groups ( Fig 4D). Bile acid metabolism is closely linked to liver size control via Fgf15 and Hippo signaling [23]. Expression of Ctgf (YAP target gene [23]) was higher in SL ( Fig 4E).

Hepatic expression markers of cellular and mitochondrial stress
We determined gene expression markers of cellular and mitochondrial stress in hepatic tissue of SL and NL mice at PN70. Hepatic expression of Ddit3 (often called Chop), a cellular stress marker, was higher in SL. Markers of endoplasmic reticulum (ER) stress (Atf4, Atf6, Gadd34, Trib3) were similar (Fig 5A). In classic ER-stress, molecular chaperones are upregulated, which was not the case in SL mice, where those were downregulated ( Fig 5B). Instead, Atf5, the mitochondrial unfolded protein response (UPR) mediator, and Hmox1, involved in protection against oxidative stress, were higher in the liver of SL mice.

Course of SL
For investigating the timeframe wherein the SL phenotype develops, we assessed body weights, body composition and food intake from weaning until PN70. For this, we used the solo- housed PN70 cohort. Bodyweight ( Fig 6A) and lean and fat mass ( Fig 6B) suggest that SL mice grow slower immediately after weaning. Food intake was similar, but the food efficiency ratio (BW gain per gram food per week) was lower in SL mice in the first week following weaning (PN21-28, Fig 6C). Relative liver weight was unimodal at PN21 and 42 ( Fig 6D). Plasma bile acid levels varied widely in the PN21 cohort (2.7-74 μM), in the PN42 cohort (0.3-176 μM), and in the PN63 cohort (0.2-191 μM; Fig 6E). Analysis of relative liver weight, plasma lipids, bile acids and hepatic fatty acyl-chain profiles using principal component analysis (PCA), showed that NL and SL mice from PN63 and PN70 appeared to form discrete clusters based on their NL/SL status. The age (PN63 versus PN70) had no apparent effect on these principal components. The PN42 cohort, which does not display bimodality with respect to wet or relative liver weight, does display 2 clusters on the PCA scatter plot, closely positioned to the NL and SL clusters. Mice from PN21 cluster together on a scatter plot, separate from NL and SL clusters (Fig 6F). Hierarchical cluster analysis of relative liver weights, plasma parameters and hepatic fatty acyl-chain profiles showed that the PN42 cohort contained mice (4/14; 29%) that  were, similarly to PN70, lower in SL versus NL mice. Liver acylcarnitines (+47%, p<0.01), in particular C5DC and C6DC (both +200%. p<0.01) were higher in cluster 2 versus cluster 1 mice at PN42 (S1 Table). At PN21, all mice had similar plasma lipid levels (S5 Fig). At PN21, all mice had similar liver acylcarnitine profiles (S5 Table). The cohorts (PN42: 29%, PN63: 53%, PN70: 35%) suggest that the SL phenotype occurs in approximately one-third of C57BL/6JOlaHsd mice fed semisynthetic diet. In a control experiment, WT C57BL/6J mice, reared on chow, were fed a low-fat semisynthetic control diet (D12450J, Research Diets Inc. USA) from PN56 until PN140 (S6 Fig). Hierarchical cluster analysis of relative liver weights, plasma lipids and total plasma bile acids indicated that among the 9 mice tested, 2 mice appeared to cluster separately from the remaining 7 mice (S6A Fig). Histological analysis revealed that these 2 mice, compared to the 7 others, had more prevalent karyocytomegaly, mitosis, single cell death and scattered moderate mononuclear infiltration (S6B- S6E Fig).

Discussion
In this study, we show spontaneous divergence in liver weight in experimentally and nutritionally identically treated C57BL/6JOlaHsd mice fed a commonly used semisynthetic low-fat diet (AIN-93G). This diet is typically not considered to induce a model of liver disease. The divergence in liver weights resulted in mice with a small liver (SL) or normal liver NL) and coincided with profound metabolic differences in terms of lipid and bile acid metabolism. Our data demonstrate that the SL phenotype resembles, to some extent, the biochemical and histological changes observed in human with chronic liver disease, what ultimately may lead to chronic liver failure. Our data indicate that the first aspects of the phenotype become notable at or immediately after weaning. Heterogeneity and heterogeneous responses to stimuli in inbred mice have been reported in the context of high-fat diet feeding [27,28], in neurobiology [29,30], and in the apoE � 3-Leiden.CETP model [24,31]. Clear parallels exist between our study and those performed in apoE � 3-L.CETP mice [24,31], which are empirically subdivided (at PN42) into phenotypical responders (R) and non-responders (NR) based on plasma lipids [24]. Upon reassessing this literature, however, it became apparent that, despite the suggestive name, R and NR mice already differ in plasma lipids on chow [24]. This illustrates that "non-responding" to a dietary challenge is preceded by an existing spontaneous phenotype [24]. NR apoE � 3-L.CETP show, similar to our SL, lower liver weight, higher liver TG, an inflammatory liver pathology, higher plasma liver enzymes and higher plasma BA, compared to R or NL mice [24]. It may thus well be that the spontaneous bimodal liver size distribution observed in our study and the (non-)responders in apoE � 3L.CETP mice have a similar etiology and, thus, that the phenotype is not facility-dependent nor strain-dependent.
Histologically and based on gene expression patterns, the livers from SL mice showed an inflammatory phenotype and mild steatosis. Marked karyocytomegaly in the livers of SL mice likely indicates polyploidy and higher rate of cell division. The karyocytomegaly was linked to hepatocellular swelling that may reflect glycogen storage or acute hydropic degeneration [32]. Despite the suggested higher hepatocyte division rate, livers of SL are smaller. Thus, higher rates of hepatocyte cell division may be counteracted by, or even compensatory for, higher hepatocyte cell death (degeneration). The presence of pigmented macrophages in the livers of SL mice, together with elevated hepatic F4/80, Tnf-α, Mcp-1 and Col1a1 expression, are suggestive of an immune response which activates fibrogenesis [33]. Mild steatosis was observed in the livers of SL mice, together with lower VLDL-TG secretion rate and higher plasma NEFA and hepatic gene expression of Cd36. Hepatic expression and protein levels of Cd36 are higher in, and thought to contribute to, non-alcoholic steatosis [34]. The higher hepatic TG levels, lower VLDL-TG secretion rates, and higher Cd36 expression levels in SL mice are suggestive of a (net) flux towards hepatic TG stores.
Supplementing DHA in an acute model of liver damage (CCL 4 ) lowers the fibrotic and inflammatory response [22]. The enzymes necessary for the conversion of 18:3ω3 to DHA (ω3) also convert ω6 species and thereby synthesize osbond acid. Higher osbond acid concentrations, and high osbond:DHA ratios have been observed in human NAFLD liver biopsies [35,36]. A high osbond:DHA ratio is considered a marker of functional DHA deficiency [35][36][37]. We noted higher concentrations of osbond acyl moieties in liver, plasma, heart and skeletal muscle of SL mice. High osbond:DHA ratios were observed in the aforementioned tissues of SL mice. As the used diet did not contain osbond nor DHA moieties, the mice rely on its endogenous synthesis from dietary linoleic (18:2ω6) and α-linolenic acid (18:3ω3), respectively. As ω6 and ω3 metabolism shares enzymes, we interpret the higher concentration of osbond acyl moieties as an endogenous attempt to synthesize (among others) DHA.
Considering that DHA can dampen fibrotic and inflammatory responses [22], we speculate that its endogenous synthesis is higher in SL mice, possibly for the purpose of generating antiinflammatory ligands, such as resolvins [38]. Therefore, the data seem to suggest that the synthesis of long-chain polyunsaturated fatty acids (LC-PUFAs) is activated secondary to the inflammatory process.
The plasma BA levels seen in SL mice were (much) higher than what can be expected in non-diseased pre-and postprandial human plasma [39]. Minor experimental differences, such as variations in fasting times (0, 4 or 9 hours fasting) cannot be held responsible for the profoundly elevated plasma BA levels [39]. Plasma lipids, in particular the triglyceride levels, may be mildly elevated in the postprandial state [40]. When comparing unfasted PN21 to (for instance) 4 h fasted PN42, plasma triglyceride levels tend to show higher levels and a higher variability in the unfasted mice (S5 Fig). Other parameters, such as wet liver weight or the relative composition of essential fatty acids (dependent on slow multi-step enzymatic elongations and desaturations [41]) are not expected to differ much between unfasted and 9 h fasted mice.
Liver fibrosis and bile duct proliferation, observed in SL, are elements observed during the onset of certain types of cholestasis [25,42]. Cholestasis is the reduction or stagnation of bile secretion and bile flow. The accumulation of bile acids (BA) in the plasma, seen in SL mice, is typical for cholestasis [43]. We did not measure bile flow, thus cannot conclusively state that cholestasis is part of the SL phenotype. Alternatively, elevated plasma BA levels can be a consequence of a deficiency in basolateral BA transporters, such as an NTCP deficiency [44,45]. We did not assess transhepatic BA fluxes, thus cannot conclusively say whether these fluxes are lower in SL mice. Expression levels of Bsep and Ntcp were lower in the livers of SL mice, but these transporters are known to have an excess transport capacity [45]. Plasma BA composition in SL mice indicated a high abundance of primary BA species. Primary BA species, such as cholate and chenodeoxycholate, are synthesized by the liver and excreted into the bile. In the gut, primary BA are partially converted by microbiota, into secondary BA species, such as deoxycholate and lithocholate. The main enteral BA uptake transporter, ASBT, does not have a strong substrate specificity for conjugated primary or secondary BA [46]. Thus, the lower relative secondary BA abundance may be indicative of either a lower microbial biotransformation capacity [47,48], or less exposure of the BA pool to the microbiota. The latter could indicate that less BA are excreted into the bile, or that the time BA spend inside the lumen (exposed to the microbiota) is shorter in SL mice. Under physiological conditions, the terminal ileum reabsorbs BA and thereby releases FGF15 protein into the portal circulation which exerts a negative feedback signal to the liver for BA synthesis [49,50]. BA and FGF15 downregulate de novo BA synthesis via Fxr and Shp and via FGFR4 [49,50], respectively. The lower Shp and higher Cyp7a1 expression levels suggest that FXR signaling is lower in SL mice [50]. The (much) higher relative and absolute tauro-beta-muricholic acid (Tβ-MCA) plasma levels in SL mice (Fig 4) may inhibit FXR signaling [50]. Tβ-MCA is a powerful (gut microbiota-sensitive) FXR antagonist in mice [50]. Considering the high plasma BA levels, the low abundance of secondary BA species, the high hepatic expression of BA synthesis genes, and the low expression of hepatic BA transporters, it is tempting to speculate that the transhepatic (and enterohepatic) BA flux is lower in SL mice.
Liver size is, at least in part, regulated by FGF15-Hippo signaling along the gut-liver axis [23]. FGF15 suppresses YAP signaling by activating Hippo signaling [23]. When YAP signaling is not suppressed, it upregulates genes necessary for bile duct and hepatocyte proliferation [23]. Our data suggest YAP signaling is activated in the livers of SL mice (Fig 4E), suggesting that the low liver weight in SL mice is likely not caused by a suppression of YAP signaling. Instead, activated YAP signaling may be in accordance with the higher rates of hepatocyte cell division seen in the livers of SL mice.
High levels of acylcarnitines have been described in liver biopsies from NASH, but not from NAFLD patients, which has been linked to mitochondrial dysfunction [51]. We noted substantially higher acylcarnitine species in the livers of SL mice, as well as higher hepatic expression of Pgc1a (the master regulator of peroxisomal and mitochondrial biogenesis) and Cpt1a (part of the carnitine shuttle, rate-limiting factor in β-oxidation). It is tempting to speculate that the inflammatory process exerts a relatively high energy demand, for example generated via mitochondrial β-oxidation. The apparent accumulation of acylcarnitine species, however, may indicate that β-oxidation or the tricarboxylic acid (TCA) cycle is not operating optimally in the livers of SL mice [51]. It could be that β-oxidation outpaces the TCA cycle, upon which incompletely oxidized acyl-carnitine intermediates can accumulate [52]. Alternatively, higher levels of acetylcarnitine may reflect higher peroxisomal β-oxidation rate [53]. Upon impairment of the carnitine shuttle, fatty acids are directed to the peroxisomes for oxidation [53]. The higher hepatic acylcarnitines levels coincided with higher expression markers of the mitochondrial unfolded protein response (Atf5 and Ddit3, higher expressed in the livers of SL mice) [54], which suggest mitochondrial dysfunction [22,55]. A potential contributor to mitochondrial dysfunction may be the absence of Nicotinamide Nucleotide Transhydrogenase (Nnt), which generates the antioxidant compounds glutathione and thioredoxin in the mitochondria [19]. Nnt is defective in C57BL/6J [19], but it appeared intact in NL and SL C57BL/6JOlaHsd mice (S4 Fig), and is therefore not likely to contribute to the observed SL phenotype.
We describe a spontaneous, but likely pathological liver phenotype in C57BL/6JOlaHsd mice, bred and reared on a semisynthetic control diet (AIN-93G). Rodent diets are typically available in 'growth' (AIN-93G) and 'maintenance' (AIN-93M) formulation [4]. The "G" variant was formulated to be suitable for growth, pregnancy and lactation [4]. The G and M formulations differ in composition, namely higher protein (20% versus 14%), fat (7% versus 4%) and minerals (mainly calcium carbonate and ferric citrate) in the former, at the expense of corn starch and maltodextrin [4]. Of note, the AIN-93 macro and micro-nutrient composition formulations [4] serve as reference formulations for other semisynthetic diets including ubiquitously used (semisynthetic) high-fat diets and their low-fat controls. We cannot conclude that the choice of the diet is an important factor for the development of the SL phenotype. We observed the SL and NL phenotype in independent cohorts fed AIN-93G from 2 different vendors. Experimental and nutritional conditions were identical for the mice within each cohort. In our view, this makes an environmental or nutritional origin of the SL phenotype unlikely. Though, we cannot exclude that our experimental conditions amplified an otherwise unremarkable genetic heterogeneity. Mice fed AIN-93G compared to a non-purified (chow) control diet, for unknown reasons, have a lower liver weight [56]. Breeders were obtained from a commercial specific-pathogen-free stock, and tested negative for major bacterial, viral and parasitic etiologic agents. We consider it unlikely that the SL phenotype is caused by an infection; we did observe the two phenotypes also within individual cages, housing different mice. The SL and NL phenotype appeared randomly distributed across (surrogate) dams. Due to randomization at PN2, we do not know whether SL/NL is linked to intrauterine conditions. Due to observation of the two phenotypes both in pair-housed and in solo-housed mice, differences in housing conditions or in social hierarchy do not seem to be a contributing factor. Through phenotyping, we were able to identify (an early form of) the SL phenotype at PN42. Such a distinction could not yet be made at weaning (PN21). However, the absence of a phenotype (at PN21) does not fully rule out that the cause for the phenotype is not present. The immediate post-weaning difference in food efficiency ratio and lean mass growth suggests that the cause of the SL phenotype is actually present at weaning. We speculate that the SL phenotype is caused by (subtle) early-life events or variations in C57BL/6JOlaHsd's (epi)genome. As the occurrence rates varied between cohorts; i.e. PN42 (~29%), PN63 (~53%) and PN70 (~36%), it is tempting to speculate that (an) uncontrolled variable(s) potentiate(s) the cause.
The biochemical differences between SL and NL mice resemble, to some extent, (the development of) liver cirrhosis or an early stage of liver failure in humans. Our histological analyses in mice, however, did not (yet) provide clear evidence of liver cirrhosis at PN63 (Fig 2). Future extended longitudinal studies in mice could demonstrate whether or not the SL phenotype evolves to liver cirrhosis and ultimately to end-stage liver failure and early mortality. As liver function appeared impaired, at least in regard to VLDL secretion (Fig 1), we speculate that SL mice may cope less effectively with (models of human) diseases and with the conditions that are associated with ageing. It was beyond the scope of the present study to completely assess whether the same phenomenon also occurs in other mouse (sub-) strains and/or under other nutritional conditions. We observed the SL phenotype in C57BL/6JOlaHsd mice fed semisynthetic AIN-93G. A similar phenotype appears to occur in C57BL/6J mice fed semisynthetic low-fat control diet (S6 Fig). A similar phenomenon appears to have been described in apoE � 3-L.CETP mice (C57BL/6J background strain) fed chow [24]. In a large phenotyping study comprising 44 inbred mouse strains (Paigen1, The Jackson Laboratory), challenged to a high fat (semisynthetic) lithogenic diet for 8 weeks, plasma bile acids showed large variability between strains: 0.9 to 350 μM but also within strains [57]. In Paigen1, male C57BL/6J liver weight distribution was bimodal (p<0.05), whereas other strains did not show this either due to low statistical power or absence of bimodality per se [57]. The cohorts PN63 and PN70 were highly similar in the assessed parameters (Fig 6F, S2 Fig), although these cohorts originated from separate breeders, and had been studied during different years (2017 and 2018, respectively). This suggests that, by keeping experimental conditions similar, the occurrence of the SL phenotype and its metabolic consequences have been highly reproducible and that it is unlikely that incidental environmental (stress) factors are a major contributing factor to our conclusion. These data suggest that the SL phenotype is likely not specific to our commercial animal supplier, to the 'JOlaHsd' C57BL/6 sub-strain, or to incidental environmental (stress) factors at our animal facility.
Theoretically, the cause of the SL phenotype could still be due to specific environmental conditions present at our facility, which could then affect each new imported cohort. Singlecenter experiments, like the ones described here, are more vulnerable to biases and methodological pitfalls compared to multi-center experiments [58]. The ApoE � 3L.CETP mice results obtained in another facility [24] suggest, however, but do not prove, that the current observation are not unique for our facility. In addition, why such an environmental factor would only affect some mice of a cohort, but not others, would then still need to be resolved. It is therefore not (yet) possible to generalize our observations and characterisations to all C57BL/6JOlaHsd mice or indeed to the ubiquitously used C57BL/6J sub-strain.
We identified markers for the identification of the SL phenotype that could potentially be used at postnatal day 42 and later ages. Plasma triglycerides, total cholesterol, bile acid levels or composition, hepatic FAME osbond:DHA ratio, C5DC or C6DC acylcarnitine species, or Col1a1 gene expression (or other fibrogenesis markers) could serve as potential markers for the (early identification) of the SL phenotype. We were not able to pinpoint any one cause for the SL phenotype and this would certainly have further strengthened our study. Our data ( Fig  6) suggest that the cause of the SL phenotype is present at weaning (PN21). We were not able to rule out a genetic origin. If the SL phenotype is, in fact, caused by a genetic variation or defect, then the identification of such a locus by means of genotyping would aid in the early (even pre-experiment) exclusion of those mice from experiments. Societal discussions surrounding the usage of animals in medical research have emphasized the importance of the '3Rs': replacement, reduction, and refinement. The presence of (congenital) defects, which cause high levels of phenotypic heterogeneity, in otherwise highly homogenous mouse strains is potentially detrimental to the quality of preclinical studies [30,59]. Reducing phenotypic heterogeneity within cohorts of C57BL/6 mice, or at least the identification of distinct subgroups within cohorts [24,59], would contribute to the overall quality of animal research and to the more efficient usage of mice, time and resources.
In summary, our data show that approximately one-third of C57BL/6JOlaHsd mice fed semisynthetic diet develop spontaneous liver disease. This correlates with low liver weight, low VLDL secretion, high plasma bile acids, liver steatosis, fibrosis and inflammation, and mitochondrial dysfunction. We advise researchers who use this strain (and possibly other C57BL/6 sub-strains, including C57BL/6J) to be aware of this spontaneous phenotype that may interfere with the interpretation of results. Our data suggest that lower plasma triglyceride and cholesterol concentrations (among others) may function as surrogate parameters for early identification of mice with the apparent SL phenotype. Heatmap containing the normalized relative liver weight, plasma parameters (lipids, total protein, liver enzymes, bile acids) and hepatic fatty acyl-chain profiles. Each column contains data from one mouse. Each row represents a discrete biometric or biochemical parameter. Each square in the heatmap represents a normalized value (parameter value divided by the average of that parameter in all animals). Missing values are shown as an X in grey. Red color indicates the values is higher than the average, whereas blue color indicates the value is lower than the average of all mice for that particular parameter. Hierarchical clusters were computed using the unweighted pair group method with arithmetic mean (UPGMA) on the Gower's similarity coefficient for mixed data. Data represent the pair-housed PN21 (n = 12), PN42 (n = 14), PN63 (n = 15) and PN70 (n = 14) cohorts. Cophenetic correlation coefficient = 0.842. Mice dissected at PN42 were split into clusters 1 (NL) and 2 (SL) based on principal component analysis (Fig 6F). Data represent the pair-housed PN21 (n = 12 total), PN42 (n = 10 cluster 1, n = 4 cluster 2), and PN63 (n = 7 NL, n = 7 SL) cohorts. Exact two-sided Mann-Whitney U test � : p<0.05, �� p<0.01, ��� p<0.001. (TIF) S6 Fig. Classical hierarchical cluster analysis on relative liver weight and assessed plasma lipid parameters of a pair-housed WT C57BL/6J cohort. WT C57BL/6J mice, reared on chow, were fed a low-fat semisynthetic control diet (D12450J) from PN56 until PN140. The mice were sacrificed at PN140. Heatmap containing the normalized relative liver weight and plasma lipids and total bile acids. Each column contains data from one mouse. Each row represents a discrete biometric or biochemical parameter. Each square in the heatmap represents a normalized value (parameter value divided by the average of that parameter in all animals). Red color indicates the values is higher than the average, whereas blue color indicates the value is lower than the average of all mice for that particular parameter. Hierarchical clusters were computed using the unweighted pair group method with arithmetic mean (UPGMA) on the Gower's similarity coefficient for mixed data (A). Data represents the pair-housed WT C57BL/ 6J cohort. Cophenetic correlation coefficient = 0.897. Histological staining using hematoxylin and eosin ( C57BL/6JOlaHsd mice were bred inhouse and nests were culled to 4 males and 2 female pups at PN2. Male pups were weaned at PN21 and either pair-housed (entire PN42 cohort, entire PN63 cohort, pair-housed PN70 cohort) or solo-housed (solo-housed PN70 cohort only). Mice were sacrificed at PN21, PN42, PN63 or PN70. In a subset of the PN70 cohort, upon intraperitoneal injection of the lipoprotein lipase inhibitor poloxamer-407, retro-orbital blood was drawn at 0, 1, 2, 3 & 5 h for determining the Very-low density lipoprotein-triglyceride (VLDL-TG) secretion rate. PN21: n = 12, PN42: n = 14, PN63: n = 15, PN70, pair-housed: n = 14, PN70, solo-housed: n = 8, PN70, pair-housed, VLDL-TG experiment: n = 14, PN70, pair-housed, VLDL-TG experiment: n = 8. (TIF) S8 Fig. Plasma samples at PN70, made from cardiac puncture blood samples. Each tube contained exactly 400 μl of plasma and was not diluted. The left sample (above the '5') was obtained from a mouse which was later labeled "SL", whereas the right sample (above the '6') was obtained from a mouse which was later labeled "NL". (TIF) S1 Raw image. (TIF) S1 Table. Liver acylcarnitine species in SL and NL mice at PN42. Mice dissected at PN42 were split into clusters 1 (NL) and 2 (SL) based on principal component analysis ( Fig 6F)