Thyroid hormone receptor beta (THRB) dependent regulation of diurnal hepatic lipid metabolism in adult male mice

Introduction
Thyroid hormones (THs)—thyroxine (T4) and its main biologically active derivative, 3,3’,5-triiodothyronine (T3)—regulate systemic energy homeostasis through transcriptional programs in metabolic tissues. TH action in target tissues mainly stems from T3 binding to nuclear hormone receptors alpha and beta (THRA and THRB—the latter being the predominantly expressed form in liver). Upon T3 binding, transcriptional co-factors are recruited and target gene transcription is initiated1,2,3,4. In the liver, THs stimulate lipogenesis and triglyceride storage, lipid breakdown by boosting β-oxidation, and ATP generation by targeting β-oxidation-derived acetyl-CoA to the TCA cycle. THs further modulate cholesterol metabolism by augmenting its biosynthesis, serum uptake, and conversion into bile acids. Overall, an elevated liver TH state results in reduced fatty acid and cholesterol levels, thus making it an attractive pharmacological target for the treatment of metabolic dysfunction associated steatotic liver disease (MASLD), previously known as non-alcoholic fatty liver disease (NAFLD)5,6.
Temporal regulation of physiological functions is vital for an organism’s adaptation to cyclic environmental changes. Animals have developed a so-called circadian timing system capable of detecting and predicting external daytime and adjusting physiological processes accordingly. At the molecular level, the circadian system is based on a set of interlocked transcriptional and translational feedback loops that drive rhythmic expression of tissue-specific clock target genes7,8.
Many endocrine signals, such as cortisol or melatonin, exhibit marked circadian secretion patterns9,10. However, to which extent THs are subject to circadian control is still debated. Previous studies show that blood T3 and T4 levels show no or only weak circadian variations. At the same time, however, substantial circadian alterations in systemic energy state, and in hepatic transcriptome regulation—including TH target genes involved in fatty acid and cholesterol metabolism—have been observed. Interestingly, changes in TH state have only minor effects on the regulation of the hepatic diurnal clock machinery11.
In this study, we investigated the role of THRB as the primary liver TH receptor isoform in regulating diurnal transcriptome and lipidome rhythms in this organ. Our findings indicate an interplay between the diurnal clock and hepatic TH action. Loss of THRB led to a combination of baseline and phasing effects in the transcriptional signature of metabolic processes indicative of a pro-steatotic state. In line with this, lack of THRB significantly impacted the diurnal lipidome of the liver with a net-increase in fatty acid and cholesterol levels. In vitro agonization of THRB in steatotic hepatocytes confirmed this role of THRB action as therapeutic target for MASLD.
Results
Loss of THRB has modest effects on liver transcriptome state
To evaluate the consequence of impaired THRB signaling on liver physiology, we used a knockout model in which the DNA binding activity of THRB is impaired12,13. Adult male control (THRB+/−, herein THRBWT) and Thrb knockout mice (THRBKO) were kept under a standard 24-hour light:dark cycle and received a regular chow diet and water ad libitum. Mice were culled at 4-hour intervals for 24 h, and liver samples were subjected to RNA sequencing. To achieve a time-of-day independent comparison of the transcriptome state between both genotypes, all sampling points were pooled (Fig. 1A, upper panel). A total of 99 differently expressed genes (DEGs; up n = 50; down n = 49) were identified in THRBKO livers, including the TH modulator Dio1 (down) and the lipoprotein receptor Lrp2 (up) (Fig. 1B, Supplementary Data 1). Compared against DEGs from previous studies in mice with systemically low (methimazole plus potassium perchlorate treatment; MMI14) or high (0.5 mg/mL supplementation of T311) TH activity, the extent of regulation in THRBKO mice was relatively modest, suggesting THRB-independent effects of TH action on liver transcription regardless of statistical thresholds (Fig. 1C). Importantly, the same conclusions were reached when considering the factor time when assessing DEGs in either THRBKO, MMI, or T3 datasets. In this refined analysis, a total of 306 DEGs were identified in THRBKO mice. Interestingly, approximately half of these transcripts were exclusively identified in THRBKO livers, while the other half was present in livers from MMI or T3 mice (Supplementary Fig. 1A, B; Supplementary Data 1).

