A preliminary transcriptome analysis suggests a transitory effect of vitamin D on mitochondrial function in obese young Finnish subjects.

OBJECTIVE: The effect of vitamin D at the transcriptome level is poorly understood, and furthermore, it is unclear if it differs between obese and normal-weight subjects. The objective of the study was to explore the transcriptome effects of vitamin D supplementation. DESIGN AND METHODS: We analysed peripheral blood gene expression using GlobinLock oligonucleotides followed by RNA sequencing in individuals participating in a 12-week randomised double-blinded placebo-controlled vitamin D intervention study. The study involved 18 obese and 18 normal-weight subjects (of which 20 males) with mean (±s.d.) age 20.4 (±2.5) years and BMIs 36 (±10) and 23 (±4) kg/m2, respectively. The supplemental daily vitamin D dose was 50 µg (2000 IU). Data were available at baseline, 6- and 12-week time points and comparisons were performed between the vitamin D and placebo groups separately in obese and normal-weight subjects.
RESULTS: Significant transcriptomic changes were observed at 6 weeks, and only in the obese subjects: 1724 genes were significantly upregulated and 186 genes were downregulated in the vitamin D group compared with placebo. Further analyses showed several enriched gene categories connected to mitochondrial function and metabolism, and the most significantly enriched pathway was related to oxidative phosphorylation (adjusted P value 3.08 × 10-14). Taken together, our data suggest an effect of vitamin D supplementation on mitochondrial function in obese subjects.
CONCLUSIONS: Vitamin D supplementation affects gene expression in obese, but not in normal-weight subjects. The altered genes are enriched in pathways related to mitochondrial function. The present study increases the understanding of the effects of vitamin D at the transcriptome level.

Introduction

Obesity and vitamin D deficiency represent frequent nutrition concerns worldwide and recent evidence links these together (1, 2, 3). Moreover, obesity and excess adipose tissue gain are suggested as a cause of vitamin D deficiency (3), defined as circulating 25-hydroxy vitamin D (25-OHD) concentration below 50 nmol/L (4). The underlying mechanism may involve a dilution effect of adipose tissue on circulating hormone concentration, and higher consumption of vitamin D by inflammatory processes in adipose tissue (5).

It is crucial to understand the metabolism and targeted effects of vitamin D and especially how these are affected in obesity. To date, we know that vitamin D supplementation increases circulating 25-OHD in obese individuals, although the effect is less pronounced than in normal-weight individuals. In contrast, the effect on biologically active, free 25-OHD is similar in the two groups (6). On the other hand, weight loss also has a beneficial effect on vitamin D status (7, 8, 9). However, studies comparing vitamin D supplementation and its effects on obese and normal-weight subjects are scarce, and well-planned studies are needed to reveal the effect on predefined outcomes. Defining the sufficiency of vitamin D should not be done solely based on the concentration of 25-OHD, but should be supported by additional information on the downstream effects of vitamin D. These mechanisms are likely to involve the genomic actions of vitamin D binding to vitamin D receptor (VDR), which is found in over 40 tissues (10) and modulates the expression of a large number of genes (11). Previous studies have mainly focused on the actions of vitamin D in specific cell lines as reviewed by Bouillon & colleagues (11), while the global transcriptomic effects of vitamin D are poorly understood.

Vitamin D deficiency is linked to a number of adverse health consequences (5, 12). For some of these, including cardiovascular disease, type 2 diabetes mellitus and metabolic syndrome, obesity is also a recognised contributor. It is essential to understand these risk factors and related molecular responses more in depth, and whether they are independent or interdependent. Recent applications of functional genomic approaches allow the systematic assessment of molecular responses in randomised trials, which unveil mechanisms of action in humans. This information will contribute to deeper understanding of the aetiology of the conditions, increase the efficacy of clinical trials and be valuable in the future prevention of chronic diseases.

We expanded a previously described 12-week randomised double-blinded placebo-controlled vitamin D intervention study involving both obese and normal-weight subjects (6). The objective of the present study was to examine the effect of vitamin D supplementation on gene expression levels in white blood cells separately in obese and normal-weight subjects, at three time points. Our findings yield insight into transcriptome changes connected to mitochondrial function unique to obese subjects at the 6-week time point.

