Vitamin D treatment induces in vitro and ex vivo transcriptomic changes indicating anti‐tumor effects

Vitamin D deficiency is associated with risk of several common cancers, including colorectal cancer (CRC). Here we have utilized patient derived epithelial organoids (ex vivo) and CRC cell lines (in vitro) to show that calcitriol (1,25OHD) increased the expression of the CRC tumor suppressor gene, CDH1, at both the transcript and protein level. Whole genome expression analysis demonstrated significant differential expression of a further six genes after 1,25OHD treatment, including genes with established links to carcinogenesis GADD45, EFTUD1 and KIAA1199. Furthermore, gene ontologies relevant to carcinogenesis were enriched by 1,25OHD treatment (e.g., ‘regulation of Wnt signaling pathway’, ‘regulation of cell death’), with common enriched processes across in vitro and ex vivo cultures including ‘negative regulation of cell proliferation’, ‘regulation of cell migration’ and ‘regulation of cell differentiation’. Our results identify genes and pathways that are modifiable by calcitriol that have links to CRC tumorigenesis. Hence the findings provide potential mechanism to the epidemiological and clinical trial data indicating a causal association between vitamin D and CRC. We suggest there is strong rationale for further well‐designed trials of vitamin D supplementation as a novel CRC chemopreventive and chemotherapeutic agent.


| INTRODUCTION
Vitamin D deficiency is associated with risk of several common cancers, with the strongest available evidence supporting a link between vitamin D and colorectal cancer. 1,2 Observational data implicating vitamin D deficiency in CRC etiology or survival are limited by potential bias as environmental risk factors associated with CRC are also associated with vitamin D status (i.e., co-causality; e.g., physical activity), while the colorectal cancer or its treatment may impair vitamin D status (i.e., reverse causation). Indeed, a definitive causal relationship remains unproven, 3 with existing mendelian randomization (MR) studies hampered by weak genetic instrumental variables [3][4][5] and large trials (e.g., VITAL Trial 6 and ViDA study 7 ) failing to demonstrate an effect of vitamin D supplementation on CRC risk. However, a recent randomized-control trial (RCT) reported an association between supplementation, vitamin D receptor genotype and risk of colorectal adenoma, supporting the premise that the beneficial effect may be causal. 8 Meanwhile, vitamin D-related genetic variation has been shown to influence the association between 25-OHD level and CRC survival, [9][10][11] and a recent meta-analysis of RCT data strongly supports a causal effect for vitamin D supplementation on CRC mortality. 12,13 A variety of pre-clinical experimental studies have provided some understanding of the biological functions of vitamin D in relation to cancer initiation or progression. In vitro studies using CRC cell lines have demonstrated the influence of vitamin D on relevant cellular processes including proliferation and apoptosis. 14 Meanwhile, differential gene expression is reported in CRC and adenoma tissue compared to normal colorectal tissue, 15,16 with genes involved in metabolism, transcription and translation and cellular processes commonly altered 17 and vitamin D is known to influence gene expression through activation of the ligand-activated transcription factor VDR. [18][19][20] Palmer et al., report induction of CDH1 in response to 1-alpha,25(OH)2D(3) treatment of SW480 (CRC) cell lines, 21 providing possible mechanism via influences on E-cadherin and Wnt/β-catenin signaling pathway. Yu et al., report impacts of calcitriol on cell migration, apoptosis, autophagy, and epithelial-mesenchymal transition, sensitizing CRC cells to ionizing radiation. 22 Other reports demonstrate similarly broad influences on gene expression with vitamin D impacting cancer cell growth in vitro. 23,24 The limitations of comparative cell line studies in the investigation of gene expression are widely acknowledged; namely that accumulated cellular genomic and karyotypic abnormalities may significantly influence results from experimental studies. It is not surprising, therefore that epithelial organoids are now emerging as a non-aberrant ex vivo model system for the investigation of tumor initiation. The colonic epithelial organoid culture is a long-term culture system derived from intestinal crypt stem cells that aims to maintain basic crypt physiology and display the hallmarks of intestinal epithelium in terms of architecture, cell type composition, and self-renewal dynamics. A small number of studies have now demonstrated gene expression changes in patient-derived organoids treated with vitamin D, supporting the use of this as a model for vitamin D intervention studies. [25][26][27] Here, we investigate the effects of vitamin D on in vitro and ex vivo epithelial cell gene expression and assess for links to anti-tumor effects as a putative mechanism of action of vitamin D, as a chemopreventive or chemotherapeutic agent in CRC.