A Experimental design and modes of DEG analysis. B Volcano plot of global differentially expressed genes (i.e., DEGs irrespective of sampling time; padj < 0.05). Red and blue represent significantly UP- and DOWN-regulated genes, respectively. C Comparison of global DEG numbers between livers from T3– and MMI-treated and THRB-deficient animals stratified for significance threshold. N = 3 independent samples per time point and condition for THRBWT and THRBKO mice. The mouse picture described in A was obtained from Biorender.
Loss of THRB alters diurnal phase of liver biological processes and transcript rhythms
To get more insight into the effects of THRB deletion on liver transcription, we stratified our data for sampling times and analyzed changes in diurnal rhythm parameters (Fig. 1A, lower panel). A combination of complementary algorithms (JTK cycle, Metacycle, DryR, and CircaN) was used to assess diurnal rhythmicity in transcript regulation, enhancing reliability and mitigating specific limitations inherent to each individual method15. Genes identified in at least one method were considered as rhythmically expressed. This approach identified 9746 and 9271 genes showing significant diurnal regulation in THRBWT and THRBKO liver samples, respectively (Fig. 2A; Supplementary Data 2). Notably, peak phase distributions of all exclusive rhythmic genes for each group of diurnal genes were consistent across genotypes, with the lowest activity in the first half of the dark phase (ZT12-18; Fig. 2B).

A Rhythmic genes after diurnal profiling (using JTK_cycle, Metacycle, DyrR, and CircaN) are shown in a Venn diagram for THRBWT (grey) and THRBKO livers (blue). Acrophase distribution of exclusive rhythmic transcripts in livers of THRBWT and THRBKO livers (B) and phase difference for robustly rhythmic genes between both genotypes (C). D THRB-dependent acrophase shifts of rhythmic transcripts of different physiological processes are depicted. ****p < 0.0001. N = 3 independent samples per time point and condition.
For those genes with significant rhythmicity in both genotypes, we observed a slight but significant phase delay of 0.4 hours in THRBKO compared to THRBWT livers (Fig. 2B, C). Using the phase set enrichment analysis (PSEA) algorithm16, we further assessed the impact of THRB on the timing of specific biological processes. Based on the overall median of the genes enriched in a biological process, an overall activity phase for a biological process was estimated. Similar processes were manually combined and averaged into broader biological categories to simplify the data. This approach yielded a general delay of up to 3.5 h in expression rhythms associated with metabolism while immune system-related processes showed a marked advance of around 6 h in THRBKO mice (Fig. 2D; Supplementary Data 2).
Loss of THRB alters diurnal rhythms of fatty acid and cholesterol metabolism
For a more detailed analysis of diurnal deregulation in THRB-deficient livers, we directly compared transcript rhythm parameters between both genotypes. CircaCompare17 analysis identified 1446, 314, and 369 genes with changes in mesor (baseline expression), amplitude (time-of-day dependent variation), and acrophase (time of peak expression), respectively, termed as differentially rhythmic genes (DRGs, Fig. 3A; Supplementary Fig. 1; Supplementary Data 3). Changes in mesor or amplitude were determined for genes rhythmic in at least one, while phase difference was determined only for genes rhythmic in both conditions11. Loss of THRB affected liver transcriptome rhythms predominantly at the mesor level (764 mesor up and 682 mesor down) with slightly lower numbers for phase (98 advanced and 271 delayed genes) and amplitude (207 and 107 genes with lower and higher amplitude, respectively) (Fig. 3A). Enrichment analysis of DRG classes for biological processes revealed a complex regulation at the transcriptome level in the absence of THRB. For instance, DRGs with reduced and increased mesor enriched for mRNA catabolic processes and phosphorylation, respectively. Conversely, DRGs with reduced amplitudes enriched for AKT pathways (PI3K pathway and AKT signaling), protein ubiquitination, and response to glucose. DRGs with increased amplitude, on the other hand, were enriched for transcriptional regulation, proteolysis, and protein kinase activity (Fig. 3B; Supplementary Data 3). Most strikingly, diurnal profiling allowed us to identify a dual regulation at the mesor level for transcripts associated with fatty acid and cholesterol regulation. Enrichment analysis of DRGs with increased mesor effect enriched for cholesterol and bile acid metabolism (Acnat2, Abcb11, Akr1d1, Baat, Cyp7a1, Npc1, Slc27a2, Hmgcs2, Vldlr) as well as fatty acid metabolism and beta-oxidation (Abcd2, Acox2, Acot3, Acot12, Cpt1a, Ehhadh, Slc27a2). Conversely, mesor down DRGs enriched for acyl-coA biosynthesis, hydrolysis, and fatty acid transport (Acsm1, Acsm2, Acaa1a, Abhd6, Acox3, Fabp1, Fabp2, Slc27a3, Slc27a4) (Fig. 3B). These analyses predicted that DRGs with increased mesor favor higher fatty acid synthesis and oxidation (Acox2, Acnat2, Acad11, Acot3, Acot12, Cpt1a, Ehhadh, Scd2, Slc27a2, Hsd17b4, Hpgd) whereas DRGs with reduced mesor point to a decreased fatty acid breakdown and transport (Acsm1, Acsm2, Acot6, Acox3, Abhd6, Fabp1, Fabp2, Mboat1, Slc27a3, Slc27a4). Genes pertaining to cholesterol synthesis and bile acid secretion (Baat, Akr1d1, Cyp7a1, Hmgcs2, Npc1, Vldlr) were identified only in DRGs with increased mesor, which suggests increased cholesterol metabolization and secretion as bile acids (Fig. 3C, D). Diurnal profiling of pathway associated genes confirmed the mesor effects with higher effects towards the end of the dark phase (ZT18-24) (Fig. 3E; Supplementary Data 3). In summary, detailed diurnal profiling predicted that THRBKO livers would tend towards increased fatty acid and cholesterol levels.