Materials and methods

Materials

We analysed the material collected for a 12-week vitamin D intervention trial in obese and normal-weight subjects, which has been described in detail elsewhere (6). In short, the study was conducted at the Children’s Hospital, Helsinki University Central Hospital, Finland. All participants were of Caucasian origin. Subjects with childhood-onset obesity (in the text referred to as obese) were identified from the Hospital’s patient registry and those fulfilling the following criteria were eligible to participate: (i) weight-for-height ratio exceeding 60% (corresponding to age and sex-specific BMI ≥99th percentile) before the age of 7 years, according to Finnish growth standards (13) and persistence of severe obesity for at least 3 years in childhood and (ii) follow-up during childhood at the Children’s Hospital, Helsinki University Hospital due to severe obesity. Among obese subjects, known endocrine and genetic disorders underlying obesity were excluded. We recruited normal-weight controls with a similar age range (in the text referred to as normal-weight) from the national population register and excluded subjects with childhood-onset obesity (weight-for-height ratio above 40% before the age of 10 years). Other exclusion criteria included sunny holidays in the 2 months preceding the trial. Habitual supplement use was not an exclusion criterion in the present study but subjects were advised to maintain their current dietary habits during the study. All study participants and guardians of minors gave informed written consent before joining the study. The Research Ethics Committee of the Hospital District of Helsinki and Uusimaa approved the study (295/13/03/03/2012), and it is registered at ClinicalTrials.gov (NCT02549326).

Intervention was double-blinded, obese and normal-weight subjects were randomised into two groups, receiving either placebo or cholecalciferol 50 µg (2000 IU) daily (Minisun®, Verman, Kerava, Finland). The study protocol included three study visits: baseline (1st time point in the text), at 6 weeks (2nd time point) and at 12 weeks (3rd time point). During each visit, a fasting blood sample was collected between 8 and 10 h in two PAXgene Blood RNA Tubes (2.5 mL, PreAnalytiX GmbH, Hembrechtikon, Switzerland). Two hours after sampling, the PAXgene tubes were transferred to a −80°C freezer for storage and later analysis.

Methods

RNA samples

Total RNA from the PaxGene RNA tube was extracted using the PaxGene RNA kit (PreAnalytiX GmbH, Hembrechtikon, Switzerland) according to the manufacturer’s instructions. Genomic DNA was removed from the whole-blood samples using the DNA-free kit (Ambion), leaving the total RNA fraction. RNA integrity numbers (RINs) were assessed on an Agilent Bioanalyzer 2100 at the Biomedicum Functional Genomics Unit (FuGU), Helsinki, Finland. RNA concentrations were measured using a Thermo Fisher Qubit 2.0 Fluorometer.

Library prep

An 80 ng sample of whole-blood RNA was treated using GlobinLock oligonucleotides (14) to mask globin transcripts before cDNA synthesis and significantly reduce the proportion of globin α and β reads in the sequencing data. This was followed by RNA sequencing using a modified version of the single-cell tagged reverse transcription (STRT) method (15) described in detail elsewhere (16).

Sequencing

The Illumina-compatible library pool was sequenced on four lanes each using the Illumina HiSeq2000 instrument and Illumina TruSeq v3 (Illumina) 60-bp single-read protocol. Sequencing was carried out at the BEA core facility at Karolinska Institutet, Huddinge, Sweden. Sequence data were converted to fastq files using Casava 1.8.2 (Illumina), and quality control performed using the STRTprep pipeline available at https://github.com/shka/STRTprep (16).

Gene ontology/pathway analyses