| Organoid culture, histopathological assessment and treatment
The culture of colonic organoids from human tissue has been described previously. 29 All participants provided informed written consent, and research was approved by local research ethics committees (13/SS/0248) and National Health Service management (2014/0058). Human colonic mucosa was removed during resectional colorectal surgery or via rectal biopsy using a rigid sigmoidoscope and rectal biopsy forceps. After washing, mucosa was incubated in crypt chelating solution (1× PBS, 45 mmol/L sucrose, 55 mmol/L D-sorbitol, 500 μmol/L DLdithiothreitol, 5  In brief, organoids were spun down at 500 g for 5 mins, fixed with 10% neutral buffered formalin overnight and embedded in 2% agarose prior to processing and paraffin embedding. For immunohistochemistry analysis, 5 µm sections were cut and stained by hematoxylin and eosin to assess morphology. Immunostaining was performed for βcatenin (1:200; BD Transduction), Ki67 (1:300; Abcam) and Cytokeratin 19 (1:500; Abcam) with the chromogenic signal was detected using DAB (#SK-4100, Vector Labs) and slides imaged using Nanozoomer (Hamamatsu). RNA Scope™ for Lgr5, an intestinal stem cell marker, was performed according to the manufacturers' protocol (ACD Newark, CA, USA). Organoid cultures were treated with calcitriol or ethanol control as per cell-line experiments, with time-course and dose-response experiments performed.

| Immunoblotting
Cell line protein extraction and western blotting for the initial dose-response experiment was performed as previously described 30 using a rabbit monoclonal E-cadherin antibody at a concentration of 1:100 at 4°C overnight (Cell signalling, Danvers, Massachusetts, USA, 3195). The secondary antibody used was goat anti-rabbit antibody tagged with the enzyme horseradish peroxidase (at 1:1000) at room temperature for an hour and antigen-antibody complexes were visualized with chemiluminescence. βactin was used as a loading control and was detected by a mouse monoclonal antibody at 1:4000 (Sigma), with secondary goat anti-mouse antibody at 1:1000 (Sigma). This was incubated at 37°C for 30 min and 95°C for 5 min.

| Quantitative real time polymerase chain reaction (qRT-PCR)
cDNA samples were diluted to produce a working stock (1:20 dilution for cDNA from cell lines, 1:10 from organoids). qRT-PCR was undertaken using TaqMan® Gene Expression Assays (Applied Biosystems; CDH1, Hs01023894_m1 and Hs01023895_m1; CYP24A1, Hs00167999_m1). For each Taqman® assay, a master mix was made up using Taqman® gene expression master mix (containing AmpliTaq Gold® DNA Polymerase, Uracil-DNA glycosylase, dNTPs (with dUTP), ROX™ Passive Reference, and optimized buffer components), the Taqman expression assay and RNase-free DNase-free water with stock cDNA. Standard curves were derived from undiluted stock solutions of cDNA (HCT116 colorectal cancer cell line, untreated or treated with calcitriol) to determine the linearity and amplification efficiency of each assay. Dilution of stock cDNA was initially 1:1, with a further six 1:4 serial dilutions. Water was also included as a negative control. Threshold cycles (Ct) were plotted against the logged cDNA quantity, and the coefficient of correlation obtained for the fitted calibration curves (R 2 ) were calculated to assess assay performance.
All qRT-PCR reactions were performed in triplicate, with amplification and detection of the product undertaken using the ABI PRISM® 7900HT Sequence Detection System (thermocycler; 50°C-2 min, 95°C-10 min, 95°C-15 s, 60°C-1 min repeated 40 times) and collected with SDS v2.3 software, using the FAM detector. Analysis of results included calculation against the standard curve and normalization to a reference gene derived via an online reference gene finder testing BestKeeper 31 and NormFinder. 32 Standard deviation between triplicates was calculated as a measure of methodological quality.

| Whole transcriptome profiling using HumanHT-12 v4.0 Expression BeadChip Arrays
Whole transcriptome profiling was performed on triplicate ETOH and calcitriol (50 nM) treated samples of HCT116, SW480 and LS174T cell lines and organoid samples using DNA microarray. RNA was quantified by NanoDrop (NanoDrop Technologies, Wilmington, DE, USA) and yield and integrity was assessed using the Agilent 2100 Bioanalyzer. The RNA yield and integrity of samples to be submitted for microarray expression/RNAseq analysis was assessed using the 2100 Bioanalyzer®. For HT12 analysis, total RNA was converted to double-stranded cDNA, followed by in vitro transcription amplification to generate labelled cRNA (Illumina TotalPrep™ RNA Amplification Kit). Gene expression profiling was undertaken using the HumanHT-12 v4.0 Expression BeadChip Arrays (Illumina) and IScan NO660 scanner, providing coverage of 47 231 transcripts and >31 000 annotated genes derived from the National Centre for Biotechnology Information Reference Sequence RefSeq Release 38 (November 7, 2009). Microarray data was exported from BeadStudio (Illumina) and processed in R using the limma package. In brief, the steps involved were background correction using negative controls, quantile normalization to remove technical variation and finally log2 transformation. Control probes and probes that were not expressed in at least three arrays to a detection-value of ≥5% were excluded. A standard approach to batch correction using ComBat was performed to control for batch effects. 33

| RNA sequencing and analysis
Whole-genome transcriptomic patterns were analyzed on total RNA. Cell line RNA samples were subjected to 150 bp paired-end total RNA-seq (155M reads) in a single batch. RNA integrity and yield was quantified using the 2100 Bioanalyzer®. RNA samples were submitted to the Edinburgh Genomics sequencing facility, where QC, ribosomal-depletion, strand-aware library preparation and Illumina adapter ligation was performed. Ribosomal RNA was depleted using the New England Biolabs NEBNext rRNA Depletion Kit according to the manufacturer's protocol. Samples were sequenced on the Illumina HiSeq 2500 platform in 'rapid mode' with 150 bp pairedend reads. Transcript indexing and quantification from RNA-seq reads was performed using Salmon v1.1.0. 34 The Salmon index was generated using the reference dataset GRCh38.primary_assembly.fa and the Ensembl 96 release, (ftp://ftp.ensem bl.org/pub/relea se-96/fasta/ homo_ sapiens). In brief, the Ensembl (cDNA and ncRNA) and reference FASTA were concatenated, the 'decoys.txt' was prepared from GRCh38 and the Salmon index run on the concatenated FASTA file in the R statistical environment (v4.0.0). Paired reads 1 and read 2 fastq files respectively were then concatenated for each sample, and Salmon run using auto-detect strandedness and 'validateMappings' flag.
Salmon transcript level quantitation files were combined and summarized to gene-level counts using txiimport 35 and biomaRt. 36,37 Transcript level counts were aggregated into gene levels counts and imported into R using the tximport package. A Differential Gene Expression object (DGEList) was created using tximport generated gene-level counts and sample metadata for analysis using edgeR (3.30.0) and limma (3.44.1). Genes were filtered using the filterByExpr function in the edgeR package (3.30.0) to include genes with at least 10 reads per gene. The calcNormFactors function was used to calculate Trimmed mean of M-values (TMM) between-sample normalization to scale the library size. The limma 'voom' function was used to transform normalized count data to log2-counts per million (logCPM), estimate the meanvariance relationship and used to compute appropriate observation-level weights.

| Genotype at the vitamin D pathway loci
Cell line genomic DNA was genotyped for rs2282679 known to be associated with 25-OHD level. 38 Three functionally relevant polymorphisms in the VDR gene which are established as functionally relevant, impacting VDR function, DNA binding and calcitriol binding (rs1544410, rs10735810, rs7975232) were also genotyped by using the Infinium Omni5 BeadChip (Illumina, San Diego, CA, USA).

| Mass-spectrometry
Cell line samples were prepared for mass-spectrometry analysis using a Guanidinium lysis followed by LysC and Trypsin digest as described. 39 Briefly, the cells were lysed in 50 µl 6M GuHCl, 100 mM Tris pH 8.5 with 1 mg/ml Chloracetamide and 1.5 mg/ml TCEP.
The cells were probe sonicated until lysed then heated to 95ºC for 5 min. and digested with 1/200 (Enzyme to Substrate) LysC (Wako Chemicals, Japan) for 4 h at 37 ºC. Samples were diluted with 250 µl containing 1/200 (Enzyme to Substrate) porcine Trypsin (Promega, UK) and the digest was continued overnight. Peptides samples were acidified with TFA to 0.5% (final concentration) and cleared by centrifugation for 10 min on a benchtop centrifuge. Samples were desalted using Stage-tips as described. 40 Mass spectrometry was performed using a Lumos Fusion (Thermo) mass spectrometer coupled to a RLS-nano uHPLC (Thermo). Peptides were separated by a 140 min linear gradient from 5% to 30% acetonitrile, 0.05% acetic acid. Proteins were identified and quantified the MaxQuant software suite 41 using label-free quantification and searching against the Human Uniprot database.

| Statistical analysis
All statistical analysis was undertaken in R. 42 Investigation of differential gene expression in response to calcitriol treatment was undertaken using the lmFit and eBayes functions within the 'limma' package. 43 Differentially expressed gene-sets were assessed for enrichment across expression analysis methodologies and between in vitro and ex vivo datasets using the geneSetTest function in R. Ranked gene lists were assessed for functional relevance using the 'GOrilla', Gene Ontology enRIchment anaLysis and visuaLizAtion tool. 44 Process ontologies were validated by Gene Set Enrichment Analysis using the gseGO function within the 'clusterprofiler' package in R. 45

and E-cadherin expression in colorectal cancer cell lines
We optimized treatment conditions, using SW480 cells with CYP24A1 and CDH1 expression as readouts. CYP24A1, which encodes 24-hydroxylase a key enzyme in vitamin D metabolism, was chosen as positive control, while CDH1, a tumor suppressor gene which encodes E-cadherin, was assessed given reported induction in SW480 cells by calcitriol. 21 Significant induction of both CYP24A1 (fold-change [FC] = 4.45 95% CI 3.00-5.89) and CDH1 (FC = 3.60 95% CI 2.17-5.03) was observed with 50 nM calcitriol treatment over 16-h, adopted hereafter as standard treatment conditions (Figures S1 and 1).
Next, we tested for CDH1 induction in a panel of cell lines. The 10 cell lines treated represent CRC as a whole, given the broad tumors of origin (see Supplementary  Methods) and are known to differ in their mutation profile e.g., microsatellite instability, CpG island methylator phenotype status and molecular subtype. 46 Significant induction of CDH1 by calcitriol treatment was observed in SW480, HCT116 and LS174T cell lines, with induction of E-Cadherin seen in SW480 and HCT116 cells but not LS174T (Figure 1). No association between cell genotype at vitamin D pathway SNPs (rs2282679, rs1544410, rs10735810, rs7975232) and either baseline CDH1, or induction with calcitriol was observed (p > .05).

| Calcitriol influences expression of 30 genes in CRC cell lines
We chose three cell lines for HT12 microarray analysis based on CDH1 induction on qRT-PCR (SW480, HCT116, LS174T). Analysis confirmed the induction of both CDH1 and CYP24A1 by 50 nM calcitriol when results were pooled for the three cell lines (CDH1 logFC = 0.80, FDR = 0.0001, CYP24A1 logFC = 3.70 FDR = 0.02; Figure S2), demonstrating strong correlation with qRT-PCR expression data: (CDH1 r = 0.91, p < 2.2e−16; CYP24A1 r = 0.88, p < 2.2e−16). In total, 30 genes were differentially expressed after calcitriol treatment (FDR < 0.05, Table S1). When considered separately, there were 22, 972 and 417 significantly differentially F I G U R E 1 (A) Replicate expression fold-change of CYP24A1 and CDH1 in SW480 cells with 16 h 500 pM-50 nM. Cells were treated with 500 pM or 50 nM calcitriol or equivalent volume of ETOH for 16 h and harvested for RNA. Up to11 biological replicate experiments were performed. qRT-PCR was performed to assess expression of CYP24A1 and CDH1 normalized to ACTB and average fold-change expression relative to ETOH control was charted. Student's t-test was performed on fold-changes per treatment. *p < .05, **p < .01, ***p < .001. Error bars represent standard error of mean (SEM). (B) Induction of CDH1 by calcitriol in panel of CRC cell lines. 10 cell lines were treated in triplicate with 50 nM calcitriol or equivalent volume of ETOH for 16 h and harvested for RNA. qRT-PCR was performed to assess expression of CDH1 normalized to ACTB and average fold-change in expression relative to ETOH control was charted. Error bars represent SEM of biological replicate samples. Student's t-test was performed on fold-changes per treatment. *p < .05, **p < .01, ***p < .001. (C) Western blot of E cadherin expression in ETOH treated cell lines. Cells were treated for 16 h ETOH and harvested for protein. Western blotting was performed with a commercial E-cadherin antibody (Cell Signalling, #3195), with βactin (Sigma, #A1978) as loading control. (D) Western blot of E cadherin expression with 50 nM calcitriol treatment. Cells were treated for 16 h with 50 nM calcitriol or ETOH and harvested for protein. Western blotting was performed with a commercial E-cadherin antibody (Cell Signalling, #3195), with βactin (Sigma, #A1978) as loading control expressed genes HCT116, LS174T and SW480 cell samples respectively. We sought technical replication of our results using RNA-sequencing. In total, 59 named genes were differentially expressed (Table S2), with CDH1 the fourth hit (logFC-0.85, FDR = 0.002, Figure S3) and strong correlation between expression values for differentially expressed genes between HT12/RNAseq (r = 0.85, p = 1.8e−08, Figure S4). Significant enrichment of the 30 differentially expressed genes on HT12 in the RNA-seq dataset was observed (2.39e−10, Figure S5).

| GO term enrichment indicates dysregulation of anti-tumor effects with calcitriol
Functional annotation of the HT12 gene list identified 63 significantly enriched ontologies, with the top 10 and additional relevant terms listed in Table 1. Results from GOrilla analysis were replicated using GSEA, which identified 131 significantly enriched terms and confirmed enrichment of 20 common processes seen with GOrilla, as depicted in Figure S6. A number of additional terms identified by GSEA had links to anti-tumor processes, e.g., 'regulation of Wnt signaling pathway', 'columnar/cuboidal epithelial cell differentiation' and 'cell adhesion' (Table S3). GOrilla and GSEA analysis were repeated using the ranked list from RNA-seq analysis. 127 significantly enriched terms were identified using GOrilla (Table S4), with 21 terms in common with analysis of the HT12 dataset and notable terms including 'regulation of cell migration' and 'negative regulation of cell adhesion'. 65 terms were enriched when using GSEA (Table S5) Note: HCT116, LS174T and SW480 cell lines were treated in triplicate with 50 nM calcitriol or ETOH and gene expression in harvested RNA assessed using HT12 microarray. Differential expression was calculated using Student's t-test and genes ranked by significance entered into the GO enrichment software GOrilla. FDR q-value' is the correction of the above p-value for multiple testing using the Benjamini and Hochberg method. 47 For the ith term (ranked according to p-value) the FDR q-value is (p-value * number of GO terms)/i. Enrichment = (the number of genes in the intersection/number of genes in the top of the input list)/(total number of genes associated with a specific GO term/total number of genes). Additional relevant processes not ranking in top 10 are given (rank).

T A B L E 1
Gene ontology analysis of differentially expressed genes in calcitriol treated cell line samples n = 3889, SW480 n = 3845). Increased expression of 1,25OHD 24-hydroxylase, encoded by CYP24A1, was seen across all cell lines (p < .0005) and the top ten significantly altered proteins in each cell line are given in Table S6. On pooling cell line data, no significantly altered proteins were identified after adjustment for multiple testing, yet enrichment analysis of this ranked list using REACTOME V57 48 revealed 216 significant entities (FDR < 0.05), with a number of pathways linked to tumorigenesis, e.g., 'Cell cycle' and 'Apoptosis' (Table S7 and Figure S7).

CDH1 and E-cadherin expression
We used patient-derived colorectal organoids to validate the effect of calcitriol on epithelial cell gene expression and processes relevant to tumorigenesis (Figures 2 and  S8). Time-course and dose-response experiment demonstrated induction of CYP24A1 and CDH1 after 50-100 nM calcitriol treatment, validating in vitro observations of both calcitriol metabolism and CDH1 induction. Calcitriol treatment (100 nM, 24 h) of 6 further patient-derived organoid cultures ( Figures S1 and 2), consistently induced both CYP24A1 and CDH1 expression.

| Calcitriol influences expression of 111 genes in patient-derived organoids
RNA from paired calcitriol treated/untreated non-FAP organoids (patient ID B, C & D) was submitted for microarray analysis. This identified 111 differentially expressed genes (top hit CD14, logFC = 4.91, FDR 4.47e−05, Table S8). The differentially expressed gene-set from cell line data was significantly enriched (p = 3.97e−08) in the organoid dataset indicating conservation of gene expression responses to calcitriol across these two model systems, with six individual genes differentially expressed (FDR < 0.05) in both cell line and organoid samples (  (Tables S9 and  S10). GO enrichment (GOrilla) revealed 328 significantly enriched terms, with many relevant to tumorigenesis, including 'cell migration', 'regulation of cell adhesion', regulation of cell death' and 'regulation of cell differentiation' (Table S11). Seventeen GO terms were enriched (FDR p < .05) in both cell line and organoid samples, including many relevant to tumorigenesis (e.g., 'regulation of cell proliferation,' 'regulation of cell differentiation'; selected relevant terms Figure 3). Replication of enrichment analysis using GSEA confirmed enrichment of similar terms (e.g., 'epithelial cell differentiation', 'apoptotic signalling pathway', 'regulation of cell death') and identified 38 confirmed common terms with in vitro dataset.

| DISCUSSION
We report for the first-time, common genes and pathways modulated by calcitriol treatment across both in vitro (cell lines) and ex vivo (organoids) models including CDH1 and several other genes linked to tumorigenesis. The enrichment of pathways relating to cell differentiation, proliferation, migration and death in addition to immune processes support epidemiological data of the role of vitamin D as a chemopreventive agent. These data provide strong rationale for further well-designed trials of vitamin D supplementation for CRC prevention. The CDH1 gene, codes for E-cadherin, a glycoprotein required for effective cell-cell adhesion. E-cadherin also regulates βcatenin signalling in the Wnt pathway, the aberrant activation of which has been shown to be an initiating event in the development of CRC. 49 The loss of E-cadherin expression is also central to the epithelialto-mesenchymal transition, a key step in the metastatic progression of human cancers, and reduced CDH1 expression through mutation or promoter methylation is seen in several epithelial cancers, resulting in increased invasiveness and a reduced survival in vivo. Inactivating germline mutations are reported in familial aggregations of gastric, colonic and breast cancers. Here we show induction of CDH1 in several established CRC cell lines and in organoids derived from normal mucosa sampled from various sites along the colorectum. Our study expands on work by Palmer et al., who reported an induction of E-Cadherin in SW480 colorectal cancer cells treated with 1,25OHD. 21 Here, we identify additional CRC cell lines with disparate molecular characteristics and mutation profile in which CDH1 is induced by 1,25OHD, providing a potential mechanism to underlie risk/survival associations between vitamin D status/supplementation and colorectal cancer.
Transcriptomic analysis identifies a number of other genes which are modulated by calcitriol treatment. These genes provide targets for validation and future mechanistic investigation. Within the common gene-set both SULT1C2 50 and CYP24A1 are known to be induced by calcitriol providing confidence that calcitriol is activating the vitamin D pathway in both culture models. Of particular interest is the consistent upregulation of GADD45A by calcitriol. The Growth Arrest and DNA Damageinducible 45 (GADD45) proteins have been implicated in regulation of many cellular functions including DNA repair, cell cycle control, senescence and genotoxic stress. Furthermore, aberrations in the GADD45 pathway and dysregulation of GADD45 expression is associated with the initiation and progression of malignancies and inhibition of GADD45 through promoter methylation or NF-κB activation, and is considered a critical step in cancer development. Meanwhile, induction of GADD45 expression is reported to be an essential step for mediating anti-cancer activity of multiple chemotherapeutic drugs. 51 Induction of GADD45 has previously been reported in prostate, 52 ovarian, 53 squamous 54 and colorectal cancer 55 cell lines supporting our current findings. Similarly, EFTUD1 is also known to be upregulated by calcitriol in breast cancer and breast tumor explants and elevated EFTUD1 predicts prolonged breast cancer survival. 56 Finally, KIAA1199, which is downregulated by calcitriol in our study, is a known oncogene which plays a role in proliferation, apoptosis, invasion and migration of various tumor cells. 57 It is overexpressed in colorectal cancer with higher levels associated with increased tumor invasion depth, stage, and poor prognosis in CRC. To the best of our knowledge this is the first study to report suppression of intestinal epithelial cell KIAA1199 expression with calcitriol, providing an exciting target for further study.
Our ex vivo dataset adds to a small number of studies demonstrating gene expression changes in patient-derived organoids treated with vitamin D. 25-27 F I G U R E 3 (A) Enrichment plot of differentially expressed genes from in vitro dataset and enrichment in ex vivo dataset. Enrichment plot y axis charts relative enrichment. Each black bar within the bar code represents a different candidate gene which are ranked from left to right by increasing T statistic from a moderated T test. The curve (or worm) above the barcode shows the relative local enrichment of the bars in each part of the plot. The dotted horizontal line indicates neutral enrichment; the worm above the dotted line shows enrichment while the worm below the dotted line shows depletion. (B) Venn diagram showing intersect of differentially expressed genes with calcitriol treatment of cell lines and organoids. Three cell lines (HCT116, LS174T, SW480) were treated in triplicate with 50 nM calcitriol or equivalent volume of ETOH for 16 h and harvested for RNA. HT12 microarray used to assess gene expression and paired results pooled to assess for in vitro responses to calcitriol. Human intestinal organoids were derived from normal mucosa from three patients and treated with calcitriol for 24 h. Organoids from six wells were pooled for each RNA sample, and HT12 microarray used to assess gene expression and paired results pooled to assess for in vitro responses to calcitriol. (C) Heatmap of common genes differentially expressed after in vitro and ex vivo calcitriol treatment. (D) Processes relevant to tumorigenesis significantly enriched after both in vitro and ex vivo calcitriol treatment. Treatment and gene expression analysis performed as described for (A). GOrilla enrichment analysis performed using ranked gene lists and compared between the in vitro and ex vivo datasets Fernández-Barral et al., demonstrated upregulation of stemness-related genes and inhibition of cell proliferation in patient-derived colonic organoids, with differential expression of >2000 genes. This gene-set was significantly expressed after calcitriol treatment of our organoids (p = 2.49e−93; Figure S9). The intersect of common differentially expressed genes was 66, including GADD45A, with complete concordance in direction of modulation (i.e., up-or down-regulated). Enrichment of our recently reported gene-set associated with 25-OHD level in human rectal mucosa was also seen in our organoid data, providing further validation of our results 58 (p = .0075, Figure S9).
Pathway analysis in the current study confirms implied links with carcinogenesis at individual gene level. Such terms as 'cell migration', 'regulation of cell adhesion', 'regulation of cell death', 'regulation of cell differentiation' and 'cell cycle' support numerous previous in vitro studies demonstrating vitamin D-induced growth arrest and apoptosis of colorectal cancer cells, modulation of the Wnt signalling pathway, DNA repair and immunomodulation, 59 thus lending support to a causal relationship between vitamin D and cancer. Several clinical studies have replicated these in vitro findings in vivo-Protiva et al., reported upregulation of genes involved in inflammation, immune response, extracellular matrix, and cell adhesion in response to 1,25OHD, 60 while a randomized study of vitamin D3 supplementation reported a reduction in tumor promoting inflammation biomarkers, decreased oxidative DNA damage, increased cell differentiation and apoptosis and modification of the APC/β-catenin pathway in the normal human colorectal epithelium. [61][62][63][64][65] Thus, our current findings confirm further possible mechanisms to underlie the observed association between vitamin D deficiency and increased CRC risk.
Whilst the complementary in vitro/ex vivo design of our study has many positive attributes, we acknowledge several limitations. First, we acknowledge the limitations of comparative cell line studies including the risk that accumulated cellular genomic and karyotypic abnormalities may significantly influence results from experimental studies. However, the replication of several common genes and pathways in our ex vivo studies, i.e., a non-aberrant mode system, provides confidence in our results. Next, for statistical power, we pooled the different cell lines and organoid for transcriptomic profiling/sequencing which ignores potential differential responses between cultures and cell lines. However, even despite this, and the known heterogeneity in gene expression between cell lines or colorectal epithelial cell sampling site, 66 relevant differential gene expression and pathway enrichment is observed, indicating conservation of responses to vitamin D irrespective of cell line differences or sampling site, thus supporting the broad potential of vitamin D as valuable adjunctive chemotherapeutic agent in CRC. We acknowledge that the smaller number of differentially expressed genes in the in vitro studies likely reflects the different characteristics of the cell lines studied, confirmed by the marked difference in significantly differentially expressed genes between the three cell lines. Further detailed investigation of the relevance of these differences in respect of baseline cell line characteristics and mechanistic studies to further define the importance of potential vitamin D target genes identified here are required. Finally, we acknowledge that RNA-sequencing is now the predominant method of transcriptome gene expression analysis. HT12 analysis was used as our primary analysis platform for the exploratory analyses providing data on >31 000 annotated genes as the most practical and cost-effective method of transcriptomic profiling. We replicated the induction of CDH1 at a technical level (qRT-PCR, HT12, RNA-seq), while the modulation of other genes of interest (GADD45, KIAA1199, GADD45A) was conserved across both microarray and RNA-seq platforms. Other groups have published reports of good high levels of reproducibility between sequencing and microarray platforms. 67 Even if correlation at a transcript/probe level was not observed, RNA-seq and microarray can be seen to be complementary. 68 As such, the differences in platform and methodology that we have utilized, with reassuring agreement in differential expression and gene-set enrichment with RNAseq/HT12 data provides confidence that these are true effects.
In conclusion, we show that calcitriol induces induction of CDH1, GADD45 and EFTUD1 in both CRC cell lines and patient-derived epithelial cell organoids, with suppression of KIAA1199 also observed; genes with recognized links to CRC tumorigenesis. Furthermore, several pathways modulated by calcitriol treatment are intimately linked to tumorigenesis. These findings support epidemiological data and provide further strong rationale for well-designed trials of vitamin D supplementation as a novel CRC chemopreventive and chemotherapeutic agent.