A Diurnal profiling of liver transcripts (using CircaCompare) was performed, and numbers of differentially rhythmic genes (DRGs) are represented as UpSet plot. B Enrichment analysis of the top-5 processes for each DRG class (amplitude, phase, mesor) are shown. C, D Manual curation of lipid metabolism-associated DRGs for mesor was performed and shown as heat maps. Chord plots show gene names and associated biological processes. E Genes associated with the biological processes were averaged and plotted across 24 h and diurnal analysis was performed using CircaCompare. BAF bile acid formation, CHO cholesterol metabolism, FAO fatty acid oxidation, FAS fatty acid synthesis, LPC lysophosphatidylcholine, PI3K phosphoinositide 3-kinase, PI and PS phosphatidyl and phosphatidylserine metabolism. N = 3 independent samples per time point and condition. R = presence of rhythmicity. NR = absence of rhythmicity.
Loss of THRB alters liver lipidome rhythms
To test this prediction, we performed diurnal untargeted lipidomic profiling on liver samples from both genotypes. In total, 285 distinct lipid species were detected of which 101 and 107 were significantly rhythmic in THRBWT and THRBKO, respectively (Fig. 4A). Assessment of lipid classes showed main genotype differences for ceramides (more rhythmic species in THRBKO), phosphatidylethanolamines (more), and diacyl- and triacylglycerides (less rhythmic species in THRBKO) (Fig. 4B). No overall phase difference between robustly rhythmic lipids was identified between the genotypes (mean = 0.15 h, p = 0.60, Fig. 4C; Supplementary Data 4).