The DAVID database (v6.8 at https://david.ncifcrf.gov/home.jsp, accessed 22 October 2017) was used to identify enriched gene ontology (GO) categories (GO BP_FAT, GO MF_FAT and GO_CC_FAT) and KEGG (17) pathways within the set of differentially expressed genes. The full list of significantly upregulated or downregulated genes was used in all pathway analyses.

qPCR validation

Five genes, with robust expression and significantly different expression levels in obese subjects in placebo and vitamin D group at the 6-week time point, were selected for technical validation by reverse transcriptase qPCR analysis. We included three obese subjects from the vitamin D group and three from the placebo group at the second time point and used the same RNA that was used for the initial whole-transcriptome analysis. cDNA was transcribed from 1 μg of total RNA using the QuantiTect Reverse Transcription Kit cDNA synthesis kit (Qiagen) according to the manufacturer’s protocol. qPCR assays were performed in quadruplicates with using the CXF96 Real-Time system (Bio-Rad Laboratories) and using the TaqMan ATP6V1F, COX5A, NDUFA2, ACTB and TBP gene expression assays Hs00855096_g1, Hs00362067_m1, Hs04187282_g1, Hs99999903_m1 and Hs00427620_m1, respectively (Applied Biosystems). ACTB and TBP were used as reference genes for data normalization. Threshold cycle (Ct) values were determined using CFX Manager Software (Bio-Rad Laboratories). Relative expression was calculated using the comparative Ct or 2ΔΔCt method (18, 19).

Serum biochemical markers

Laboratory methods concerning measurements of total 25-OHD, free 25-OHD, insulin, adiponectin and high-sensitive C-reactive protein (hs-CRP) are described in detail elsewhere (6, 20). In brief, serum total 25-OHD was measured with an automated IDS-iSYS analyser (IDS Ltd., Boldon, UK) with coefficient of variation (CV) <5% and compared against in-house liquid chromatography in tandem with mass spectrometry (20), and free 25-OHD was assessed using a two-step immunoassay (Future Diagnostics BV, Wijchen, Netherlands) with CV <8.5%. Plasma glucose was analysed by spectrophotometric a hexokinase and glucose-6-phosphate dehydrogenase assay (Gluko-quant glucose/hexokinase, Roche Diagnostics) with a Hitachi Modular automatic analyser (Indianapolis, USA). Serum insulin was measured using a time-resolved immunofluorometric assay (Perkin Elmer Life Sciences) with a detection limit of 0.5 mU/L and an interassay-CV <4%. Serum adiponectin was determined using a Human Total Adiponectin/Acrp30 Quantikine ELISA Kit and serum leptin with Human Leptin R Quantikine ELISA Kit (R&D Systems) with intra- and inter-assay CV of <12%. hs-CRP was determined by an immunoturbidimetric assay on a Roche-automated clinical chemistry analysed at the central laboratory of the Helsinki University Central Hospital.

Statistical methods

Baseline characteristics between obese and normal-weight groups were compared separately with t-tests, and in case of categorical variables, using a Fisher’s exact test. Changes in total and free 25-OHD, insulin, hs-CRP, adiponectin and the proportion of white blood cells between three time points were compared using the repeated-measures ANOVA. We tested changes over time and responses as area under the curve between groups of interest. Differences were considered significant at P < 0.05. All statistical analyses were conducted using the IBM SPSS program for Windows, version 22 (IBM).

The STRTprep pipeline uses SAMstrt (21) to identify differentially expressed genes; the differences were considered significant at fluctuation P values <0.05 and difference q values <0.05. The fluctuation is a test of alternative hypothesis that the biological CV in the samples is larger than the square of the technical CV; the technical CV is estimated by spike-in RNAs (16). The difference between the sample groups is a test of the alternative hypothesis that the distribution of expression levels in one group is stochastically greater than the other (21). The fluctuation P values were corrected by the Benjamini and Hochberg method, and difference q values were corrected using the Storey and Tibshirani method.

The significance of pathway and GO category enrichments was assessed using q values, which are P values adjusted using the Benjamini–Hochberg method for multiple hypothesis testing.

Results

Clinical characteristics of the intervention cohort

Baseline characteristics are shown in Table 1 for each of the four groups: obese + placebo (O_P), obese with vitamin D (O_VD), normal-weight + placebo (NW_P) and normal-weight + vitamin D (NW_VD). Serum biochemical markers and anthropometric measurement did not differ between obese groups, except for hs-CRP, which was higher in O_P compared with O_VD. Correspondingly, the characteristics did not differ between normal-weight groups. At baseline, the fraction of subjects with vitamin D deficiency (25-OHD ≤50 nmol/L) was 3/9, 5/7, 1/8 and 3/8 in the O_P, O_VD, NW_P and NW_VD groups, respectively. The response to vitamin D supplementation was observed in both the total and free 25-OHD concentrations (P < 0.001): The obese subjects had a lower response in total 25-OHD but a similar free 25-OHD when comparing the O_VD and NW_VD groups (Fig. 1). There were no changes in insulin (repeated-measures ANOVA: P = 0.475), hs-CRP (P = 0.062) or adiponectin (P = 0.263) levels during the 12 weeks and these did not differ between the four groups (ANOVA P = 0.278, P = 0.753 and P = 0.516, for insulin, hs-CRP and adiponectin, respectively).

Table 1

Cohort characteristics with mean (s.d.), unless otherwise noted.

	O_P	O_VD	NW_P	NW_VD	P value_O	P value_NW
N	9	7	8	8
Maleb	7/9	4/7	4/8	5/8	0.596a	0.608a
Age (years)	18.9 (2.4)	19.8 (3.2)	21.9 (2.2)	21.1 (2.4)	0.533	0.502
BMI (kg/m2)	38.8 (13.3)	34.4 (6.6)	22.0 (2.8)	25.2 (5.3)	0.403	0.168
Waist (cm)	117.8 (29.6)	104.7 (19.3)	74.0 (7.0)	83.2 (14.0)	0.306	0.127
S-25-OHD (nmol/L)	55.8 (13.6)	44.9 (13.0)	66.0 (23.7)	51.8 (17.4)	0.128	0.197
Free 25-OHD (pg/mL)	3.4 (1.3)	2.6 (1.1)	5.0 (1.9)	4.1 (1.7)	0.219	0.323
Free/total 25-OHD × 102	6.20 (2.20)	5.80 (1.51)	7.78 (2.33)	7.84 (1.27)	0.682	0.945
Insulin (mU/L)	13.8 (7.8)	15.8 (9.2)	5.4 (2.7)	9.7 (5.9)	0.647	0.091
hs-CRP (mg/L)	7.55 (7.2)	0.87 (0.5)	0.79 (0.6)	0.68 (0.8)	0.024	0.773
Adiponectin (ng/mL)	3760 (2327)	5180 (2689)	10221 (7316)	10615 (6216)	0.288	0.909
Gc genotype 1/1b	5/9	5/7	8/8	6/8	0.633a	0.467a
White blood cell count, cells × 109/L	7.8 (2.4)	6.8 (1.1)	4.4 (1.2)	5.4 (1.9)	0.341	0.231
Neutrophils	4.2 (1.3)	3.2 (1.1)	2.1 (0.9)	2.6 (1.1)	0.125	0.289
Lymphocytes	2.7 (1.1)	2.7 (0.8)	1.7 (0.4)	2.1 (0.7)	0.971	0.213
Monocytes	0.55 (0.17)	0.62 (0.14)	0.39 (0.08)	0.49 (0.17)	0.398	0.163
Eosinophils	0.25 (0.17)	0.25 (0.21)	0.25 (0.21)	0.19 (0.13)	0.992	0.524
Basophils	0.04 (0.03)	0.02 (0.01)	0.03 (0.02)	0.04 (0.03)	0.200	0.949

aFisher’s exact test; bfraction.

NW_P, normal-weight + placebo; NW_VD, normal-weight + vitamin D; O_P, obese + placebo; O_VD, obese + vitamin D; P value_O, P value for the comparison of obese groups; P value_NW, P value for the comparison of normal-weight groups.

Figure 1

Change in serum total (A) and free (B) 25-OHD concentrations and their adjusted ratio (C) in the four groups: obese + placebo (O_P), obese + vitamin D (O_VD), normal-weight + placebo (NW_P), and normal-weight + vitamin D (NW_VD) with mean and SEM during the 12-week study.

RNA samples and sequencing

The RNA samples were of high quality (RIN >8, in all samples) and similar in all groups. Supplementary Table 1 (see section on supplementary data given at the end of this article) shows the results of the sequencing analysis in both normal-weight and obese individuals, at three time points, including number of reads per sample and proportion of sequencing reads aligning to the 5´ends of coding genes. Five samples were removed from further analysis (two due to low RNA quality, three due to low yield of RNA ERCC spike-in). Each sample had on average >6 million sequence raw-reads, and on average, >5 millions of those mapped to the human genome (hg19). No significant differences between the quality of the libraries, between different treatments or time points were seen (data not shown).

Differential expression analysis

A comparison of the gene expression in the vitamin D and placebo groups was performed separately in the obese and normal-weight subjects at the three time points. Principal component analysis (PCA) plots for each comparison are shown in Supplementary Fig. 1 and the full dataset of normalised expression values for each subject and time point, and each detected gene, can be found in Supplementary Dataset 1. In normal-weight subjects, the differences in gene expression between the vitamin D and placebo groups were not significant at any time point.

Significant fold-changes were observed only at the second time point in the obese individuals (Table 2 and Supplementary Table 2), between the placebo (O_P) and vitamin D (O_VD) groups, and thus, we proceeded further only with these samples. Figure 2 shows a hierarchical clustering of the samples from obese individuals at the second time point and a heatmap of the differentially expressed genes (1724 upregulated genes and 186 downregulated, in vitamin D supplemented compared with placebo receiving individuals). Individuals receiving vitamin D and placebo cluster separately, except for one individual in placebo group. The study was underpowered to look at changes in gene expression within subjects (data not shown). Still, we provide pairwise analyses of the top three upregulated and downregulated genes in Supplementary Figs 2 and 3 to support our group-wise comparisons.

Table 2

List of top 50 up- and downregulated genes in the obese group with vitamin D (see Supplementary Table 2 for full list).

Upregulated			Downregulated
Gene	DiffexpScore	q value	Gene	DiffexpScore	q value
TAPBP	24.00	0.01	HBQ1	−16.15	0.01
XAB2	23.95	0.01	TUBGCP3	−14.10	0.01
CHERP	23.90	0.01	KANK2	−14.00	0.01
ETHE1	23.90	0.01	FUNDC2	−13.00	0.01
ITGB1BP1	23.75	0.01	AK1	−11.75	0.01
DECR1	23.50	0.01	SLC25A37	−11.60	0.01
EFTUD2	23.35	0.01	FAM210B	−11.30	0.01
RASGRP4	23.25	0.01	MTRNR2L9	−11.30	0.01
CFLAR	23.15	0.01	EPB42	−11.20	0.01
ABI3	23.10	0.01	UBB	−10.85	0.01
FAM127A	23.05	0.01	HEMGN	−10.70	0.01
MAZ	22.95	0.01	GCAT	−10.05	0.01
SIRT7	22.85	0.01	BIRC2	−10.00	0.01
PSAP	22.75	0.01	TMEM176B	−9.80	0.01
IFI27L2	22.70	0.01	HMBS	−9.20	0.01
TMEM161A	22.70	0.01	HBB	−8.90	0.01
CXCR3	22.65	0.01	MICALL2	−8.85	0.01
VIM	22.65	0.01	AHSP	−8.40	0.01
MYO1G	22.60	0.01	OPTN	−8.35	0.01
SLC25A3	22.60	0.01	RNF10	−8.30	0.01
CALM1	22.55	0.01	HBG2	−8.25	0.01
SERPING1	22.55	0.01	MYL4	−8.20	0.01
C12orf75	22.50	0.01	GMPR	−8.15	0.01
CRTC2	22.45	0.01	KLC3	−8.10	0.01
ZFP36L2	22.45	0.01	TMEM176A	−7.95	0.01
GSTP1	22.30	0.01	HBD	−7.90	0.01
SEPT9	22.30	0.01	GSTM1	−7.70	0.01
PEF1	22.25	0.01	ISCA1	−7.60	0.01
RHBDF2	22.25	0.01	CA1	−7.15	0.01
PRAM1	22.20	0.01	LYPD2	−6.80	0.01
ANXA2	22.15	0.01	SLC22A16	−6.80	0.01
CHI3L1	22.15	0.01	GSTM3	−6.75	0.01
DHCR7	22.15	0.01	HAGH	−6.75	0.01
MRPL43	22.15	0.01	MTRNR2L2	−6.75	0.01
TP53I13	22.15	0.01	PPP1R3B	−6.45	0.01
HMGA1	22.10	0.01	PNP	−6.40	0.01
SQSTM1	22.05	0.01	YBX1	−6.35	0.01
CD7	22.00	0.01	C7orf73	−6.25	0.01
CD74	22.00	0.01	PRDX6	−6.15	0.01
MGST3	22.00	0.01	DHX16	−6.00	0.01
PLD3	21.95	0.01	ADIPOR1	−5.95	0.01
EMD	21.90	0.01	SLC2A1	−5.95	0.01
TRAPPC2L	21.90	0.01	RUNX2	−5.90	0.01
CHTF8	21.85	0.01	STRADB	−5.80	0.01
ARSA	21.80	0.01	TESC	−5.80	0.01
ATP5G3	21.80	0.01	FAM104A	−5.75	0.01
SSH3	21.80	0.01	C6orf48	−5.70	0.01
ARHGAP1	21.75	0.01	EIF1AY	−5.70	0.01
ATF5	21.75	0.01	RNF175	−5.60	0.01
PDHB	21.70	0.01	ACKR1	−5.55	0.01

Figure 2

Hierarchical clustering of samples based on the expression of differentially expressed genes (vitamin D treatment compared with placebo). Red for higher and blue for lower expression.

The qPCR results of ATP6V1F, COX5A and NDUFA2 gene expression were in line with the results of the STRT RNA-seq data (Supplementary Fig. 4). The expression of ATP6V1F (Independent samples t-test: P = 0.003) and COX5A (P = 0.002) differed between the vitamin D (O_VD) and placebo (O_P) groups: a significant upregulation was observed with vitamin D compared with the placebo group at the 6-week time point in the obese ones. There were no significant differences between groups in the expression of NDUFA2 (P = 0.29) most likely due to low sample size.

Pathways

In order to understand the potential mechanism underlying vitamin D response at the transcriptome level, we performed GO category and KEGG pathway analysis of the differentially expressed genes. Supplementary Dataset 2 shows the significant or top 25 GO categories (BP, MF and CC) that were enriched in the list of genes upregulated upon vitamin D supplementation in the obese individuals. The enriched GO-BP categories (representing biological processes) include the immune response and nucleoside metabolic processes. The cellular component categories (GO-CC) highlight the mitochondrial membrane. Lastly, the enriched GO-MF (molecular function) categories involve RNA binding, cadherin binding and transmembrane ion transport. Supplementary Dataset 2 shows the GO categories that were enriched in downregulated genes.

Figure 3 shows the top 10 upregulated KEGG categories in the vitamin D-supplemented obese group (O_VD) at the 6-week time point. Supplementary Dataset 3 shows the full list (n = 14) of significant pathways. The most significantly enriched category was oxidative phosphorylation (hsa:00190) with an adjusted P value of 3.08 × 10−14. Figure 4 summarises all significantly differentially expressed genes in this pathway, highlighted with a red asterisk. The other top enriched pathways include Alzheimer’s disease, Huntington’s disease and Parkinson’s disease. While the names of these pathways imply a connection to a specific disease, looking further into the genes involved instead shows that they show considerable overlap with the oxidative phosphorylation pathway (Supplementary Fig. 5), with 51.8% of the differentially expressed genes in these pathways being shared by all four pathways. Thus, it is clear that the four pathways largely overlap, despite their names, and the common factor is genes involved in mitochondrial function. A second group of enriched KEGG pathways is related to immune responses including phagosome (hsa04145) and antigen processing and presentation (hsa04612) pathways. Only one KEGG pathway, hsa03010:ribosome, was significantly enriched in the set of downregulated genes (adjusted P value 6.54 × 10−7).

Figure 3

The ten most significantly enriched KEGG pathways in genes upregulated by vitamin D treatment compared to placebo. The number of upregulated genes within each pathway is shown.

Figure 4

Schematic of the major components of the KEGG oxidative phosphorylation pathway, from the DAVID bioinformatics portal (https://david.ncifcrf.gov/). Genes significantly upregulated at the second time point in obese subjects receiving vitamin D are marked with a red star.

Discussion

In the present study, we applied the recently published GlobinLock technology (14) to minimise globin RNA contamination when performing transcriptome analysis on RNA from circulating white blood cells. Our results do not allow the study of the tissue-specific actions of vitamin D, for example, in adipose tissue, but they give an overview of the global transcriptomic effects of vitamin D. Our study shows for the first time that vitamin D supplementation has an impact on gene expression in circulating white blood cells of obese individuals, an effect which was not observed in normal-weight subjects. A part of the affected genes are enriched in pathways related to immune responses, which is in line with previous findings (10, 11, 22). The strongest transcriptomic effects, however, are those connected to mitochondrial function.

Oxidative pathway

An activation of the oxidative phosphorylation pathway by vitamin D supplementation was observed in obese subject in our study. To our knowledge, there is no comparable study available. Oxidative phosphorylation is a part of the energy-producing pathway within the mitochondrial inner membrane and exerts its actions on cellular energy production. Several studies have reported reduced mitochondrial number, size and activity in skeletal muscle and adipose tissue in patients with obesity and type 2 diabetes (23, 24, 25). In fact, mitochondrial dysfunction is suggested to precede the development of insulin resistance (26, 27). In a recent metabolomics study, the metabolic signature differed between obese children with and without insulin resistance, and the enriched and the most altered pathways included the urea cycle, alanine metabolism and the glucose–alanine cycle, also connected to mitochondrial function (28).

Using a systematic approach, we detected significant upregulation of altogether 53 of the approx. 130 genes in the KEGG oxidative phosphorylation pathway in obese subjects receiving vitamin D supplementation. A link between vitamin D and mitochondria function has been reported previously in experimental models (29, 30, 31, 32, 33, 34), this was recently summarised by Ricca et al. (35): 1,25(OH)2D limits excessive mitochondrial respiratory activity and ROS production and turns cell metabolism toward biosynthesis and proliferation. A number of these effects are highly cell type specific. Our finding that the uncoupling protein 2, UCP2, is upregulated upon vitamin D supplementation is not in line with previous observations (30, 31, 32). This may be due to differences in timing, discordance in phenotypes or for example in the severity of the vitamin D deficiency.

Obesity-specific effects

Interestingly, the transcriptomic effect of vitamin D was seen only in obese individuals. A plausible explanation for this is that the obese are more likely to have dysfunction of mitochondria and in the oxidative phosphorylation than normal-weight individuals. If we use fasting insulin concentration as a proxy for metabolic health, we notice that obese individuals in our cohort had higher median values and more variation in insulin concentrations compared with controls. Abnormally high fasting insulin concentrations (≥12 mU/L) (36) were also mostly observed in the obese subjects. Furthermore, it has been shown that mitochondrial dysfunction decreases adiponectin secretion in adipocytes, resulting in lower glucose uptake in muscle (37) and adipocytes (38). In line with this, we observed lower adiponectin concentrations in the obese groups compared with the normal-weight groups. At baseline, the frequency of vitamin D deficiency was more common in obese than in normal-weight subjects, which is an alternative explanation to our finding. Hence, the vitamin D supplementation might have a stronger corrective effect in those with vitamin D deficiency. Limited sample size precluded addressing the confounders in more detail. On the other hand, the ratio of free-to-total 25-OHD, reflecting availability of vitamin D at the tissue level, differed at the second time point between the obese groups. This might also be the explanation for the observed transcriptomic effects.

Implication for further studies or for other phenotypes

Our findings on transient effects are novel, but these need to be replicated. However, they potentially have several implications for future research not only in obesity but also in public health. Mounting evidence links vitamin D deficiency to the development of type 2 diabetes (39, 40), but the exact mechanism has remained inconclusive at least in humans. It is widely accepted that impaired mitochondrial function and insulin resistance are related. Thus, vitamin D supplementation may have far-reaching consequences on insulin resistance, which is supported by some animal (41) and human studies (42), although not consistently (43). Correspondingly, vitamin D deficiency has been linked to Parkinson’s disease (44, 45) and Alzheimer’s disease (44) but only a single study has linked it to Huntington’s disease (46). Our study supports these associations, presumably through the shared mitochondrial mechanism.

Timing

The most pronounced effect of vitamin D on the transcriptome is seen at 6 weeks after the initiation of the vitamin D supplementation, which is a novel finding. The effect was already less clear again at the 12-week time point. This might be related to the vitamin D turnover rate in the body (or saturation level), which depends on the dose of vitamin D, baseline 25-OHD concentration and polymorphism of vitamin D-related genes (47, 48, 49). Similarly, we and others have demonstrated with multiple doses of vitamin D that the effect on circulating 25-OHD plateaus after 6 weeks (47, 50). The catabolizing enzyme CYP24A1 is known to initiate the degradation pathway of 25-OHD (10), but the expression of CYP24A1 was not detected in the present study, likely due to the rather tissue-specific expression of this protein. In a Norwegian study, long-term vitamin D supplementation induced only minor changes in whole-blood gene expression in 47 overweight senior citizens with prediabetes (51). The global transcriptome was analysed at a single 50-month time point. In a subgroup analysis of the low versus high 25-OHD group, differences in the expression of 198 mRNAs was observed. The affected genes were enriched within pathways involving oxidative stress response, apoptosis signalling and gonadotropin-releasing hormone receptor. These authors suggested that short-term effects of vitamin D supplementation, which is what we look at in our current study, are likely to have a more prominent effect on global transcriptome than a long-term study.

On the contrary, a more recent, randomised, placebo-controlled trial, supplementing older English men with high doses of vitamin D for 12 months detected no significant effect on gene expression or circulating cytokine concentrations (52). It is likely that they missed the early transcriptomic changes, since the sampling was performed only at baseline and 12 months.

Strengths and limitations

Our study has several limitations, including the descriptive nature of the study and due to removal of samples failing quality control, the groups were not equal-sized, which resulted in a loss of statistical power. The total amount of RNA in the sequencing was the same in each sample, and since the proportion of cell types did not differ between the groups and did not change with time within the groups, our findings are not likely to be confounded by cell type distributions. Reassuringly, we did not observe significant transcriptome differences between the obese groups at baseline. A more general limitation is that the transcriptome of blood cells serves merely as a proxy of an endocrine phenotype. The mRNA transcriptome furthermore represents only a snapshot of the ongoing processes. Despite the intervention trial, the study could not prove causality, and the possibility of indirect effects mediated by calcium fluxes remains unsolved (31).

The current study can be regarded as a pilot study or descriptive study on the transcriptomic effects of vitamin D supplementation due to lack of functional outcomes. Future work should ideally be expanded to larger cohorts and to include analyses of multiple tissues and to verify these findings at the protein level and link them with endocrine and metabolic parameters.

Conclusions

Here, we demonstrate that vitamin D supplementation affects gene expression in young obese subjects supplemented with a daily dose of vitamin D (50 μg) and the affected genes are enriched within pathways related to mitochondrial function. The present study increases our understanding of the effects of vitamin D at the transcriptome level but is not conclusive and should be regarded as pilot data for future studies. Our results may be valuable for considering the effects of vitamin D supplementation in obese individuals and when evaluating the adequacy and safety of the supplementation dose.

Declaration of interest

E E, M P, S K, O M and H V have nothing to declare. K K and J K are co-inventors of GlobinLock international patent application.

Funding

The authors also acknowledge support from the Academy of Finland, Foundation for Pediatric Research, Päivikki and Sakari Sohlberg Foundation, Folkhälsan Research Foundation, the Swedish Research Council, Swedish Childhood Cancer Foundation and the Sigrid Jusélius Foundation. The computations were performed on resources provided by SNIC through the Uppsala Multidisciplinary Center for Advanced Computational Science (UPPMAX) under project b2014069. They thank the Karolinska Bioinformatics and Expression Analysis (BEA) core facility, Huddinge, Sweden.