A Rhythmic lipids after diurnal profiling (using JTK_cycle, Metacycle, and CircaN) are shown in a Venn diagram for THRBWT (grey) and THRBKO livers (blue). B Distribution of lipid classes identified as rhythmic in THRBWT (left) or THRBKO livers. C Rose plot depicting acrophase differences of robustly rhythmic lipids between genotypes. D UpSet plot representing numbers of lipids with specific rhythm parameter changes (mesor, amplitude, phase) identified by CircaCompare. E Examples of lipids with alterations in specific rhythm parameters. N = 3 independent samples per time point and condition.
Detailed comparison of lipid rhythms revealed marked changes mostly at the mesor level (Fig. 4D; Supplementary Data 4). Lipid classes such as phosphatidylcholine (PC 42:2, 31:1, 18:0_22:5, 20:0_20:4, 22:0_20:4), plasmalogen PC (38:2), phosphatidylinositol (PI, 18:0_22:6, 20:0_20:4) showed mesor down effects in THRBKO. Increased mesor, in turn, was observed for cholesterol esters and sulfates (CE 18:1, ST 28:1;O;S, ST 29:1;O;S), triacylglycerols (TG 58:9, TG 58:8, TG 53:4, TG 58:7, TG 52:5, TG 55:5, TG 51:4), diacylglycerols (DG 18:1_22:6, DG 18.2_20.2, DG 18.2_20.4, DG 18.1_18.2, DG 18.0_18.2, DG 16.0_18.2), phosphatidylglycerols (PG 18:2_22:6, PG 18:1_22:6, PG 18:1_18:1), free fatty acids (18:2), ceramides (Cer 42:3;2O; Cer 18:1;2O/23:0), and phosphatidylcholines (PC 18:1_20:3, PC 18:1_18:2). Increased rhythm amplitudes were identified for phosphatidylglycerol (PG 16:0_18:2) while a phase advance of ca. 4 h (associated with a mesor increase) was found for two diacylglycerol entities, DG 18:0_18:2 and DG 18:2_20:4, with similar tendencies (p-value < 0.1) for TG 60.4 and DG 36.4 (Fig. 4E; Supplementary Fig. 2 and 3; Supplementary Data 4).
In line with transcriptome findings, diurnal lipidomic profiling suggested slightly increased lipid contents (TAG, DAG, free fatty acids, cholesterol) in THRBKO mice. It further identified additional lipid classes, such as PCs, PIs, PGs, and ceramides whose mesor levels were also affected by the absence of THRB. Our findings suggest that loss of THRB results in changes in metabolic rhythms yielding a higher susceptibility to steatogenic stimuli.
Activation of THRB reduces steatosis in hepatocytes in vitro
To directly assess the role of hepatic THRB action in the regulation of steatosis, we treated AML12 hepatocytes with palmitate with or without pharmacological activation of THRB. Palmitate treatment for 48 hours resulted in consistent levels of steatosis (Fig. 5A, B). Using a highly sensitive flow cytometry approach, we concurrently measured cell death, triglyceride, and cholesterol content at the individual cell level. Palmitate treatment increased cell mortality, which was significantly reversed by co-stimulation with T3 or the THRB agonist resmetirom. Remarkably, resmetirom reduced cell death even below baseline levels of the BSA-treated group, underscoring a potent hepatoprotective effect mediated by THRB activation (Fig. 5B, C, left panels). As anticipated, palmitate increased, whereas T3 and resmetirom rescued triglyceride and cholesterol levels when compared to BSA controls (Fig. 5B, C, center and right panels). The reduction by T3 and resmetirom fully reverted lipid content whereas cholesterol decrease was only partial (Fig. 5C). Collectively, our findings suggest a protective role of THRB in lipid metabolism under steatotic conditions.

A Representative micrographs of AML-12 hepatocyte triglyceride staining (AdipoRed) under control (left) and with palmitate (right). Scale bars depict 100 µm. B Gating strategy representing viability, triglyceride and cholesterol signal. C Percentage of dead cells and mediant intensity fluorescence (MIF) of lipid -and cholesterol-positive cells. N = 5–12 independent samples per condition.
THRB-independent TH-regulated genes are associated with immune function
Our initial findings suggested that loss of THRB has only moderate effects on overall transcriptome regulation compared to conditions with altered systemic TH state (Fig. 1). THRB is recognized as the primary thyroid hormone receptor in the liver, which is largely attributed to its high prevalence in hepatocytes of the liver parenchyma18. While THRA is also expressed in the liver, it has primarily been associated with other cell types19. Single-cell RNA sequencing data from the Liver Cell Atlas20 supported these findings, indicating that Thrb expression is mainly found in hepatocytes and Kupffer cells, while Thra is primarily expressed in Kupffer, T, and endothelial cells (Supplementary Fig. 4).
Our data suggested that many genes exhibit increased expression under high-TH conditions while being downregulated under low-TH conditions11,14. Interestingly, 277 of these putative TH target genes did not show significant alterations in the absence of THRB (Fig. 6A). Enrichment analysis for these TH-regulated but THRB-insensitive genes highlighted several biological processes, such as immune system, inflammation, and lipid and cholesterol metabolism (Fig. 6B, Supplementary Data 4). Targeted pathway analysis confirmed this notion. Interestingly, at the pathway level immune system-associated transcripts showed phase changes in line with previous PSEA data from THRBKO mice (Figs. 6C, 2D). These changes suggest an overall THRB-independent regulation of associated pathways by TH—possibly via THRA—which would be consistent with the predominant expression of THRA in non-hepatocyte cells.

A UpSet plot comparing liver global DEGs down-regulated under low-TH (MMI), up-regulated under high-TH (T3), and unchanged in THRB. B Enrichment analysis of putative THRB-independent TH targets (n = 277). C Overall pathway expression profiles of the THRB-independent TH-responsive genes in the immune system, lipid and cholesterol metabolism in different groups. Data from low (MMI) and high (T3) levels were obtained from GSE199998. N = 3–4 independent samples per time point and condition.
Discussion
We show that loss of THRB results in subtle alterations in liver transcriptome regulation. Detailed diurnal profile analysis unveiled alterations in the expression profiles of metabolism-associated and other genes. Particularly, in the absence of THRB, immune-related genes were predominantly phase-advanced. In contrast, genes involved in metabolic functions, such as lipid metabolism, were phase-delayed. Numerous genes associated with lipid metabolism pathways were significantly altered with regard to rhythm parameters, especially mesor. Diurnal lipidomic profiling supported these transcriptomic insights and underscored an elevation in TAGs, DAGs, and cholesterol. Cell-based assays confirm that THRB activation mitigates lipid and cholesterol build-up, underscoring its protective role against steatosis. Leveraging existing datasets, we discerned a subset of genes regulated by TH independent of THRB that are predominantly implicated in immunological and metabolic processes.
The hepatic transcriptomic response to low TH levels is much smaller compared to high-TH conditions14, and this reduction is even more pronounced in THRBKO mice. A previous microarray study identified that 201 out of 1900 genes as responsive to T3 treatment, of which 60% show a THRB-dependent effect21. The pronounced hepatic response to high TH versus the diminished response to low TH state—whether induced pharmacologically or genetically (THRB knockout)—further substantiates the notion of low hepatic TH action under physiological conditions.
We uncovered more fine-grained THRB deletion effects when considering sampling time in our analysis. Diurnal profiling yielded distinct changes in transcript regulation under THRB-deficient conditions. Such an approach has previously been shown to enhance the detection of transcriptomic alterations in response to a high-fat diet15. Nevertheless, the effects of THRB deletion were overall less pronounced compared to those of altered TH conditions. Focusing on lipid metabolism, we observed bidirectional mesor level alterations: genes driving fatty acid oxidation and cholesterol metabolism were upregulated, while transcripts involved in fatty acid transport/metabolism were downregulated.
Diurnal lipidomics validated these findings, revealing accompanying changes following THRB deletion at the metabolite level. Overall, we identified 285 lipid species, with 35% exhibiting significant rhythmicity across both genotypes. This exceeds previous findings by Adamovich et al.22 which reported 159 lipid entities with 17% rhythmicity. Variations in experimental protocols, such as sample processing and lipid identification, along with the methodologies for assessing rhythmicity, likely account for these differences. In a recent study, Sprenger et al.23 showed that half of the liver lipidome is rhythmic, employing ANOVA for temporal analysis due to their high-resolution data and irregular sampling intervals. Direct comparisons with our results are challenging because of methodological differences. However, a consistent rhythmicity observed in TAGs across studies is notable.
33 lipids (11%) showed significant mesor variations between genotypes, and a considerable increase in cholesterol, TAGs, DAGs, and FAs in THRBKO livers was identified. We further observed significant alterations in PCs, PGs, and PIs, with most changes presenting as increased mesor levels. Notably, a specific PG species (PG 16:0_18:2) exhibited a gain of rhythmicity in THRBKO mice, while some PC species showed mesor alterations (i.e., increase or decrease), thus representing a mixed effect. Similarly, specific PG and PI species were distinctly impacted in THRBKO livers, displaying elevated or diminished mesor levels, respectively. Together, our findings suggest that the absence of THRB has a more pronounced impact on the liver lipidome than the transcriptome. This metabolic reprogramming lends further support to THRB’s protective role against steatosis. Lower PC/PE ratios and PC levels have been associated with metabolic syndrome-associated steatohepatitis development24,25,26,27, and insulin resistance has been implicated in MASLD development28 through changes in DAGs and ceramides29,30,31. THRB KO livers show increased ceramide and TAG levels, which could predispose to MASLD development.
The ablation of THRB resulted in elevated hepatic lipid and cholesterol levels, underscoring its significant role in lipid regulation5,32. Pharmacologically activating THRB has consistently been shown to decrease lipid and cholesterol levels18,33,34,35. Contrastingly, receptor-level manipulations yielded variable results. For example, mice with a dominant-negative THRB mutation (TRBPV) display a THRB protective role32, while a recent THRB knockout study did not show a preventive effect on metabolic dysfunction-associated steatohepatitis (MASH) development. In contrast, THRBKO mice were further protected against fibrosis compared to wild-type controls36. These discrepancies could be due to differences in experimental approaches, such as the type of knockout and temperature conditions during experiments. Importantly, thermoneutrality-kept THRBKO mice show normalized T3 and T4 levels, whereas mice kept under non-thermoneutral conditions are hyperthyroid36. It has been suggested that the availability of the thyroid hormone ligand is more critical than receptor availability in MASH. This hypothesis aligns with findings that increased hepatic TH production via upregulated DIO1 has a protective effect against MASH37. Our in vitro experiments utilizing remestirom to assess THRB action revealed that THRB agonization reduces lipid and cholesterol levels, mimicking the hepatic action of T3.
THRB deletion leads to significant disruptions in lipid metabolism diurnal rhythms, without concurrent changes in the expression of core diurnal clock genes. This uncoupling of clock and clock outputs is reiterated in the liver under varied TH conditions, as previously shown11,14, suggesting that TH may influence diurnal physiology downstream of or independently from the diurnal machinery. This concept remains to be further elucidated, e.g., in clock mutant models. Moreover, our study has identified several DRGs that are directly or indirectly affected by THRB. However, bulk RNAseq data cannot distinguish the main cell types that may contribute to this effect. Combining our diurnal dataset with (non-temporal) single-cell RNAseq data could provide interesting insights and hypotheses for further experimentation. For example, genes with a mesor UP effect, such as Cyp7a1, that encode a crucial enzyme that catalyzes the first step in the conversion of cholesterol into bile acids, and Npc1, an intracellular cholesterol transporter, are predominantly expressed in Thrb-positive hepatocytes. On the other hand, Baat, which encodes an enzyme involved in the conjugation of bile acids, is expressed in both Thrb-positive and Thrb-negative hepatocytes (Supplementary Fig. 4).
Our research is not without its limitations. The use of systemic THRB knockout mice introduces known non-hepatic variables, including increased serum T3 levels12. Increased serum T3 could activate THRA in the liver (and other organs) and contribute to the observed effects in THRBKO mice. Additional studies with hepatocyte-specific THRBKO mice and mice with deletions in the clock gene machinery are required to clarify how hepatocyte THRB regulation interacts with local clocks to shape diurnal transcriptome and lipidome regulation. Our differential rhythm analysis was based on an uncorrected p-value, which may result in higher false positive rates. To compensate for this, a stricter p-value cut-off (<0.01) was used to establish rhythmicity, which was subsequently used for differential rhythm analysis. Importantly, our classification of liver as a low-TH action organ is mainly based on gene expression signatures. Direct quantification of intra-hepatic TH levels is required to substantiate this concept and enzymatic profiling under different steatosis conditions to better assess the physiological consequences.
Methods
Mouse model and experimental conditions
THRB knockout (for both isoforms) mice and controls (THRB+/−), which are phenotypically similar to wild-type (THRB+/+)12,13 (B6.129S1-Thrbtm1Df/J Jackson # 003462) were bred on the C57Bl/6NCrl background in the Gemeinsame Tierhaltung (GTH) at the University of Lübeck. Two to three-months-old male controls (THRB+/−) and knockout (THRB−/−) were obtained from group-housed under 12 h/12 h light/dark conditions (LD, 200–400 lux) at 22 ± 2 °C with food and water provided ad libitum (normal chow, 5% fat, 1314, Altromin, Germany). Mice were culled by cervical dislocation at 4-hour intervals, tissues were kept in dry ice, and stored at −80 °C. Animal experiments were ethically approved by the Animal Health and Care Committee of the Government of Schleswig-Holstein and in line with international guidelines for the ethical use of animals. Data was reported according to the ARRIVE guidelines. For each timepoint, a total of 3 mice were used.
Transcriptome analyses
Total RNA was extracted from livers using TRIzol reagent (ThermoFisher, Germany) and columns (Zymo Research, USA) according to the manufactory’s instructions. Truseq RNAseq comprising 40 million reads per sample (12GB) was performed at Novogene (Germany). Samples underwent quality control and only those with RNA integrity number (RIN) higher than 6.5 were used. Messenger RNA was purified from total RNA using poly-T oligo-attached magnetic beads. After fragmentation, the first strand cDNA was synthesized using random hexamer primers, followed by the second strand cDNA synthesis using dTTP for non-directional library. The library was checked with Qubit and real-time PCR for quantification and bioanalyzer for size distribution detection. The clustering of the index-coded samples was performed according to the manufacturer’s instructions. After cluster generation, the library preparations were sequenced on an Illumina platform and paired-end reads were generated. Raw data (raw reads) of fastq format were firstly processed through in-house perl scripts. Clean data (clean reads) were obtained by removing reads containing adapter, reads containing poly-N and low-quality reads from raw data. Q20, Q30 and GC content the clean data were calculated.
Reference genome (GRCm38/mm10) and gene model annotation files were downloaded from genome website browser (NCBI/UCSC/ENSEMBL) directly. Paired-end clean reads were mapped to the reference genome using HISAT2 software. All the downstream analyses were based on clean data with high quality. FeatureCounts was used to count the reads numbers mapped to each gene. Genes containing a sum of reads among two groups lower than 100 were excluded from the analysis. The remaining genes were log2 transformed in DESeq2 using vst function. A total of 17,994 transcripts were considered for subsequent analyses.
Lipidomics analysis of the liver via liquid chromatography coupled to tandem mass spectrometry (LC-MS/MS)
All measurements were performed using a Dionex Ultimate 3000 RS LC-system coupled to an Orbitrap mass spectrometer (QExactive, ThermoFisher Scientific, Bremen, Germany) equipped with a heated-electrospray ionization (HESI-II) probe. All solvents were of LC-MS grade quality and were purchased from Merck (Darmstadt, Germany). Zirconium oxide beats and 500 µl of water were added to the frozen tissue (70–150 mg). Liver samples were homogenized for 3 min at speed level 8 (Bullet Blender, Next Advance, New York, USA) and centrifuged (2000 g/10 min/4 °C). The supernatant (400 µl) was transferred into a new reaction tube and diluted accordingly to reach a protein concentration between 0.6 mg/ml and 1.2 mg/ml. Protein content was determined using a Bradford assay (Sigma-Aldrich, USA).
Lipid profiling was performed as described previously38. Briefly, lipids were extracted by adding 1 ml of methanol/methyl tertiary-butyl ether/chloroform (1.33:1:1, v/v), containing butylated hydroxytoluene (100 mg/L, Sigma-Aldrich, Germany) and 2.5 ml/L of SPLASH internal standard mix (Avanti Polar Lipids, Alabaster, USA), to 50 µl of tissue homogenate. After incubation and centrifugation, supernatants were dried under vacuum and resuspended in 50 µl of methanol/isopropyl alcohol (1:1, v/v) for LC-MS/MS analysis. Extracted lipids were separated on an Accucore C30 RP column (150 × 2.1 mm, 2.6 µm) using acetonitrile/water 6:4 (v/v) as eluent A and isopropyl alcohol/acetonitrile (9:1, v/v) as eluent B; 10 mM ammonium acetate (Sigma-Aldrich, Germany) and 0.1% formic acid (Biosolve, France) were added to both eluents. Ionization and data acquisition with data-dependent MS² scans (top 15) was performed as described previously38,39. Lipid species were identified according to exact mass, isotopic pattern, retention time and fragmentation ions using Compound Discoverer 3.3 (Thermo Fisher Scientific) and two in silico databases40. The area under the peak was normalized to the internal standard and the protein content. An extraction blank and quality control samples were used to ensure data quality and linearity.
Differentially expressed gene (DEG) analysis
Global DEGs analysis was performed independently of sampling time. DESeq241 with an adjusted p value < 0.05 (Benjamini-Hochberg Procedure, FDR) and a log2 fold-change higher or low than zero was considered as upregulated or downregulated genes, respectively. To account for the effects of time in detecting DEGs, the factors time and genotype were included in DESEq2. The likelihood ratio test (LRT) was used to identify the influence of the genotype, considering the time factor. An FDR < 0.05 and a log2 fold-change higher or lower than zero were considered upregulated or downregulated genes, respectively. For the microarray data from low or high TH state mice (GSE199998), Two-Way ANOVA was performed in a pair-wise comparison (CON vs. T3 and MMI vs. T3). Genes with a significant effect caused by treatment (MMI or T3, FDR < 0.05) were considered significant.
Rhythm analyses
To identify diurnal (24-hour) oscillations in gene expression, we combined four independent algorithms: CircaN42, JTK cycle43, Metacycle44, and DryR45. A single list of rhythmic genes (identified as rhythmic in at least one method) was input in CircaCompare17 to compare rhythmic parameters between groups. For mesor and amplitude comparison, CircaCompare was allowed to fit a sine curve irrespective of meeting rhythmicity thresholds. Phase comparison was only performed in robustly rhythmic genes as described11. CircaN, JTK cycle, and Metacycle algorithms were run in CircaN algorithm using the default parameters. An exact period of 24 h was used. The significance of rhythmicity was set as p < 0.01 for each method. For all comparisons in CircaCompare a p value of 0.05 was used. The same pipeline was done for lipidomics samples, except DryR method was not included. For heatmap and rose plots, phase estimation was extracted from CircaCompare. Rhythmic gene visualization was made using CircaCompare algorithm or ggplot (method “lm” and formula = y ~ 1 + sin(2*pi*x/24) + cos(2*pi*x/24). Methodological constraints led us to perform differential rhythm analysis using CircaCompare pairwise, without adjusting for multiple comparisons, potentially inflating the rate of false-positive results.
Gene enrichment and phase set enrichment analysis (PSEA)
Enrichment analyses were performed using Gene Ontology (GO) annotations for Biological Processes with DAVID 6.8 software46 and cut-offs of 3 genes per process and a p-value < 0.05. REVIGO (reduction of 0.7 and default conditions) was used to reduce redundant calls in the enrichment analyses47. The phases of rhythmic genes were estimated using CircaSingle, round up to exact number, and used for PSEA analysis16. Mouse gene names were converted to human with the assistance of geneName (version 0.2.3). A minimum of 10 genes and a p value of < 0.01 were selected (against the background). Biological processes IDs were extracted using the biological processes’ names using the version 7.5.1 (http://geneontology.org/docs/download-ontology/). Processes were filtered in REVIGO (0.7 reduction). The phase of biological processes was estimated by calculating the median peak expression from all genes enriched to a biological process. Enriched processes from each genotype were merged into a single list, manually curated by assigning each biological process to a higher hierarchical biological process class. Phase differences between genotypes were determined by comparing process phases between groups.
In vitro experiments
AML-12 cells, sourced from the ATCC Biobank (catalog no. CRL-2254), were cultured in Gibco Dulbecco’s Modified Eagle Medium supplemented with Nutrient Mixture F12 (DMEM/F12, ThermoFisher). The growth medium was enriched with 1% penicillin-streptomycin, 1% insulin-transferrin-selenium (ITS, ThermoFisher), and 10% non-heat-inactivated fetal bovine serum (FBS, ThermoFisher). Additionally, 10 nM dexamethasone (Sigma-Aldrich) was included to maintain the cells, which were incubated at 37 °C in a humidified atmosphere with 5% CO2. For specific experiments, dexamethasone was excluded, and regular FBS was substituted with charcoal/dextrane-treated FBS (Hyclone) to reduce T4 and T3 concentrations by 40% and 90%, respectively, in addition to reducing endogenous hormones.
Initially, 105 AML-12 cells were seeded per well in 12-well plates. After 24 hours, the cells underwent synchronization with 200 nM dexamethasone for a duration of 2 hours. Post-synchronization, cells were washed with PBS, and drug pre-treatments were administered for 48 h. Subsequent to the drug pre-treatments, a 0.25 mM concentration of BSA-palmitate saturated fatty acid complex (5 mM, 7:1, Cayman Chemical) was introduced. At 48 hours post-palmitate addition, both cells and supernatant were collected. The collection protocol involved centrifugation and washing steps, followed by incubation with Zombie NIR dye (1:500 dilution, Biolegend), Bodipy-Cholesterol (1 µM, MedChem Express), and adipored (0.3 µL per mL, Lonza) at room temperature. Each incubation step was combined with centrifugation and washing. Cells were then fixed in 4% paraformaldehyde (PFA) for 30 minutes at room temperature, subsequently stored for fluorescence-activated cell sorting (FACS) analysis. All drugs were dissolved in dimethyl sulfoxide (DMSO), ensuring the final DMSO concentration remained below 1% in all treatment conditions. The resmetiron (MedchemExpress) concentration was determined based on known EC50 (0.2 µM). T3 concentration was 100 nM.
Samples were assessed in Cytek Aurora at the Cell Analysis Core Facility (CAnaCore facility at the University of Lübeck). Data were acquired using predefined templates with adjustments for optimal event rates and detector gains. Spectral unmixing was applied using the manufacturer’s software with single-stained controls to resolve fluorescent populations. At least 10,000 events were captured for each sample. Gating strategies were implemented using FlowJO V10 software, with each gate clearly defined in terms of fluorescence. Data was analyzed using Two-Way ANOVA followed by false discovery rate (Two-stage linear step-up procedure of Benjamini, Krieger, and Yekutieli) correction for multiple comparisons in Prisma (Version 10).
Responses